Diagonalisation quantique de Krylov basée sur des échantillons (SKQD)

Cette leçon sur la diagonalisation quantique de Krylov basée sur l'échantillonnage (SKQD) combine les méthodes expliquées dans les méthodes précédentes. Il s'agit d'un exemple unique qui s'appuie sur le cadre de modèles Qiskit :

• Étape 1 : Tracer le problème à l'aide de circuits et d'opérateurs quantiques
• Étape 2 : Optimisation pour le matériel cible
• Étape 3 : Exécution à l'aide des primitives « IBM Quantum »
• Étape 4 : Post-traitement

Une étape importante de la méthode de diagonalisation quantique basée sur l'échantillon consiste à générer des vecteurs de qualité pour le sous-espace. Dans la leçon précédente, nous avons utilisé l'ansatz LUCJ pour générer des vecteurs de sous-espace pour un hamiltonien de chimie. Dans cette leçon, nous utiliserons les états de Krylov quantiques [1], comme nous l'avons vu dans la leçon 2. Tout d'abord, nous verrons comment créer l'espace de Krylov sur un ordinateur quantique en utilisant des opérations d'évolution temporelle. Nous en tirerons ensuite un échantillon. Nous projetterons l'hamiltonien du système sur le sous-espace échantillonné et le diagonaliserons pour estimer l'énergie de l'état fondamental. L'algorithme converge de manière prouvée et efficace vers l'état fondamental, sous les hypothèses décrites dans la leçon 2.

0. L'espace de Krylov

Rappelons qu'un espace de Krylov $\mathcal{K}^r$ d'ordre $r$ est l'espace couvert par les vecteurs obtenus en multipliant les puissances supérieures d'une matrice $A$, jusqu'à $r-1$, avec un vecteur de référence $\vert v \rangle$.

Si la matrice $A$ est l'hamiltonien $H$, l'espace correspondant est appelé espace de Krylov des puissances $\mathcal{K}_P$. Dans le cas où $A$ est l'opérateur d'évolution temporelle généré par l'hamiltonien $U=e^{-iH(dt)}$, l'espace est appelé espace de Krylov unitaire $\mathcal{K}_U$. Le sous-espace de Krylov puissance ne peut pas être généré directement sur un ordinateur quantique car $H$ n'est pas un opérateur unitaire. Au lieu de cela, nous pouvons utiliser l'opérateur d'évolution temporelle $U = e^{-iH(dt)}$, dont on peut montrer qu'il offre des garanties de convergence similaires à celles de l'espace de Krylov puissance. Les puissances de $U$ deviennent alors des pas de temps différents $U^k = e^{-iH(k dt)}$ où $k = 0, 1, 2, ..., (r-1)$.

1. Mapper le problème sur des circuits quantiques et des opérateurs

Dans cette leçon, nous considérons l'hamiltonien de la chaîne antiferromagnétique XX-Z spin-1/2 avec des sites $L = 22$ avec la condition de limite périodique :

from qiskit.transpiler import CouplingMap
from qiskit_addon_utils.problem_generators import generate_xyz_hamiltonian

num_spins = 22
coupling_map = CouplingMap.from_ring(num_spins)
H_op = generate_xyz_hamiltonian(coupling_map, coupling_constants=(0.3, 0.3, 1.0))

Pour construire l'espace de Krylov, nous avons besoin de trois ingrédients principaux :

1. Choix de la dimension de Krylov ( $r$ ) et du pas de temps ( $dt$ ).
2. Un état initial (de référence) (vecteur $\vert v \rangle$ ci-dessus) avec un chevauchement polynomial avec l'état cible (au sol), lorsque l'état cible est peu dense. Cette exigence de chevauchement polynomial est la même que dans l'algorithme d'estimation de la phase quantique.
3. Opérateurs d'évolution temporelle $U^{k}=e^{-iH(k * dt)}$ ( $k = 0, 1, 2, ..., r-1$ ).

Pour une valeur choisie de $r$ (et, $dt$ ), nous créerons $r$ circuits quantiques distincts et les échantillonnerons. Chaque circuit quantique est créé en joignant la représentation du circuit quantique de l'état de référence et l'opérateur d'évolution temporelle pour une valeur $k$.

Une dimension de Krylov plus importante améliore la convergence de l'énergie estimée. Dans cette leçon, nous avons fixé la dimension à $5$ pour illustrer la tendance à la convergence.

Ref [2] a montré qu'un pas de temps suffisamment petit pour la KQD est $\pi / \vert \vert H \vert \vert$, et qu'il est préférable de sous-estimer cette valeur plutôt que de la surestimer. D'autre part, si l'on choisit $dt$ comme étant trop petit, le conditionnement du sous-espace de Krylov est moins bon, car les vecteurs de base de Krylov diffèrent moins d'un pas de temps à l'autre. En outre, bien que ce choix de $dt$ soit adéquat pour la convergence de SKQD, dans ce contexte basé sur l'échantillonnage, le choix optimal de $dt$ dans la pratique est un sujet d'étude en cours. Dans cette leçon, nous avons défini $dt = 0.15$.

Outre la dimension de Krylov et le pas de temps, nous devons définir le nombre de pas de Trotter pour l'évolution temporelle. L'utilisation d'un nombre insuffisant d'étapes conduit à des erreurs de trotterisation plus importantes, tandis qu'un nombre trop élevé d'étapes conduit à des circuits plus profonds. Dans cette leçon, nous avons fixé le nombre de pas de Trotter à $6$.

# Set parameters for quantum Krylov algorithm
krylov_dim = 5  # size of krylov subspace
dt = 0.15
num_trotter_steps = 6

Ensuite, nous devons choisir un état de référence $\vert \psi \rangle$ qui présente un certain chevauchement avec l'état fondamental. Pour cet hamiltonien, nous utilisons l'état de Neel avec une alternance de 1s et 0s $\vert ...101...010...101 \rangle$ comme état de référence.

# Prep `Neel` state as the reference state for evolution
from qiskit import QuantumCircuit

qc_state_prep = QuantumCircuit(num_spins)
for i in range(num_spins):
    if i % 2 == 0:
        qc_state_prep.x(i)

Enfin, nous devons faire correspondre l'opérateur d'évolution temporelle à un circuit quantique. Cela a été fait dans la leçon 2, mais ici nous allons utiliser des méthodes de Qiskit, en particulier une méthode appelée synthèse. Il existe différentes méthodes pour synthétiser les opérateurs mathématiques en circuits quantiques avec des portes quantiques. De nombreuses techniques de ce type sont disponibles dans le module de synthèse Qiskit. Nous utiliserons l'approche LieTrotter de synthèse [3] [4].

from qiskit.circuit import QuantumRegister
from qiskit.circuit.library import PauliEvolutionGate
from qiskit.synthesis import LieTrotter

evol_gate = PauliEvolutionGate(
    H_op, time=(dt / num_trotter_steps), synthesis=LieTrotter(reps=num_trotter_steps)
)  # `U` operator

qr = QuantumRegister(num_spins)
qc_evol = QuantumCircuit(qr)
qc_evol.append(evol_gate, qargs=qr)

circuits = []
for rep in range(krylov_dim):
    circ = qc_state_prep.copy()

    # Repeating the `U` operator to implement U^0, U^1, U^2, and so on, for power Krylov space
    for _ in range(rep):
        circ.compose(other=qc_evol, inplace=True)

    circ.measure_all()
    circuits.append(circ)

circuits[1].decompose().draw("mpl", fold=-1)

Output:

circuits[2].decompose().draw("mpl", fold=-1)

Output:

2. Optimisation pour le matériel cible

Maintenant que nous avons créé les circuits, nous pouvons les optimiser pour un matériel cible. Nous choisissons une QPU à l'échelle d'un service public.

import warnings

from qiskit import generate_preset_pass_manager
from qiskit_ibm_runtime import QiskitRuntimeService

warnings.filterwarnings("ignore")

service = QiskitRuntimeService()
# Use the least-busy backend or specify a quantum computer using the syntax commented out below.
backend = service.least_busy(operational=True, simulator=False)
# backend = service.backend("ibm_brisbane")

Nous transposons ensuite les circuits vers le backend cible à l'aide d'un gestionnaire de passe prédéfini.

pm = generate_preset_pass_manager(backend=backend, optimization_level=3)
isa_circuits = pm.run(circuits=circuits)

3. Exécuter sur le matériel cible

Après avoir optimisé les circuits pour l'exécution matérielle, nous sommes prêts à les exécuter sur le matériel cible et à collecter des échantillons pour l'estimation de l'énergie de l'état fondamental.

from qiskit_ibm_runtime import SamplerV2 as Sampler

sampler = Sampler(mode=backend)
job = sampler.run(isa_circuits, shots=100_000)  # Takes approximately 2m 58s of QPU time

counts_all = [job.result()[k].data.meas.get_counts() for k in range(krylov_dim)]

4. Résultats post-traitement

Ensuite, nous agrégeons les chiffres pour les dimensions de Krylov croissantes de manière cumulative. En utilisant les comptes cumulés, nous couvrirons des sous-espaces pour une dimension de Krylov croissante et analyserons le comportement de convergence.

from collections import Counter

counts_cumulative = []
for i in range(krylov_dim):
    counter = Counter()
    for d in counts_all[: i + 1]:
        counter.update(d)

    counts = dict(counter)
    counts_cumulative.append(counts)

Pour projeter et diagonaliser l'hamiltonien, nous utilisons les capacités de qiskit-addon-sqd. L'addon offre des fonctionnalités permettant de projeter les hamiltoniens basés sur la chaîne de Pauli sur un sous-espace et de résoudre les valeurs propres à l'aide de SciPy.

from qiskit_addon_sqd.counts import counts_to_arrays
from qiskit_addon_sqd.qubit import solve_qubit

En principe, nous pouvons filtrer les chaînes de bits présentant un motif incorrect avant de couvrir le sous-espace. Par exemple, l'état fondamental de l'hamiltonien antiferromagnétique de cette leçon présente généralement un nombre égal de spins "up" et "down", c'est-à-dire que le nombre de "1" dans la chaîne de bits doit être exactement égal à la moitié du nombre total de bits (spins) dans le système. La fonction suivante filtre les chaînes de bits dont le nombre de "1" est incorrect.

# Filters out bitstrings that do not have specified number (`num_ones`) of `1` bits.
def postselect_counts(counts, num_ones):
    filtered_counts = {}
    for bitstring, freq in counts.items():
        if bitstring.count("1") == num_ones:
            filtered_counts[bitstring] = freq

    return filtered_counts

En utilisant des chaînes de bits avec le nombre correct d'électrons ascendants/descendants, nous couvrons des sous-espaces et calculons les valeurs propres pour une dimension de Krylov croissante. En fonction de la taille du problème et des ressources classiques disponibles, il peut s'avérer nécessaire d'adopter un sous-échantillonnage (similaire à la leçon sur le SQD ) pour contrôler la dimension du sous-espace. De plus, nous pouvons appliquer la notion de récupération de configuration similaire à la leçon 4. Nous pouvons calculer l'occupation électronique par site à partir des états propres reconstruits et utiliser ces informations pour corriger les chaînes de bits comportant un nombre incorrect d'électrons ascendants/descendants. Nous laissons cet exercice aux lecteurs intéressés.

import numpy as np

num_batches = 10
rand_seed = 0
scipy_kwargs = {"k": 2, "which": "SA"}

ground_state_energies = []
for idx, counts in enumerate(counts_cumulative):
    counts = postselect_counts(counts, num_ones=num_spins // 2)
    bitstring_matrix, probs = counts_to_arrays(counts=counts)

    eigenvals, eigenstates = solve_qubit(
        bitstring_matrix, H_op, verbose=False, **scipy_kwargs
    )
    gs_en = np.min(eigenvals)
    ground_state_energies.append(gs_en)

Ensuite, nous traçons l'énergie calculée en fonction de la dimension de Krylov et la comparons à l'énergie exacte. L'énergie exacte est calculée séparément à l'aide d'une méthode classique de force brute. Nous pouvons constater que l'énergie estimée de l'état fondamental converge avec l'augmentation de la dimension de l'espace de Krylov. Bien que la dimension de Krylov de $5$ soit limitée, les résultats montrent toujours une convergence impressionnante, qui devrait s'améliorer avec une dimension de Krylov plus grande [1].

import matplotlib.pyplot as plt

exact_gs_en = -23.934184
plt.plot(
    range(1, krylov_dim + 1),
    ground_state_energies,
    color="blue",
    linestyle="-.",
    label="estimate",
)
plt.plot(
    range(1, krylov_dim + 1),
    [exact_gs_en] * krylov_dim,
    color="red",
    linestyle="-",
    label="exact",
)
plt.xticks(range(1, krylov_dim + 1), range(1, krylov_dim + 1))
plt.legend()
plt.xlabel("Krylov space dimension")
plt.ylabel("Energy")
plt.ylim([-24, -22.50])
plt.title(
    "Estimating Ground state energy with Sample-based Krylov Quantum Diagonalization"
)
plt.show()

Output:

Vérifiez votre compréhension

Lisez les questions ci-dessous, réfléchissez à vos réponses, puis cliquez sur les triangles pour trouver les solutions.

• Que pourrait-on faire pour améliorer la convergence dans le graphique ci-dessus?

  Réponse :

  Augmenter la dimension de Krylov. En général, on pourrait aussi augmenter le nombre de tirs, mais celui-ci est déjà assez élevé dans le calcul ci-dessus.

• Quels sont les principaux avantages du SKQD par rapport (a) au SQD et (b) au KQD?

  Réponse :

  Il peut y avoir d'autres réponses valables, mais les réponses complètes doivent comprendre les éléments suivants :

  (a) SKQD offre des garanties de convergence que SQD n'offre pas. En SQD, vous devez soit faire une très bonne supposition pour votre ansatz qui a un excellent chevauchement avec le support de l'état fondamental dans la base de calcul, soit introduire une composante variationnelle dans le calcul pour échantillonner une famille d'ansatz.

  (b) SKQD nécessite beaucoup moins de temps QPU, car il évite le calcul coûteux des éléments de la matrice via le test de Hadamard.

5. Résumé

• L'estimation de l'énergie de l'état fondamental par l'échantillonnage des états de base de Krylov est très bien adaptée aux modèles de treillis, y compris les systèmes de spin, les problèmes de matière condensée et les théories de jauge sur treillis. Cette approche s'adapte beaucoup mieux que l'EQV, car elle ne nécessite pas d'optimisation sur de nombreux paramètres dans un ansatz variationnel comme dans l'EQV, ou dans la SQD basée sur un ansatz heuristique (par exemple, le problème de chimie de la leçon précédente).
  • Pour réduire la profondeur des circuits, il est judicieux d'aborder les problèmes de treillis qui se prêtent à l'utilisation de matériel tolérant aux pannes.
• SKQD ne pose pas de problème de mesure quantique comme dans VQE. Il n'y a pas de groupes d'opérateurs de Pauli commutés à estimer.
• La méthode SKQD est résistante aux échantillons bruités, car il est possible d'utiliser une routine de post-sélection spécifique au problème (par exemple, en filtrant les chaînes binaires qui ne respectent pas les motifs propres au problème) ou d'accepter la surcharge liée à la diagonalisation classique (c'est-à-dire de diagonaliser dans un sous-espace plus grand) afin d'éliminer efficacement l'effet du bruit.

Références

[1] Jeffery Yu et al, "Algorithme centré sur le quantum pour la diagonalisation de Krylov basée sur l'échantillonnage" (2025). arxiv:quant-ph/2501.09702.

[2] Ethan N. Epperly, Lin Lin et Yuji Nakatsukasa. "Une théorie de la diagonalisation du sous-espace quantique". SIAM Journal on Matrix Analysis and Applications 43, 1263-1290 (2022).

[2] N. Hatano et M. Suzuki, "Finding Exponential Product Formulas of Higher Orders" (2005). arXiv:math-ph/0506007.

[4] D. Berry, G. Ahokas, R. Cleve et B. Sanders, "Efficient quantum algorithms for simulating sparse Hamiltonians" (2006). arXiv:quant-ph/0508139.