1640 lines
72 KiB
Plaintext
1640 lines
72 KiB
Plaintext
{
|
||
"cells": [
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "509f7bd9-b597-4d23-b3af-0a76a7b4d33d",
|
||
"metadata": {},
|
||
"source": [
|
||
"{/* cspell:ignore prefactors */}\n",
|
||
"\n",
|
||
"# Krylov quantum diagonalization of lattice Hamiltonians\n",
|
||
"*Usage estimate: 20 minutes on IBM Fez (NOTE: This is an estimate only. Your runtime may vary.)*"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": null,
|
||
"id": "71d42c57-d1aa-424a-8116-7a3b77bd5fd9",
|
||
"metadata": {
|
||
"tags": [
|
||
"remove-cell"
|
||
]
|
||
},
|
||
"outputs": [],
|
||
"source": [
|
||
"# This cell is hidden from users – it disables some lint rules\n",
|
||
"# ruff: noqa: E402 E722 F601"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "921c7b04-5b5d-4cfc-aba0-15a61334e619",
|
||
"metadata": {},
|
||
"source": [
|
||
"## Background\n",
|
||
"\n",
|
||
"This tutorial demonstrates how to implement the Krylov Quantum Diagonalization Algorithm (KQD) within the context of Qiskit patterns. You will first learn about the theory behind the algorithm and then see a demonstration of its execution on a QPU.\n",
|
||
"\n",
|
||
"Across disciplines, we're interested in learning ground state properties of quantum systems. Examples include understanding the fundamental nature of particles and forces, predicting and understanding the behavior of complex materials and understanding bio-chemical interactions and reactions. Because of the exponential growth of the Hilbert space and the correlation that arise in entangled systems, classical algorithm struggle to solve this problem for quantum systems of increasing size. At one end of the spectrum, existing approach that take advantage of the quantum hardware focus on variational quantum methods (e.g., [variational quantum eigen-solver](https://learning.quantum.ibm.com/tutorial/variational-quantum-eigensolver)). These techniques face challenges with current devices because of the high number of function calls required in the optimization process, which adds a large overhead in resources once advanced error mitigation techniques are introduced, thus limiting their efficacy to small systems. At the other end of the spectrum, there are fault-tolerant quantum methods with performance guarantees (e.g., [quantum phase estimation](https://arxiv.org/abs/quant-ph/0604193)) which require deep circuits that can be executed only on a fault-tolerant device. For these reasons, we introduce here a quantum algorithm based on subspace methods (as described in this [review paper](https://arxiv.org/abs/2312.00178)), the Krylov quantum diagonalization (KQD) algorithm. This algorithm performs well at large scale [\\[1\\]](#references) on existing quantum hardware, shares similar [performance guarantees](https://arxiv.org/abs/2110.07492) as phase estimation, is compatible with advanced error mitigation techniques and could provide results that are classically inaccessible."
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "5f698c82-95ca-4dc4-b9fa-d6e741e2c02c",
|
||
"metadata": {},
|
||
"source": [
|
||
"## Requirements\n",
|
||
"Before starting this tutorial, be sure you have the following installed:\n",
|
||
"\n",
|
||
"- Qiskit SDK v1.0 or later, with visualization support ( `pip install 'qiskit[visualization]'` )\n",
|
||
"- Qiskit Runtime 0.22 or later ( `pip install qiskit-ibm-runtime` )"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "c44956e4-47ab-4b0f-9d6d-553080110062",
|
||
"metadata": {},
|
||
"source": [
|
||
"## Setup"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": null,
|
||
"id": "ce7d5adc-3ef2-4654-b865-14d5141ce41a",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"import numpy as np\n",
|
||
"import scipy as sp\n",
|
||
"import matplotlib.pylab as plt\n",
|
||
"from typing import Union, List\n",
|
||
"import itertools as it\n",
|
||
"import copy\n",
|
||
"from sympy import Matrix\n",
|
||
"import warnings\n",
|
||
"\n",
|
||
"warnings.filterwarnings(\"ignore\")\n",
|
||
"\n",
|
||
"from qiskit.quantum_info import SparsePauliOp, Pauli, StabilizerState\n",
|
||
"from qiskit.circuit import Parameter\n",
|
||
"from qiskit import QuantumCircuit, QuantumRegister\n",
|
||
"from qiskit.circuit.library import PauliEvolutionGate\n",
|
||
"from qiskit.synthesis import LieTrotter\n",
|
||
"from qiskit.transpiler import Target, CouplingMap\n",
|
||
"from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager\n",
|
||
"\n",
|
||
"\n",
|
||
"from qiskit_ibm_runtime import (\n",
|
||
" QiskitRuntimeService,\n",
|
||
" EstimatorOptions,\n",
|
||
" EstimatorV2 as Estimator,\n",
|
||
")\n",
|
||
"\n",
|
||
"\n",
|
||
"def solve_regularized_gen_eig(\n",
|
||
" h: np.ndarray,\n",
|
||
" s: np.ndarray,\n",
|
||
" threshold: float,\n",
|
||
" k: int = 1,\n",
|
||
" return_dimn: bool = False,\n",
|
||
") -> Union[float, List[float]]:\n",
|
||
" \"\"\"\n",
|
||
" Method for solving the generalized eigenvalue problem with regularization\n",
|
||
"\n",
|
||
" Args:\n",
|
||
" h (numpy.ndarray):\n",
|
||
" The effective representation of the matrix in the Krylov subspace\n",
|
||
" s (numpy.ndarray):\n",
|
||
" The matrix of overlaps between vectors of the Krylov subspace\n",
|
||
" threshold (float):\n",
|
||
" Cut-off value for the eigenvalue of s\n",
|
||
" k (int):\n",
|
||
" Number of eigenvalues to return\n",
|
||
" return_dimn (bool):\n",
|
||
" Whether to return the size of the regularized subspace\n",
|
||
"\n",
|
||
" Returns:\n",
|
||
" lowest k-eigenvalue(s) that are the solution of the regularized generalized eigenvalue problem\n",
|
||
"\n",
|
||
"\n",
|
||
" \"\"\"\n",
|
||
" s_vals, s_vecs = sp.linalg.eigh(s)\n",
|
||
" s_vecs = s_vecs.T\n",
|
||
" good_vecs = np.array(\n",
|
||
" [vec for val, vec in zip(s_vals, s_vecs) if val > threshold]\n",
|
||
" )\n",
|
||
" h_reg = good_vecs.conj() @ h @ good_vecs.T\n",
|
||
" s_reg = good_vecs.conj() @ s @ good_vecs.T\n",
|
||
" if k == 1:\n",
|
||
" if return_dimn:\n",
|
||
" return sp.linalg.eigh(h_reg, s_reg)[0][0], len(good_vecs)\n",
|
||
" else:\n",
|
||
" return sp.linalg.eigh(h_reg, s_reg)[0][0]\n",
|
||
" else:\n",
|
||
" if return_dimn:\n",
|
||
" return sp.linalg.eigh(h_reg, s_reg)[0][:k], len(good_vecs)\n",
|
||
" else:\n",
|
||
" return sp.linalg.eigh(h_reg, s_reg)[0][:k]\n",
|
||
"\n",
|
||
"\n",
|
||
"def single_particle_gs(H_op, n_qubits):\n",
|
||
" \"\"\"\n",
|
||
" Find the ground state of the single particle(excitation) sector\n",
|
||
" \"\"\"\n",
|
||
" H_x = []\n",
|
||
" for p, coeff in H_op.to_list():\n",
|
||
" H_x.append(set([i for i, v in enumerate(Pauli(p).x) if v]))\n",
|
||
"\n",
|
||
" H_z = []\n",
|
||
" for p, coeff in H_op.to_list():\n",
|
||
" H_z.append(set([i for i, v in enumerate(Pauli(p).z) if v]))\n",
|
||
"\n",
|
||
" H_c = H_op.coeffs\n",
|
||
"\n",
|
||
" print(\"n_sys_qubits\", n_qubits)\n",
|
||
"\n",
|
||
" n_exc = 1\n",
|
||
" sub_dimn = int(sp.special.comb(n_qubits + 1, n_exc))\n",
|
||
" print(\"n_exc\", n_exc, \", subspace dimension\", sub_dimn)\n",
|
||
"\n",
|
||
" few_particle_H = np.zeros((sub_dimn, sub_dimn), dtype=complex)\n",
|
||
"\n",
|
||
" sparse_vecs = [\n",
|
||
" set(vec) for vec in it.combinations(range(n_qubits + 1), r=n_exc)\n",
|
||
" ] # list all of the possible sets of n_exc indices of 1s in n_exc-particle states\n",
|
||
"\n",
|
||
" m = 0\n",
|
||
" for i, i_set in enumerate(sparse_vecs):\n",
|
||
" for j, j_set in enumerate(sparse_vecs):\n",
|
||
" m += 1\n",
|
||
"\n",
|
||
" if len(i_set.symmetric_difference(j_set)) <= 2:\n",
|
||
" for p_x, p_z, coeff in zip(H_x, H_z, H_c):\n",
|
||
" if i_set.symmetric_difference(j_set) == p_x:\n",
|
||
" sgn = ((-1j) ** len(p_x.intersection(p_z))) * (\n",
|
||
" (-1) ** len(i_set.intersection(p_z))\n",
|
||
" )\n",
|
||
" else:\n",
|
||
" sgn = 0\n",
|
||
"\n",
|
||
" few_particle_H[i, j] += sgn * coeff\n",
|
||
"\n",
|
||
" gs_en = min(np.linalg.eigvalsh(few_particle_H))\n",
|
||
" print(\"single particle ground state energy: \", gs_en)\n",
|
||
" return gs_en"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "76a1c616-a749-48a2-b5fd-8beff760760b",
|
||
"metadata": {},
|
||
"source": [
|
||
"## Step 1: Map classical inputs to a quantum problem"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "a166e2a1-3003-4799-9c15-59f5e78b6ea8",
|
||
"metadata": {},
|
||
"source": [
|
||
"### The Krylov space\n",
|
||
"\n",
|
||
"The Krylov space $\\mathcal{K}^r$ of order $r$ is the space spanned by vectors obtained by multiplying higher powers of a matrix $A$, up to $r-1$, with a reference vector $\\vert v \\rangle$.\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\mathcal{K}^r = \\left\\{ \\vert v \\rangle, A \\vert v \\rangle, A^2 \\vert v \\rangle, ..., A^{r-1} \\vert v \\rangle \\right\\}\n",
|
||
"$$\n",
|
||
"\n",
|
||
"If the matrix $A$ is the Hamiltonian $H$, we'll refer to the corresponding space as the power Krylov space $\\mathcal{K}_P$. In the case where $A$ is the time-evolution operator generated by the Hamiltonian $U=e^{-iHt}$, we'll refer to the space as the unitary Krylov space $\\mathcal{K}_U$. The power Krylov subspace that we use classically cannot be generated directly on a quantum computer as $H$ is not a unitary operator. Instead, we can use the time-evolution operator $U = e^{-iHt}$ which can be shown to give similar [convergence guarantees](https://arxiv.org/abs/2110.07492) as the power method. Powers of $U$ then become different time steps $U^k = e^{-iH(kt)}$.\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\mathcal{K}_U^r = \\left\\{ \\vert \\psi \\rangle, U \\vert \\psi \\rangle, U^2 \\vert \\psi \\rangle, ..., U^{r-1} \\vert \\psi \\rangle \\right\\}\n",
|
||
"$$\n",
|
||
"\n",
|
||
"See the Appendix for a detailed derivation of how the unitary Krylov space allows to represents low-energy eigenstates accurately."
|
||
]
|
||
},
|
||
{
|
||
"attachments": {},
|
||
"cell_type": "markdown",
|
||
"id": "5573ca7e-ab16-4488-88b3-d8d1eba9e20c",
|
||
"metadata": {},
|
||
"source": [
|
||
"### Krylov quantum diagonalization algorithm\n",
|
||
"\n",
|
||
"Given an Hamiltonian $H$ that we wish to diagonalize, first we consider the corresponding unitary Krylov space $\\mathcal{K}_U$. The goal is to find a compact representation of the Hamiltonian in $\\mathcal{K}_U$, which we'll refer to as $\\tilde{H}$. The matrix elements of $\\tilde{H}$, the projection of the Hamiltonian in the Krylov space, can be calculated by calculating the following expectation values\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\tilde{H}_{mn} = \\langle \\psi_m \\vert H \\vert \\psi_n \\rangle =\n",
|
||
"$$\n",
|
||
"$$\n",
|
||
"= \\langle \\psi \\vert e^{i H t_m} H e^{-i H t_n} \\vert \\psi \\rangle\n",
|
||
"$$\n",
|
||
"$$\n",
|
||
"= \\langle \\psi \\vert e^{i H m dt} H e^{-i H n dt} \\vert \\psi \\rangle\n",
|
||
"$$\n",
|
||
"\n",
|
||
"Where $\\vert \\psi_n \\rangle = e^{-i H t_n} \\vert \\psi \\rangle$ are the vectors of the unitary Krylov space and $t_n = n dt$ are the multiples of the time step $dt$ chosen. On a quantum computer, the calculation of each matrix elements can be done with any algorithm which allows to obtain overlap between quantum states. In this tutorial we'll focus on the Hadamard test. Given that the $\\mathcal{K}_U$ has dimension $r$, the Hamiltonian projected into the subspace will have dimensions $r \\times r$. With $r$ small enough (generally $r<<100$ is sufficient to obtain convergence of estimates of eigenenergies) we can then easily diagonalize the projected Hamiltonian $\\tilde{H}$. However, we cannot directly diagonalize $\\tilde{H}$ because of the non-orthogonality of the Krylov space vectors. We'll have to measure their overlaps and construct a matrix $\\tilde{S}$\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\tilde{S}_{mn} = \\langle \\psi_m \\vert \\psi_n \\rangle\n",
|
||
"$$\n",
|
||
"\n",
|
||
"\n",
|
||
"This allows us to solve the eigenvalue problem in a non-orthogonal space (also called generalized eigenvalue problem)\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\tilde{H} \\ \\vec{c} = E \\ \\tilde{S} \\ \\vec{c}\n",
|
||
"$$\n",
|
||
"\n",
|
||
"One can then obtain estimates of the eigenvalues and eigenstates of $H$ by looking at the ones of $\\tilde{H}$. For example, the estimate of the ground state energy is obtained by taking the smallest eigenvalue $c$ and the ground state from the corresponding eigenvector $\\vec{c}$. The coefficients in $\\vec{c}$ determines the contribution of the different vectors that span $\\mathcal{K}_U$.\n",
|
||
"\n",
|
||
"\n",
|
||
"\n",
|
||
"The Figure shows a circuit representation of the modified Hadamard test, a method that is used to compute the overlap between different quantum states. For each matrix element $\\tilde{H}_{i,j}$, a Hadamard test between the state $\\vert \\psi_i \\rangle$, $\\vert \\psi_j \\rangle$ is carried out. This is highlighted in the figure by the color scheme for the matrix elements and the corresponding $\\text{Prep} \\; \\psi_i$, $\\text{Prep} \\; \\psi_j$ operations. Thus, a set of Hadamard tests for all the possible combinations of Krylov space vectors is required to compute all the matrix elements of the projected Hamiltonian $\\tilde{H}$. The top wire in the Hadamard test circuit is an ancilla qubit which is measured either in the X or Y basis, its expectation value determines the value of the overlap between the states. The bottom wire represents all the qubits of the system Hamiltonian. The $\\text{Prep} \\; \\psi_i$ operation prepares the system qubit in the state $\\vert \\psi_i \\rangle$ controlled by the state of the ancilla qubit (similarly for $\\text{Prep} \\; \\psi_j$) and the operation $P$ represents Pauli decomposition of the system Hamiltonian $H = \\sum_i P_i$. A more detailed derivation of the operations calculated by the Hadamard test is given below."
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "1a6d7a4f-6c1f-4069-93d1-f5b670645d7f",
|
||
"metadata": {},
|
||
"source": [
|
||
"#### Define Hamiltonian\n",
|
||
"\n",
|
||
"Let's consider the Heisenberg Hamiltonian for $N$ qubits on a linear chain: $H= \\sum_{i,j}^N X_i X_j + Y_i Y_j - J Z_i Z_j$"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 2,
|
||
"id": "82163249-bafd-4bc6-9741-20a4077971a7",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
"[('ZZIIIIIIIIIIIIIIIIIIIIIIIIIIII', 1), ('IZZIIIIIIIIIIIIIIIIIIIIIIIIIII', 1), ('IIZZIIIIIIIIIIIIIIIIIIIIIIIIII', 1), ('IIIZZIIIIIIIIIIIIIIIIIIIIIIIII', 1), ('IIIIZZIIIIIIIIIIIIIIIIIIIIIIII', 1), ('IIIIIZZIIIIIIIIIIIIIIIIIIIIIII', 1), ('IIIIIIZZIIIIIIIIIIIIIIIIIIIIII', 1), ('IIIIIIIZZIIIIIIIIIIIIIIIIIIIII', 1), ('IIIIIIIIZZIIIIIIIIIIIIIIIIIIII', 1), ('IIIIIIIIIZZIIIIIIIIIIIIIIIIIII', 1), ('IIIIIIIIIIZZIIIIIIIIIIIIIIIIII', 1), ('IIIIIIIIIIIZZIIIIIIIIIIIIIIIII', 1), ('IIIIIIIIIIIIZZIIIIIIIIIIIIIIII', 1), ('IIIIIIIIIIIIIZZIIIIIIIIIIIIIII', 1), ('IIIIIIIIIIIIIIZZIIIIIIIIIIIIII', 1), ('IIIIIIIIIIIIIIIZZIIIIIIIIIIIII', 1), ('IIIIIIIIIIIIIIIIZZIIIIIIIIIIII', 1), ('IIIIIIIIIIIIIIIIIZZIIIIIIIIIII', 1), ('IIIIIIIIIIIIIIIIIIZZIIIIIIIIII', 1), ('IIIIIIIIIIIIIIIIIIIZZIIIIIIIII', 1), ('IIIIIIIIIIIIIIIIIIIIZZIIIIIIII', 1), ('IIIIIIIIIIIIIIIIIIIIIZZIIIIIII', 1), ('IIIIIIIIIIIIIIIIIIIIIIZZIIIIII', 1), ('IIIIIIIIIIIIIIIIIIIIIIIZZIIIII', 1), ('IIIIIIIIIIIIIIIIIIIIIIIIZZIIII', 1), ('IIIIIIIIIIIIIIIIIIIIIIIIIZZIII', 1), ('IIIIIIIIIIIIIIIIIIIIIIIIIIZZII', 1), ('IIIIIIIIIIIIIIIIIIIIIIIIIIIZZI', 1), ('IIIIIIIIIIIIIIIIIIIIIIIIIIIIZZ', 1), ('XXIIIIIIIIIIIIIIIIIIIIIIIIIIII', 1), ('IXXIIIIIIIIIIIIIIIIIIIIIIIIIII', 1), ('IIXXIIIIIIIIIIIIIIIIIIIIIIIIII', 1), ('IIIXXIIIIIIIIIIIIIIIIIIIIIIIII', 1), ('IIIIXXIIIIIIIIIIIIIIIIIIIIIIII', 1), ('IIIIIXXIIIIIIIIIIIIIIIIIIIIIII', 1), ('IIIIIIXXIIIIIIIIIIIIIIIIIIIIII', 1), ('IIIIIIIXXIIIIIIIIIIIIIIIIIIIII', 1), ('IIIIIIIIXXIIIIIIIIIIIIIIIIIIII', 1), ('IIIIIIIIIXXIIIIIIIIIIIIIIIIIII', 1), ('IIIIIIIIIIXXIIIIIIIIIIIIIIIIII', 1), ('IIIIIIIIIIIXXIIIIIIIIIIIIIIIII', 1), ('IIIIIIIIIIIIXXIIIIIIIIIIIIIIII', 1), ('IIIIIIIIIIIIIXXIIIIIIIIIIIIIII', 1), ('IIIIIIIIIIIIIIXXIIIIIIIIIIIIII', 1), ('IIIIIIIIIIIIIIIXXIIIIIIIIIIIII', 1), ('IIIIIIIIIIIIIIIIXXIIIIIIIIIIII', 1), ('IIIIIIIIIIIIIIIIIXXIIIIIIIIIII', 1), ('IIIIIIIIIIIIIIIIIIXXIIIIIIIIII', 1), ('IIIIIIIIIIIIIIIIIIIXXIIIIIIIII', 1), ('IIIIIIIIIIIIIIIIIIIIXXIIIIIIII', 1), ('IIIIIIIIIIIIIIIIIIIIIXXIIIIIII', 1), ('IIIIIIIIIIIIIIIIIIIIIIXXIIIIII', 1), ('IIIIIIIIIIIIIIIIIIIIIIIXXIIIII', 1), ('IIIIIIIIIIIIIIIIIIIIIIIIXXIIII', 1), ('IIIIIIIIIIIIIIIIIIIIIIIIIXXIII', 1), ('IIIIIIIIIIIIIIIIIIIIIIIIIIXXII', 1), ('IIIIIIIIIIIIIIIIIIIIIIIIIIIXXI', 1), ('IIIIIIIIIIIIIIIIIIIIIIIIIIIIXX', 1), ('YYIIIIIIIIIIIIIIIIIIIIIIIIIIII', 1), ('IYYIIIIIIIIIIIIIIIIIIIIIIIIIII', 1), ('IIYYIIIIIIIIIIIIIIIIIIIIIIIIII', 1), ('IIIYYIIIIIIIIIIIIIIIIIIIIIIIII', 1), ('IIIIYYIIIIIIIIIIIIIIIIIIIIIIII', 1), ('IIIIIYYIIIIIIIIIIIIIIIIIIIIIII', 1), ('IIIIIIYYIIIIIIIIIIIIIIIIIIIIII', 1), ('IIIIIIIYYIIIIIIIIIIIIIIIIIIIII', 1), ('IIIIIIIIYYIIIIIIIIIIIIIIIIIIII', 1), ('IIIIIIIIIYYIIIIIIIIIIIIIIIIIII', 1), ('IIIIIIIIIIYYIIIIIIIIIIIIIIIIII', 1), ('IIIIIIIIIIIYYIIIIIIIIIIIIIIIII', 1), ('IIIIIIIIIIIIYYIIIIIIIIIIIIIIII', 1), ('IIIIIIIIIIIIIYYIIIIIIIIIIIIIII', 1), ('IIIIIIIIIIIIIIYYIIIIIIIIIIIIII', 1), ('IIIIIIIIIIIIIIIYYIIIIIIIIIIIII', 1), ('IIIIIIIIIIIIIIIIYYIIIIIIIIIIII', 1), ('IIIIIIIIIIIIIIIIIYYIIIIIIIIIII', 1), ('IIIIIIIIIIIIIIIIIIYYIIIIIIIIII', 1), ('IIIIIIIIIIIIIIIIIIIYYIIIIIIIII', 1), ('IIIIIIIIIIIIIIIIIIIIYYIIIIIIII', 1), ('IIIIIIIIIIIIIIIIIIIIIYYIIIIIII', 1), ('IIIIIIIIIIIIIIIIIIIIIIYYIIIIII', 1), ('IIIIIIIIIIIIIIIIIIIIIIIYYIIIII', 1), ('IIIIIIIIIIIIIIIIIIIIIIIIYYIIII', 1), ('IIIIIIIIIIIIIIIIIIIIIIIIIYYIII', 1), ('IIIIIIIIIIIIIIIIIIIIIIIIIIYYII', 1), ('IIIIIIIIIIIIIIIIIIIIIIIIIIIYYI', 1), ('IIIIIIIIIIIIIIIIIIIIIIIIIIIIYY', 1)]\n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"# Define problem Hamiltonian.\n",
|
||
"n_qubits = 30\n",
|
||
"J = 1 # coupling strength for ZZ interaction\n",
|
||
"\n",
|
||
"# Define the Hamiltonian:\n",
|
||
"H_int = [[\"I\"] * n_qubits for _ in range(3 * (n_qubits - 1))]\n",
|
||
"for i in range(n_qubits - 1):\n",
|
||
" H_int[i][i] = \"Z\"\n",
|
||
" H_int[i][i + 1] = \"Z\"\n",
|
||
"for i in range(n_qubits - 1):\n",
|
||
" H_int[n_qubits - 1 + i][i] = \"X\"\n",
|
||
" H_int[n_qubits - 1 + i][i + 1] = \"X\"\n",
|
||
"for i in range(n_qubits - 1):\n",
|
||
" H_int[2 * (n_qubits - 1) + i][i] = \"Y\"\n",
|
||
" H_int[2 * (n_qubits - 1) + i][i + 1] = \"Y\"\n",
|
||
"H_int = [\"\".join(term) for term in H_int]\n",
|
||
"H_tot = [(term, J) if term.count(\"Z\") == 2 else (term, 1) for term in H_int]\n",
|
||
"\n",
|
||
"# Get operator\n",
|
||
"H_op = SparsePauliOp.from_list(H_tot)\n",
|
||
"print(H_tot)"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "f26578ba-24ad-42f3-baad-c348d3c05699",
|
||
"metadata": {},
|
||
"source": [
|
||
"#### Set parameters for the algorithm"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "8d376c0e-fb7a-4a41-838a-ae041d5c9afa",
|
||
"metadata": {},
|
||
"source": [
|
||
"We heuristically choose a value for the time-step `dt` (based on upper bounds on the Hamiltonian norm). Ref [\\[2\\]](#references) showed that a sufficiently small timestep is $\\pi/\\vert \\vert H \\vert \\vert$, and that it is preferable up to a point to underestimate this value rather than overestimate, since overestimating can allow contributions from high-energy states to corrupt even the optimal state in the Krylov space. On the other hand, choosing $dt$ to be too small leads to worse conditioning of the Krylov subspace, since the Krylov basis vectors differ less from timestep to timestep."
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": null,
|
||
"id": "963dd2e9-f4ea-456d-b4ed-e3d05a114082",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"data": {
|
||
"text/plain": [
|
||
"0.10833078115826872"
|
||
]
|
||
},
|
||
"execution_count": 3,
|
||
"metadata": {},
|
||
"output_type": "execute_result"
|
||
}
|
||
],
|
||
"source": [
|
||
"# Get Hamiltonian restricted to single-particle states\n",
|
||
"single_particle_H = np.zeros((n_qubits, n_qubits))\n",
|
||
"for i in range(n_qubits):\n",
|
||
" for j in range(i + 1):\n",
|
||
" for p, coeff in H_op.to_list():\n",
|
||
" p_x = Pauli(p).x\n",
|
||
" p_z = Pauli(p).z\n",
|
||
" if all(\n",
|
||
" p_x[k] == ((i == k) + (j == k)) % 2 for k in range(n_qubits)\n",
|
||
" ):\n",
|
||
" sgn = (\n",
|
||
" (-1j) ** sum(p_z[k] and p_x[k] for k in range(n_qubits))\n",
|
||
" ) * ((-1) ** p_z[i])\n",
|
||
" else:\n",
|
||
" sgn = 0\n",
|
||
" single_particle_H[i, j] += sgn * coeff\n",
|
||
"for i in range(n_qubits):\n",
|
||
" for j in range(i + 1, n_qubits):\n",
|
||
" single_particle_H[i, j] = np.conj(single_particle_H[j, i])\n",
|
||
"\n",
|
||
"# Set dt according to spectral norm\n",
|
||
"dt = np.pi / np.linalg.norm(single_particle_H, ord=2)\n",
|
||
"dt"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "d2bcc2b9-ca7c-4147-8e6b-0f15208297ce",
|
||
"metadata": {},
|
||
"source": [
|
||
"And set other parameters of the algorithm. For the sake of this tutorial, we'll limit ourselves to using a Krylov space with only 5 dimensions, which is quite limiting."
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 4,
|
||
"id": "ddde76b8-446f-4cd5-bc4b-8f00e2ec726c",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"# Set parameters for quantum Krylov algorithm\n",
|
||
"krylov_dim = 5 # size of Krylov subspace\n",
|
||
"num_trotter_steps = 6\n",
|
||
"dt_circ = dt / num_trotter_steps"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "3c0dcc08-6963-4e14-a50c-21c5dd7f8ade",
|
||
"metadata": {},
|
||
"source": [
|
||
"#### State preparation\n",
|
||
"Pick a reference state $\\vert \\psi \\rangle$ that has some overlap with the ground state. For this Hamiltonian, We use the a state with an excitation in the middle qubit $\\vert 00..010...00 \\rangle$ as our reference state."
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 5,
|
||
"id": "70161afe-8ace-4642-894a-cd21ed77a3b9",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"data": {
|
||
"text/plain": [
|
||
"<Image src=\"/docs/images/tutorials/krylov-quantum-diagonalization/extracted-outputs/70161afe-8ace-4642-894a-cd21ed77a3b9-0.avif\" alt=\"Output of the previous code cell\" />"
|
||
]
|
||
},
|
||
"execution_count": 5,
|
||
"metadata": {},
|
||
"output_type": "execute_result"
|
||
}
|
||
],
|
||
"source": [
|
||
"qc_state_prep = QuantumCircuit(n_qubits)\n",
|
||
"qc_state_prep.x(int(n_qubits / 2) + 1)\n",
|
||
"qc_state_prep.draw(\"mpl\", scale=0.5)"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "1786de6b-476f-47ea-9c20-3464d3cfbbdb",
|
||
"metadata": {},
|
||
"source": [
|
||
"#### Time evolution\n",
|
||
"\n",
|
||
"We can realize the time-evolution operator generated by a given Hamiltonian: $U=e^{-iHt}$ via the [Lie-Trotter approximation](/docs/api/qiskit/qiskit.synthesis.LieTrotter)."
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 6,
|
||
"id": "a23b9e9c-5dc8-447f-8c73-b4fa01630c8f",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"data": {
|
||
"text/plain": [
|
||
"<qiskit.circuit.instructionset.InstructionSet at 0x136639060>"
|
||
]
|
||
},
|
||
"execution_count": 6,
|
||
"metadata": {},
|
||
"output_type": "execute_result"
|
||
}
|
||
],
|
||
"source": [
|
||
"t = Parameter(\"t\")\n",
|
||
"\n",
|
||
"## Create the time-evo op circuit\n",
|
||
"evol_gate = PauliEvolutionGate(\n",
|
||
" H_op, time=t, synthesis=LieTrotter(reps=num_trotter_steps)\n",
|
||
")\n",
|
||
"\n",
|
||
"qr = QuantumRegister(n_qubits)\n",
|
||
"qc_evol = QuantumCircuit(qr)\n",
|
||
"qc_evol.append(evol_gate, qargs=qr)"
|
||
]
|
||
},
|
||
{
|
||
"attachments": {},
|
||
"cell_type": "markdown",
|
||
"id": "41c2cc43-51a2-42c3-88ee-a3b7b30eabb4",
|
||
"metadata": {},
|
||
"source": [
|
||
"#### Hadamard test\n",
|
||
"\n",
|
||
"\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\begin{equation*}\n",
|
||
" |0\\rangle|0\\rangle^N \\quad\\longrightarrow\\quad \\frac{1}{\\sqrt{2}}\\Big(|0\\rangle + |1\\rangle \\Big)|0\\rangle^N \\quad\\longrightarrow\\quad \\frac{1}{\\sqrt{2}}\\Big(|0\\rangle|0\\rangle^N+|1\\rangle |\\psi_i\\rangle\\Big) \\quad\\longrightarrow\\quad \\frac{1}{\\sqrt{2}}\\Big(|0\\rangle |0\\rangle^N+|1\\rangle P |\\psi_i\\rangle\\Big) \\quad\\longrightarrow\\quad\\frac{1}{\\sqrt{2}}\\Big(|0\\rangle |\\psi_j\\rangle+|1\\rangle P|\\psi_i\\rangle\\Big)\n",
|
||
"\\end{equation*}\n",
|
||
"$$\n",
|
||
"\n",
|
||
"Where $P$ is one of the terms in the decomposition of the Hamiltonian $H=\\sum P$ and $\\text{Prep} \\; \\psi_i$, $\\text{Prep} \\; \\psi_j$ are controlled operations that prepare $|\\psi_i\\rangle$, $|\\psi_j\\rangle$ vectors of the unitary Krylov space, with $|\\psi_k\\rangle = e^{-i H k dt } \\vert \\psi \\rangle = e^{-i H k dt } U_{\\psi} \\vert 0 \\rangle^N$. To measure $X$, first apply $H$...\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\begin{equation*}\n",
|
||
" \\longrightarrow\\quad\\frac{1}{2}|0\\rangle\\Big( |\\psi_j\\rangle + P|\\psi_i\\rangle\\Big) + \\frac{1}{2}|1\\rangle\\Big(|\\psi_j\\rangle - P|\\psi_i\\rangle\\Big)\n",
|
||
"\\end{equation*}\n",
|
||
"$$\n",
|
||
"\n",
|
||
"... then measure:\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\begin{equation*}\n",
|
||
"\\begin{split}\n",
|
||
" \\Rightarrow\\quad\\langle X\\rangle &= \\frac{1}{4}\\Bigg(\\Big\\|| \\psi_j\\rangle + P|\\psi_i\\rangle \\Big\\|^2-\\Big\\||\\psi_j\\rangle - P|\\psi_i\\rangle\\Big\\|^2\\Bigg) \\\\\n",
|
||
" &= \\text{Re}\\Big[\\langle\\psi_j| P|\\psi_i\\rangle\\Big].\n",
|
||
"\\end{split}\n",
|
||
"\\end{equation*}\n",
|
||
"$$\n",
|
||
"\n",
|
||
"From the identity $|a + b\\|^2 = \\langle a + b | a + b \\rangle = \\|a\\|^2 + \\|b\\|^2 + 2\\text{Re}\\langle a | b \\rangle$. Similarly, measuring $Y$ yields\n",
|
||
"$$\n",
|
||
"\\begin{equation*}\n",
|
||
" \\langle Y\\rangle = \\text{Im}\\Big[\\langle\\psi_j| P|\\psi_i\\rangle\\Big].\n",
|
||
"\\end{equation*}\n",
|
||
"$$"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 7,
|
||
"id": "7c1efca7-7db9-43a9-bcb7-053ba274d6f6",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
"Circuit for calculating the real part of the overlap in S via Hadamard test\n"
|
||
]
|
||
},
|
||
{
|
||
"data": {
|
||
"text/plain": [
|
||
"<Image src=\"/docs/images/tutorials/krylov-quantum-diagonalization/extracted-outputs/7c1efca7-7db9-43a9-bcb7-053ba274d6f6-1.avif\" alt=\"Output of the previous code cell\" />"
|
||
]
|
||
},
|
||
"execution_count": 7,
|
||
"metadata": {},
|
||
"output_type": "execute_result"
|
||
}
|
||
],
|
||
"source": [
|
||
"## Create the time-evo op circuit\n",
|
||
"evol_gate = PauliEvolutionGate(\n",
|
||
" H_op, time=dt, synthesis=LieTrotter(reps=num_trotter_steps)\n",
|
||
")\n",
|
||
"\n",
|
||
"## Create the time-evo op dagger circuit\n",
|
||
"evol_gate_d = PauliEvolutionGate(\n",
|
||
" H_op, time=dt, synthesis=LieTrotter(reps=num_trotter_steps)\n",
|
||
")\n",
|
||
"evol_gate_d = evol_gate_d.inverse()\n",
|
||
"\n",
|
||
"# Put pieces together\n",
|
||
"qc_reg = QuantumRegister(n_qubits)\n",
|
||
"qc_temp = QuantumCircuit(qc_reg)\n",
|
||
"qc_temp.compose(qc_state_prep, inplace=True)\n",
|
||
"for _ in range(num_trotter_steps):\n",
|
||
" qc_temp.append(evol_gate, qargs=qc_reg)\n",
|
||
"for _ in range(num_trotter_steps):\n",
|
||
" qc_temp.append(evol_gate_d, qargs=qc_reg)\n",
|
||
"qc_temp.compose(qc_state_prep.inverse(), inplace=True)\n",
|
||
"\n",
|
||
"# Create controlled version of the circuit\n",
|
||
"controlled_U = qc_temp.to_gate().control(1)\n",
|
||
"\n",
|
||
"# Create hadamard test circuit for real part\n",
|
||
"qr = QuantumRegister(n_qubits + 1)\n",
|
||
"qc_real = QuantumCircuit(qr)\n",
|
||
"qc_real.h(0)\n",
|
||
"qc_real.append(controlled_U, list(range(n_qubits + 1)))\n",
|
||
"qc_real.h(0)\n",
|
||
"\n",
|
||
"print(\n",
|
||
" \"Circuit for calculating the real part of the overlap in S via Hadamard test\"\n",
|
||
")\n",
|
||
"qc_real.draw(\"mpl\", fold=-1, scale=0.5)"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "c72fd5e5-a586-498e-bbc1-0f48ecea49d5",
|
||
"metadata": {},
|
||
"source": [
|
||
"The Hadamard test circuit can be a deep circuit once we decompose to native gates (which will increase even more if we account for the topology of the device)"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 8,
|
||
"id": "06f22c05-6e1c-4540-9dac-884b3400a50e",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
"Number of layers of 2Q operations 120003\n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"print(\n",
|
||
" \"Number of layers of 2Q operations\",\n",
|
||
" qc_real.decompose(reps=2).depth(lambda x: x[0].num_qubits == 2),\n",
|
||
")"
|
||
]
|
||
},
|
||
{
|
||
"attachments": {},
|
||
"cell_type": "markdown",
|
||
"id": "8e15903d-d4a8-4fbe-8f9e-e960e686629e",
|
||
"metadata": {},
|
||
"source": [
|
||
"## Step 2: Optimize problem for quantum hardware execution\n",
|
||
"\n",
|
||
"### Efficient Hadamard test\n",
|
||
"We can optimize the deep circuits for the Hadamard test that we have obtained by introducing some approximations and relying on some assumption about the model Hamiltonian. For example, consider the following circuit for the Hadamard test:\n",
|
||
"\n",
|
||
"\n",
|
||
"\n",
|
||
"Assume we can classically calculate $E_0$, the eigenvalue of $|0\\rangle^N$ under the Hamiltonian $H$.\n",
|
||
"This is satisfied when the Hamiltonian preserves the U(1) symmetry. Although this may seem like a strong assumption, there are many cases where it is safe to assume that there is a vacuum state (in this case it maps to the $|0\\rangle^N$ state) which is unaffected by the action of the Hamiltonian. This is true for example for chemistry Hamiltonians that describe stable molecule (where the number of electrons is conserved).\n",
|
||
"Given that the gate $\\text{Prep} \\; \\psi$, prepares the desired reference state $\\ket{psi} = \\text{Prep} \\; \\psi \\ket{0} = e^{-i H 0 dt} U_{\\psi} \\ket{0}$, e.g., to prepare the HF state for chemistry $\\text{Prep} \\; \\psi$ would be a product of single-qubit NOTs, so controlled-$\\text{Prep} \\; \\psi$ is just a product of CNOTs.\n",
|
||
"Then the circuit above implements the following state prior to measurement:\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\begin{equation}\n",
|
||
"\\begin{split}\n",
|
||
" \\ket{0} \\ket{0}^N\\xrightarrow{H}&\\frac{1}{\\sqrt{2}}\n",
|
||
" \\left(\n",
|
||
" \\ket{0}\\ket{0}^N+ \\ket{1} \\ket{0}^N\n",
|
||
" \\right)\\\\\n",
|
||
" \\xrightarrow{\\text{1-ctrl-init}}&\\frac{1}{\\sqrt{2}}\\left(|0\\rangle|0\\rangle^N+|1\\rangle|\\psi\\rangle\\right)\\\\\n",
|
||
" \\xrightarrow{U}&\\frac{1}{\\sqrt{2}}\\left(e^{i\\phi}\\ket{0}\\ket{0}^N+\\ket{1} U\\ket{\\psi}\\right)\\\\\n",
|
||
" \\xrightarrow{\\text{0-ctrl-init}}&\\frac{1}{\\sqrt{2}}\n",
|
||
" \\left(\n",
|
||
" e^{i\\phi}\\ket{0} \\ket{\\psi}\n",
|
||
" +\\ket{1} U\\ket{\\psi}\n",
|
||
" \\right)\\\\\n",
|
||
" =&\\frac{1}{2}\n",
|
||
" \\left(\n",
|
||
" \\ket{+}\\left(e^{i\\phi}\\ket{\\psi}+U\\ket{\\psi}\\right)\n",
|
||
" +\\ket{-}\\left(e^{i\\phi}\\ket{\\psi}-U\\ket{\\psi}\\right)\n",
|
||
" \\right)\\\\\n",
|
||
" =&\\frac{1}{2}\n",
|
||
" \\left(\n",
|
||
" \\ket{+i}\\left(e^{i\\phi}\\ket{\\psi}-iU\\ket{\\psi}\\right)\n",
|
||
" +\\ket{-i}\\left(e^{i\\phi}\\ket{\\psi}+iU\\ket{\\psi}\\right)\n",
|
||
" \\right)\n",
|
||
"\\end{split}\n",
|
||
"\\end{equation}\n",
|
||
"$$\n",
|
||
"\n",
|
||
"where we have used the classical simulable phase shift $ U\\ket{0}^N = e^{i\\phi}\\ket{0}^N$ in the third line. Therefore the expectation values are obtained as\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\begin{equation}\n",
|
||
"\\begin{split}\n",
|
||
" \\langle X\\otimes P\\rangle&=\\frac{1}{4}\n",
|
||
" \\Big(\n",
|
||
" \\left(e^{-i\\phi}\\bra{\\psi}+\\bra{\\psi}U^\\dagger\\right)P\\left(e^{i\\phi}\\ket{\\psi}+U\\ket{\\psi}\\right)\n",
|
||
" \\\\\n",
|
||
" &\\qquad-\\left(e^{-i\\phi}\\bra{\\psi}-\\bra{\\psi}U^\\dagger\\right)P\\left(e^{i\\phi}\\ket{\\psi}-U\\ket{\\psi}\\right)\n",
|
||
" \\Big)\\\\\n",
|
||
" &=\\text{Re}\\left[e^{-i\\phi}\\bra{\\psi}PU\\ket{\\psi}\\right],\n",
|
||
"\\end{split}\n",
|
||
"\\end{equation}\n",
|
||
"$$\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\begin{equation}\n",
|
||
"\\begin{split}\n",
|
||
" \\langle Y\\otimes P\\rangle&=\\frac{1}{4}\n",
|
||
" \\Big(\n",
|
||
" \\left(e^{-i\\phi}\\bra{\\psi}+i\\bra{\\psi}U^\\dagger\\right)P\\left(e^{i\\phi}\\ket{\\psi}-iU\\ket{\\psi}\\right)\n",
|
||
" \\\\\n",
|
||
" &\\qquad-\\left(e^{-i\\phi}\\bra{\\psi}-i\\bra{\\psi}U^\\dagger\\right)P\\left(e^{i\\phi}\\ket{\\psi}+iU\\ket{\\psi}\\right)\n",
|
||
" \\Big)\\\\\n",
|
||
" &=\\text{Im}\\left[e^{-i\\phi}\\bra{\\psi}PU\\ket{\\psi}\\right].\n",
|
||
"\\end{split}\n",
|
||
"\\end{equation}\n",
|
||
"$$\n",
|
||
"\n",
|
||
"Using these assumptions we were able to write the expectation values of operators of interest with fewer controlled operations. In fact, we only need to implement the controlled state preparation $\\text{Prep} \\; \\psi$ and not controlled time evolutions. Reframing our calculation as above will allow us to greatly reduce the depth of the resulting circuits."
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "2e29da5f-b5f5-4b9e-96fc-3fa7ab698398",
|
||
"metadata": {},
|
||
"source": [
|
||
"### Decompose time-evolution operator with Trotter decomposition\n",
|
||
"Instead of implementing the time-evolution operator exactly we can use the Trotter decomposition to implement an approximation of it. Repeating several times a certain order Trotter decomposition gives us further reduction of the error introduced from the approximation. In the following, we directly build the Trotter implementation in the most efficient way for the interaction graph of the Hamiltonian we are considering (nearest neighbor interactions only). In practice we insert Pauli rotations $R_{xx}$, $R_{yy}$, $R_{zz}$ with a parametrized angle $t$ which correspond to the approximate implementation of $e^{-i (XX + YY + ZZ) t}$. Given the difference in definition of the Pauli rotations and the time-evolution that we are trying to implement, we'll have to use the parameter $2*dt$ to achieve a time-evolution of $dt$. Furthermore, we reverse the order of the operations for odd number of repetitions of the Trotter steps, which is functionally equivalent but allows for synthesizing adjacent operations in a single $SU(2)$ unitary. This gives a much shallower circuit than what is obtained using the generic `PauliEvolutionGate()` functionality."
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 7,
|
||
"id": "267716dc-fa23-41bd-abe4-6d4e0499a0f4",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"data": {
|
||
"text/plain": [
|
||
"<Image src=\"/docs/images/tutorials/krylov-quantum-diagonalization/extracted-outputs/267716dc-fa23-41bd-abe4-6d4e0499a0f4-0.avif\" alt=\"Output of the previous code cell\" />"
|
||
]
|
||
},
|
||
"execution_count": 7,
|
||
"metadata": {},
|
||
"output_type": "execute_result"
|
||
}
|
||
],
|
||
"source": [
|
||
"t = Parameter(\"t\")\n",
|
||
"\n",
|
||
"# Create instruction for rotation about XX+YY-ZZ:\n",
|
||
"Rxyz_circ = QuantumCircuit(2)\n",
|
||
"Rxyz_circ.rxx(t, 0, 1)\n",
|
||
"Rxyz_circ.ryy(t, 0, 1)\n",
|
||
"Rxyz_circ.rzz(t, 0, 1)\n",
|
||
"Rxyz_instr = Rxyz_circ.to_instruction(label=\"RXX+YY+ZZ\")\n",
|
||
"\n",
|
||
"interaction_list = [\n",
|
||
" [[i, i + 1] for i in range(0, n_qubits - 1, 2)],\n",
|
||
" [[i, i + 1] for i in range(1, n_qubits - 1, 2)],\n",
|
||
"] # linear chain\n",
|
||
"\n",
|
||
"qr = QuantumRegister(n_qubits)\n",
|
||
"trotter_step_circ = QuantumCircuit(qr)\n",
|
||
"for i, color in enumerate(interaction_list):\n",
|
||
" for interaction in color:\n",
|
||
" trotter_step_circ.append(Rxyz_instr, interaction)\n",
|
||
" if i < len(interaction_list) - 1:\n",
|
||
" trotter_step_circ.barrier()\n",
|
||
"reverse_trotter_step_circ = trotter_step_circ.reverse_ops()\n",
|
||
"\n",
|
||
"qc_evol = QuantumCircuit(qr)\n",
|
||
"for step in range(num_trotter_steps):\n",
|
||
" if step % 2 == 0:\n",
|
||
" qc_evol = qc_evol.compose(trotter_step_circ)\n",
|
||
" else:\n",
|
||
" qc_evol = qc_evol.compose(reverse_trotter_step_circ)\n",
|
||
"\n",
|
||
"qc_evol.decompose().draw(\"mpl\", fold=-1, scale=0.5)"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "0a9c3a4d-a678-41e4-b9dd-995fe34341fb",
|
||
"metadata": {},
|
||
"source": [
|
||
"### Use an optimized circuit for state preparation"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 8,
|
||
"id": "70411715-eed3-4cf5-961d-06a6f1e04efc",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"data": {
|
||
"text/plain": [
|
||
"<Image src=\"/docs/images/tutorials/krylov-quantum-diagonalization/extracted-outputs/70411715-eed3-4cf5-961d-06a6f1e04efc-0.avif\" alt=\"Output of the previous code cell\" />"
|
||
]
|
||
},
|
||
"execution_count": 8,
|
||
"metadata": {},
|
||
"output_type": "execute_result"
|
||
}
|
||
],
|
||
"source": [
|
||
"control = 0\n",
|
||
"excitation = int(n_qubits / 2) + 1\n",
|
||
"controlled_state_prep = QuantumCircuit(n_qubits + 1)\n",
|
||
"controlled_state_prep.cx(control, excitation)\n",
|
||
"controlled_state_prep.draw(\"mpl\", fold=-1, scale=0.5)"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "2421ff11-d97a-488a-a382-a652df8c94d6",
|
||
"metadata": {},
|
||
"source": [
|
||
"### Template circuits for calculating matrix elements of $\\tilde{S}$ and $\\tilde{H}$ via Hadamard test\n",
|
||
"The only difference between the circuits used in the Hadamard test will be the phase in the time-evolution operator and the observables measured. Therefore we can prepare a template circuit which represent the generic circuit for the Hadamard test, with placeholders for the gates that depend on the time-evolution operator."
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 9,
|
||
"id": "2f102112-4ddc-41ea-999c-db5863bc77ac",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"# Parameters for the template circuits\n",
|
||
"parameters = []\n",
|
||
"for idx in range(1, krylov_dim):\n",
|
||
" parameters.append(2 * dt_circ * (idx))"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 10,
|
||
"id": "33ec7c29-904e-4445-a654-405214349a4d",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"data": {
|
||
"text/plain": [
|
||
"<Image src=\"/docs/images/tutorials/krylov-quantum-diagonalization/extracted-outputs/33ec7c29-904e-4445-a654-405214349a4d-0.avif\" alt=\"Output of the previous code cell\" />"
|
||
]
|
||
},
|
||
"execution_count": 10,
|
||
"metadata": {},
|
||
"output_type": "execute_result"
|
||
}
|
||
],
|
||
"source": [
|
||
"# Create modified hadamard test circuit\n",
|
||
"qr = QuantumRegister(n_qubits + 1)\n",
|
||
"qc = QuantumCircuit(qr)\n",
|
||
"qc.h(0)\n",
|
||
"qc.compose(controlled_state_prep, list(range(n_qubits + 1)), inplace=True)\n",
|
||
"qc.barrier()\n",
|
||
"qc.compose(qc_evol, list(range(1, n_qubits + 1)), inplace=True)\n",
|
||
"qc.barrier()\n",
|
||
"qc.x(0)\n",
|
||
"qc.compose(\n",
|
||
" controlled_state_prep.inverse(), list(range(n_qubits + 1)), inplace=True\n",
|
||
")\n",
|
||
"qc.x(0)\n",
|
||
"\n",
|
||
"qc.decompose().draw(\"mpl\", fold=-1)"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 13,
|
||
"id": "157356c9-06bb-411c-87ef-cd1d2f73be6f",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
"The optimized circuit has 2Q gates depth: 74\n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"print(\n",
|
||
" \"The optimized circuit has 2Q gates depth: \",\n",
|
||
" qc.decompose().decompose().depth(lambda x: x[0].num_qubits == 2),\n",
|
||
")"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "5a59bbc4-72c7-4a55-9a9e-57bb7b469021",
|
||
"metadata": {},
|
||
"source": [
|
||
"We have considerably reduced the depth of the Hadamard test with a combination of Trotter approximation and uncontrolled unitaries"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "ca78bfe3-684f-4d3a-b8c2-3d4e55c9ec30",
|
||
"metadata": {},
|
||
"source": [
|
||
"## Step 3: Execute using Qiskit primitives"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "a082f023-d752-40e4-b693-c4d0da5ba102",
|
||
"metadata": {},
|
||
"source": [
|
||
"Instantiate the backend and set runtime parameters"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": null,
|
||
"id": "0d90e4df-e262-4852-a811-ff6a1d3232ae",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"service = QiskitRuntimeService()\n",
|
||
"backend = service.least_busy(\n",
|
||
" operational=True, simulator=False, min_num_qubits=127\n",
|
||
")"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "2b3f7227-e6ef-4a11-93d6-d3696054b9cc",
|
||
"metadata": {},
|
||
"source": [
|
||
"### Transpiling to a QPU\n",
|
||
"First, let's pick subsets of the coupling map with \"good\" performing qubits (where \"good\" is pretty arbitraty here, we mostly want to avoid really poor performing qubits) and create a new target for transpilation"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": null,
|
||
"id": "e123fda1-6454-4893-9612-2b8591e8cfb9",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"target = backend.target\n",
|
||
"cmap = target.build_coupling_map(filter_idle_qubits=True)\n",
|
||
"cmap_list = list(cmap.get_edges())\n",
|
||
"\n",
|
||
"cust_cmap_list = copy.deepcopy(cmap_list)\n",
|
||
"for q in range(target.num_qubits):\n",
|
||
" meas_err = target[\"measure\"][(q,)].error\n",
|
||
" t2 = target.qubit_properties[q].t2 * 1e6\n",
|
||
" if meas_err > 0.05 or t2 < 30:\n",
|
||
" # print(q)\n",
|
||
" for q_pair in cmap_list:\n",
|
||
" if q in q_pair:\n",
|
||
" try:\n",
|
||
" cust_cmap_list.remove(q_pair)\n",
|
||
" except:\n",
|
||
" continue\n",
|
||
"\n",
|
||
"for q in cmap_list:\n",
|
||
" twoq_gate_err = target[\"cz\"][q].error\n",
|
||
" if twoq_gate_err > 0.015:\n",
|
||
" for q_pair in cmap_list:\n",
|
||
" if q == q_pair:\n",
|
||
" try:\n",
|
||
" cust_cmap_list.remove(q)\n",
|
||
" except:\n",
|
||
" continue\n",
|
||
"\n",
|
||
"\n",
|
||
"cust_cmap = CouplingMap(cust_cmap_list)\n",
|
||
"cust_target = Target.from_configuration(\n",
|
||
" basis_gates=backend.configuration().basis_gates,\n",
|
||
" coupling_map=cust_cmap,\n",
|
||
")"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "96f79e35-a265-4de4-b98c-f9c84f58f0d5",
|
||
"metadata": {},
|
||
"source": [
|
||
"Then transpile the virtual circuit to the best physical layout in this new target"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": null,
|
||
"id": "62109c95-79fe-4076-b591-9bd920fd51f4",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"basis_gates = list(target.operation_names)\n",
|
||
"pm = generate_preset_pass_manager(\n",
|
||
" optimization_level=3,\n",
|
||
" target=cust_target,\n",
|
||
" basis_gates=basis_gates,\n",
|
||
")\n",
|
||
"\n",
|
||
"qc_trans = pm.run(qc)\n",
|
||
"\n",
|
||
"print(\"depth\", qc_trans.depth(lambda x: x[0].num_qubits == 2))\n",
|
||
"print(\"num 2q ops\", qc_trans.count_ops())\n",
|
||
"print(\n",
|
||
" \"physical qubits\",\n",
|
||
" sorted(\n",
|
||
" [\n",
|
||
" idx\n",
|
||
" for idx, qb in qc_trans.layout.initial_layout.get_physical_bits().items()\n",
|
||
" if qb._register.name != \"ancilla\"\n",
|
||
" ]\n",
|
||
" ),\n",
|
||
")"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "9b9b0170-b96a-47b9-b5ab-5112d02338a3",
|
||
"metadata": {},
|
||
"source": [
|
||
"### Create PUBs for execution with Estimator"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": null,
|
||
"id": "bcfd6693-d1fb-44e4-9b06-c70e6766e877",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"# Define observables to measure for S\n",
|
||
"observable_S_real = \"I\" * (n_qubits) + \"X\"\n",
|
||
"observable_S_imag = \"I\" * (n_qubits) + \"Y\"\n",
|
||
"\n",
|
||
"observable_op_real = SparsePauliOp(\n",
|
||
" observable_S_real\n",
|
||
") # define a sparse pauli operator for the observable\n",
|
||
"observable_op_imag = SparsePauliOp(observable_S_imag)\n",
|
||
"\n",
|
||
"layout = qc_trans.layout # get layout of transpiled circuit\n",
|
||
"observable_op_real = observable_op_real.apply_layout(\n",
|
||
" layout\n",
|
||
") # apply physical layout to the observable\n",
|
||
"observable_op_imag = observable_op_imag.apply_layout(layout)\n",
|
||
"observable_S_real = (\n",
|
||
" observable_op_real.paulis.to_labels()\n",
|
||
") # get the label of the physical observable\n",
|
||
"observable_S_imag = observable_op_imag.paulis.to_labels()\n",
|
||
"\n",
|
||
"observables_S = [[observable_S_real], [observable_S_imag]]\n",
|
||
"\n",
|
||
"\n",
|
||
"# Define observables to measure for H\n",
|
||
"# Hamiltonian terms to measure\n",
|
||
"observable_list = []\n",
|
||
"for pauli, coeff in zip(H_op.paulis, H_op.coeffs):\n",
|
||
" # print(pauli)\n",
|
||
" observable_H_real = pauli[::-1].to_label() + \"X\"\n",
|
||
" observable_H_imag = pauli[::-1].to_label() + \"Y\"\n",
|
||
" observable_list.append([observable_H_real])\n",
|
||
" observable_list.append([observable_H_imag])\n",
|
||
"\n",
|
||
"layout = qc_trans.layout\n",
|
||
"\n",
|
||
"observable_trans_list = []\n",
|
||
"for observable in observable_list:\n",
|
||
" observable_op = SparsePauliOp(observable)\n",
|
||
" observable_op = observable_op.apply_layout(layout)\n",
|
||
" observable_trans_list.append([observable_op.paulis.to_labels()])\n",
|
||
"\n",
|
||
"observables_H = observable_trans_list\n",
|
||
"\n",
|
||
"\n",
|
||
"# Define a sweep over parameter values\n",
|
||
"params = np.vstack(parameters).T\n",
|
||
"\n",
|
||
"\n",
|
||
"# Estimate the expectation value for all combinations of\n",
|
||
"# observables and parameter values, where the pub result will have\n",
|
||
"# shape (# observables, # parameter values).\n",
|
||
"pub = (qc_trans, observables_S + observables_H, params)"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "8bd27e84-d673-4326-9004-f3be5a65aed3",
|
||
"metadata": {},
|
||
"source": [
|
||
"### Run circuits\n",
|
||
"Circuits for $t=0$ are classically calculable"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": null,
|
||
"id": "768cd443-2ed8-4ba2-ba9d-98491b36fa54",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
"(25+0j)\n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"qc_cliff = qc.assign_parameters({t: 0})\n",
|
||
"\n",
|
||
"\n",
|
||
"# Get expectation values from experiment\n",
|
||
"S_expval_real = StabilizerState(qc_cliff).expectation_value(\n",
|
||
" Pauli(\"I\" * (n_qubits) + \"X\")\n",
|
||
")\n",
|
||
"S_expval_imag = StabilizerState(qc_cliff).expectation_value(\n",
|
||
" Pauli(\"I\" * (n_qubits) + \"Y\")\n",
|
||
")\n",
|
||
"\n",
|
||
"# Get expectation values\n",
|
||
"S_expval = S_expval_real + 1j * S_expval_imag\n",
|
||
"\n",
|
||
"H_expval = 0\n",
|
||
"for obs_idx, (pauli, coeff) in enumerate(zip(H_op.paulis, H_op.coeffs)):\n",
|
||
" # Get expectation values from experiment\n",
|
||
" expval_real = StabilizerState(qc_cliff).expectation_value(\n",
|
||
" Pauli(pauli[::-1].to_label() + \"X\")\n",
|
||
" )\n",
|
||
" expval_imag = StabilizerState(qc_cliff).expectation_value(\n",
|
||
" Pauli(pauli[::-1].to_label() + \"Y\")\n",
|
||
" )\n",
|
||
" expval = expval_real + 1j * expval_imag\n",
|
||
"\n",
|
||
" # Fill-in matrix elements\n",
|
||
" H_expval += coeff * expval\n",
|
||
"\n",
|
||
"\n",
|
||
"print(H_expval)"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "794edf4c-d539-4864-865d-52a88a450b5c",
|
||
"metadata": {},
|
||
"source": [
|
||
"Execute circuits for $S$ and $\\tilde{H}$ with the Estimator"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 20,
|
||
"id": "9059dc9e-9203-4572-92c4-76e9c85dfed9",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"# Experiment options\n",
|
||
"num_randomizations = 300\n",
|
||
"num_randomizations_learning = 30\n",
|
||
"shots_per_randomization = 100\n",
|
||
"noise_factors = [1, 1.2, 1.4]\n",
|
||
"learning_pair_depths = [0, 4, 24, 48]\n",
|
||
"\n",
|
||
"\n",
|
||
"experimental_opts = {}\n",
|
||
"experimental_opts[\"execution_path\"] = \"gen3-turbo\"\n",
|
||
"experimental_opts[\"resilience\"] = {\n",
|
||
" \"measure_mitigation\": True,\n",
|
||
" \"measure_noise_learning\": {\n",
|
||
" \"num_randomizations\": num_randomizations_learning,\n",
|
||
" \"shots_per_randomization\": shots_per_randomization,\n",
|
||
" },\n",
|
||
" \"zne_mitigation\": True,\n",
|
||
" \"zne\": {\"noise_factors\": noise_factors},\n",
|
||
" \"layer_noise_learning\": {\n",
|
||
" \"max_layers_to_learn\": 10,\n",
|
||
" \"layer_pair_depths\": learning_pair_depths,\n",
|
||
" \"shots_per_randomization\": shots_per_randomization,\n",
|
||
" \"num_randomizations\": num_randomizations_learning,\n",
|
||
" },\n",
|
||
" \"zne\": {\n",
|
||
" \"amplifier\": \"pea\",\n",
|
||
" \"return_all_extrapolated\": True,\n",
|
||
" \"return_unextrapolated\": True,\n",
|
||
" \"extrapolated_noise_factors\": [0] + noise_factors,\n",
|
||
" },\n",
|
||
"}\n",
|
||
"experimental_opts[\"twirling\"] = {\n",
|
||
" \"num_randomizations\": num_randomizations,\n",
|
||
" \"shots_per_randomization\": shots_per_randomization,\n",
|
||
" \"strategy\": \"all\",\n",
|
||
" # 'strategy':'active-accum'\n",
|
||
"}\n",
|
||
"\n",
|
||
"options = EstimatorOptions(experimental=experimental_opts)\n",
|
||
"estimator = Estimator(mode=backend, options=options)\n",
|
||
"\n",
|
||
"\n",
|
||
"job = estimator.run([pub])"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "9524e609-b457-42e9-9de7-8dbd4464ac74",
|
||
"metadata": {},
|
||
"source": [
|
||
"## Step 4: Post-process and return result in desired classical format"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 46,
|
||
"id": "30d72032-02e4-40e7-af28-cadda2f7f3bd",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"results = job.result()[0]"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "8149d431-f069-4567-a0be-0cc4ed8d523b",
|
||
"metadata": {},
|
||
"source": [
|
||
"### Calculate Effective Hamiltonian and Overlap matrices\n",
|
||
"First calculate the phase accumulated by the $\\vert 0 \\rangle$ state during the uncontrolled time evolution"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 47,
|
||
"id": "34ba8b05-d52c-4b0f-8d28-acaf7a5ddffc",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"prefactors = [\n",
|
||
" np.exp(-1j * sum([c for p, c in H_op.to_list() if \"Z\" in p]) * i * dt)\n",
|
||
" for i in range(1, krylov_dim)\n",
|
||
"]"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "36875c8b-36c8-464b-9f0e-7a79e17d5fef",
|
||
"metadata": {},
|
||
"source": [
|
||
"Once we have the results of the circuit executions we can post-process the data to calculate the matrix elements of $S$"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 48,
|
||
"id": "16dd2534-8d5f-40f7-a938-5d519475fd8d",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"# Assemble S, the overlap matrix of dimension D:\n",
|
||
"S_first_row = np.zeros(krylov_dim, dtype=complex)\n",
|
||
"S_first_row[0] = 1 + 0j\n",
|
||
"\n",
|
||
"# Add in ancilla-only measurements:\n",
|
||
"for i in range(krylov_dim - 1):\n",
|
||
" # Get expectation values from experiment\n",
|
||
" expval_real = results.data.evs[0][0][\n",
|
||
" i\n",
|
||
" ] # automatic extrapolated evs if ZNE is used\n",
|
||
" expval_imag = results.data.evs[1][0][\n",
|
||
" i\n",
|
||
" ] # automatic extrapolated evs if ZNE is used\n",
|
||
"\n",
|
||
" # Get expectation values\n",
|
||
" expval = expval_real + 1j * expval_imag\n",
|
||
" S_first_row[i + 1] += prefactors[i] * expval\n",
|
||
"\n",
|
||
"S_first_row_list = S_first_row.tolist() # for saving purposes\n",
|
||
"\n",
|
||
"\n",
|
||
"S_circ = np.zeros((krylov_dim, krylov_dim), dtype=complex)\n",
|
||
"\n",
|
||
"# Distribute entries from first row across matrix:\n",
|
||
"for i, j in it.product(range(krylov_dim), repeat=2):\n",
|
||
" if i >= j:\n",
|
||
" S_circ[j, i] = S_first_row[i - j]\n",
|
||
" else:\n",
|
||
" S_circ[j, i] = np.conj(S_first_row[j - i])"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 49,
|
||
"id": "01c6563d-87ca-4bd6-b487-dcdaece2d8c2",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"data": {
|
||
"text/latex": [
|
||
"$\\displaystyle \\left[\\begin{matrix}1.0 & -0.723052998582984 - 0.345085413575966 i & 0.467051960502366 + 0.516197865254034 i & -0.180546747798251 - 0.492624093654174 i & 0.0012070853532697 + 0.312052218182462 i\\\\-0.723052998582984 + 0.345085413575966 i & 1.0 & -0.723052998582984 - 0.345085413575966 i & 0.467051960502366 + 0.516197865254034 i & -0.180546747798251 - 0.492624093654174 i\\\\0.467051960502366 - 0.516197865254034 i & -0.723052998582984 + 0.345085413575966 i & 1.0 & -0.723052998582984 - 0.345085413575966 i & 0.467051960502366 + 0.516197865254034 i\\\\-0.180546747798251 + 0.492624093654174 i & 0.467051960502366 - 0.516197865254034 i & -0.723052998582984 + 0.345085413575966 i & 1.0 & -0.723052998582984 - 0.345085413575966 i\\\\0.0012070853532697 - 0.312052218182462 i & -0.180546747798251 + 0.492624093654174 i & 0.467051960502366 - 0.516197865254034 i & -0.723052998582984 + 0.345085413575966 i & 1.0\\end{matrix}\\right]$"
|
||
],
|
||
"text/plain": [
|
||
"Matrix([\n",
|
||
"[ 1.0, -0.723052998582984 - 0.345085413575966*I, 0.467051960502366 + 0.516197865254034*I, -0.180546747798251 - 0.492624093654174*I, 0.0012070853532697 + 0.312052218182462*I],\n",
|
||
"[-0.723052998582984 + 0.345085413575966*I, 1.0, -0.723052998582984 - 0.345085413575966*I, 0.467051960502366 + 0.516197865254034*I, -0.180546747798251 - 0.492624093654174*I],\n",
|
||
"[ 0.467051960502366 - 0.516197865254034*I, -0.723052998582984 + 0.345085413575966*I, 1.0, -0.723052998582984 - 0.345085413575966*I, 0.467051960502366 + 0.516197865254034*I],\n",
|
||
"[-0.180546747798251 + 0.492624093654174*I, 0.467051960502366 - 0.516197865254034*I, -0.723052998582984 + 0.345085413575966*I, 1.0, -0.723052998582984 - 0.345085413575966*I],\n",
|
||
"[0.0012070853532697 - 0.312052218182462*I, -0.180546747798251 + 0.492624093654174*I, 0.467051960502366 - 0.516197865254034*I, -0.723052998582984 + 0.345085413575966*I, 1.0]])"
|
||
]
|
||
},
|
||
"execution_count": 49,
|
||
"metadata": {},
|
||
"output_type": "execute_result"
|
||
}
|
||
],
|
||
"source": [
|
||
"Matrix(S_circ)"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "a34c72c0-2984-47fb-8579-089ea795ef39",
|
||
"metadata": {},
|
||
"source": [
|
||
"And the matrix elements of $\\tilde{H}$"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 50,
|
||
"id": "9cde2419-8a29-4d7a-b5b1-9ef2d4c126d2",
|
||
"metadata": {},
|
||
"outputs": [],
|
||
"source": [
|
||
"# Assemble S, the overlap matrix of dimension D:\n",
|
||
"H_first_row = np.zeros(krylov_dim, dtype=complex)\n",
|
||
"H_first_row[0] = H_expval\n",
|
||
"\n",
|
||
"for obs_idx, (pauli, coeff) in enumerate(zip(H_op.paulis, H_op.coeffs)):\n",
|
||
" # Add in ancilla-only measurements:\n",
|
||
" for i in range(krylov_dim - 1):\n",
|
||
" # Get expectation values from experiment\n",
|
||
" expval_real = results.data.evs[2 + 2 * obs_idx][0][\n",
|
||
" i\n",
|
||
" ] # automatic extrapolated evs if ZNE is used\n",
|
||
" expval_imag = results.data.evs[2 + 2 * obs_idx + 1][0][\n",
|
||
" i\n",
|
||
" ] # automatic extrapolated evs if ZNE is used\n",
|
||
"\n",
|
||
" # Get expectation values\n",
|
||
" expval = expval_real + 1j * expval_imag\n",
|
||
" H_first_row[i + 1] += prefactors[i] * coeff * expval\n",
|
||
"\n",
|
||
"H_first_row_list = H_first_row.tolist()\n",
|
||
"\n",
|
||
"H_eff_circ = np.zeros((krylov_dim, krylov_dim), dtype=complex)\n",
|
||
"\n",
|
||
"# Distribute entries from first row across matrix:\n",
|
||
"for i, j in it.product(range(krylov_dim), repeat=2):\n",
|
||
" if i >= j:\n",
|
||
" H_eff_circ[j, i] = H_first_row[i - j]\n",
|
||
" else:\n",
|
||
" H_eff_circ[j, i] = np.conj(H_first_row[j - i])"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 51,
|
||
"id": "5800d892-5554-4a92-bfe0-3750ef0a3ed7",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"data": {
|
||
"text/latex": [
|
||
"$\\displaystyle \\left[\\begin{matrix}25.0 & -14.2437089383409 - 6.50486277982165 i & 10.2857217968584 + 9.0431912203186 i & -5.15587257589417 - 8.88280836036843 i & 1.98818301405581 + 5.8897614762563 i\\\\-14.2437089383409 + 6.50486277982165 i & 25.0 & -14.2437089383409 - 6.50486277982165 i & 10.2857217968584 + 9.0431912203186 i & -5.15587257589417 - 8.88280836036843 i\\\\10.2857217968584 - 9.0431912203186 i & -14.2437089383409 + 6.50486277982165 i & 25.0 & -14.2437089383409 - 6.50486277982165 i & 10.2857217968584 + 9.0431912203186 i\\\\-5.15587257589417 + 8.88280836036843 i & 10.2857217968584 - 9.0431912203186 i & -14.2437089383409 + 6.50486277982165 i & 25.0 & -14.2437089383409 - 6.50486277982165 i\\\\1.98818301405581 - 5.8897614762563 i & -5.15587257589417 + 8.88280836036843 i & 10.2857217968584 - 9.0431912203186 i & -14.2437089383409 + 6.50486277982165 i & 25.0\\end{matrix}\\right]$"
|
||
],
|
||
"text/plain": [
|
||
"Matrix([\n",
|
||
"[ 25.0, -14.2437089383409 - 6.50486277982165*I, 10.2857217968584 + 9.0431912203186*I, -5.15587257589417 - 8.88280836036843*I, 1.98818301405581 + 5.8897614762563*I],\n",
|
||
"[-14.2437089383409 + 6.50486277982165*I, 25.0, -14.2437089383409 - 6.50486277982165*I, 10.2857217968584 + 9.0431912203186*I, -5.15587257589417 - 8.88280836036843*I],\n",
|
||
"[ 10.2857217968584 - 9.0431912203186*I, -14.2437089383409 + 6.50486277982165*I, 25.0, -14.2437089383409 - 6.50486277982165*I, 10.2857217968584 + 9.0431912203186*I],\n",
|
||
"[-5.15587257589417 + 8.88280836036843*I, 10.2857217968584 - 9.0431912203186*I, -14.2437089383409 + 6.50486277982165*I, 25.0, -14.2437089383409 - 6.50486277982165*I],\n",
|
||
"[ 1.98818301405581 - 5.8897614762563*I, -5.15587257589417 + 8.88280836036843*I, 10.2857217968584 - 9.0431912203186*I, -14.2437089383409 + 6.50486277982165*I, 25.0]])"
|
||
]
|
||
},
|
||
"execution_count": 51,
|
||
"metadata": {},
|
||
"output_type": "execute_result"
|
||
}
|
||
],
|
||
"source": [
|
||
"Matrix(H_eff_circ)"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "896553fa-3c28-44cc-89ad-5031d2f8d35f",
|
||
"metadata": {},
|
||
"source": [
|
||
"Finally, we can solve the generalized eigenvalue problem for $\\tilde{H}$:\n",
|
||
"\n",
|
||
"$$\\tilde{H} \\vec{c} = c S \\vec{c}$$\n",
|
||
"\n",
|
||
"and get an estimate of the ground state energy $c_{min}$"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 58,
|
||
"id": "8b997d15-ee25-40eb-80e8-1f054a298ee9",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
"The estimated ground state energy is: 25.0\n",
|
||
"The estimated ground state energy is: 22.572154819954875\n",
|
||
"The estimated ground state energy is: 21.691509219286587\n",
|
||
"The estimated ground state energy is: 21.23882298756386\n",
|
||
"The estimated ground state energy is: 20.965499325470294\n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"gnd_en_circ_est_list = []\n",
|
||
"for d in range(1, krylov_dim + 1):\n",
|
||
" # Solve generalized eigenvalue problem for different size of the Krylov space\n",
|
||
" gnd_en_circ_est = solve_regularized_gen_eig(\n",
|
||
" H_eff_circ[:d, :d], S_circ[:d, :d], threshold=9e-1\n",
|
||
" )\n",
|
||
" gnd_en_circ_est_list.append(gnd_en_circ_est)\n",
|
||
" print(\"The estimated ground state energy is: \", gnd_en_circ_est)"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "b39943f8-b75f-4bd4-b246-833ea4f9b710",
|
||
"metadata": {},
|
||
"source": [
|
||
"For a single-particle sector, we can efficiently calculate the ground state of this sector of the Hamiltonian classically"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 59,
|
||
"id": "fcfe07e5-99e4-4276-a4d6-6f1c7e28c5f2",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"name": "stdout",
|
||
"output_type": "stream",
|
||
"text": [
|
||
"n_sys_qubits 30\n",
|
||
"n_exc 1 , subspace dimension 31\n",
|
||
"single particle ground state energy: 21.021912418526906\n"
|
||
]
|
||
}
|
||
],
|
||
"source": [
|
||
"gs_en = single_particle_gs(H_op, n_qubits)"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "code",
|
||
"execution_count": 60,
|
||
"id": "4bc52594-0376-497f-8a61-0949415a1fe0",
|
||
"metadata": {},
|
||
"outputs": [
|
||
{
|
||
"data": {
|
||
"text/plain": [
|
||
"<Image src=\"/docs/images/tutorials/krylov-quantum-diagonalization/extracted-outputs/4bc52594-0376-497f-8a61-0949415a1fe0-0.avif\" alt=\"Output of the previous code cell\" />"
|
||
]
|
||
},
|
||
"metadata": {},
|
||
"output_type": "display_data"
|
||
}
|
||
],
|
||
"source": [
|
||
"plt.plot(\n",
|
||
" range(1, krylov_dim + 1),\n",
|
||
" gnd_en_circ_est_list,\n",
|
||
" color=\"blue\",\n",
|
||
" linestyle=\"-.\",\n",
|
||
" label=\"KQD estimate\",\n",
|
||
")\n",
|
||
"plt.plot(\n",
|
||
" range(1, krylov_dim + 1),\n",
|
||
" [gs_en] * krylov_dim,\n",
|
||
" color=\"red\",\n",
|
||
" linestyle=\"-\",\n",
|
||
" label=\"exact\",\n",
|
||
")\n",
|
||
"plt.xticks(range(1, krylov_dim + 1), range(1, krylov_dim + 1))\n",
|
||
"plt.legend()\n",
|
||
"plt.xlabel(\"Krylov space dimension\")\n",
|
||
"plt.ylabel(\"Energy\")\n",
|
||
"plt.title(\n",
|
||
" \"Estimating Ground state energy with Krylov Quantum Diagonalization\"\n",
|
||
")\n",
|
||
"plt.show()"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "2ef2fac9-b4b9-4457-a551-c98999a88710",
|
||
"metadata": {},
|
||
"source": [
|
||
"## Appendix: Krylov subspace from real time-evolutions\n",
|
||
"\n",
|
||
"The unitary Krylov space is defined as\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\mathcal{K}_U(H, |\\psi\\rangle) = \\text{span}\\left\\{ |\\psi\\rangle, e^{-iH\\,dt} |\\psi\\rangle, \\dots, e^{-irH\\,dt} |\\psi\\rangle \\right\\}\n",
|
||
"$$\n",
|
||
"\n",
|
||
"for some timestep $dt$ that we will determine later. Temporarily assume $r$ is even: then define $d=r/2$. Notice that when we project the Hamiltonian into the Krylov space above, it is indistinguishable from the Krylov space\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\mathcal{K}_U(H, |\\psi\\rangle) = \\text{span}\\left\\{ e^{i\\,d\\,H\\,dt}|\\psi\\rangle, e^{i(d-1)H\\,dt} |\\psi\\rangle, \\dots, e^{-i(d-1)H\\,dt} |\\psi\\rangle, e^{-i\\,d\\,H\\,dt} |\\psi\\rangle \\right\\},\n",
|
||
"$$\n",
|
||
"\n",
|
||
"i.e., where all the time-evolutions are shifted backward by $d$ timesteps.\n",
|
||
"The reason it is indistinguishable is because the matrix elements\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\tilde{H}_{j,k} = \\langle\\psi|e^{i\\,j\\,H\\,dt}He^{-i\\,k\\,H\\,dt}|\\psi\\rangle=\\langle\\psi|He^{i(j-k)H\\,dt}|\\psi\\rangle\n",
|
||
"$$\n",
|
||
"\n",
|
||
"are invariant under overall shifts of the evolution time, since the time-evolutions commute with the Hamiltonian. For odd $r$, we can use the analysis for $r-1$.\n",
|
||
"\n",
|
||
"We want to show that somewhere in this Krylov space, there is guaranteed to be a low-energy state. We do so by way of the following result, which is derived from Theorem 3.1 in [\\[3\\]](#references):\n",
|
||
"\n",
|
||
"**Claim 1:** there exists a function $f$ such that for energies $E$ in the spectral range of the Hamiltonian (i.e., between the ground state energy and the maximum energy)...\n",
|
||
"1. $f(E_0)=1$\n",
|
||
"2. $|f(E)|\\le2\\left(1 + \\delta\\right)^{-d}$ for all values of $E$ that lie $\\ge\\delta$ away from $E_0$, i.e., it is exponentially suppressed\n",
|
||
"3. $f(E)$ is a linear combination of $e^{ijE\\,dt}$ for $j=-d,-d+1,...,d-1,d$\n",
|
||
"\n",
|
||
"We give a proof below, but that can be safely skipped unless one wants to understand the full, rigorous argument. For now we focus on the implications of the above claim. By property 3 above, we can see that the shifted Krylov space above contains the state $f(H)|\\psi\\rangle$. This is our low-energy state. To see why, write $|\\psi\\rangle$ in the energy eigenbasis:\n",
|
||
"\n",
|
||
"$$\n",
|
||
"|\\psi\\rangle = \\sum_{k=0}^{N}\\gamma_k|E_k\\rangle,\n",
|
||
"$$\n",
|
||
"\n",
|
||
"where $|E_k\\rangle$ is the kth energy eigenstate and $\\gamma_k$ is its amplitude in the initial state $|\\psi\\rangle$. Expressed in terms of this, $f(H)|\\psi\\rangle$ is given by\n",
|
||
"\n",
|
||
"$$\n",
|
||
"f(H)|\\psi\\rangle = \\sum_{k=0}^{N}\\gamma_kf(E_k)|E_k\\rangle,\n",
|
||
"$$\n",
|
||
"\n",
|
||
"using the fact that we can replace $H$ by $E_k$ when it acts on the eigenstate $|E_k\\rangle$. The energy error of this state is therefore\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\text{energy error} = \\frac{\\langle\\psi|f(H)(H-E_0)f(H)|\\psi\\rangle}{\\langle\\psi|f(H)^2|\\psi\\rangle}\n",
|
||
"$$\n",
|
||
"\n",
|
||
"$$\n",
|
||
"= \\frac{\\sum_{k=0}^{N}|\\gamma_k|^2f(E_k)^2(E_k-E_0)}{\\sum_{k=0}^{N}|\\gamma_k|^2f(E_k)^2}.\n",
|
||
"$$\n",
|
||
"\n",
|
||
"To turn this into an upper bound that is easier to understand, we first separate the sum in the numerator into terms with $E_k-E_0\\le\\delta$ and terms with $E_k-E_0>\\delta$:\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\text{energy error} = \\frac{\\sum_{E_k\\le E_0+\\delta}|\\gamma_k|^2f(E_k)^2(E_k-E_0)}{\\sum_{k=0}^{N}|\\gamma_k|^2f(E_k)^2} + \\frac{\\sum_{E_k> E_0+\\delta}|\\gamma_k|^2f(E_k)^2(E_k-E_0)}{\\sum_{k=0}^{N}|\\gamma_k|^2f(E_k)^2}.\n",
|
||
"$$\n",
|
||
"\n",
|
||
"We can upper bound the first term by $\\delta$,\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\frac{\\sum_{E_k\\le E_0+\\delta}|\\gamma_k|^2f(E_k)^2(E_k-E_0)}{\\sum_{k=0}^{N}|\\gamma_k|^2f(E_k)^2} < \\frac{\\delta\\sum_{E_k\\le E_0+\\delta}|\\gamma_k|^2f(E_k)^2}{\\sum_{k=0}^{N}|\\gamma_k|^2f(E_k)^2} \\le \\delta,\n",
|
||
"$$\n",
|
||
"\n",
|
||
"where the first step follows because $E_k-E_0\\le\\delta$ for every $E_k$ in the sum, and the second step follows because the sum in the numerator is a subset of the sum in the denominator. For the second term, first we lower bound the denominator by $|\\gamma_0|^2$, since $f(E_0)^2=1$: adding everything back together, this gives\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\text{energy error} \\le \\delta + \\frac{1}{|\\gamma_0|^2}\\sum_{E_k>E_0+\\delta}|\\gamma_k|^2f(E_k)^2(E_k-E_0).\n",
|
||
"$$\n",
|
||
"\n",
|
||
"To simplify what is left, notice that for all these $E_k$, by the definition of $f$ we know that $f(E_k)^2 \\le 4\\left(1 + \\delta\\right)^{-2d}$. Additionally upper bounding $E_k-E_0<2\\|H\\|$ and upper bounding $\\sum_{E_k>E_0+\\delta}|\\gamma_k|^2<1$ gives\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\text{energy error} \\le \\delta + \\frac{8}{|\\gamma_0|^2}\\|H\\|\\left(1 + \\delta\\right)^{-2d}.\n",
|
||
"$$\n",
|
||
"\n",
|
||
"This holds for any $\\delta>0$, so if we set $\\delta$ equal to our goal error, then the error bound above converges towards that exponentially with the Krylov dimension $2d=r$. Also note that if $\\delta<E_1-E_0$ then the $\\delta$ term actually goes away entirely in the above bound.\n",
|
||
"\n",
|
||
"To complete the argument, we first note that the above is just the energy error of the particular state $f(H)|\\psi\\rangle$, rather than the energy error of the lowest energy state in the Krylov space. However, by the (Rayleigh-Ritz) variational principle, the energy error of the lowest energy state in the Krylov space is upper bounded by the energy error of any state in the Krylov space, so the above is also an upper bound on the energy error of the lowest energy state, i.e., the output of the Krylov quantum diagonalization algorithm.\n",
|
||
"\n",
|
||
"A similar analysis as the above can be carried out that additionally accounts for noise and the thresholding procedure discussed in the notebook. See [\\[2\\]](#references) and [\\[4\\]](#references) for this analysis.\n",
|
||
"\n",
|
||
"\n",
|
||
"## Appendix: proof of Claim 1\n",
|
||
"The following is mostly derived from [\\[3\\]](#references), Theorem 3.1: let $ 0 < a < b $ and let $ \\Pi^*_d $ be the space of residual polynomials (polynomials whose value at 0 is 1) of degree at most $ d $. The solution to\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\beta(a, b, d) = \\min_{p \\in \\Pi^*_d} \\max_{x \\in [a, b]} |p(x)| \\quad\n",
|
||
"$$\n",
|
||
"\n",
|
||
"is\n",
|
||
"\n",
|
||
"$$\n",
|
||
"p^*(x) = \\frac{T_d\\left(\\frac{b + a - 2x}{b - a}\\right)}{T_d\\left(\\frac{b + a}{b - a}\\right)}, \\quad\n",
|
||
"$$\n",
|
||
"\n",
|
||
"and the corresponding minimal value is\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\beta(a, b, d) = T_d^{-1}\\left(\\frac{b + a}{b - a}\\right).\n",
|
||
"$$\n",
|
||
"\n",
|
||
"We want to convert this into a function that can be expressed naturally in terms of complex exponentials, because those are the real time-evolutions that generate the quantum Krylov space.\n",
|
||
"To do so, it is convenient to introduce the following transformation of energies within the spectral range of the Hamiltonian to numbers in the range $[0,1]$: define\n",
|
||
"\n",
|
||
"$$\n",
|
||
"g(E) = \\frac{1-\\cos\\big((E-E_0)dt\\big)}{2},\n",
|
||
"$$\n",
|
||
"\n",
|
||
"where $dt$ is a timestep such that $-\\pi < E_0dt < E_\\text{max}dt < \\pi$.\n",
|
||
"Notice that $g(E_0)=0$ and $g(E)$ grows as $E$ moves away from $E_0$.\n",
|
||
"\n",
|
||
"Now using the polynomial $p^*(x)$ with the parameters a, b, d set to $a = g(E_0 + \\delta)$, $b = 1$, and d = int(r/2), we define the function:\n",
|
||
"\n",
|
||
"$$\n",
|
||
"f(E) = p^* \\left( g(E) \\right) = \\frac{T_d\\left(1 + 2\\frac{\\cos\\big((E-E_0)dt\\big) - \\cos\\big(\\delta\\,dt\\big)}{1 +\\cos\\big(\\delta\\,dt\\big)}\\right)}{T_d\\left(1 + 2\\frac{1-\\cos\\big(\\delta\\,dt\\big)}{1 + \\cos\\big(\\delta\\,dt\\big)}\\right)}\n",
|
||
"$$\n",
|
||
"\n",
|
||
"where $E_0$ is the ground state energy. We can see by inserting $\\cos(x)=\\frac{e^{ix}+e^{-ix}}{2}$ that $f(E)$ is a trigonometric polynomial of degree $d$, i.e., a linear combination of $e^{ijE\\,dt}$ for $j=-d,-d+1,...,d-1,d$. Furthermore, from the definition of $p^*(x)$ above we have that $f(E_0)=p(0)=1$ and for any $E$ in the spectral range such that $\\vert E-E_0 \\vert > \\delta$ we have\n",
|
||
"\n",
|
||
"$$\n",
|
||
"|f(E)| \\le \\beta(a, b, d) = T_d^{-1}\\left(1 + 2\\frac{1-\\cos\\big(\\delta\\,dt\\big)}{1 + \\cos\\big(\\delta\\,dt\\big)}\\right)\n",
|
||
"$$\n",
|
||
"\n",
|
||
"$$\n",
|
||
"\\leq 2\\left(1 + \\delta\\right)^{-d} = 2\\left(1 + \\delta\\right)^{-\\lfloor k/2\\rfloor}.\n",
|
||
"$$"
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "ccd6c21d-ed1e-4894-9cbb-8d3a9d91afe2",
|
||
"metadata": {},
|
||
"source": [
|
||
"## References\n",
|
||
"\n",
|
||
"[1] N. Yoshioka, M. Amico, W. Kirby et al. \"Diagonalization of large many-body Hamiltonians on a quantum processor\". [arXiv:2407.14431](https://arxiv.org/abs/2407.14431)\n",
|
||
"\n",
|
||
"[2] Ethan N. Epperly, Lin Lin, and Yuji Nakatsukasa. \"A theory of quantum subspace diagonalization\". SIAM Journal on Matrix Analysis and Applications 43, 1263–1290 (2022).\n",
|
||
"\n",
|
||
"[3] Å. Björck. \"Numerical methods in matrix computations\". Texts in Applied Mathematics. Springer International Publishing. (2014).\n",
|
||
"\n",
|
||
"[4] William Kirby. \"Analysis of quantum Krylov algorithms with errors\". Quantum 8, 1457 (2024)."
|
||
]
|
||
},
|
||
{
|
||
"cell_type": "markdown",
|
||
"id": "74948cc7-041f-412c-ab16-2554bf164061",
|
||
"metadata": {},
|
||
"source": [
|
||
"## Tutorial survey\n",
|
||
"\n",
|
||
"Please take one minute to provide feedback on this tutorial. Your insights will help us improve our content offerings and user experience.\n",
|
||
"\n",
|
||
"[Link to survey](https://your.feedback.ibm.com/jfe/form/SV_82nennpKIjjD8rQ)"
|
||
]
|
||
}
|
||
],
|
||
"metadata": {
|
||
"description": "Implement the Krylov Quantum Diagonalization Algorithm (KQD) within the context of Qiskit patterns.",
|
||
"kernelspec": {
|
||
"display_name": "Python 3",
|
||
"language": "python",
|
||
"name": "python3"
|
||
},
|
||
"language_info": {
|
||
"codemirror_mode": {
|
||
"name": "ipython",
|
||
"version": 3
|
||
},
|
||
"file_extension": ".py",
|
||
"mimetype": "text/x-python",
|
||
"name": "python",
|
||
"nbconvert_exporter": "python",
|
||
"pygments_lexer": "ipython3",
|
||
"version": "3"
|
||
},
|
||
"platform": "cloud",
|
||
"title": "Krylov quantum diagonalization of lattice Hamiltonians\n"
|
||
},
|
||
"nbformat": 4,
|
||
"nbformat_minor": 4
|
||
}
|