Estimation de l'énergie de l'état fondamental de la chaîne de Heisenberg avec VQE

Estimation d'utilisation : 37 minutes sur un processeur Heron (REMARQUE : il s'agit d'une estimation uniquement. Ton temps d'exécution peut varier.)

Contexte

Ce tutoriel montre comment construire, déployer et exécuter un workflow de développement appelé patron Qiskit pour simuler une chaîne de Heisenberg et estimer son énergie de l'état fondamental en utilisant l'optimiseur SPSA.

Prérequis

Avant de commencer ce tutoriel, assure-toi que les éléments suivants sont installés :

• Qiskit SDK v2.0 ou version ultérieure, avec le support de visualisation
• Qiskit Runtime v0.44 ou version ultérieure (pip install qiskit-ibm-runtime)

Configuration

# Added by doQumentation — required packages for this notebook
!pip install -q matplotlib numpy qiskit qiskit-ibm-runtime

import numpy as np
import matplotlib.pyplot as plt
from typing import Sequence

from qiskit import QuantumCircuit
from qiskit.quantum_info import SparsePauliOp
from qiskit.primitives import BaseEstimatorV2
from qiskit.circuit.library import XGate
from qiskit.circuit.library import efficient_su2
from qiskit.transpiler import PassManager
from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager
from qiskit.transpiler.passes.scheduling import (
    ALAPScheduleAnalysis,
    PadDynamicalDecoupling,
)
from qiskit_ibm_runtime import QiskitRuntimeService, Session, EstimatorV2

def visualize_results(results):
    plt.plot(results["cost_history"], lw=2)
    plt.xlabel("Number of function evaluations")
    plt.ylabel("Energy")
    plt.show()

Étape 1 : Transposer les entrées classiques en un problème quantique

• Entrée : Nombre de spins
• Sortie : Ansatz et hamiltonien modélisant la chaîne de Heisenberg

Construis un ansatz et un hamiltonien qui modélisent une chaîne de Heisenberg à 10 spins. Tout d'abord, nous importons quelques paquets génériques et créons quelques fonctions utilitaires.

num_spins = 10
ansatz = efficient_su2(num_qubits=num_spins, reps=2)

service = QiskitRuntimeService(
    channel="ibm_cloud",
    token="<YOUR_API_TOKEN>",  # Replace with your actual API token
    instance="<YOUR_INSTANCE_NAME>",  # Replace with your instance name if needed
)
backend = service.least_busy(
    operational=True, min_num_qubits=num_spins, simulator=False
)

coupling = backend.target.build_coupling_map()
reduced_coupling = coupling.reduce(list(range(num_spins)))

edge_list = reduced_coupling.graph.edge_list()
ham_list = []

for edge in edge_list:
    ham_list.append(("ZZ", edge, 0.5))
    ham_list.append(("YY", edge, 0.5))
    ham_list.append(("XX", edge, 0.5))

for qubit in reduced_coupling.physical_qubits:
    ham_list.append(("Z", [qubit], np.random.random() * 2 - 1))

hamiltonian = SparsePauliOp.from_sparse_list(ham_list, num_qubits=num_spins)

ansatz.draw("mpl", style="iqp")

Étape 2 : Optimiser le problème pour l'exécution sur du matériel quantique

• Entrée : Circuit abstrait, observable
• Sortie : Circuit et observable cibles, optimisés pour le QPU sélectionné

Utilise la fonction generate_preset_pass_manager de Qiskit pour générer automatiquement une routine d'optimisation pour notre circuit par rapport au QPU sélectionné. Nous choisissons optimization_level=3, qui fournit le niveau d'optimisation le plus élevé des gestionnaires de passes prédéfinis. Nous incluons également les passes de planification ALAPScheduleAnalysis et PadDynamicalDecoupling pour supprimer les erreurs de décohérence.

target = backend.target
pm = generate_preset_pass_manager(optimization_level=3, target=target)
pm.scheduling = PassManager(
    [
        ALAPScheduleAnalysis(durations=target.durations()),
        PadDynamicalDecoupling(
            durations=target.durations(),
            dd_sequence=[XGate(), XGate()],
            pulse_alignment=target.pulse_alignment,
        ),
    ]
)
isa_ansatz = pm.run(ansatz)
isa_observable = hamiltonian.apply_layout(isa_ansatz.layout)
isa_ansatz.draw("mpl", scale=0.6, style="iqp", fold=-1, idle_wires=False)

Étape 3 : Exécuter à l'aide des primitives Qiskit

• Entrée : Circuit et observable cibles
• Sortie : Résultats de l'optimisation

Minimise l'énergie estimée de l'état fondamental du système en optimisant les paramètres du circuit. Utilise la primitive Estimator de Qiskit Runtime pour évaluer la fonction de coût pendant l'optimisation.

Puisque nous avons optimisé le circuit pour le backend à l'étape 2, nous pouvons éviter la transpilation sur le serveur Runtime en définissant skip_transpilation=True et en passant le circuit optimisé. Pour cette démonstration, nous exécuterons sur un QPU en utilisant les primitives qiskit-ibm-runtime. Pour exécuter avec les primitives basées sur le vecteur d'état de qiskit, remplace le bloc de code utilisant les primitives Qiskit Runtime par le bloc commenté.

Dans ce tutoriel, nous utilisons l'approximation stochastique par perturbation simultanée (SPSA), qui est un optimiseur basé sur le gradient. Ci-après, nous donnons une brève introduction et fournissons le code pour implémenter SPSA avec Qiskit v2.0.

Présentation de SPSA

L'approximation stochastique par perturbation simultanée (SPSA) [1] est un algorithme d'optimisation qui approxime l'intégralité du vecteur gradient en utilisant seulement deux appels de fonction à chaque itération. Soit $f:\mathbb{R}^p\rightarrow \mathbb{R}$ la fonction de coût avec $p$ paramètres à optimiser, et $x_i\in \mathbb{R}^p$ le vecteur de paramètres à la $i^{\text{ième}}$ étape de l'itération. Pour calculer le gradient, un vecteur aléatoire $\Delta_i$ de taille $p$ est créé, où chaque élément $\Delta_{ij}$, $\forall$ $j\in \{1,2,...,p\}$, est échantillonné uniformément dans $\{-1, 1\}$. Ensuite, chaque élément du vecteur aléatoire $\Delta_i$ est multiplié par une petite valeur $c_i$ pour créer une perturbation aléatoire. Le gradient est alors estimé comme

$[\nabla f(x_i)]_j \approx \frac{f(x_i + c_i \Delta_i) - f(x_i - c_i \Delta_i)}{2c_i\Delta_{ij}}.$

Intuitivement, puisqu'une perturbation aléatoire est appliquée lors de l'estimation du gradient, on s'attend à ce que de petites déviations dans les valeurs exactes de $f$ dues au bruit puissent être tolérées et prises en compte. En fait, SPSA est particulièrement connu pour sa robustesse au bruit, et ne nécessite que deux appels au matériel par itération. Il est donc l'un des optimiseurs très préférés pour l'implémentation d'algorithmes variationnels.

Dans ce tutoriel, les hyperparamètres pour la $i^{\text{ième}}$ itération, $a_i$ et $c_i$, sont calculés comme $a_i = a/(A + i + 1)^\alpha$ et $c_i = c / (i+1)^\gamma$, où les valeurs constantes sont $A = 30$, $\alpha = 0,9$, $a = 0,3$, $c = 0,1$, et $\gamma = 0,4$. Ces valeurs sont tirées de [2]. Un réglage approprié des hyperparamètres est nécessaire pour obtenir de bonnes performances avec SPSA.

def spsa(
    fun, x0, args=(), A=30, alpha=0.9, a=0.3, c=0.1, gamma=0.4, maxiter=100
):
    nparams = len(x0)
    x = np.copy(x0)

    for i in range(maxiter):
        a_i = a / (A + i + 1) ** alpha
        c_i = c / (i + 1) ** gamma
        delta_i = np.random.choice([-1, 1], nparams)

        # two hardware calls
        eval_1 = fun(x + c_i * delta_i, *args)
        eval_2 = fun(x - c_i * delta_i, *args)

        # compute the gradient and update the parameters
        grad = (eval_1 - eval_2) / (2 * c_i) * np.reciprocal(delta_i)
        x = x - a_i * grad

    return x

def cost_func(
    params: Sequence,
    ansatz: QuantumCircuit,
    hamiltonian: SparsePauliOp,
    estimator: BaseEstimatorV2,
    cost_history_dict: dict,
) -> float:
    """Ground state energy evaluation."""
    energy = (
        estimator.run([(ansatz, hamiltonian, [params])]).result()[0].data.evs
    )

    cost_history_dict["iters"] += 1
    cost_history_dict["prev_vector"] = list(params)
    cost_history_dict["cost_history"].append(float(energy[0]))

    print(
        f"Fx Iters. done: {cost_history_dict['iters']} [Current cost: {round(energy[0], 5)}]",
        end="\r",
    )

    return energy

def solve(x0, isa_ansatz, isa_observable, maxiter=150):
    cost_history_dict = {
        "prev_vector": None,
        "iters": 0,
        "cost_history": [],
        "y_min": None,
    }

    # Evaluate the problem using a QPU via Qiskit IBM Runtime
    with Session(backend=backend) as session:
        estimator = EstimatorV2(mode=session)
        estimator.skip_transpilation = True
        x_opt = spsa(
            cost_func,
            x0=x0,
            args=(isa_ansatz, isa_observable, estimator, cost_history_dict),
            maxiter=maxiter,
        )

        y_min = cost_func(
            x_opt, isa_ansatz, isa_observable, estimator, cost_history_dict
        )

    return y_min, cost_history_dict

np.random.seed(42)
num_params = ansatz.num_parameters
params = 2 * np.pi * np.random.random(num_params)

Ici nous définissons maxiter = 50. Remarque que chaque itération nécessite deux appels à la fonction pour calculer le gradient, donc le nombre total d'appels de fonction sera $2 \times \text{maxiter}$. Le paramètre maxiter peut être augmenté à n'importe quelle valeur supérieure pour une meilleure estimation de l'énergie.

maxiter = 50
spsa_min, spsa_history = solve(
    params, isa_ansatz, isa_observable, maxiter=maxiter
)

Fx Iters. done: 101 [Current cost: -2.19621]

Étape 4 : Post-traiter et renvoyer le résultat dans le format classique souhaité

• Entrée : Estimations de l'énergie de l'état fondamental pendant l'optimisation
• Sortie : Énergie estimée de l'état fondamental

print(f"Estimated ground state energy: {spsa_min}")

Estimated ground state energy: [-2.19621239]

results = {
    "spsa": spsa_history,
}

visualize_results(spsa_history)

Références

[1] Spall, J. C. (2002). Implementation of the simultaneous perturbation algorithm for stochastic optimization. IEEE Transactions on Aerospace and Electronic Systems, 34(3), 817-823.

[2] Sahin, M. Emre, et al. (2025). Qiskit Machine Learning: an open-source library for quantum machine learning tasks at scale on quantum hardware and classical simulators. arXiv:2505.17756.

Prochaines étapes

Recommandations

Si ce travail t'a intéressé, tu pourrais être intéressé par le matériel suivant :

• Trouver l'énergie de l'état fondamental d'un hamiltonien sparse en suivant le tutoriel SQD
• Suivre le cours Conception d'algorithmes variationnels sur IBM Learning