PEC avec TREX

La Probabilistic error cancellation (PEC) est une technique d'atténuation du bruit quantique dans le calcul de valeurs moyennes d'observables. Elle te permet de calculer une estimation non biaisée des valeurs moyennes. Cette technique consiste à exécuter un ensemble de circuits bruités dans lesquels on injecte l'inverse d'un modèle de bruit appris. Le modèle de bruit est ensuite approximé par un canal de Pauli grâce au Pauli twirling (également connu sous le nom de randomized compiling). De tels modèles de bruit peuvent être appris efficacement et atténués avec la PEC (voir Probabilistic error cancellation with sparse Pauli-Lindblad models on noisy quantum processors pour plus de détails).

Twirled Readout Error eXtinction (TREX) est une technique qui atténue les erreurs de mesure et qui est utile pour estimer les valeurs moyennes d'observables de Pauli. Elle repose sur le twirling de mesure (injection aléatoire d'une Gate X avant la mesure et inversion du résultat mesuré), qui diagonalise la matrice de transfert des erreurs de lecture et la rend facile à inverser. Un terme de remise à l'échelle est appris à partir du canal de bruit diagonal en mesurant des circuits de calibration twirlés initialisés dans l'état zéro.

Ce guide montre comment utiliser la PEC pour atténuer le bruit des Gates et la combiner avec TREX pour atténuer aussi le bruit de mesure lors de l'estimation des valeurs moyennes d'observables.

Le flux de travail est le suivant :

1. Spécifier un circuit et une observable.
2. Calculer les valeurs moyennes exactes avec le simulateur de vecteur d'état de Qiskit.
3. Calculer les valeurs moyennes bruitées et non atténuées.
4. Apprendre le bruit des couches d'intrication avec qiskit-noise-learning.
5. Préparer un programme PEC + TREX et y ajouter un circuit de calibration TREX.
6. Comparer les valeurs moyennes atténuées.

Spécifier un circuit et une observable

Prépare un circuit miroir d'une évolution temporelle trotterisée d'un modèle d'Ising linéaire en champ transverse sur huit qubits. Applique le miroir au circuit pour annuler l'unitaire du circuit et former un circuit identité. Cela te permet de calculer facilement la valeur moyenne idéale de chaque observable. De plus, le miroir double la profondeur du circuit, ce qui double à son tour l'effet du bruit.

Utilise l'angle $\theta= \frac{\pi}{2}$ pour que toutes les Gates du circuit appartiennent au groupe de Clifford, ce qui permet une simulation plus rapide avec le simulateur stabilizer.

# Added by doQumentation — installs the packages this notebook needs if they are missing
import importlib.util

_needed = {"matplotlib": "matplotlib", "numpy": "numpy", "qiskit": "qiskit", "qiskit_aer": "qiskit-aer", "qiskit_ibm_runtime": "qiskit-ibm-runtime", "qiskit_mitigation": "qiskit-mitigation", "qiskit_noise_learning": "qiskit-noise-learning", "samplomatic": "samplomatic"}
_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 QuantumCircuit
from qiskit.quantum_info import SparsePauliOp

warnings.filterwarnings("ignore")

def construct_ising_circuit(
    num_qubits: int,
    num_trotter_steps: int,
    rx_angle: float,
    barrier: bool = True,
) -> QuantumCircuit:
    circuit = QuantumCircuit(num_qubits)

    for _step in range(num_trotter_steps):
        circuit.rx(rx_angle, range(num_qubits))
        if barrier:
            circuit.barrier()
        for first_qubit in (1, 2):
            for idx in range(first_qubit, num_qubits, 2):
                # equivalent to Rzz(-pi/2):
                circuit.sdg([idx - 1, idx])
                circuit.cz(idx - 1, idx)
        if barrier:
            circuit.barrier()

    return circuit

def mirror_circuit(circuit: QuantumCircuit) -> QuantumCircuit:
    mirror_circ = circuit.copy_empty_like()
    mirror_circ.compose(circuit, inplace=True)
    mirror_circ.barrier()
    mirror_circ.compose(circuit.inverse(), inplace=True)
    mirror_circ.measure_active()
    return mirror_circ

# define observables
num_qubits = 8
observables = [
    SparsePauliOp("I" * ((num_qubits - 1) - q) + "Z" + "I" * q)
    for q in range(num_qubits)
]

# Instantiate circuit
num_trotter_steps = 5
rx_angle = np.pi / 2

circuit = construct_ising_circuit(
    num_qubits, num_trotter_steps, rx_angle, barrier=False
)
mirrored_circuit = mirror_circuit(circuit)
print("Visualization of 1 mirrored Trotter step")
mirror_circuit(
    construct_ising_circuit(num_qubits, 1, rx_angle, barrier=False)
).draw("mpl", scale=0.7)

Visualization of 1 mirrored Trotter step

Calculer les valeurs moyennes bruitées et 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 place le circuit dans des boîtes selon les options de boxing données. Définis l'option de boxing measure_annotations sur "change_basis" pour faire pivoter les mesures vers les bases requises, et ajoute add_tags: "unique_box" afin que le simulateur puisse associer un modèle de bruit aux boîtes du circuit. Le programme est soumis par défaut au Pauli twirling, ce qui, sur le matériel, transforme le bruit cohérent en bruit de Pauli supposé par 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 la limite de shots par circuit.

Si tu as accès à un QPU, tu peux exécuter le programme avec un Executor construit à partir d'un Backend du QiskitRuntimeService. Cet exemple utilise un AerSimulator pour exécuter le programme localement. Son option de simulateur layer_noise_model injecte un modèle de bruit de Pauli-Lindblad aléatoire sur chaque couche de Gates d'intrication, généré avec random_pauli_lindblad_noise de qiskit-mitigation avec une erreur par Gate en couche (EPLG) de 0.003, typique des QPU Heron actuels.

import matplotlib.pyplot as plt
from qiskit.circuit import CircuitInstruction
from qiskit.quantum_info import PauliLindbladMap, pauli_basis
from qiskit.transpiler import PassManager
from qiskit_aer import AerSimulator
from qiskit_ibm_runtime import Executor
from qiskit_mitigation import MitigationTask, find_combined_unique_layers
from qiskit_noise_learning.aer_executor.insert_noise_pass import (
    InsertNoisePass,
)
from samplomatic import Tag
from samplomatic.utils import get_annotation

def generate_random_pauli_lindblad_noise(
    layers,
    num_qubits: int,
    eplg: float = 3e-3,
    readout_noise_scale: float = 3e-2,
    seed: int | None = None,
) -> tuple[
    list[tuple[CircuitInstruction, PauliLindbladMap]], 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 not contain_measurement:
            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))
    # generate a readout noise model
    rng = np.random.default_rng(seed=seed)
    rates = rng.random(num_qubits) * readout_noise_scale
    readout_noise_model = PauliLindbladMap.from_sparse_list(
        [("X", [q], rate) for q, rate in enumerate(rates)],
        num_qubits=num_qubits,
    )
    return layer_noise_model, readout_noise_model

def add_readout_noise(program, readout_noise_model):
    """Injects readout noise to the circuits in a quantum program."""
    readout_noise_map = {}
    for item in program.items:
        readout_tag = item.circuit[-1].label[
            -4:
        ]  # extract the readout layer tag from the circuit
        readout_noise_map[readout_tag] = readout_noise_model
        item.circuit = PassManager(
            [InsertNoisePass(noise_dict=readout_noise_map, noise_after=False)]
        ).run(item.circuit)

num_randomizations = 300
shots_per_randomization = 64

boxing_options = {"add_tags": "unique_box"}

# find the unique layers to generate the layered noise for
box_layers = find_combined_unique_layers(
    [mirrored_circuit], custom_boxing_options=boxing_options
)
layer_noise_model, readout_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(method="stabilizer"),
    options={"simulator": {"layer_noise_model": layer_noise_model}},
)

task = MitigationTask()
program = task.prepare(
    circuit=mirrored_circuit,
    observables=observables,
    parameters=None,
    num_randomizations=num_randomizations,
    shots_per_randomization=shots_per_randomization,
    custom_boxing_options=boxing_options,
)
# For the simulation - we inject readout noise to the circuits of the program
add_readout_noise(program, readout_noise_model)

program._semantic_role = "simulation"  # flag the program for the simulator
unmitigated_res = task.postprocess(noisy_executor.run(program).result())

x_vals = np.array(range(num_qubits))
fig, ax = plt.subplots(figsize=(7, 4))
ax.plot(x_vals, [1] * len(x_vals), "k-", alpha=0.5, label="Exact")
ax.errorbar(
    x_vals,
    unmitigated_res.data.evs,
    unmitigated_res.data.stds,
    fmt="o",
    label="Unmitigated",
)

x_tick_labels = [r"$Z_{ {i} }$" for i in range(num_qubits)]
ax.set_xlabel(r"$\theta/\pi$")
ax.set_ylabel("Expectation value")
ax.legend()
ax.set_xticks(x_vals, x_tick_labels)
plt.tight_layout()
plt.show()

Apprendre le bruit des couches d'intrication

La PEC commence par apprendre le modèle de bruit twirlé de chaque couche de Gates d'intrication du circuit. Si tu supposes un modèle de bruit markovien, la profondeur du circuit est supposée ne pas affecter le canal de bruit d'une couche d'intrication donnée. Par conséquent, seul le modèle de bruit de chaque couche unique de Gates d'intrication doit être appris. Le circuit utilisé dans ce guide peut être divisé en deux couches uniques. Concevoir des circuits avec un ensemble minimal de couches d'intrication uniques est important si tu veux apprendre un modèle de bruit, car les couches uniques doivent être apprises indépendamment. Si ton circuit contient trop de couches uniques, l'apprentissage du bruit peut devenir irréalisable en raison du coût d'échantillonnage de l'apprentissage de nombreuses couches.

qiskit-noise-learning est utilisé pour apprendre le bruit. Le module d'apprentissage répète chaque couche fragment_depths fois dans ses circuits d'apprentissage, applique le twirling à chaque circuit num_randomizations fois et ajuste la décroissance des fidélités mesurées à un modèle de Pauli-Lindblad. Ces paramètres peuvent affecter la précision du modèle appris. Plus le modèle de bruit appris est précis, plus l'annulation du bruit sera précise. Le code suivant dessine le modèle de bruit appris avec les deux couches fusionnées en une seule carte. Chaque qubit est un diagramme circulaire de ses taux X, Y et Z de poids 1 additionnés sur les deux couches, et chaque arête indique les taux de poids 2 appris. Par exemple, un taux de 3e-4 signifie que l'erreur se produit environ trois fois sur 10 000 utilisations de la couche.

L'atténuation du bruit de lecture est effectuée en post-traitement avec TREX, contrairement au modèle de bruit des Gates, qui est requis pour la préparation du programme afin d'injecter l'anti-bruit correspondant au modèle de bruit appris. Par conséquent, le bruit de lecture peut être exécuté et appris soit lors de l'apprentissage du bruit des Gates, soit lors de l'exécution du programme atténué. Apprendre le bruit de lecture en même temps que le bruit des Gates peut être bénéfique pour apprendre un modèle de bruit auto-cohérent (voir par exemple Disambiguating Pauli noise in quantum computers). Dans les cas où le bruit des Gates et le bruit de lecture sont séparés dans le modèle de bruit, il est recommandé d'apprendre le bruit de lecture lors de l'exécution du programme atténué afin d'éviter une dérive du bruit entre le moment où il a été appris et celui où il est exécuté.

Dans ce guide, le bruit de lecture est appris pendant l'exécution du programme atténué grâce au circuit de calibration TREX.

from qiskit.providers.fake_provider import GenericBackendV2
from qiskit.quantum_info import PauliList
from qiskit_ibm_runtime.results import LayerError, PauliLindbladError
from qiskit_ibm_runtime.visualization.embeddings import Embedding
from qiskit_mitigation import PEC
from qiskit_noise_learning.aer_executor import AerExecutor
from qiskit_noise_learning.noise_learner import LearningOptions, NoiseLearner
from samplomatic import InjectNoise

line_coupling = [(i, i + 1) for i in range(num_qubits - 1)]
backend = GenericBackendV2(
    num_qubits=num_qubits,
    basis_gates=["id", "rz", "sx", "x", "cz"],
    coupling_map=line_coupling,
)

task_pec = PEC()
# first learn only the layers PEC will inject noise into, i.e. not the measurement layer
pec_layers = [
    layer
    for layer in task_pec.find_unique_layers(mirrored_circuit, boxing_options)
    if get_annotation(layer.operation, InjectNoise)
]
layers = task_pec.find_unique_layers(mirrored_circuit, boxing_options)

learning_options = LearningOptions(
    fragment_depths=[2, 4, 8, 16],
    num_randomizations=250,
    shots_per_randomizations=16,
)

# qiskit-noise-learning is not compatible with the Executor in local mode,
# and requires the simulator present in the package for running simulation
noise_dict = {
    get_annotation(layer.operation, Tag).ref: noise_layer
    for layer, noise_layer in layer_noise_model
}
noise_learner_simulator = AerExecutor(
    qasm_simulator=AerSimulator(method="stabilizer"), noise_dict=noise_dict
)
learner = NoiseLearner(
    backend, options=learning_options, executor=noise_learner_simulator
)
learned_noise = (
    learner.run(pec_layers).result().to_dict()
)  # keyed by InjectNoise reference

# Draw the learned model on the line, laid out as a line with the qubits in order.
# Composing the two layers' maps sums their rates, giving one map over all the edges.
embedding = Embedding(
    [(0, 0), (0, 2), (0, 4), (2, 4), (4, 4), (4, 2), (4, 0), (2, 0)],
    line_coupling,
)
first, second = learned_noise.values()
learned_model = first.compose(second).simplify()
error = PauliLindbladError(
    PauliList(
        learned_model.get_qubit_sparse_pauli_list_copy().to_pauli_list()
    ),
    learned_model.rates,
)
fig = LayerError(circuit, list(range(num_qubits)), error).draw_map(
    embedding, radius=0.4, height=460, width=620
)
fig.update_traces(
    marker={
        "cmin": 0,
        "cmax": max(learned_model.rates),
        "colorbar": {"tickformat": ".1e", "title": "rate"},
    },
    selector={"mode": "markers"},
)
fig.show()

Préparer un programme combiné PEC et TREX

Ici, nous préparons à la fois une tâche TREX et une tâche PEC, et nous ajoutons la tâche TREX comme variable à la tâche PEC. Cela garantit que toutes les propriétés requises pour l'atténuation TREX sont ajoutées à l'élément atténué par PEC (comme le twirling de mesure). Le bruit injecté est échantillonné une fois par randomisation ; l'erreur statistique de la PEC est donc contrôlée par le nombre de randomisations plutôt que par les shots par randomisation. Ainsi, de nombreuses randomisations sont utilisées - chacune avec peu de shots - en gardant le nombre total de shots identique à celui du programme non atténué.

Le circuit de calibration TREX est ajouté au programme afin d'apprendre le bruit de lecture pendant l'exécution du programme. Le circuit de calibration TREX est très peu profond, ce qui rend des randomisations supplémentaires relativement peu coûteuses et permet d'apprendre le bruit de lecture plus précisément.

from qiskit_mitigation import TREX

trex = TREX()

program_pec = task_pec.prepare(
    circuit=mirrored_circuit,
    observables=observables,
    parameters=None,
    num_randomizations=1200,
    shots_per_randomization=16,
    custom_boxing_options=boxing_options,
    trex=trex,
    noise_maps=learned_noise,  # PEC option
    noise_gain=0,  # PEC option
    max_sampling_overhead=100,  # PEC option
)
program_pec = trex.prepare(
    num_randomizations=1500,
    quantum_program=program_pec,
    custom_boxing_options=boxing_options,
)

# For the simulation - we inject readout noise to the circuits of the program
add_readout_noise(program_pec, readout_noise_model)

Examiner le surcoût en shots, exécuter le programme et post-traiter

La PEC renvoie une estimation non biaisée de la valeur moyenne, au prix d'un grand surcoût en shots. Le surcoût en shots requis peut être approximé en calculant le facteur gamma du circuit, qui reflète l'effet du modèle de bruit appris sur le circuit. Par défaut, le nombre de randomisations de la tâche PEC est mis à l'échelle par la quantité approximative de surcoût en shots nécessaire pour obtenir une estimation non biaisée (~$\gamma^2$). Pour éviter cette mise à l'échelle automatique, tu peux définir le paramètre scale_randomizations_by_gamma sur False. Le surcoût en shots est plafonné par le paramètre max_sampling_overhead.

Tu devrais examiner le nombre de shots requis par la PEC pour atteindre une estimation non biaisée avant l'exécution afin de décider si cette quantité est abordable, car il croît exponentiellement avec les taux de bruit et la profondeur du circuit.

print(f"The circuit calculated gamma is: {task_pec.gamma}")
print(f"The circuit shot overhead is ~{task_pec.gamma**2}")

The circuit calculated gamma is: 1.2774846896089953
The circuit shot overhead is ~1.631967132185391

program_pec._semantic_role = (
    "simulation"  # flag the program for the simulator
)
job_pec_results = noisy_executor.run(program_pec).result()

# post-processing
pec_res = task_pec.postprocess(job_pec_results)
# In order to post process the results without TREX, we call `postprocess` again, but instead of the default options we
# pass only the PEC item result and the `measure_noise_data` as `None`
pec_without_trex_res = task_pec.postprocess(
    job_pec_results[task_pec._program_item_index], measure_noise_data=None
)

Comparer les valeurs moyennes atténuées

Enfin, trace la valeur moyenne de chaque observable. La PEC + TREX retrouve les valeurs moyennes exactes dans ses barres d'erreur, tandis que le résultat non atténué est biaisé bien en dessous.

x = np.array(range(num_qubits))
x_tick_labels = [r"$Z_{ {i} }$" for i in range(num_qubits)]

fig, ax = plt.subplots(figsize=(8, 5))
ax.errorbar(
    x,
    unmitigated_res.data.evs,
    unmitigated_res.data.stds,
    fmt="o",
    label="Twirl only",
)
ax.errorbar(
    x,
    pec_without_trex_res.data.evs,
    pec_without_trex_res.data.stds,
    fmt="o",
    label="PEC only",
)
ax.errorbar(
    x,
    pec_res.data.evs,
    pec_res.data.stds,
    fmt="o",
    label="PEC + TREX",
)

ax.plot(x, [1] * len(x), "--", color="black", alpha=0.5, label="Exact")

tick_labels = [f"$Z_{ {i} }$" for i in range(num_qubits)]
ax.set_xticks(x, tick_labels)

ax.set_ylabel("Expectation value")
ax.legend()
ax.set_title("PEC+TREX Mitigated Expectation Values")
plt.tight_layout()
plt.show()