Extrapolation à bruit nul avec repliement de portes et extrapolation personnalisée
L'extrapolation à bruit nul (ZNE) est une technique d'atténuation du bruit quantique dans les calculs de valeurs d'espérance d'observables. Elle consiste à exécuter un circuit quantique à des échelles de bruit croissantes, à calculer les valeurs d'espérance pour ces échelles de bruit, puis à extrapoler pour estimer une valeur d'espérance à bruit nul. Les étapes d'amplification du bruit et d'extrapolation peuvent être implémentées de différentes manières. Ce guide amplifie le bruit par « repliement numérique de portes », ce qui signifie que les portes à deux qubits sont remplacées par des séquences équivalentes de la porte et de son inverse. Par exemple, remplacer une unitaire par donnerait un facteur d'amplification du bruit de 3. Tu peux aussi implémenter l'amplification par amplification probabiliste d'erreur (PEA), qui nécessite une étape supplémentaire d'apprentissage d'un modèle de bruit et est traitée dans le guide sur l'amplification probabiliste d'erreur. Ce guide extrapole à l'aide de quelques stratégies d'ajustement intégrées (linéaire et exponentielle), et montre aussi comment implémenter ta propre fonction d'ajustement personnalisée.
Garde à l'esprit que, si la ZNE peut être efficace pour atténuer le bruit quantique avec un coût d'échantillonnage inférieur à celui d'autres méthodes, elle ne donne aucune borne sur le biais.
Le flux de travail est le suivant :
- Spécifier un circuit et une observable
- Calculer les valeurs d'espérance exactes avec le simulateur de vecteur d'état de Qiskit
- Calculer les valeurs d'espérance bruitées, non atténuées
- Préparer un programme ZNE avec repliement numérique de portes et échantillonner des circuits d'échelles de bruit variées
- Extrapoler vers la limite à bruit nul avec des fonctions d'ajustement intégrées et personnalisées
- Comparer les valeurs d'espérance atténuées
Spécifier un circuit et une observable
Calcule l'aimantation moyenne par site d'un modèle d'Ising à champ transverse (TFIM) trotterisé sur huit qubits en topologie en anneau. Tu vas balayer plusieurs valeurs de l'intensité du champ transverse et tracer l'aimantation en fonction de l'intensité du champ.
# Added by doQumentation — installs the packages this notebook needs if they are missing
import importlib.util
_needed = {"IPython": "IPython", "matplotlib": "matplotlib", "numpy": "numpy", "qiskit": "qiskit", "qiskit_aer": "qiskit-aer", "qiskit_ibm_runtime": "qiskit-ibm-runtime", "qiskit_mitigation": "qiskit-mitigation", "scipy": "scipy"}
_missing = [pip for module, pip in _needed.items()
if importlib.util.find_spec(module) is None]
# One at a time, so a package that fails to install does not block the others
for _pip in _missing:
%pip install -q {_pip}
if not _missing:
print("\u2713 All required packages are installed")
import warnings
import numpy as np
from qiskit.circuit import Parameter, QuantumCircuit
from qiskit.quantum_info import SparsePauliOp
warnings.filterwarnings("ignore")
def trotter_circuit(num_qubits, theta, layers, num_steps):
"""Trotterized transverse-field Ising evolution: RX on every site, then a CZ layer per coupling layer."""
circuit = QuantumCircuit(num_qubits)
for _ in range(num_steps):
circuit.rx(theta, range(num_qubits))
for layer in layers:
circuit.barrier()
for edge in layer:
circuit.cz(*edge)
circuit.sdg(edge)
circuit.barrier()
return circuit
num_qubits = 8
trotter_steps = 15
ring_coupling = [(i, (i + 1) % num_qubits) for i in range(num_qubits)]
layers = [
ring_coupling[0::2],
ring_coupling[1::2],
] # two disjoint layers of CZ gates
theta = Parameter("theta")
circuit = trotter_circuit(num_qubits, theta, layers, trotter_steps)
# Parameter values: theta in [0, 0.25 pi] - values for the transverse-field sweep
data_points = 6
parameter_values = np.linspace(0, 0.25 * np.pi, data_points).reshape(
(data_points, 1)
)
# Define observables - Z on each site
observables = [
SparsePauliOp("I" * (num_qubits - q - 1) + "Z" + "I" * q)
for q in range(num_qubits)
]
print("Visualization of 1 Trotter step")
trotter_circuit(num_qubits, theta, layers, 1).draw("mpl", scale=0.7)
Visualization of 1 Trotter step
Calculer les valeurs d'espérance exactes avec le simulateur de vecteur d'état de Qiskit
Utilise la classe Statevector de Qiskit pour calculer la valeur d'espérance exacte de sur chaque site pour chaque valeur de .
import matplotlib.pyplot as plt
from qiskit.quantum_info import Statevector
exact_evs = np.array(
[
[
Statevector(circuit.assign_parameters([value]))
.expectation_value(obs)
.real
for value in parameter_values.ravel()
]
for obs in observables
]
) # shape (observables, parameter values)
x_vals = parameter_values.ravel() / np.pi
fig, ax = plt.subplots(figsize=(7, 4))
ax.plot(x_vals, exact_evs.mean(axis=0), "ko-", label="Exact")
ax.set_xlabel(r"$\theta/\pi$")
ax.set_ylabel(r"$\langle Z \rangle_{\text{avg}}$")
ax.legend()
plt.tight_layout()
plt.show()
Calculer les valeurs d'espérance bruitées, non atténuées
Un QuantumProgram contient un circuit quantique et un samplex dans un format adapté à la primitive Executor. MitigationTask.prepare construit un QuantumProgram à partir du circuit, des observables et des valeurs de paramètres, et encadre le circuit dans des boîtes avec les options de boxing données. Tu définis l'option de boxing measure_annotations sur "change_basis", qui fait pivoter les mesures vers les bases requises, et tu ajoutes add_tags: "unique_box" pour que le simulateur puisse appliquer un modèle de bruit aux boîtes du circuit. Le programme est par défaut soumis à un twirling de Pauli, qui sur du matériel transforme le bruit cohérent en bruit de Pauli, celui que supposent les méthodes d'atténuation. num_randomizations définit combien d'instances twirlées de chaque circuit sont échantillonnées et shots_per_randomization définit combien de shots reçoit chaque instance ; leur produit est donc le budget de shots par circuit.
Un utilisateur ayant accès à un QPU exécuterait le programme avec un Executor construit à partir d'un backend du QiskitRuntimeService. Ici, tu vas passer un AerSimulator à l'Executor, qui exécute le programme localement, et utiliser son option de simulateur layer_noise_model pour injecter un modèle de bruit de Pauli-Lindblad aléatoire sur chaque couche de portes d'intrication, généré avec random_pauli_lindblad_noise de qiskit-mitigation avec une erreur par porte en couche (EPLG) de 0.003, typique des QPU Heron actuels.
from qiskit.circuit import CircuitInstruction
from qiskit.quantum_info import PauliLindbladMap, pauli_basis
from qiskit_aer import AerSimulator
from qiskit_ibm_runtime import Executor
from qiskit_mitigation import MitigationTask, find_combined_unique_layers
def generate_random_pauli_lindblad_noise(
layers,
num_qubits: int,
eplg: float = 3e-3,
seed: int | None = None,
) -> list[tuple[CircuitInstruction, PauliLindbladMap]]:
"""Generate random Pauli-Lindblad noise over the full Pauli basis."""
layer_noise_model = []
for layer in layers:
layer_instructions = layer.operation.params[0]
contain_measurement = any(
inst.name.startswith("measure") for inst in layer_instructions
)
if contain_measurement: # do not add measurement noise
layer_noise = PauliLindbladMap.identity(num_qubits)
else:
edges = [
(inst.qubits[0]._index, inst.qubits[1]._index)
for inst in layer_instructions
] # extract the edges of the layer
basis_paulis = [p for p in pauli_basis(2) if np.sum(p.x + p.z)]
rng = np.random.default_rng(seed=seed)
num_generators = len(basis_paulis)
rates = rng.random(len(basis_paulis)) * eplg / num_generators
layer_noise = PauliLindbladMap.from_sparse_list(
[
(pauli.to_label(), edge, rate)
for pauli, rate in zip(basis_paulis, rates, strict=True)
for edge in edges
],
num_qubits=num_qubits,
)
layer_noise_model.append((layer, layer_noise))
return layer_noise_model
num_randomizations = 50
shots_per_randomization = 64
boxing_options = {
"measure_annotations": "change_basis",
"add_tags": "unique_box",
}
# find the unique layers to generate the layered noise for
box_layers = find_combined_unique_layers(
[circuit], custom_boxing_options=boxing_options
)
layer_noise_model = generate_random_pauli_lindblad_noise(
box_layers, num_qubits, seed=42
)
# a local noisy simulator in place of a QPU backend
noisy_executor = Executor(
AerSimulator(),
options={"simulator": {"layer_noise_model": layer_noise_model}},
)
task = MitigationTask()
program = task.prepare(
circuit=circuit,
observables=observables,
parameters=parameter_values,
num_randomizations=num_randomizations,
shots_per_randomization=shots_per_randomization,
custom_boxing_options=boxing_options,
)
program._semantic_role = "simulation" # flag the program for the simulator
unmitigated_res = task.postprocess(noisy_executor.run(program).result())
fig, ax = plt.subplots(figsize=(7, 4))
ax.plot(x_vals, exact_evs.mean(axis=0), "k-", label="Exact")
ax.errorbar(
x_vals,
np.mean(unmitigated_res.data.evs, axis=0),
np.mean(unmitigated_res.data.stds, axis=0),
fmt="o-",
label="Unmitigated",
)
ax.set_xlabel(r"$\theta/\pi$")
ax.set_ylabel(r"$\langle Z \rangle_{\text{avg}}$")
ax.legend()
plt.tight_layout()
plt.show()
Préparer un programme ZNE avec repliement numérique de portes et échantillonner des circuits d'échelles de bruit variées
Prépare un programme ZNE contenant des circuits d'échelles de bruit croissantes. Chaque facteur de bruit représente le multiple de portes d'intrication que contient le circuit (un facteur de bruit de 3.0 contient trois fois plus de portes d'intrication que le circuit d'origine). Le circuit de chaque facteur de bruit est twirlé avec le même nombre de randomisations et de shots que le programme non atténué.
from qiskit_mitigation import GateFolding
task_gf = GateFolding()
program_zne = task_gf.prepare(
circuit=circuit,
observables=observables,
parameters=parameter_values,
num_randomizations=num_randomizations,
shots_per_randomization=shots_per_randomization,
custom_boxing_options=boxing_options,
noise_factors=[1, 3, 5], # ZNE option
)
Dessine la première étape de Trotter du modèle twirlé avec le facteur de bruit 1 et avec le facteur de bruit 3 pour visualiser l'amplification. Pour le facteur de bruit 3.0, chaque couche du circuit est remplacée par la couche, son inverse, puis de nouveau la couche.
import IPython.display
def first_trotter_step(template, noise_factor):
"""The template up to the closing barrier of its first Trotter step (two layers, each folded ``noise_factor`` times)."""
toy = template.copy_empty_like()
for instruction in template.data:
toy.append(instruction)
if (
instruction.operation.name == "barrier"
and instruction.operation.label.startswith(
f"R{2 * noise_factor - 1}@"
)
):
break
return toy
IPython.display.display(
first_trotter_step(program_zne.items[0].circuit, 1).draw(
"mpl", scale=0.7, plot_barriers=False
)
)
first_trotter_step(program_zne.items[1].circuit, 3).draw(
"mpl", fold=-1, plot_barriers=False
)
Enfin, exécute le programme sur le simulateur bruité. Ses résultats alimentent l'extrapolation de l'étape suivante.
program_zne._semantic_role = (
"simulation" # flag the program for the simulator
)
job_zne_results = noisy_executor.run(program_zne).result()
Extrapoler vers la limite à bruit nul avec des fonctions d'ajustement intégrées et personnalisées
Maintenant que tu as des données pour chaque échelle de bruit, tu peux ajuster une fonction aux données et extrapoler jusqu'à la valeur à bruit nul pour estimer la valeur d'espérance de l'observable sans bruit. Il existe de nombreuses façons d'ajuster des données, et le bon choix de fonction dépend de la nature des données échantillonnées. GateFolding.postprocess fournit plusieurs fonctions d'extrapolation intégrées, dont linear, exponential et double_exponential. La méthode postprocess prend une liste de méthodes d'extrapolation, les évalue dans l'ordre et renvoie la première valeur physiquement valide et de variance finie. Cela signifie que l'ordre dans lequel tu spécifies les extrapolateurs a un grand impact sur celui qui est choisi. Soumettre fallback comme extrapolateur signifie que si aucune réponse valide n'est trouvée, la valeur d'espérance du plus petit facteur de bruit sera choisie et renvoyée. Si fallback n'est pas soumis et qu'aucun extrapolateur ne donne de valeur physiquement valide, le résultat avec le plus petit std est renvoyé. Sous bruit de Pauli, une observable décroît exponentiellement avec le facteur de bruit ; avec 30 couches (15 étapes de Trotter, 2 couches par étape) de portes d'intrication, l'extrapolateur exponentiel est donc le choix naturel ; tu peux aussi extrapoler linéairement, ce qui est un choix courant lorsque la décroissance est faible.
Pour utiliser une fonction d'ajustement personnalisée, passe l'extrapolateur "custom" avec un tuple custom_fit contenant la fonction d'ajustement, une liste avec une estimation initiale pour chacun de ses paramètres, et des bornes. Ceux-ci sont transmis à scipy.optimize.curve_fit. Le post-traitement peut être répété avec différents extrapolateurs sans ré-exécuter le programme. Ici, tu utilises une fonction d'ajustement personnalisée en racine carrée. Pour visualiser l'extrapolation, trace l'aimantation moyenne par site calculée à chaque facteur de bruit pour , les trois ajustements, et la moyenne sur les sites des valeurs d'espérance extrapolées au point de bruit nul. Les trois points de données semblent presque linéaires, et pourtant l'extrapolation linéaire tombe nettement en dessous de la valeur exacte. Comme prévu, l'ajustement exponentiel est celui qui estime le mieux la valeur d'espérance sans bruit.
from scipy.optimize import curve_fit
zne_mitigated_exp = task_gf.postprocess(
job_zne_results, extrapolator=["exponential", "linear"]
)
zne_mitigated_lin = task_gf.postprocess(
job_zne_results, extrapolator=["linear"]
)
def sqrt_fit(x, a, b):
return a * np.sqrt(x) + b
custom_fit = (sqrt_fit, [1.0, 1.0], (-np.inf, np.inf))
zne_mitigated_custom = task_gf.postprocess(
job_zne_results,
extrapolator=["custom", "exponential", "linear"],
custom_fit=custom_fit,
)
# The result holds the expectation value at each noise factor and the extrapolated value, but
# not the fit parameters, so we refit the same models here to draw the curves.
noise_factors = np.array([1, 3, 5])
angle = 3 # theta = 0.15 pi, the angle with the largest unmitigated error
amplifications = zne_mitigated_exp.data.evs_noise_factors[:, angle, :].mean(
axis=0
)
x_fit = np.linspace(0, 5, 51)
fits = {
"exponential": (
lambda x, a, b: a * np.exp(b * x),
zne_mitigated_exp,
"black",
),
"linear": (lambda x, a, b: a * x + b, zne_mitigated_lin, "blue"),
"square-root": (sqrt_fit, zne_mitigated_custom, "green"),
}
plt.plot(
noise_factors, amplifications, "o", color="red", label="amplifications"
)
for name, (fit, res, color) in fits.items():
params, _ = curve_fit(fit, noise_factors, amplifications, p0=[1.0, -0.1])
plt.plot(
x_fit, fit(x_fit, *params), "--", color=color, label=f"{name} fit"
)
plt.plot(
[0],
res.data.evs[:, angle].mean(),
"*",
color=color,
markersize=15,
label=f"{name} extrapolation",
)
plt.plot(
[0],
exact_evs[:, angle].mean(),
"_",
color="red",
markersize=15,
markeredgewidth=3,
label="exact",
)
plt.xlabel("noise factor")
plt.ylabel(r"$\langle Z \rangle_{\text{avg}}$")
plt.legend()
plt.show()
Comparer les valeurs d'espérance atténuées
Trace l'aimantation moyenne par site de chaque méthode sur le balayage des angles de rotation d'Ising. L'extrapolation exponentielle retrouve l'aimantation exacte dans ses barres d'erreur, tandis que l'extrapolation linéaire supprime une grande partie du biais mais sous-estime systématiquement, car la vraie décroissance est exponentielle. L'ajustement en racine carrée surestime fortement.
results = {
"Unmitigated": unmitigated_res,
"ZNE exponential": zne_mitigated_exp,
"ZNE linear": zne_mitigated_lin,
"ZNE square-root": zne_mitigated_custom,
}
fig, ax = plt.subplots(figsize=(7, 4))
ax.plot(x_vals, exact_evs.mean(axis=0), "k-", label="Exact")
for label, res in results.items():
ax.errorbar(
x_vals,
np.mean(res.data.evs, axis=0),
np.mean(res.data.stds, axis=0),
fmt="o-",
label=label,
)
ax.set_xlabel(r"$\theta/\pi$")
ax.set_ylabel(r"$\langle Z \rangle_{\text{avg}}$")
ax.legend()
ax.set_title(
f"Error-Mitigated Magnetization [{num_qubits}-qubit, {trotter_steps} Trotter steps]"
)
plt.tight_layout()
plt.show()
