Premiers pas avec les formules multi-produits (MPF)

Versions des packages

Le code de cette page a été développé avec les prérequis suivants. Nous recommandons d'utiliser ces versions ou des versions plus récentes.

qiskit[all]~=2.3.0
qiskit-addon-utils~=0.3.0
qiskit-addon-mpf~=0.3.0
scipy~=1.16.3

Ce guide montre comment utiliser le package qiskit-addon-mpf, en prenant l'évolution temporelle d'un modèle d'Ising comme exemple. Avec ce package, tu peux construire une formule multi-produits (MPF) qui permet d'atteindre une erreur de Trotter plus faible sur les mesures d'observables. Les outils fournis te permettent de déterminer les poids d'une MPF choisie, qui peuvent ensuite être utilisés pour recombiner les valeurs d'espérance estimées à partir de plusieurs circuits d'évolution temporelle, chacun avec un nombre différent de pas de Trotter.

Commence par considérer le Hamiltonien d'un modèle d'Ising à 10 sites :

$H_{\text{Ising}} = \sum_{i=1}^9 J_{i,(i+1)}Z_iZ_{i+1} + \sum_{i=1}^{10} h_i X_i$

où $J_{i,(i+1)}$ est la force de couplage et $h_i$ est l'intensité du champ magnétique externe. Pour formuler le problème, l'observable à mesurer sera la magnétisation totale du système :

$\langle M \rangle = \sum_{i=1}^{10} \langle Z_i \rangle.$

L'extrait de code ci-dessous prépare le Hamiltonien de la chaîne d'Ising à l'aide du package qiskit-addon-utils, et définit l'observable qui sera mesurée.

# Added by doQumentation — required packages for this notebook
!pip install -q numpy qiskit qiskit-addon-mpf qiskit-addon-utils scipy

from qiskit.transpiler import CouplingMap
from qiskit.synthesis import SuzukiTrotter
from qiskit.quantum_info import SparsePauliOp
from qiskit.primitives import StatevectorEstimator
from qiskit.providers.fake_provider import GenericBackendV2
from qiskit.transpiler import generate_preset_pass_manager
from qiskit_addon_utils.problem_generators import (
    generate_xyz_hamiltonian,
    generate_time_evolution_circuit,
)
from qiskit_addon_mpf.costs import (
    setup_exact_problem,
    setup_sum_of_squares_problem,
)
from qiskit_addon_mpf.static import setup_static_lse

from scipy.linalg import expm
import numpy as np

# Generate some coupling map to use for this example
coupling_map = CouplingMap.from_line(10, bidirectional=False)

# Get a qubit operator describing the Ising field model
hamiltonian = generate_xyz_hamiltonian(
    coupling_map,
    coupling_constants=(0.0, 0.0, 1.0),
    ext_magnetic_field=(0.4, 0.0, 0.0),
)
print(f"Hamiltonian:\n {hamiltonian}\n")

L = coupling_map.size()
observable = SparsePauliOp.from_sparse_list(
    [("Z", [i], 1 / L / 2) for i in range(L)], num_qubits=L
)
print(f"Observable:\n {observable}")

Hamiltonian:
 SparsePauliOp(['IIIIIIIZZI', 'IIIIIZZIII', 'IIIZZIIIII', 'IZZIIIIIII', 'IIIIIIIIZZ', 'IIIIIIZZII', 'IIIIZZIIII', 'IIZZIIIIII', 'ZZIIIIIIII', 'IIIIIIIIIX', 'IIIIIIIIXI', 'IIIIIIIXII', 'IIIIIIXIII', 'IIIIIXIIII', 'IIIIXIIIII', 'IIIXIIIIII', 'IIXIIIIIII', 'IXIIIIIIII', 'XIIIIIIIII'],
              coeffs=[1. +0.j, 1. +0.j, 1. +0.j, 1. +0.j, 1. +0.j, 1. +0.j, 1. +0.j, 1. +0.j,
 1. +0.j, 0.4+0.j, 0.4+0.j, 0.4+0.j, 0.4+0.j, 0.4+0.j, 0.4+0.j, 0.4+0.j,
 0.4+0.j, 0.4+0.j, 0.4+0.j])

Observable:
 SparsePauliOp(['IIIIIIIIIZ', 'IIIIIIIIZI', 'IIIIIIIZII', 'IIIIIIZIII', 'IIIIIZIIII', 'IIIIZIIIII', 'IIIZIIIIII', 'IIZIIIIIII', 'IZIIIIIIII', 'ZIIIIIIIII'],
              coeffs=[0.05+0.j, 0.05+0.j, 0.05+0.j, 0.05+0.j, 0.05+0.j, 0.05+0.j, 0.05+0.j,
 0.05+0.j, 0.05+0.j, 0.05+0.j])

Ensuite, tu prépares la MPF. Le premier choix à faire est de déterminer si les coefficients seront statiques (indépendants du temps) ou dynamiques ; ce tutoriel utilise une MPF statique. Le choix suivant porte sur l'ensemble des valeurs $k_j$ . Cela détermine le nombre de termes dans la MPF, ainsi que le nombre de pas de Trotter que chaque terme utilise pour simuler l'évolution temporelle. Le choix des valeurs $k_j$ est heuristique, tu dois donc trouver ton propre ensemble de « bonnes » valeurs $k_j$ . Des lignes directrices pour trouver un bon ensemble de valeurs sont disponibles à la fin de la page de démarrage.

Une fois les valeurs $k_j$ déterminées, tu peux préparer le système d'équations $Ax=b$ à résoudre. La matrice $A$ est également déterminée par la formule produit à utiliser. Les choix portent ici sur son ordre, fixé à $2$ dans cet exemple, et sur le caractère symétrique ou non de la formule produit, fixé à True dans cet exemple. L'extrait de code ci-dessous sélectionne un temps total d'évolution du système, les valeurs $k_j$ à utiliser, et l'ensemble des équations à résoudre à l'aide de la méthode qiskit_addon_mpf.static.setup_static_lse.

time = 8.0
trotter_steps = (8, 12, 19)

lse = setup_static_lse(trotter_steps, order=2, symmetric=True)
print(lse)

LSE(A=array([[1.00000000e+00, 1.00000000e+00, 1.00000000e+00],
       [1.56250000e-02, 6.94444444e-03, 2.77008310e-03],
       [2.44140625e-04, 4.82253086e-05, 7.67336039e-06]]), b=array([1., 0., 0.]))

Une fois le système d'équations linéaires instancié, il peut être résolu soit de façon exacte, soit via un modèle approché utilisant une somme des carrés (ou la norme de Frobenius pour les coefficients dynamiques ; voir la référence API pour plus d'informations). Le recours à un modèle approché intervient généralement lorsque la norme des coefficients pour l'ensemble de valeurs $k_j$ choisi est jugée trop élevée et qu'il n'est pas possible de choisir un autre ensemble de valeurs $k_j$ . Ce guide présente les deux approches afin de comparer les résultats.

model_exact, coeffs_exact = setup_exact_problem(lse)
model_approx, coeffs_approx = setup_sum_of_squares_problem(
    lse, max_l1_norm=3.0
)
model_exact.solve()
model_approx.solve()
print(f"Exact solution: {coeffs_exact.value}")
print(f"Approximate solution: {coeffs_approx.value}")

Exact solution: [ 0.17239057 -1.19447005  2.02207947]
Approximate solution: [-0.40454257  0.57553173  0.8290123 ]

remarque

L'objet LSE possède également une méthode LSE.solve() qui résout le système d'équations de façon exacte. Si setup_exact_problem() est utilisé dans ce guide, c'est pour illustrer l'interface fournie par les autres méthodes approchées.

Configurer et exécuter les circuits de Trotter

Maintenant que les coefficients $x_j$ ont été obtenus, la dernière étape consiste à générer les circuits d'évolution temporelle pour l'ordre et l'ensemble de pas $k_j$ choisis pour la MPF. Le package qiskit-addon-utils peut accélérer ce processus.

circuits = []
for k in trotter_steps:
    circ = generate_time_evolution_circuit(
        hamiltonian,
        synthesis=SuzukiTrotter(order=2, reps=k),
        time=time,
    )
    circuits.append(circ)

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

Output of the previous code cell

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

Output of the previous code cell

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

Output of the previous code cell

Une fois ces circuits construits, tu peux les transpiler et les exécuter sur un QPU. Dans cet exemple, nous utilisons simplement l'un des simulateurs sans bruit pour illustrer la réduction de l'erreur de Trotter.

backend = GenericBackendV2(num_qubits=10)
transpiler = generate_preset_pass_manager(
    optimization_level=2, backend=backend
)

transpiled_circuits = [transpiler.run(circ) for circ in circuits]

estimator = StatevectorEstimator()
job = estimator.run([(circ, observable) for circ in transpiled_circuits])
result = job.result()

mpf_evs = [res.data.evs for res in result]
print(mpf_evs)

[array(0.23799162), array(0.35754312), array(0.38649906)]

Reconstruire les résultats

Maintenant que les circuits ont été exécutés, la reconstruction des résultats est assez simple. Comme mentionné sur la page de présentation de la MPF, notre observable est reconstruite par la somme pondérée :

$\langle M \rangle = \sum_j x_j \langle M_j \rangle.$

où $x_j$ sont les coefficients trouvés et $\langle M_j \rangle$ est l'estimation de l'observable $\sum_i \langle Z_i \rangle$ pour chaque circuit. On peut ensuite comparer les résultats obtenus avec la valeur exacte à l'aide du package scipy.linalg.

exp_H = expm(-1j * time * hamiltonian.to_matrix())
initial_state = np.zeros(exp_H.shape[0])
initial_state[0] = 1.0

time_evolved_state = exp_H @ initial_state
exact_obs = (
    time_evolved_state.conj() @ observable.to_matrix() @ time_evolved_state
)

# Print out the different observable measurements
print(f"Exact value: {exact_obs.real}")
print(f"PF with 19 steps: {mpf_evs[-1]}")
print(f"MPF using exact solution: {mpf_evs @ coeffs_exact.value}")
print(f"MPF using approximate solution: {mpf_evs @ coeffs_approx.value}")

Exact value: 0.4006024248789992
PF with 19 steps: 0.3864990619977402
MPF using exact solution: 0.3954847855979902
MPF using approximate solution: 0.4299121425348959

On peut constater que la MPF a réduit l'erreur de Trotter par rapport à celle obtenue avec une simple formule produit à $k_j=19$ pas. Cependant, le modèle approché a produit une valeur d'espérance moins précise que le modèle exact. Cela illustre l'importance d'utiliser des critères de convergence stricts pour le modèle approché et de trouver un « bon » ensemble de valeurs $k_j$ .