{ "cells": [ { "cell_type": "markdown", "id": "ee388de5-2507-48ff-97a6-6777698d6256", "metadata": {}, "source": [ "---\n", "title: \"Probenbasierte Krylov-Quantendiagonalisierung eines fermionischen Gittermodells\"\n", "description: \"Verwenden Sie den samplebasierten Quanten-Diagonalisierungsalgorithmus, um das Anderson-Modell mit einer einzelnen Verunreinigung unter Verwendung von rauschbehafteter Quantenhardware zu simulieren.\"\n", "---\n", "\n", "{/* cspell:ignore fontdict fontsize nocc SQKD DMRG textrm varepsilon vecs pqrs ijkl */}\n", "\n", "\n", "\n", "# Probenbasierte Krylov-Quantendiagonalisierung eines fermionischen Gittermodells\n", "\n", "*Geschätzte Nutzungsdauer: Neun Sekunden auf einem Heron r2 Prozessor (HINWEIS: Dies ist nur eine Schätzung. Ihre Laufzeit kann variieren.)*\n", "\n" ] }, { "cell_type": "markdown", "id": "a1b2c3d4-e5f6-7890-abcd-ef1234567890", "metadata": {}, "source": [ "\n", "\n", "## Lernziele\n", "\n", "* Wie man das [SQD-Qiskit-Add-on](https://github.com/Qiskit/qiskit-addon-sqd) verwendet, um die Grundzustandsenergie eines Gittermodells anhand von Bitstrings zu approximieren, die von einer Quantenverarbeitungseinheit (QPU) abgetastet wurden.\n", "* Wie man mit [ffsim](https://github.com/qiskit-community/ffsim) Zeitentwicklungsschaltungen für die Fermionensimulation erstellt.\n", "* Wie man Signale aus mehreren Schaltkreisen für die Weiterverarbeitung mit dem algorithmischen Verfahren der samplebasierten Krylov-Diagonalisierung (SKQD) zusammenführt.\n", "\n", "\n", "\n", "## Voraussetzungen\n", "\n", "* [Quantendiagonalisierung eines Chemie-Hamiltonians auf der Grundlage von Stichproben](/docs/tutorials/sample-based-quantum-diagonalization)\n", "* [Krylov-Quanten-Diagonalisierung von Gitter-Hamilton-Operatoren](/docs/tutorials/krylov-quantum-diagonalization)\n", "* [Qiskit primitives](/docs/guides/primitives)\n", "\n" ] }, { "cell_type": "markdown", "id": "dc5cc74e-06bf-45ac-a69b-81778138e08f", "metadata": {}, "source": [ "\n", "\n", "## Hintergrund\n", "\n", "Dieses Tutorial zeigt, wie man die stichprobenbasierte Quantendiagonalisierung (SQD) verwendet, um die Grundzustandsenergie eines fermionischen Gittermodells abzuschätzen. Konkret untersuchen wir das eindimensionale Single-Impurity-Anderson-Modell (SIAM), das zur Beschreibung von in Metallen eingebetteten magnetischen Verunreinigungen verwendet wird.\n", "\n", "Dieses Tutorial folgt einem ähnlichen Arbeitsablauf wie das verwandte Tutorial [Probenbasierte Quantendiagonalisierung eines Chemie-Hamiltonians](/docs/tutorials/sample-based-quantum-diagonalization). Ein wesentlicher Unterschied besteht jedoch darin, wie die Quantenschaltungen aufgebaut sind. Das andere Tutorial verwendet einen heuristischen Variationsansatz, der für chemische Hamiltonianer mit potenziell Millionen von Wechselwirkungstermen interessant ist. Andererseits werden in diesem Tutorium Schaltungen verwendet, die die Zeitentwicklung durch den Hamiltonian annähern. Solche Schaltkreise können tief sein, was diesen Ansatz für Anwendungen auf Gittermodelle besser geeignet macht. Die von diesen Schaltkreisen vorbereiteten Zustandsvektoren bilden die Grundlage für einen [Krylov-Unterraum](https://en.wikipedia.org/wiki/Krylov_subspace), so dass der Algorithmus unter geeigneten Annahmen nachweislich und effizient gegen den Grundzustand konvergiert.\n", "\n", "Der in diesem Tutorium verwendete Ansatz kann als eine Kombination aus den Techniken der SQD und der [Krylov-Quanten-Diagonalisierung (KQD](https://arxiv.org/abs/2407.14431) ) betrachtet werden. Der kombinierte Ansatz wird manchmal als stichprobenbasierte Krylov-Quantendiagonalisierung (SQKD) bezeichnet. Siehe [Krylov-Quanten-Diagonalisierung von Gitter-Hamiltonianern](/docs/tutorials/krylov-quantum-diagonalization) für eine Einführung in die KQD-Methode.\n", "\n", "Dieses Tutorial basiert auf der Arbeit [\"Quantum-Centric Algorithm for Sample-Based Krylov Diagonalization\"](https://arxiv.org/abs/2501.09702), die für weitere Details herangezogen werden kann.\n", "\n", "\n", "\n", "### Anderson-Modell mit einer einzigen Verunreinigung (SIAM)\n", "\n", "Der eindimensionale SIAM-Hamiltonian ist eine Summe aus drei Termen:\n", "\n", "$$\n", "H = H_{\\textrm{imp}}+ H_\\textrm{bath} + H_\\textrm{hyb},\n", "$$\n", "\n", "wo\n", "\n", "$$\n", "\\begin{align*}\n", " H_\\textrm{imp} &= \\varepsilon \\left( \\hat{n}_{d\\uparrow} + \\hat{n}_{d\\downarrow} \\right) + U \\hat{n}_{d\\uparrow}\\hat{n}_{d\\downarrow}, \\\\\n", " H_\\textrm{bath} &= -t \\sum_{\\substack{\\mathbf{j} = 0\\\\ \\sigma\\in \\{\\uparrow, \\downarrow\\}}}^{L-1} \\left(\\hat{c}^\\dagger_{\\mathbf{j}, \\sigma}\\hat{c}_{\\mathbf{j}+1, \\sigma} + \\hat{c}^\\dagger_{\\mathbf{j}+1, \\sigma}\\hat{c}_{\\mathbf{j}, \\sigma} \\right), \\\\\n", " H_\\textrm{hyb} &= V\\sum_{\\sigma \\in \\{\\uparrow, \\downarrow \\}} \\left(\\hat{d}^\\dagger_\\sigma \\hat{c}_{0, \\sigma} + \\hat{c}^\\dagger_{0, \\sigma} \\hat{d}_{\\sigma} \\right).\n", "\\end{align*}\n", "$$\n", "\n", "Hier sind $c^\\dagger_{\\mathbf{j},\\sigma}/c_{\\mathbf{j},\\sigma}$ die fermionischen Erzeugungs-/Annihilationsoperatoren für die $\\mathbf{j}^{\\textrm{th}}$ Badstelle mit dem Spin $\\sigma$, $\\hat{d}^\\dagger_{\\sigma}/\\hat{d}_{\\sigma}$ sind Erzeugungs-/Annihilationsoperatoren für den Störstellenmodus und $\\hat{n}_{d\\sigma} = \\hat{d}^\\dagger_{\\sigma} \\hat{d}_{\\sigma}$. $t$, $U$ und $V$ sind reelle Zahlen, die die Sprung-, Vor-Ort- und Hybridisierungswechselwirkungen beschreiben, und $\\varepsilon$ ist eine reelle Zahl, die das chemische Potenzial angibt.\n", "\n", "Es ist zu beachten, dass der Hamiltonian eine spezifische Instanz des allgemeinen Wechselwirkungs-Elektronen-Hamiltonians ist,\n", "\n", "$$\n", "\\begin{align*}\n", " H &= \\sum_{\\substack{p, q \\\\ \\sigma}} h_{pq} \\hat{a}^\\dagger_{p\\sigma} \\hat{a}_{q\\sigma} + \\sum_{\\substack{p, q, r, s \\\\ \\sigma \\tau}} \\frac{h_{pqrs}}{2} \\hat{a}^\\dagger_{p\\sigma} \\hat{a}^\\dagger_{q\\tau} \\hat{a}_{s\\tau} \\hat{a}_{r\\sigma} \\\\\n", " &= H_1 + H_2,\n", "\\end{align*}\n", "$$\n", "\n", "wobei $H_1$ aus Einkörpertermen besteht, die quadratisch in den fermionischen Erzeugungs- und Annihilationsoperatoren sind, und $H_2$ aus Zweikörpertermen besteht, die quartisch sind. Für die SIAM,\n", "\n", "$$\n", "H_2 = U \\hat{n}_{d\\uparrow}\\hat{n}_{d\\downarrow}\n", "$$\n", "\n", "und $H_1$ enthält den Rest der Terme im Hamiltonian. Um den Hamiltonian programmatisch darstellen zu können, speichern wir die Matrix $h_{pq}$ und den Tensor $h_{pqrs}$.\n", "\n", "\n", "\n", "### Positions- und Impulsbasen\n", "\n", "Aufgrund der annähernden Translationssymmetrie in $H_\\textrm{bath}$ erwarten wir nicht, dass der Grundzustand in der Positionsbasis (der Orbitalbasis, in der der Hamiltonian oben angegeben ist) spärlich ist. Die Leistung von SQD ist nur dann gewährleistet, wenn der Grundzustand spärlich ist, d. h. nur eine kleine Anzahl von Basiszuständen eine signifikante Bedeutung hat. Um die Sparsamkeit des Grundzustands zu verbessern, führen wir die Simulation in der Orbitalbasis durch, in der $H_\\textrm{bath}$ diagonal ist. Wir nennen diese Basis die *Impulsbasis*. Da $H_\\textrm{bath}$ ein quadratischer fermionischer Hamiltonian ist, kann er durch eine Orbitalrotation effizient diagonalisiert werden.\n", "\n", "\n", "\n", "### Ungefähre zeitliche Entwicklung durch den Hamilton-Operator\n", "\n", "Zur Annäherung der Zeitentwicklung durch den Hamiltonian verwenden wir eine Trotter-Suzuki-Zerlegung zweiter Ordnung,\n", "\n", "$$\n", " e^{-i \\Delta t H} \\approx e^{-i\\frac{\\Delta t}{2} H_2} e^{-i\\Delta t H_1} e^{-i\\frac{\\Delta t}{2} H_2}.\n", "$$\n", "\n", "Unter der [Jordan-Wigner-Transformation](https://en.wikipedia.org/wiki/Jordan%E2%80%93Wigner_transformation) läuft die Zeitentwicklung durch $H_2$ auf ein einziges [CPhase-Gate](/docs/api/qiskit/qiskit.circuit.library.CPhaseGate) zwischen den Spin-up- und Spin-down-Orbitalen an der Verunreinigungsstelle hinaus. Da es sich bei $H_1$ um einen quadratischen fermionischen Hamiltonian handelt, kommt die Zeitentwicklung über $H_1$ einer Orbitalrotation gleich.\n", "\n", "Die Krylov-Basiszustände $\\{ |\\psi_k\\rangle \\}_{k=0}^{D-1}$, wobei $D$ die Dimension des Krylov-Unterraums ist, werden durch wiederholte Anwendung eines einzigen Trotter-Schritts gebildet, so dass\n", "\n", "$$\n", " |\\psi_k\\rangle \\approx \\left[e^{-i\\frac{\\Delta t}{2} H_2} e^{-i\\Delta t H_1} e^{-i\\frac{\\Delta t}{2} H_2} \\right]^k\\ket{\\psi_0}.\n", "$$\n", "\n", "Im folgenden SQD-basierten Arbeitsablauf werden wir aus diesem Satz von Schaltkreisen Stichproben nehmen und den kombinierten Satz von Bitstrings mit SQD nachbearbeiten. Diese Herangehensweise steht im Gegensatz zu der im verwandten Tutorial [Sample-based quantum diagonalization of a chemistry Hamiltonian](/docs/tutorials/sample-based-quantum-diagonalization) verwendeten, bei der die Proben aus einem einzigen heuristischen Variationskreis gezogen wurden.\n", "\n" ] }, { "cell_type": "markdown", "id": "df66d697-102c-4a2c-80f3-f67fdda05573", "metadata": {}, "source": [ "\n", "\n", "## Anforderungen\n", "\n", "Vergewissern Sie sich vor Beginn dieses Tutorials, dass Sie die folgenden Komponenten installiert haben:\n", "\n", "* Qiskit SDK v1.0 oder später, mit [Visual](/docs/api/qiskit/visualization) isierungsunterstützung\n", "* Qiskit Runtime v0.22 oder später (`pip install qiskit-ibm-runtime`)\n", "* SQD Qiskit-Add-on v0.11 oder höher (`pip install qiskit-addon-sqd`)\n", "* ffsim v0.0.72 oder höher (`pip install ffsim`)\n", "\n" ] }, { "cell_type": "markdown", "id": "c3d4e5f6-a7b8-9012-cdef-123456789012", "metadata": {}, "source": [ "\n", "\n", "## Beispiel für einen Simulator im kleinen Maßstab\n", "\n" ] }, { "cell_type": "markdown", "id": "8540487a-8033-49c2-9f30-022336105f64", "metadata": {}, "source": [ "\n", "\n", "### Schritt 1: Problem auf einen Quantenschaltkreis abbilden\n", "\n" ] }, { "cell_type": "markdown", "id": "4e6652bd-c97d-4f2f-98a7-d44857805cf1", "metadata": {}, "source": [ "Zunächst erzeugen wir den SIAM-Hamiltonschen in der Positionsbasis. Der Hamiltonian wird durch die Matrix $h_{pq}$ und den Tensor $h_{pqrs}$ dargestellt. Dann drehen wir ihn in die Impulsbasis um. In der Positionsbasis setzen wir die Verunreinigung an die erste Stelle. Wenn wir jedoch zur Impulsbasis wechseln, verschieben wir die Verunreinigung an einen zentralen Ort, um Wechselwirkungen mit anderen Orbitalen zu erleichtern.\n", "\n" ] }, { "cell_type": "code", "execution_count": 1, "id": "b5cb9c28-4721-4141-8665-96885038e210", "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "import pyscf.fci\n", "\n", "\n", "def siam_hamiltonian(\n", " norb: int,\n", " hopping: float,\n", " onsite: float,\n", " hybridization: float,\n", " chemical_potential: float,\n", ") -> tuple[np.ndarray, np.ndarray]:\n", " \"\"\"Hamiltonian for the single-impurity Anderson model.\"\"\"\n", " # Place the impurity on the first site\n", " impurity_orb = 0\n", "\n", " # One body matrix elements in the \"position\" basis\n", " h1e = np.zeros((norb, norb))\n", " np.fill_diagonal(h1e[:, 1:], -hopping)\n", " np.fill_diagonal(h1e[1:, :], -hopping)\n", " h1e[impurity_orb, impurity_orb + 1] = -hybridization\n", " h1e[impurity_orb + 1, impurity_orb] = -hybridization\n", " h1e[impurity_orb, impurity_orb] = chemical_potential\n", "\n", " # Two body matrix elements in the \"position\" basis\n", " h2e = np.zeros((norb, norb, norb, norb))\n", " h2e[impurity_orb, impurity_orb, impurity_orb, impurity_orb] = onsite\n", "\n", " return h1e, h2e\n", "\n", "\n", "def momentum_basis(norb: int) -> np.ndarray:\n", " \"\"\"Get the orbital rotation to change from the position to the momentum basis.\"\"\"\n", " n_bath = norb - 1\n", "\n", " # Orbital rotation that diagonalizes the bath (non-interacting system)\n", " hopping_matrix = np.zeros((n_bath, n_bath))\n", " np.fill_diagonal(hopping_matrix[:, 1:], -1)\n", " np.fill_diagonal(hopping_matrix[1:, :], -1)\n", " _, vecs = np.linalg.eigh(hopping_matrix)\n", "\n", " # Expand to include impurity\n", " orbital_rotation = np.zeros((norb, norb))\n", " # Impurity is on the first site\n", " orbital_rotation[0, 0] = 1\n", " orbital_rotation[1:, 1:] = vecs\n", "\n", " # Move the impurity to the center\n", " new_index = n_bath // 2\n", " perm = np.r_[1 : (new_index + 1), 0, (new_index + 1) : norb]\n", " orbital_rotation = orbital_rotation[:, perm]\n", "\n", " return orbital_rotation\n", "\n", "\n", "def rotated(\n", " h1e: np.ndarray, h2e: np.ndarray, orbital_rotation: np.ndarray\n", ") -> tuple[np.ndarray, np.ndarray]:\n", " \"\"\"Rotate the orbital basis of a Hamiltonian.\"\"\"\n", " h1e_rotated = np.einsum(\n", " \"ab,Aa,Bb->AB\",\n", " h1e,\n", " orbital_rotation,\n", " orbital_rotation.conj(),\n", " optimize=\"greedy\",\n", " )\n", " h2e_rotated = np.einsum(\n", " \"abcd,Aa,Bb,Cc,Dd->ABCD\",\n", " h2e,\n", " orbital_rotation,\n", " orbital_rotation.conj(),\n", " orbital_rotation,\n", " orbital_rotation.conj(),\n", " optimize=\"greedy\",\n", " )\n", " return h1e_rotated, h2e_rotated\n", "\n", "\n", "# Total number of spatial orbitals, including the bath sites and the impurity\n", "# This should be an even number\n", "norb = 8\n", "\n", "# System is half-filled\n", "nelec = (norb // 2, norb // 2)\n", "# One orbital is the impurity, the rest are bath sites\n", "n_bath = norb - 1\n", "\n", "# Hamiltonian parameters\n", "hybridization = 1.0\n", "hopping = 1.0\n", "onsite = 10.0\n", "chemical_potential = -0.5 * onsite\n", "\n", "# Generate Hamiltonian in position basis\n", "h1e, h2e = siam_hamiltonian(\n", " norb=norb,\n", " hopping=hopping,\n", " onsite=onsite,\n", " hybridization=hybridization,\n", " chemical_potential=chemical_potential,\n", ")\n", "\n", "# Rotate to momentum basis\n", "orbital_rotation = momentum_basis(norb)\n", "h1e_momentum, h2e_momentum = rotated(h1e, h2e, orbital_rotation.T.conj())\n", "# In the momentum basis, the impurity is placed in the center\n", "impurity_index = n_bath // 2\n", "\n", "# Use PySCF to compute the exact ground state energy\n", "reference_energy, _ = pyscf.fci.direct_spin1.kernel(h1e, h2e, norb, nelec)" ] }, { "cell_type": "markdown", "id": "e888edf2-7865-41a8-be57-6ddb72dd0cc7", "metadata": {}, "source": [ "Als Nächstes erzeugen wir die Schaltkreise, um die Krylov-Basiszustände zu erzeugen.\n", "Für jede Spinspezies ist der Ausgangszustand $\\ket{\\psi_0}$ durch die Überlagerung aller möglichen Anregungen der drei dem Fermi-Niveau am nächsten liegenden Elektronen in die 4 nächstgelegenen leeren Moden, ausgehend vom Zustand $|00\\cdots 0011 \\cdots 11\\rangle$, gegeben und wird durch die Anwendung von sieben [XXPlusYYGates](/docs/api/qiskit/qiskit.circuit.library.XXPlusYYGate) realisiert.\n", "Die zeitlich veränderten Zustände werden durch aufeinanderfolgende Anwendungen eines Trotter-Schrittes zweiter Ordnung erzeugt.\n", "\n", "Eine ausführlichere Beschreibung dieses Modells und des Entwurfs der Schaltkreise finden Sie unter [\"Quantum-Centric Algorithm for Sample-Based Krylov Diagonalization\"](https://arxiv.org/abs/2501.09702).\n", "\n" ] }, { "cell_type": "code", "execution_count": 2, "id": "0f729f86-1814-4d5a-ae65-3f9614e103b3", "metadata": {}, "outputs": [], "source": [ "from typing import Sequence\n", "\n", "import ffsim\n", "import scipy\n", "from qiskit import QuantumCircuit, QuantumRegister\n", "from qiskit.circuit import CircuitInstruction, Qubit\n", "from qiskit.circuit.library import CPhaseGate, XGate, XXPlusYYGate\n", "\n", "\n", "def prepare_initial_state(qubits: Sequence[Qubit], norb: int, nocc: int):\n", " \"\"\"Prepare initial state.\"\"\"\n", " assert norb >= 8\n", " x_gate = XGate()\n", " rot = XXPlusYYGate(0.5 * np.pi, -0.5 * np.pi)\n", " for i in range(nocc):\n", " yield CircuitInstruction(x_gate, [qubits[i]])\n", " yield CircuitInstruction(x_gate, [qubits[norb + i]])\n", " for i in range(3):\n", " for j in range(nocc - i - 1, nocc + i, 2):\n", " yield CircuitInstruction(rot, [qubits[j], qubits[j + 1]])\n", " yield CircuitInstruction(\n", " rot, [qubits[norb + j], qubits[norb + j + 1]]\n", " )\n", " yield CircuitInstruction(rot, [qubits[j + 1], qubits[j + 2]])\n", " yield CircuitInstruction(\n", " rot, [qubits[norb + j + 1], qubits[norb + j + 2]]\n", " )\n", "\n", "\n", "def trotter_step(\n", " qubits: Sequence[Qubit],\n", " time_step: float,\n", " one_body_evolution: np.ndarray,\n", " h2e: np.ndarray,\n", " impurity_index: int,\n", " norb: int,\n", "):\n", " \"\"\"A Trotter step.\"\"\"\n", " # Assume the two-body interaction is just the on-site interaction of the impurity\n", " onsite = h2e[\n", " impurity_index, impurity_index, impurity_index, impurity_index\n", " ]\n", " # Two-body evolution for half the time\n", " yield CircuitInstruction(\n", " CPhaseGate(-0.5 * time_step * onsite),\n", " [qubits[impurity_index], qubits[norb + impurity_index]],\n", " )\n", " # One-body evolution for the full time\n", " yield CircuitInstruction(\n", " ffsim.qiskit.OrbitalRotationJW(norb, one_body_evolution), qubits\n", " )\n", " # Two-body evolution for half the time\n", " yield CircuitInstruction(\n", " CPhaseGate(-0.5 * time_step * onsite),\n", " [qubits[impurity_index], qubits[norb + impurity_index]],\n", " )\n", "\n", "\n", "# Time step\n", "time_step = 0.2\n", "# Number of Krylov basis states\n", "krylov_dim = 8\n", "\n", "# Initialize circuit\n", "qubits = QuantumRegister(2 * norb, name=\"q\")\n", "circuit = QuantumCircuit(qubits)\n", "\n", "# Generate initial state\n", "for instruction in prepare_initial_state(qubits, norb=norb, nocc=norb // 2):\n", " circuit.append(instruction)\n", "circuit.measure_all()\n", "\n", "# Create list of circuits, starting with the initial state circuit\n", "circuits = [circuit.copy()]\n", "\n", "# Add time evolution circuits to the list\n", "one_body_evolution = scipy.linalg.expm(-1j * time_step * h1e_momentum)\n", "for i in range(krylov_dim - 1):\n", " # Remove measurements\n", " circuit.remove_final_measurements()\n", " # Append another Trotter step\n", " for instruction in trotter_step(\n", " qubits,\n", " time_step,\n", " one_body_evolution,\n", " h2e_momentum,\n", " impurity_index,\n", " norb,\n", " ):\n", " circuit.append(instruction)\n", " # Measure qubits\n", " circuit.measure_all()\n", " # Add a copy of the circuit to the list\n", " circuits.append(circuit.copy())" ] }, { "cell_type": "code", "execution_count": 3, "id": "9f2cc4d4-ecac-457a-bcae-558319668e1f", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "\"Output" ] }, "execution_count": 3, "metadata": {}, "output_type": "execute_result" } ], "source": [ "circuits[0].draw(\"mpl\", scale=0.4, fold=-1)" ] }, { "cell_type": "code", "execution_count": 4, "id": "827976ec-4815-4707-80b1-e13fb2fef309", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "\"Output" ] }, "execution_count": 4, "metadata": {}, "output_type": "execute_result" } ], "source": [ "circuits[-1].draw(\"mpl\", scale=0.4, fold=-1)" ] }, { "cell_type": "markdown", "id": "b8ca6be4-61d9-47be-8099-8712c7ecc774", "metadata": {}, "source": [ "\n", "\n", "### Schritt 2: Optimierung des Problems für die Quantenausführung\n", "\n" ] }, { "cell_type": "markdown", "id": "a3304e1b-9c7c-4212-8744-d1c62292eced", "metadata": {}, "source": [ "Als Nächstes optimieren wir die Schaltung für eine bestimmte Hardware. Zunächst erstellen wir ein generisches Backend mit einer bestimmten Anzahl von Qubits und einem Gatter-Satz, in den sich die Zeitentwicklungsschaltungen auf natürliche Weise zerlegen lassen.\n", "\n" ] }, { "cell_type": "code", "execution_count": 5, "id": "2d2fdbff-1e22-45af-a2eb-c334e4328c59", "metadata": {}, "outputs": [], "source": [ "from qiskit.providers.fake_provider import GenericBackendV2\n", "\n", "backend = GenericBackendV2(\n", " 2 * norb, basis_gates=[\"cp\", \"xx_plus_yy\", \"p\", \"x\"]\n", ")" ] }, { "cell_type": "markdown", "id": "ce3e6a52-7b99-49b6-8294-b005efa59cfc", "metadata": {}, "source": [ "Jetzt verwenden wir Qiskit, um die Schaltkreise in das Ziel-Backend zu transpilieren.\n", "\n" ] }, { "cell_type": "code", "execution_count": 6, "id": "c8643533-9fec-40bf-a307-da8839b1e444", "metadata": {}, "outputs": [], "source": [ "from qiskit.transpiler import generate_preset_pass_manager\n", "\n", "pass_manager = generate_preset_pass_manager(\n", " optimization_level=3, backend=backend\n", ")\n", "isa_circuits = pass_manager.run(circuits)" ] }, { "cell_type": "markdown", "id": "6cfd3eea-e2d9-40a5-a449-1d3d790a5f2d", "metadata": {}, "source": [ "\n", "\n", "### Schritt 3: Führen Sie den Befehl „ Qiskit primitives “ aus\n", "\n" ] }, { "cell_type": "markdown", "id": "ad48e3e6-1013-45a1-942c-63766fed5819", "metadata": {}, "source": [ "Nach der Optimierung der Schaltungen für die Hardware-Ausführung können wir sie auf der Zielhardware ausführen und Proben für die Abschätzung der Grundzustandsenergie sammeln. Nachdem wir mit dem Sampler-Primitiv Bitstrings aus jedem Schaltkreis gesampelt haben, fassen wir alle Ergebnisse in einem einzigen Zählwörterbuch zusammen und zeichnen die 20 am häufigsten gesampelten Bitstrings auf.\n", "\n" ] }, { "cell_type": "code", "execution_count": null, "id": "80eee553-60d6-4258-88ab-d8d120418c36", "metadata": {}, "outputs": [], "source": [ "from qiskit.visualization import plot_histogram\n", "from qiskit.primitives import StatevectorSampler\n", "\n", "# Sample from the circuits\n", "sampler = StatevectorSampler()\n", "job = sampler.run(isa_circuits, shots=500)" ] }, { "cell_type": "code", "execution_count": 8, "id": "10af4663-7375-4b50-bae6-9f3d5106457b", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "\"Output" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "from qiskit.primitives import BitArray\n", "\n", "# Combine the shots from the individual Trotter circuits\n", "bit_array = BitArray.concatenate_shots(\n", " [result.data.meas for result in job.result()]\n", ")\n", "\n", "plot_histogram(bit_array.get_counts(), number_to_keep=20)" ] }, { "cell_type": "markdown", "id": "2aa74455-d16b-4ac3-a354-54a79d5c5759", "metadata": {}, "source": [ "\n", "\n", "### Schritt 4: Nachbearbeitung und Rückgabe des Ergebnisses im gewünschten klassischen Format\n", "\n" ] }, { "cell_type": "markdown", "id": "d3f7713c-a4e6-407f-94de-0e5cfeb1134c", "metadata": {}, "source": [ "Nun führen wir den SQD-Algorithmus mit der Funktion `diagonalize_fermionic_hamiltonian` aus. Erläuterungen zu den Argumenten für diese Funktion finden Sie in der [API-Dokumentation](https://qiskit.github.io/qiskit-addon-sqd/apidocs/qiskit_addon_sqd.fermion.html#qiskit_addon_sqd.fermion.diagonalize_fermionic_hamiltonian).\n", "\n" ] }, { "cell_type": "code", "execution_count": 9, "id": "7609d1e1-e8ef-48e1-a965-97927f403163", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Iteration 1\n", "\tSubsample 0\n", "\t\tEnergy: -13.4222953188441\n", "\t\tSubspace dimension: 529\n", "\tSubsample 1\n", "\t\tEnergy: -13.42237556285828\n", "\t\tSubspace dimension: 784\n", "\tSubsample 2\n", "\t\tEnergy: -13.422045397387413\n", "\t\tSubspace dimension: 529\n", "Iteration 2\n", "\tSubsample 0\n", "\t\tEnergy: -13.422379583305478\n", "\t\tSubspace dimension: 900\n", "\tSubsample 1\n", "\t\tEnergy: -13.422376197704326\n", "\t\tSubspace dimension: 841\n", "\tSubsample 2\n", "\t\tEnergy: -13.422421162849295\n", "\t\tSubspace dimension: 1089\n", "Iteration 3\n", "\tSubsample 0\n", "\t\tEnergy: -13.422421164670345\n", "\t\tSubspace dimension: 1156\n", "\tSubsample 1\n", "\t\tEnergy: -13.422421492737689\n", "\t\tSubspace dimension: 1156\n", "\tSubsample 2\n", "\t\tEnergy: -13.422421205869572\n", "\t\tSubspace dimension: 1156\n", "Iteration 4\n", "\tSubsample 0\n", "\t\tEnergy: -13.422421494558726\n", "\t\tSubspace dimension: 1225\n", "\tSubsample 1\n", "\t\tEnergy: -13.422421492737689\n", "\t\tSubspace dimension: 1156\n", "\tSubsample 2\n", "\t\tEnergy: -13.422421492737689\n", "\t\tSubspace dimension: 1156\n" ] } ], "source": [ "from qiskit_addon_sqd.fermion import (\n", " SCIResult,\n", " diagonalize_fermionic_hamiltonian,\n", ")\n", "\n", "# List to capture intermediate results\n", "result_history = []\n", "\n", "\n", "def callback(results: list[SCIResult]):\n", " result_history.append(results)\n", " iteration = len(result_history)\n", " print(f\"Iteration {iteration}\")\n", " for i, result in enumerate(results):\n", " print(f\"\\tSubsample {i}\")\n", " print(f\"\\t\\tEnergy: {result.energy}\")\n", " print(\n", " f\"\\t\\tSubspace dimension: {np.prod(result.sci_state.amplitudes.shape)}\"\n", " )\n", "\n", "\n", "rng = np.random.default_rng(24)\n", "result = diagonalize_fermionic_hamiltonian(\n", " h1e_momentum,\n", " h2e_momentum,\n", " bit_array,\n", " samples_per_batch=100,\n", " norb=norb,\n", " nelec=nelec,\n", " num_batches=3,\n", " max_iterations=5,\n", " symmetrize_spin=True,\n", " callback=callback,\n", " seed=rng,\n", ")" ] }, { "cell_type": "markdown", "id": "3dee9c61-fc42-48e9-8888-af8fc831cd5c", "metadata": {}, "source": [ "Die folgende Codezelle stellt die Ergebnisse grafisch dar. Die erste Grafik zeigt die berechnete Energie als Funktion der Anzahl der Iterationen zur Konfigurationswiederherstellung, und die zweite Grafik zeigt die durchschnittliche Besetzung jedes räumlichen Orbitals nach der letzten Iteration. Da es sich um ein so kleines Problem handelt, kommen wir bereits mit der ersten Iteration sehr nahe an den exakten Energiewert heran (beachte die Skalierung der y-Achse).\n", "\n" ] }, { "cell_type": "code", "execution_count": 10, "id": "b6879566-8bf5-4c28-bfb6-b2686692e3d3", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Reference energy: -13.42249\n", "SQD energy: -13.42242\n", "Absolute error: 0.00007\n" ] }, { "data": { "text/plain": [ "\"Output" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "import matplotlib.pyplot as plt\n", "\n", "min_es = [\n", " min(result, key=lambda res: res.energy).energy\n", " for result in result_history\n", "]\n", "min_id, min_e = min(enumerate(min_es), key=lambda x: x[1])\n", "\n", "# Data for energies plot\n", "x1 = range(len(result_history))\n", "\n", "# Data for avg spatial orbital occupancy\n", "y2 = np.sum(result.orbital_occupancies, axis=0)\n", "x2 = range(len(y2))\n", "\n", "fig, axs = plt.subplots(1, 2, figsize=(12, 6))\n", "\n", "# Plot energies\n", "axs[0].plot(x1, min_es, label=\"energy\", marker=\"o\")\n", "axs[0].set_xticks(x1)\n", "axs[0].set_xticklabels(x1)\n", "axs[0].axhline(\n", " y=reference_energy,\n", " color=\"#BF5700\",\n", " linestyle=\"--\",\n", " label=\"reference energy\",\n", ")\n", "axs[0].set_title(\"Approximated Ground State Energy vs SQD Iterations\")\n", "axs[0].set_xlabel(\"Iteration Index\", fontdict={\"fontsize\": 12})\n", "axs[0].set_ylabel(\"Energy\", fontdict={\"fontsize\": 12})\n", "axs[0].legend()\n", "\n", "# Plot orbital occupancy\n", "axs[1].bar(x2, y2, width=0.8)\n", "axs[1].set_xticks(x2)\n", "axs[1].set_xticklabels(x2)\n", "axs[1].set_title(\"Avg Occupancy per Spatial Orbital\")\n", "axs[1].set_xlabel(\"Orbital Index\", fontdict={\"fontsize\": 12})\n", "axs[1].set_ylabel(\"Avg Occupancy\", fontdict={\"fontsize\": 12})\n", "\n", "print(f\"Reference energy: {reference_energy:.5f}\")\n", "print(f\"SQD energy: {min_e:.5f}\")\n", "print(f\"Absolute error: {abs(min_e - reference_energy):.5f}\")\n", "plt.tight_layout()\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "0c9f4976-d770-426a-822e-e6756c2cfbe7", "metadata": {}, "source": [ "\n", "\n", "### Überprüfen Sie den Energieverbrauch\n", "\n", "Die von SQD zurückgegebene Energie ist garantiert eine Obergrenze der tatsächlichen Grundzustandsenergie. Der Wert der Energie lässt sich überprüfen, da SQD auch die Koeffizienten des Zustandsvektors zurückgibt, der den Grundzustand approximiert. Sie können die Energie aus dem Zustandsvektor anhand seiner reduzierten Dichtematrizen für ein und zwei Teilchen berechnen, wie in der folgenden Codezelle gezeigt.\n", "\n" ] }, { "cell_type": "code", "execution_count": 11, "id": "e2b9de72-61cf-49d3-a1b5-f043e4b16956", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Recomputed energy: -13.42242\n" ] } ], "source": [ "rdm1 = result.sci_state.rdm(rank=1, spin_summed=True)\n", "rdm2 = result.sci_state.rdm(rank=2, spin_summed=True)\n", "\n", "energy = np.sum(h1e_momentum * rdm1) + 0.5 * np.sum(h2e_momentum * rdm2)\n", "\n", "print(f\"Recomputed energy: {energy:.5f}\")" ] }, { "cell_type": "markdown", "id": "5221928a-79ff-4f54-90b4-1fe4d7739aae", "metadata": {}, "source": [ "\n", "\n", "## Beispiel für groß angelegte Hardware\n", "\n", "Nun führen wir ein umfangreicheres Beispiel auf einer echten QPU aus.\n", "Als Referenzenergie verwenden wir die Ergebnisse einer separat durchgeführten [DMRG](https://en.wikipedia.org/wiki/Density_matrix_renormalization_group) -Berechnung.\n", "\n" ] }, { "cell_type": "code", "execution_count": null, "id": "933037d8-847e-4986-80da-5ac8d677b2ff", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Using backend ibm_boston\n", "Iteration 1\n", "\tSubsample 0\n", "\t\tEnergy: -28.63965951544449\n", "\t\tSubspace dimension: 9801\n", "\tSubsample 1\n", "\t\tEnergy: -28.625588929202006\n", "\t\tSubspace dimension: 9409\n", "\tSubsample 2\n", "\t\tEnergy: -28.647371834135498\n", "\t\tSubspace dimension: 8281\n", "Iteration 2\n", "\tSubsample 0\n", "\t\tEnergy: -28.67213260849567\n", "\t\tSubspace dimension: 29584\n", "\tSubsample 1\n", "\t\tEnergy: -28.670340686158816\n", "\t\tSubspace dimension: 27225\n", "\tSubsample 2\n", "\t\tEnergy: -28.669976379525988\n", "\t\tSubspace dimension: 31329\n", "Iteration 3\n", "\tSubsample 0\n", "\t\tEnergy: -28.68622875601382\n", "\t\tSubspace dimension: 36100\n", "\tSubsample 1\n", "\t\tEnergy: -28.698569623143126\n", "\t\tSubspace dimension: 34225\n", "\tSubsample 2\n", "\t\tEnergy: -28.694848533971882\n", "\t\tSubspace dimension: 33856\n", "Iteration 4\n", "\tSubsample 0\n", "\t\tEnergy: -28.69883392844593\n", "\t\tSubspace dimension: 42025\n", "\tSubsample 1\n", "\t\tEnergy: -28.701289495200996\n", "\t\tSubspace dimension: 38025\n", "\tSubsample 2\n", "\t\tEnergy: -28.699319594978245\n", "\t\tSubspace dimension: 45369\n", "Iteration 5\n", "\tSubsample 0\n", "\t\tEnergy: -28.701936886834154\n", "\t\tSubspace dimension: 51076\n", "\tSubsample 1\n", "\t\tEnergy: -28.702468711812013\n", "\t\tSubspace dimension: 53824\n", "\tSubsample 2\n", "\t\tEnergy: -28.702298147575938\n", "\t\tSubspace dimension: 52900\n", "Reference energy: -28.70660\n", "SQD energy: -28.70247\n", "Absolute error: 0.00413\n" ] }, { "data": { "text/plain": [ "\"Output" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "from qiskit_ibm_runtime import SamplerV2 as Sampler\n", "from qiskit_ibm_runtime import QiskitRuntimeService\n", "\n", "# Model parameters\n", "norb = 20\n", "nelec = (norb // 2, norb // 2)\n", "n_bath = norb - 1\n", "hybridization = 1.0\n", "hopping = 1.0\n", "onsite = 10.0\n", "chemical_potential = -0.5 * onsite\n", "\n", "# Generate Hamiltonian and orbital rotation\n", "h1e, h2e = siam_hamiltonian(\n", " norb=norb,\n", " hopping=hopping,\n", " onsite=onsite,\n", " hybridization=hybridization,\n", " chemical_potential=chemical_potential,\n", ")\n", "orbital_rotation = momentum_basis(norb)\n", "h1e_momentum, h2e_momentum = rotated(h1e, h2e, orbital_rotation.T.conj())\n", "impurity_index = n_bath // 2\n", "\n", "# Set reference energy to DMRG value computed separately\n", "reference_energy = -28.70659686\n", "\n", "# Algorithm parameters\n", "time_step = 0.2\n", "krylov_dim = 8\n", "\n", "# Construct circuits\n", "qubits = QuantumRegister(2 * norb, name=\"q\")\n", "circuit = QuantumCircuit(qubits)\n", "for instruction in prepare_initial_state(qubits, norb=norb, nocc=norb // 2):\n", " circuit.append(instruction)\n", "circuit.measure_all()\n", "circuits = [circuit.copy()]\n", "one_body_evolution = scipy.linalg.expm(-1j * time_step * h1e_momentum)\n", "for i in range(krylov_dim - 1):\n", " circuit.remove_final_measurements()\n", " for instruction in trotter_step(\n", " qubits,\n", " time_step,\n", " one_body_evolution,\n", " h2e_momentum,\n", " impurity_index,\n", " norb,\n", " ):\n", " circuit.append(instruction)\n", " circuit.measure_all()\n", " circuits.append(circuit.copy())\n", "\n", "# Initialize hardware backend\n", "service = QiskitRuntimeService()\n", "backend = service.least_busy(\n", " operational=True, simulator=False, min_num_qubits=127\n", ")\n", "print(f\"Using backend {backend.name}\")\n", "\n", "# Transpile to backend\n", "pass_manager = generate_preset_pass_manager(\n", " optimization_level=3, backend=backend\n", ")\n", "isa_circuits = pass_manager.run(circuits)\n", "\n", "# Sample from the circuits\n", "sampler = Sampler(backend)\n", "sampler.options.environment.job_tags = [\"TUT_SKQD\"]\n", "job = sampler.run(isa_circuits, shots=500)\n", "\n", "# Combine the shots from the individual Trotter circuits\n", "bit_array = BitArray.concatenate_shots(\n", " [result.data.meas for result in job.result()]\n", ")\n", "\n", "# Run configuration recovery and diagonalization\n", "result_history = []\n", "\n", "\n", "def callback(results: list[SCIResult]):\n", " result_history.append(results)\n", " iteration = len(result_history)\n", " print(f\"Iteration {iteration}\")\n", " for i, result in enumerate(results):\n", " print(f\"\\tSubsample {i}\")\n", " print(f\"\\t\\tEnergy: {result.energy}\")\n", " print(\n", " f\"\\t\\tSubspace dimension: {np.prod(result.sci_state.amplitudes.shape)}\"\n", " )\n", "\n", "\n", "rng = np.random.default_rng(24)\n", "result = diagonalize_fermionic_hamiltonian(\n", " h1e_momentum,\n", " h2e_momentum,\n", " bit_array,\n", " samples_per_batch=100,\n", " norb=norb,\n", " nelec=nelec,\n", " num_batches=3,\n", " max_iterations=5,\n", " symmetrize_spin=True,\n", " callback=callback,\n", " seed=rng,\n", ")\n", "\n", "\n", "# Plot results\n", "min_es = [\n", " min(result, key=lambda res: res.energy).energy\n", " for result in result_history\n", "]\n", "min_id, min_e = min(enumerate(min_es), key=lambda x: x[1])\n", "x1 = range(len(result_history))\n", "y2 = np.sum(result.orbital_occupancies, axis=0)\n", "x2 = range(len(y2))\n", "fig, axs = plt.subplots(1, 2, figsize=(12, 6))\n", "axs[0].plot(x1, min_es, label=\"energy\", marker=\"o\")\n", "axs[0].set_xticks(x1)\n", "axs[0].set_xticklabels(x1)\n", "axs[0].axhline(\n", " y=reference_energy,\n", " color=\"#BF5700\",\n", " linestyle=\"--\",\n", " label=\"reference energy\",\n", ")\n", "axs[0].set_title(\"Approximated Ground State Energy vs SQD Iterations\")\n", "axs[0].set_xlabel(\"Iteration Index\", fontdict={\"fontsize\": 12})\n", "axs[0].set_ylabel(\"Energy\", fontdict={\"fontsize\": 12})\n", "axs[0].legend()\n", "axs[1].bar(x2, y2, width=0.8)\n", "axs[1].set_xticks(x2)\n", "axs[1].set_xticklabels(x2)\n", "axs[1].set_title(\"Avg Occupancy per Spatial Orbital\")\n", "axs[1].set_xlabel(\"Orbital Index\", fontdict={\"fontsize\": 12})\n", "axs[1].set_ylabel(\"Avg Occupancy\", fontdict={\"fontsize\": 12})\n", "print(f\"Reference energy: {reference_energy:.5f}\")\n", "print(f\"SQD energy: {min_e:.5f}\")\n", "print(f\"Absolute error: {abs(min_e - reference_energy):.5f}\")\n", "plt.tight_layout()\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "482ebea3-84b8-471b-bddc-b282e23162ad", "metadata": {}, "source": [ "\n", "\n", "## Nächste Schritte\n", "\n", "\n", " Wenn Sie diese Arbeit interessant fanden, könnten Sie auch an folgenden Materialien interessiert sein:\n", "\n", " * [Beispielbasierte Quanten-Diagonalisierung eines chemischen Hamilton-Operators](/docs/tutorials/sample-based-quantum-diagonalization) – ein dazugehöriges Tutorial, das anstelle von Trotter-Schaltkreisen einen heuristischen Variationsansatz verwendet\n", " * [Krylov-Quanten-Diagonalisierung von Gitter-Hamiltonianen](/docs/tutorials/krylov-quantum-diagonalization) – eine Einführung in die KQD-Methode\n", " * [SQD-Addon-API-Dokumentation](https://qiskit.github.io/qiskit-addon-sqd/apidocs/qiskit_addon_sqd.fermion.html#qiskit_addon_sqd.fermion.diagonalize_fermionic_hamiltonian) – Referenz zur Funktion `diagonalize_fermionic_hamiltonian`\n", " * [*„Quantum-Centric Algorithm for Sample-Based Krylov Diagonalization*](https://arxiv.org/abs/2501.09702) “ – der Artikel, auf dem dieses Tutorial basiert\n", "\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": 9 }, "nbformat": 4, "nbformat_minor": 4 }