{ "cells": [ { "cell_type": "markdown", "id": "509f7bd9-b597-4d23-b3af-0a76a7b4d33d", "metadata": {}, "source": [ "---\n", "title: \"Diagonalizzazione quantistica di Krylov degli hamiltoniani reticolari\"\n", "description: \"Implementare l'algoritmo di diagonalizzazione quantistica di Krylov (KQD) nel contesto dei modelli Qiskit.\"\n", "---\n", "\n", "{/* cspell:ignore prefactors */}\n", "\n", "\n", "\n", "# Diagonalizzazione quantistica di Krylov degli hamiltoniani reticolari\n", "\n", "*Stima di utilizzo: 20 minuti su un Heron r2 (NOTA: questa è solo una stima. Il tempo di esecuzione potrebbe variare)*\n", "\n" ] }, { "cell_type": "markdown", "id": "921c7b04-5b5d-4cfc-aba0-15a61334e619", "metadata": {}, "source": [ "\n", "\n", "## Sfondo\n", "\n", "Questa esercitazione mostra come implementare l'algoritmo di diagonalizzazione quantistica di Krylov (KQD) nel contesto dei modelli Qiskit. Verrà dapprima illustrata la teoria alla base dell'algoritmo e poi verrà mostrata una dimostrazione della sua esecuzione su una QPU.\n", "\n", "In tutte le discipline, siamo interessati a conoscere le proprietà dello stato fondamentale dei sistemi quantistici. Tra gli esempi vi sono la comprensione della natura fondamentale delle particelle e delle forze, la previsione e la comprensione del comportamento di materiali complessi e la comprensione delle interazioni e delle reazioni biochimiche. A causa della crescita esponenziale dello spazio di Hilbert e delle correlazioni che sorgono nei sistemi entangled, gli algoritmi classici faticano a risolvere questo problema per sistemi quantistici di dimensioni crescenti. A un'estremità dello spettro c'è l'approccio esistente che sfrutta l'hardware quantistico concentrandosi sui metodi quantistici variazionali (ad esempio, l' [eigensolver quantistico variazionale](/docs/tutorials/spin-chain-vqe) ). Queste tecniche si scontrano con i dispositivi attuali a causa dell'elevato numero di chiamate di funzione richieste nel processo di ottimizzazione, che aggiungono un grande overhead di risorse una volta introdotte le tecniche avanzate di mitigazione degli errori, limitando così la loro efficacia ai sistemi di piccole dimensioni. All'altra estremità dello spettro, ci sono i metodi quantistici a tolleranza di errore con garanzie di prestazioni (ad esempio, la [stima quantistica della fase](https://arxiv.org/abs/quant-ph/0604193) ), che richiedono circuiti profondi che possono essere eseguiti solo su un dispositivo a tolleranza di errore. Per questi motivi, introduciamo qui un algoritmo quantistico basato su metodi subspaziali (come descritto in questo [articolo](https://arxiv.org/abs/2312.00178) ), l'algoritmo di diagonalizzazione quantistica di Krylov (KQD). Questo algoritmo si comporta bene su larga scala [\\[1\\]](#references) sull'hardware quantistico esistente, condivide [garanzie di prestazione](https://arxiv.org/abs/2110.07492) simili a quelle della stima di fase, è compatibile con tecniche avanzate di mitigazione degli errori e potrebbe fornire risultati classicamente inaccessibili.\n", "\n" ] }, { "cell_type": "markdown", "id": "5f698c82-95ca-4dc4-b9fa-d6e741e2c02c", "metadata": {}, "source": [ "\n", "\n", "## Requisiti\n", "\n", "Prima di iniziare questa esercitazione, assicuratevi di aver installato quanto segue:\n", "\n", "* Qiskit SDK v2.0 o versioni successive, con supporto [alla visualizzazione](/docs/api/qiskit/visualization)\n", "* Qiskit Runtime v0.22 o successivamente ( `pip install qiskit-ibm-runtime` )\n", "\n" ] }, { "cell_type": "markdown", "id": "c44956e4-47ab-4b0f-9d6d-553080110062", "metadata": {}, "source": [ "\n", "\n", "## Configura\n", "\n" ] }, { "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, IfElseOp\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", " 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\n", " 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", " # list all of the possible sets of n_exc indices of 1s in\n", " # n_exc-particle states\n", " sparse_vecs = [\n", " set(vec) for vec in it.combinations(range(n_qubits + 1), r=n_exc)\n", " ]\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": [ "\n", "\n", "## Fase 1: mappare gli input classici su un problema quantistico\n", "\n" ] }, { "cell_type": "markdown", "id": "a166e2a1-3003-4799-9c15-59f5e78b6ea8", "metadata": {}, "source": [ "\n", "\n", "### Lo spazio di Krylov\n", "\n", "Lo spazio di Krylov $\\mathcal{K}^r$ di ordine $r$ è lo spazio spaziato dai vettori ottenuti moltiplicando le potenze superiori di una matrice $A$, fino a $r-1$, con un vettore di riferimento $\\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", "Se la matrice $A$ è l'hamiltoniana $H$, ci riferiremo allo spazio corrispondente come spazio di Krylov di potenza $\\mathcal{K}_P$. Nel caso in cui $A$ sia l'operatore di evoluzione temporale generato dall'hamiltoniana $U=e^{-iHt}$, ci riferiremo allo spazio come spazio di Krylov unitario $\\mathcal{K}_U$. Il sottospazio di Krylov di potenza che utilizziamo classicamente non può essere generato direttamente su un computer quantistico, poiché $H$ non è un operatore unitario. Si può invece utilizzare l'operatore di evoluzione temporale $U = e^{-iHt}$, che può dare [garanzie di convergenza](https://arxiv.org/abs/2110.07492) simili a quelle del metodo delle potenze. Le potenze di $U$ diventano quindi passi temporali diversi $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", "Si veda l'Appendice per una derivazione dettagliata di come lo spazio di Krylov unitario permetta di rappresentare accuratamente gli autostati a bassa energia.\n", "\n" ] }, { "attachments": {}, "cell_type": "markdown", "id": "5573ca7e-ab16-4488-88b3-d8d1eba9e20c", "metadata": {}, "source": [ "\n", "\n", "### Algoritmo di diagonalizzazione quantistica di Krylov\n", "\n", "Data un'hamiltoniana $H$ che vogliamo diagonalizzare, consideriamo innanzitutto il corrispondente spazio di Krylov unitario $\\mathcal{K}_U$. L'obiettivo è trovare una rappresentazione compatta dell'hamiltoniana in $\\mathcal{K}_U$, che chiameremo $\\tilde{H}$. Gli elementi della matrice di $\\tilde{H}$, la proiezione dell'hamiltoniana nello spazio di Krylov, possono essere calcolati calcolando i seguenti valori di aspettativa\n", "\n", "$$\n", "\\tilde{H}_{mn} = \\langle \\psi_m \\vert H \\vert \\psi_n \\rangle =\n", "$$\n", "\n", "$$\n", "= \\langle \\psi \\vert e^{i H t_m} H e^{-i H t_n} \\vert \\psi \\rangle\n", "$$\n", "\n", "$$\n", "= \\langle \\psi \\vert e^{i H m dt} H e^{-i H n dt} \\vert \\psi \\rangle\n", "$$\n", "\n", "Dove $\\vert \\psi_n \\rangle = e^{-i H t_n} \\vert \\psi \\rangle$ sono i vettori dello spazio di Krylov unitario e $t_n = n dt$ sono i multipli del passo temporale $dt$ scelti. Su un computer quantistico, il calcolo degli elementi di ciascuna matrice può essere effettuato con qualsiasi algoritmo che consenta di ottenere la sovrapposizione tra stati quantistici. Questa esercitazione si concentra sul test di Hadamard. Dato che $\\mathcal{K}_U$ ha dimensione $r$, l'hamiltoniana proiettata nel sottospazio avrà dimensione $r \\times r$. Con $r$ sufficientemente piccola (in genere $r<<100$ è sufficiente per ottenere la convergenza delle stime delle autovalori) possiamo quindi facilmente diagonalizzare l'hamiltoniana proiettata $\\tilde{H}$. Tuttavia, non possiamo diagonalizzare direttamente $\\tilde{H}$ a causa della non ortogonalità dei vettori dello spazio di Krylov. Dovremo misurare le loro sovrapposizioni e costruire una matrice $\\tilde{S}$\n", "\n", "$$\n", "\\tilde{S}_{mn} = \\langle \\psi_m \\vert \\psi_n \\rangle\n", "$$\n", "\n", "Questo ci permette di risolvere il problema degli autovalori in uno spazio non ortogonale (detto anche problema degli autovalori generalizzato)\n", "\n", "$$\n", "\\tilde{H} \\ \\vec{c} = E \\ \\tilde{S} \\ \\vec{c}\n", "$$\n", "\n", "Si possono quindi ottenere stime degli autovalori e degli autostati di $H$ osservando quelli di $\\tilde{H}$. Ad esempio, la stima dell'energia dello stato fondamentale si ottiene prendendo l'autovalore più piccolo $c$ e lo stato fondamentale dal corrispondente autovettore $\\vec{c}$. I coefficienti in $\\vec{c}$ determinano il contributo dei diversi vettori che coprono $\\mathcal{K}_U$.\n", "\n", "![fig1.png](https://quantum.cloud.ibm.com/docs/images/tutorials/krylov-subspace-diagonalization/fc662b76-8ad7-4a6c-8c49-5f08c125aee8.avif)\n", "\n", "La figura mostra una rappresentazione circuitale del test di Hadamard modificato, un metodo utilizzato per calcolare la sovrapposizione tra diversi stati quantistici. Per ogni elemento della matrice $\\tilde{H}_{i,j}$, si esegue un test di Hadamard tra gli stati $\\vert \\psi_i \\rangle$, $\\vert \\psi_j \\rangle$. Ciò è evidenziato nella figura dallo schema di colori degli elementi della matrice e dalle corrispondenti operazioni $\\text{Prep} \\; \\psi_i$, $\\text{Prep} \\; \\psi_j$. Pertanto, per calcolare tutti gli elementi della matrice dell'hamiltoniana proiettata $\\tilde{H}$ è necessario un insieme di test di Hadamard per tutte le possibili combinazioni di vettori dello spazio di Krylov. Il filo superiore del circuito del test di Hadamard è un qubit ancilla che viene misurato nella base X o Y; il suo valore di aspettativa determina il valore della sovrapposizione tra gli stati. Il filo inferiore rappresenta tutti i qubit dell'hamiltoniana del sistema. L'operazione $\\text{Prep} \\; \\psi_i$ prepara il qubit del sistema nello stato $\\vert \\psi_i \\rangle$ controllato dallo stato del qubit ancilla (analogamente per $\\text{Prep} \\; \\psi_j$ ) e l'operazione $P$ rappresenta la decomposizione Pauli dell'hamiltoniana del sistema $H = \\sum_i P_i$. Di seguito viene fornita una derivazione più dettagliata delle operazioni calcolate dal test di Hadamard.\n", "\n" ] }, { "cell_type": "markdown", "id": "1a6d7a4f-6c1f-4069-93d1-f5b670645d7f", "metadata": {}, "source": [ "\n", "\n", "#### Definire l'Hamiltoniano\n", "\n", "Consideriamo l'hamiltoniana di Heisenberg per $N$ qubit su una catena lineare: $H= \\sum_{i,j}^N X_i X_j + Y_i Y_j - J Z_i Z_j$\n", "\n" ] }, { "cell_type": "code", "execution_count": 3, "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": [ "\n", "\n", "#### Imposta i parametri per l'algoritmo\n", "\n" ] }, { "cell_type": "markdown", "id": "8d376c0e-fb7a-4a41-838a-ae041d5c9afa", "metadata": {}, "source": [ "Scegliamo euristicamente un valore per il passo temporale `dt` (basato sui limiti superiori della norma hamiltoniana). Ref [\\[2\\]](#references) ha dimostrato che un passo temporale sufficientemente piccolo è $\\pi/\\vert \\vert H \\vert \\vert$, e che è preferibile fino a un certo punto sottostimare questo valore piuttosto che sovrastimarlo, poiché una sovrastima può permettere ai contributi degli stati ad alta energia di corrompere anche lo stato ottimale nello spazio di Krylov. D'altra parte, la scelta di $dt$ troppo piccola porta a un condizionamento peggiore del sottospazio di Krylov, poiché i vettori base di Krylov differiscono meno da un passo all'altro.\n", "\n" ] }, { "cell_type": "code", "execution_count": 4, "id": "963dd2e9-f4ea-456d-b4ed-e3d05a114082", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "np.float64(0.10833078115826875)" ] }, "execution_count": 4, "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": [ "E impostare altri parametri dell'algoritmo. Ai fini di questa esercitazione, ci limiteremo a utilizzare uno spazio di Krylov con solo cinque dimensioni, che è abbastanza limitante.\n", "\n" ] }, { "cell_type": "code", "execution_count": 5, "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": [ "\n", "\n", "#### Preparazione dello Stato\n", "\n", "Scegliere uno stato di riferimento $\\vert \\psi \\rangle$ che abbia una certa sovrapposizione con lo stato fondamentale. Per questa hamiltoniana, utilizziamo come stato di riferimento lo stato a con un'eccitazione nel qubit centrale $\\vert 00..010...00 \\rangle$.\n", "\n" ] }, { "cell_type": "code", "execution_count": 6, "id": "70161afe-8ace-4642-894a-cd21ed77a3b9", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "\"Output" ] }, "execution_count": 6, "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": [ "\n", "\n", "#### Evoluzione temporale\n", "\n", "Possiamo realizzare l'operatore di evoluzione temporale generato da una data hamiltoniana: $U=e^{-iHt}$ tramite l' [approssimazione di Lie-Trotter](/docs/api/qiskit/qiskit.synthesis.LieTrotter).\n", "\n" ] }, { "cell_type": "code", "execution_count": 7, "id": "a23b9e9c-5dc8-447f-8c73-b4fa01630c8f", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "" ] }, "execution_count": 7, "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": [ "\n", "\n", "#### test di Hadamard\n", "\n", "![fig2.png](https://quantum.cloud.ibm.com/docs/images/tutorials/krylov-subspace-diagonalization/c5263851-6067-4ca2-8e0c-a835631cdc7f.avif)\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", "Dove $P$ è uno dei termini della decomposizione dell'hamiltoniana $H=\\sum P$ e $\\text{Prep} \\; \\psi_i$, $\\text{Prep} \\; \\psi_j$ sono operazioni controllate che preparano $|\\psi_i\\rangle$, $|\\psi_j\\rangle$ vettori dello spazio di Krylov unitario, con $|\\psi_k\\rangle = e^{-i H k dt } \\vert \\psi \\rangle = e^{-i H k dt } U_{\\psi} \\vert 0 \\rangle^N$. Per misurare $X$, applicare prima $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", "... poi misurare:\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", "Dall'identità $|a + b\\|^2 = \\langle a + b | a + b \\rangle = \\|a\\|^2 + \\|b\\|^2 + 2\\text{Re}\\langle a | b \\rangle$. Analogamente, misurando $Y$ si ottiene\n", "\n", "$$\n", "\\begin{equation*}\n", " \\langle Y\\rangle = \\text{Im}\\Big[\\langle\\psi_j| P|\\psi_i\\rangle\\Big].\n", "\\end{equation*}\n", "$$\n", "\n" ] }, { "cell_type": "code", "execution_count": 8, "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": [ "\"Output" ] }, "execution_count": 8, "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": [ "Il circuito di test Hadamard può essere un circuito profondo una volta decomposto in porte native (il che aumenterà ancora di più se teniamo conto della topologia del dispositivo)\n", "\n" ] }, { "cell_type": "code", "execution_count": 9, "id": "06f22c05-6e1c-4540-9dac-884b3400a50e", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Number of layers of 2Q operations 112753\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": [ "\n", "\n", "## Fase 2: Ottimizzazione del problema per l'esecuzione su hardware quantistico\n", "\n", "\n", "\n", "### Test di Hadamard efficiente\n", "\n", "Possiamo ottimizzare i circuiti profondi per il test di Hadamard che abbiamo ottenuto introducendo alcune approssimazioni e facendo affidamento su alcune ipotesi sull'hamiltoniana del modello. Ad esempio, si consideri il seguente circuito per il test di Hadamard:\n", "\n", "![fig3.png](https://quantum.cloud.ibm.com/docs/images/tutorials/krylov-subspace-diagonalization/35b13797-5a46-486c-b50e-97c205cc9747.avif)\n", "\n", "Supponiamo di poter calcolare classicamente $E_0$, l'autovalore di $|0\\rangle^N$ sotto l'hamiltoniana $H$. Questo è soddisfatto quando l'hamiltoniana preserva la simmetria U(1). Anche se questa può sembrare un'ipotesi forte, in molti casi è lecito supporre che esista uno stato di vuoto (che in questo caso corrisponde allo stato $|0\\rangle^N$ ) che non è influenzato dall'azione dell'hamiltoniana. Questo è vero, ad esempio, per gli hamiltoniani della chimica che descrivono molecole stabili (dove il numero di elettroni è conservato).\n", "Dato che il gate $\\text{Prep} \\; \\psi$, prepara lo stato di riferimento desiderato $\\ket{psi} = \\text{Prep} \\; \\psi \\ket{0} = e^{-i H 0 dt} U_{\\psi} \\ket{0}$, ad esempio, preparare lo stato HF per la chimica $\\text{Prep} \\; \\psi$ sarebbe un prodotto di NOT a singolo-qubit, quindi controllato- $\\text{Prep} \\; \\psi$ è solo un prodotto di CNOT.\n", "Il circuito di cui sopra implementa il seguente stato prima della misurazione:\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", "dove abbiamo utilizzato il classico spostamento di fase simulabile $ U\\ket{0}^N = e^{i\\phi}\\ket{0}^N$ nella terza riga. Pertanto i valori di aspettativa si ottengono come\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", "Utilizzando queste ipotesi, siamo stati in grado di scrivere i valori di aspettativa degli operatori di interesse con un minor numero di operazioni controllate. In effetti, dobbiamo implementare solo la preparazione controllata degli stati $\\text{Prep} \\; \\psi$ e non le evoluzioni temporali controllate. Riformulare il calcolo come sopra ci permetterà di ridurre notevolmente la profondità dei circuiti risultanti.\n", "\n" ] }, { "cell_type": "markdown", "id": "2e29da5f-b5f5-4b9e-96fc-3fa7ab698398", "metadata": {}, "source": [ "\n", "\n", "### Scomporre l'operatore di evoluzione temporale con la decomposizione di Trotter\n", "\n", "Invece di implementare esattamente l'operatore di evoluzione temporale, possiamo utilizzare la decomposizione di Trotter per implementarne un'approssimazione. Ripetendo più volte una decomposizione di Trotter di un certo ordine si ottiene un'ulteriore riduzione dell'errore introdotto dall'approssimazione. Di seguito, costruiamo direttamente l'implementazione di Trotter nel modo più efficiente per il grafo delle interazioni dell'hamiltoniana che stiamo considerando (solo interazioni di prossimità). In pratica inseriamo le rotazioni Pauli $R_{xx}$, $R_{yy}$, $R_{zz}$ con un angolo parametrizzato $t$ che corrispondono all'implementazione approssimata di $e^{-i (XX + YY + ZZ) t}$. Data la differenza di definizione delle rotazioni Pauli e dell'evoluzione temporale che stiamo cercando di implementare, dovremo utilizzare il parametro $2*dt$ per ottenere un'evoluzione temporale di $dt$. Inoltre, invertiamo l'ordine delle operazioni per un numero dispari di ripetizioni dei passi di Trotter, il che è funzionalmente equivalente ma permette di sintetizzare operazioni adiacenti in un unico $SU(2)$ unitario. In questo modo si ottiene un circuito molto meno profondo di quello che si ottiene utilizzando la funzionalità generica `PauliEvolutionGate()` .\n", "\n" ] }, { "cell_type": "code", "execution_count": 10, "id": "267716dc-fa23-41bd-abe4-6d4e0499a0f4", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "\"Output" ] }, "execution_count": 10, "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": [ "\n", "\n", "### Utilizzare un circuito ottimizzato per la preparazione dello stato\n", "\n" ] }, { "cell_type": "code", "execution_count": 11, "id": "70411715-eed3-4cf5-961d-06a6f1e04efc", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "\"Output" ] }, "execution_count": 11, "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": [ "\n", "\n", "### Circuiti modello per il calcolo degli elementi di matrice di $\\tilde{S}$ e $\\tilde{H}$ tramite il test di Hadamard\n", "\n", "L'unica differenza tra i circuiti utilizzati nel test di Hadamard sarà la fase dell'operatore di evoluzione temporale e le osservabili misurate. Possiamo quindi preparare un circuito modello che rappresenti il circuito generico per il test di Hadamard, con segnaposti per le porte che dipendono dall'operatore di evoluzione temporale.\n", "\n" ] }, { "cell_type": "code", "execution_count": 12, "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": 13, "id": "33ec7c29-904e-4445-a654-405214349a4d", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "\"Output" ] }, "execution_count": 13, "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": 14, "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": [ "Abbiamo ridotto considerevolmente la profondità del test di Hadamard con una combinazione di approssimazione di Trotter e unitarietà non controllata\n", "\n" ] }, { "cell_type": "markdown", "id": "ca78bfe3-684f-4d3a-b8c2-3d4e55c9ec30", "metadata": {}, "source": [ "\n", "\n", "## Passaggio 3: eseguire utilizzando Qiskit primitives\n", "\n" ] }, { "cell_type": "markdown", "id": "a082f023-d752-40e4-b693-c4d0da5ba102", "metadata": {}, "source": [ "Istanziare il backend e impostare i parametri di runtime\n", "\n" ] }, { "cell_type": "code", "execution_count": null, "id": "0d90e4df-e262-4852-a811-ff6a1d3232ae", "metadata": {}, "outputs": [], "source": [ "service = QiskitRuntimeService()\n", "backend = service.least_busy(operational=True, simulator=False)\n", "if (\n", " \"if_else\" not in backend.target.operation_names\n", "): # Needed as \"op_name\" could be \"if_else\"\n", " backend.target.add_instruction(IfElseOp, name=\"if_else\")\n", "print(backend.name)" ] }, { "cell_type": "markdown", "id": "2b3f7227-e6ef-4a11-93d6-d3696054b9cc", "metadata": {}, "source": [ "\n", "\n", "### Transpiling su una QPU\n", "\n", "Per prima cosa, selezioniamo sottoinsiemi della mappa di accoppiamento con qubit dalle prestazioni \"buone\" (dove \"buone\" è piuttosto arbitrario in questo caso, vogliamo principalmente evitare qubit dalle prestazioni davvero scadenti) e creiamo un nuovo obiettivo per la transpilazione\n", "\n" ] }, { "cell_type": "code", "execution_count": 32, "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.02 or t2 < 100:\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", " op_name = list(target.operation_names_for_qargs(q))[0]\n", " twoq_gate_err = target[f\"{op_name}\"][q].error\n", " if twoq_gate_err > 0.005:\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": [ "Quindi trasporre il circuito virtuale nel miglior layout fisico in questo nuovo target\n", "\n" ] }, { "cell_type": "code", "execution_count": 36, "id": "62109c95-79fe-4076-b591-9bd920fd51f4", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "depth 52\n", "num 2q ops OrderedDict([('rz', 2058), ('sx', 1703), ('cz', 728), ('x', 84), ('barrier', 8)])\n", "physical qubits [91, 92, 93, 94, 95, 98, 99, 108, 109, 110, 111, 113, 114, 115, 119, 127, 132, 133, 134, 135, 137, 139, 147, 148, 149, 150, 151, 152, 153, 154, 155]\n" ] } ], "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": [ "\n", "\n", "### Creare PUB per l'esecuzione con Estimator\n", "\n" ] }, { "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": [ "\n", "\n", "### Esegui circuiti\n", "\n", "I circuiti per $t=0$ sono calcolabili in modo classico\n", "\n" ] }, { "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": [ "Esegui i circuiti per $S$ e $\\tilde{H}$ con Estimator\n", "\n" ] }, { "cell_type": "code", "execution_count": null, "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[\"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", " \"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", "}\n", "\n", "estimator = Estimator(mode=backend, options=experimental_opts)\n", "\n", "\n", "job = estimator.run([pub])" ] }, { "cell_type": "markdown", "id": "9524e609-b457-42e9-9de7-8dbd4464ac74", "metadata": {}, "source": [ "\n", "\n", "## Fase 4: Post-elaborazione e restituzione del risultato nel formato classico desiderato\n", "\n" ] }, { "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": [ "\n", "\n", "### Calcolare le matrici Hamiltoniane efficaci e di sovrapposizione\n", "\n", "Calcolare innanzitutto la fase accumulata dallo stato $\\vert 0 \\rangle$ durante l'evoluzione temporale non controllata\n", "\n" ] }, { "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": [ "Una volta ottenuti i risultati delle esecuzioni dei circuiti, possiamo post-elaborare i dati per calcolare gli elementi della matrice di $S$\n", "\n" ] }, { "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": [ "$$\n", "\\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]\n", "$$" ], "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": [ "E gli elementi della matrice di $\\tilde{H}$\n", "\n" ] }, { "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": [ "$$\n", "\\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]\n", "$$" ], "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": [ "Infine, possiamo risolvere il problema generalizzato degli autovalori per $\\tilde{H}$ :\n", "\n", "$\\tilde{H} \\vec{c} = c S \\vec{c}$\n", "\n", "e ottenere una stima dell'energia dello stato fondamentale $c_{min}$\n", "\n" ] }, { "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": [ "Per un settore a una sola particella, possiamo calcolare in modo efficiente lo stato fondamentale di questo settore dell'hamiltoniana in modo classico\n", "\n" ] }, { "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": [ "\"Output" ] }, "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": [ "\n", "\n", "## Appendice: sottospazio di Krylov dalle evoluzioni in tempo reale\n", "\n", "Lo spazio di Krylov unitario è definito come\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", "per qualche passo temporale $dt$ che determineremo in seguito. Assumiamo temporaneamente che $r$ sia pari: allora definiamo $d=r/2$. Si noti che quando proiettiamo l'Hamiltoniana nello spazio di Krylov di cui sopra, essa è indistinguibile dallo spazio di Krylov\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", "cioè dove tutte le evoluzioni temporali sono spostate all'indietro di $d$ passi temporali.\n", "Il motivo per cui è indistinguibile è che gli elementi della matrice\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", "sono invarianti sotto spostamenti complessivi del tempo di evoluzione, poiché le evoluzioni temporali commutano con l'hamiltoniana. Per $r$ dispari, possiamo utilizzare l'analisi per $r-1$.\n", "\n", "Vogliamo dimostrare che da qualche parte in questo spazio di Krylov è garantita l'esistenza di uno stato a bassa energia. Lo facciamo attraverso il seguente risultato, derivato dal teorema 3.1 di [\\[3\\]](#references) :\n", "\n", "**Asserzione 1:** esiste una funzione $f$ tale che per le energie $E$ nell'intervallo spettrale dell'hamiltoniana (cioè tra l'energia dello stato fondamentale e l'energia massima)...\n", "\n", "1. $f(E_0)=1$\n", "2. $|f(E)|\\le2\\left(1 + \\delta\\right)^{-d}$ per tutti i valori di $E$ che si trovano a $\\ge\\delta$ distanza da $E_0$, cioè è soppresso esponenzialmente\n", "3. $f(E)$ è una combinazione lineare di $e^{ijE\\,dt}$ per $j=-d,-d+1,...,d-1,d$\n", "\n", "Di seguito viene fornita una dimostrazione, che può essere tranquillamente saltata a meno che non si voglia comprendere l'argomentazione completa e rigorosa. Per ora ci concentriamo sulle implicazioni di questa affermazione. In base alla proprietà 3, possiamo vedere che lo spazio di Krylov spostato sopra contiene lo stato $f(H)|\\psi\\rangle$. Questo è il nostro stato a bassa energia. Per capire perché, scriviamo $|\\psi\\rangle$ nell'autobase dell'energia:\n", "\n", "$$\n", "|\\psi\\rangle = \\sum_{k=0}^{N}\\gamma_k|E_k\\rangle,\n", "$$\n", "\n", "dove $|E_k\\rangle$ è l'autostato di energia k e $\\gamma_k$ è la sua ampiezza nello stato iniziale $|\\psi\\rangle$. Espresso in termini di questo, $f(H)|\\psi\\rangle$ è dato da\n", "\n", "$$\n", "f(H)|\\psi\\rangle = \\sum_{k=0}^{N}\\gamma_kf(E_k)|E_k\\rangle,\n", "$$\n", "\n", "utilizzando il fatto che possiamo sostituire $H$ con $E_k$ quando agisce sull'autostato $|E_k\\rangle$. L'errore energetico di questo stato è quindi\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", "Per trasformare questo risultato in un limite superiore di più facile comprensione, separiamo prima la somma nel numeratore in termini con $E_k-E_0\\le\\delta$ e termini con $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", "Il primo termine può essere limitato in alto da $\\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", "dove il primo passo segue perché $E_k-E_0\\le\\delta$ per ogni $E_k$ nella somma, e il secondo passo segue perché la somma nel numeratore è un sottoinsieme della somma nel denominatore. Per il secondo termine, prima di tutto abbassiamo il denominatore di $|\\gamma_0|^2$, dato che $f(E_0)^2=1$ : sommando il tutto, si ottiene\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", "Per semplificare ciò che rimane, si noti che per tutti questi $E_k$, dalla definizione di $f$ sappiamo che $f(E_k)^2 \\le 4\\left(1 + \\delta\\right)^{-2d}$. Inoltre, l'upper bounding di $E_k-E_0<2\\|H\\|$ e l'upper bounding di $\\sum_{E_k>E_0+\\delta}|\\gamma_k|^2<1$ dà\n", "\n", "$$\n", "\\text{energy error} \\le \\delta + \\frac{8}{|\\gamma_0|^2}\\|H\\|\\left(1 + \\delta\\right)^{-2d}.\n", "$$\n", "\n", "Questo vale per qualsiasi $\\delta>0$, quindi se impostiamo $\\delta$ uguale all'errore che ci prefiggiamo, il limite dell'errore di cui sopra converge verso tale valore in modo esponenziale con la dimensione di Krylov $2d=r$. Si noti inoltre che se $\\delta\n", "\n", "## Appendice: prova della rivendicazione 1\n", "\n", "Quanto segue deriva in gran parte da [\\[3\\]](#references), Teorema 3.1: Sia $0 < a < b$ e sia $\\Pi^*_d$ lo spazio dei polinomi residui (polinomi il cui valore a 0 è 1) di grado massimo $d$. La soluzione di\n", "\n", "$$\n", "\\beta(a, b, d) = \\min_{p \\in \\Pi^*_d} \\max_{x \\in [a, b]} |p(x)| \\quad\n", "$$\n", "\n", "è\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", "e il corrispondente valore minimo è\n", "\n", "$$\n", "\\beta(a, b, d) = T_d^{-1}\\left(\\frac{b + a}{b - a}\\right).\n", "$$\n", "\n", "Vogliamo convertire questa funzione in una funzione che possa essere espressa naturalmente in termini di esponenziali complessi, perché queste sono le reali evoluzioni temporali che generano lo spazio di Krylov quantistico.\n", "Per fare ciò, è conveniente introdurre la seguente trasformazione delle energie all'interno dell'intervallo spettrale dell'hamiltoniano in numeri nell'intervallo $[0,1]$ : definire\n", "\n", "$$\n", "g(E) = \\frac{1-\\cos\\big((E-E_0)dt\\big)}{2},\n", "$$\n", "\n", "dove $dt$ è un passo temporale tale che $-\\pi < E_0dt < E_\\text{max}dt < \\pi$. Si noti che $g(E_0)=0$ e $g(E)$ crescono quando $E$ si allontana da $E_0$.\n", "\n", "Ora, utilizzando il polinomio $p^*(x)$ con i parametri a, b, d impostati a $a = g(E_0 + \\delta)$, $b = 1$, e d = int( r/2 ), definiamo la funzione:\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", "dove $E_0$ è l'energia dello stato fondamentale. Inserendo $\\cos(x)=\\frac{e^{ix}+e^{-ix}}{2}$ si vede che $f(E)$ è un polinomio trigonometrico di grado $d$, cioè una combinazione lineare di $e^{ijE\\,dt}$ per $j=-d,-d+1,...,d-1,d$. Inoltre, dalla definizione di $p^*(x)$ abbiamo che $f(E_0)=p(0)=1$ e per qualsiasi $E$ nell'intervallo spettrale tale che $\\vert E-E_0 \\vert > \\delta$ abbiamo\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", "$$\n", "\n" ] }, { "cell_type": "markdown", "id": "ccd6c21d-ed1e-4894-9cbb-8d3a9d91afe2", "metadata": {}, "source": [ "\n", "\n", "## Riferimenti\n", "\n", "\\[1] N. Yoshioka, M. Amico, W. Kirby et al. \"Diagonalizzazione di grandi hamiltoniani a molti corpi su un processore quantistico\". [arXiv:2407.14431](https://arxiv.org/abs/2407.14431)\n", "\n", "\\[2] Ethan N. Epperly, Lin Lin e Yuji Nakatsukasa. \"Una teoria della diagonalizzazione dei sottospazi quantistici\". SIAM Journal on Matrix Analysis and Applications 43, 1263-1290 (2022).\n", "\n", "\\[3] Å. Björck. \"Metodi numerici nel calcolo delle matrici\". Testi di matematica applicata. Springer International Publishing. (2014).\n", "\n", "\\[4] William Kirby. \"Analisi di algoritmi di Krylov quantistici con errori\". Quantum 8, 1457 (2024).\n", "\n" ] }, { "cell_type": "markdown", "id": "74948cc7-041f-412c-ab16-2554bf164061", "metadata": {}, "source": [ "\n", "\n", "## Sondaggio tutorial\n", "\n", "Vi invitiamo a partecipare a questo breve sondaggio per fornire un feedback su questa esercitazione. Le vostre indicazioni ci aiuteranno a migliorare la nostra offerta di contenuti e l'esperienza degli utenti.\n", "\n", "[Collegamento al sondaggio](https://your.feedback.ibm.com/jfe/form/SV_82nennpKIjjD8rQ)\n", "\n" ] }, { "cell_type": "markdown", "metadata": {}, "id": "a1b8767d", "source": "© IBM Corp., 2017-2026" } ], "metadata": { "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" }, "hours": 1.5, "qpuSeconds": 1200 }, "nbformat": 4, "nbformat_minor": 4 }