PEC עם TREX

ביטול שגיאות הסתברותי (PEC) היא טכניקה למיתון רעש קוונטי בחישובי ערכי תצפית של אובזרבבלים. היא מאפשרת לך לחשב הערכה בלתי מוטה של ערכי התצפית. הטכניקה פועלת על ידי הרצת אנסמבל של מעגלים רועשים שמוזרק בהם ההופכי של מודל רעש שנלמד. לאחר מכן מודל הרעש מקורב לערוץ Pauli באמצעות Pauli twirling (שידוע גם בשם randomized compiling). אפשר ללמוד מודלי רעש כאלה ביעילות, ואפשר למתן אותם באמצעות PEC (לפרטים נוספים ראה Probabilistic error cancellation with sparse Pauli-Lindblad models on noisy quantum processors).

Twirled Readout Error eXtinction (TREX) היא טכניקה שממתנת שגיאות מדידה ושימושית להערכת ערכי תצפית של אובזרבבלי Pauli. היא מסתמכת על twirling של המדידה (הזרקה אקראית של שער X לפני המדידה והיפוך התוצאה שנמדדה), שמלכסן את מטריצת המעבר של שגיאת הקריאה והופך אותה לקלה להיפוך. איבר הסקלה מחדש נלמד מערוץ הרעש האלכסוני באמצעות מדידת מעגלי כיול מסובבים שאותחלו במצב האפס.

המדריך הזה מראה איך להשתמש ב-PEC כדי למתן רעש בשערים, ולשלב אותו עם TREX כדי למתן גם רעש מדידה בעת הערכת ערכי תצפית של אובזרבבלים.

זרימת העבודה היא כזו:

1. מגדירים מעגל ואובזרבבל.
2. מחשבים ערכי תצפית מדויקים באמצעות סימולטור וקטור המצב של Qiskit.
3. מחשבים ערכי תצפית רועשים ללא מיתון.
4. לומדים את הרעש של שכבות השזירה באמצעות qiskit-noise-learning.
5. מכינים תוכנית PEC + TREX ומוסיפים אליה מעגל כיול של TREX.
6. משווים את ערכי התצפית אחרי המיתון.

הגדרת מעגל ואובזרבבל

הכן מעגל מראה של התפתחות בזמן בקירוב Trotter של מודל Ising ליניארי עם שדה רוחבי על שמונה Qubits. שקף את המעגל כדי לבטל את היוניטרי של המעגל וליצור מעגל זהות. כך אפשר לחשב בקלות את ערך התצפית האידיאלי של כל אובזרבבל. בנוסף, השיקוף מכפיל את עומק המעגל, ובכך מכפיל את השפעת הרעש.

השתמש בזווית $\theta= \frac{\pi}{2}$ כדי שכל השערים במעגל יישארו בחבורת Clifford, מה שמאפשר סימולציה מהירה יותר באמצעות הסימולטור 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

חישוב ערכי תצפית רועשים ללא מיתון

QuantumProgram מכילה מעגל קוונטי ו-samplex בפורמט שמתאים לפרימיטיב Executor. MitigationTask.prepare בונה QuantumProgram מהמעגל, מהאובזרבבלים ומערכי הפרמטרים, ועוטף את המעגל בקופסאות (boxes) לפי אפשרויות ה-boxing שניתנו. הגדר את אפשרות ה-boxing ‏measure_annotations ל-"change_basis" כדי לסובב את המדידות לבסיסים הנדרשים, והוסף add_tags: "unique_box" כדי שהסימולטור יוכל למפות מודל רעש על הקופסאות של המעגל. התוכנית עוברת Pauli twirling כברירת מחדל, מה שבחומרה הופך רעש קוהרנטי לרעש Pauli שאליו שיטות המיתון מניחות. num_randomizations קובע כמה מופעים מסובבים של כל מעגל נדגמים, ו-shots_per_randomization קובע כמה shots מקבל כל מופע, כך שהמכפלה שלהם היא מגבלת ה-shots לכל מעגל.

אם יש לך גישה ל-QPU, אפשר להריץ את התוכנית באמצעות Executor שנבנה מ-Backend של QiskitRuntimeService. הדוגמה הזו משתמשת ב-AerSimulator כדי להריץ את התוכנית באופן מקומי. אפשרות הסימולטור layer_noise_model שלו מזריקה מודל רעש Pauli-Lindblad אקראי בכל שכבה של שערי שזירה, שנוצר באמצעות random_pauli_lindblad_noise מתוך qiskit-mitigation, עם שגיאה לכל שער בשכבה (EPLG) של 0.003, אופיינית ל-QPU-ים הנוכחיים מסוג Heron.

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()

לימוד הרעש של שכבות השזירה

PEC מתחיל בלימוד מודל הרעש המסובב (twirled) של כל שכבה של שערי שזירה במעגל. אם מניחים מודל רעש מרקובי, מניחים שעומק המעגל לא משפיע על ערוץ הרעש של שכבת שזירה נתונה. לכן צריך ללמוד רק את מודל הרעש של כל שכבה ייחודית של שערי שזירה. אפשר לפצל את המעגל שבו משתמשים במדריך הזה לשתי שכבות ייחודיות. תכנון מעגלים עם קבוצה מינימלית של שכבות שזירה ייחודיות חשוב אם רוצים ללמוד מודל רעש, כי את השכבות הייחודיות צריך ללמוד בנפרד. אם יש במעגל יותר מדי שכבות ייחודיות, לימוד הרעש עלול להפוך לבלתי אפשרי בגלל עלות הדגימה של לימוד שכבות רבות.

qiskit-noise-learning משמש ללימוד הרעש. הלומד חוזר על כל שכבה fragment_depths פעמים במעגלי הלמידה שלו, מסובב (twirl) כל מעגל num_randomizations פעמים, ומתאים את דעיכת הנאמנויות (fidelities) הנמדדות למודל Pauli-Lindblad. הפרמטרים האלה יכולים להשפיע על הדיוק של המודל שנלמד. ככל שמודל הרעש שנלמד מדויק יותר, כך ביטול הרעש יהיה מדויק יותר. הקוד הבא משרטט את מודל הרעש שנלמד, כששתי השכבות ממוזגות למפה אחת. כל Qubit הוא עוגה של שיעורי X, Y ו-Z ממשקל 1 שסוכמו על פני שתי השכבות, וכל צלע מראה את השיעורים ממשקל 2 שנלמדו. למשל, שיעור של 3e-4 אומר שהשגיאה מתרחשת בערך שלוש פעמים מכל 10,000 שימושים בשכבה.

מיתון רעש הקריאה נעשה בעיבוד מקדים באמצעות TREX, בניגוד למודל רעש השערים, שנדרש להכנת התוכנית כדי להזריק את האנטי-רעש שמתאים למודל הרעש שנלמד. לכן אפשר להריץ וללמוד את רעש הקריאה בזמן לימוד רעש השערים, או בזמן הרצת התוכנית הממותנת. לימוד רעש הקריאה יחד עם רעש השערים עשוי להועיל ללימוד מודל רעש עקבי בפני עצמו (ראה למשל Disambiguating Pauli noise in quantum computers). במקרים שבהם רעש השערים ורעש הקריאה מופרדים במודל הרעש, מומלץ ללמוד את רעש הקריאה בזמן הרצת התוכנית הממותנת, כדי להימנע מסחיפה ברעש בין הרגע שבו הוא נלמד לרגע שבו התוכנית הורצה.

במדריך הזה, רעש הקריאה נלמד בזמן הרצת התוכנית הממותנת באמצעות מעגל הכיול של 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()

הכנת תוכנית משולבת של PEC ו-TREX

כאן אנחנו מכינים גם משימת TREX וגם משימת PEC, ומוסיפים את משימת ה-TREX כמשתנה למשימת ה-PEC. זה מבטיח שכל התכונות הנדרשות למיתון TREX מתווספות לפריט הממותן של PEC (כמו twirling של המדידה). הרעש המוזרק נדגם פעם אחת לכל randomization, ולכן השגיאה הסטטיסטית של PEC נקבעת לפי מספר ה-randomizations ולא לפי מספר ה-shots לכל randomization. לכן משתמשים בהרבה randomizations - כל אחד עם מעט shots - ומשאירים את המספר הכולל של ה-shots זהה לזה של התוכנית ללא מיתון.

מעגל הכיול של TREX מתווסף לתוכנית כדי ללמוד את רעש הקריאה בזמן הרצת התוכנית. מעגל הכיול של TREX רדוד מאוד, מה שהופך randomizations נוספים עבורו לזולים יחסית ומאפשר ללמוד את רעש הקריאה בצורה מדויקת יותר.

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)

בדיקת תקורת ה-shots, הרצת התוכנית ועיבוד מקדים

PEC מחזיר הערכה בלתי מוטה של ערך התצפית, במחיר של תקורת shots גדולה. אפשר לקרב את תקורת ה-shots הנדרשת על ידי חישוב גורם הגאמה של המעגל, שמשקף את ההשפעה של מודל הרעש שנלמד על המעגל. כברירת מחדל, מספר ה-randomizations של משימת ה-PEC מוגדל בכמות המשוערת של תקורת ה-shots שנדרשת כדי לקבל הערכה בלתי מוטה (~$\gamma^2$). כדי להימנע מההגדלה האוטומטית הזו, אפשר להגדיר את הפרמטר scale_randomizations_by_gamma ל-False. תקורת ה-shots מוגבלת על ידי הפרמטר max_sampling_overhead.

כדאי לך לבדוק לפני ההרצה את מספר ה-shots ש-PEC דורש כדי להגיע להערכה בלתי מוטה, ולהחליט אם הכמות הזו משתלמת, כי היא גדלה אקספוננציאלית עם שיעורי הרעש ועומק המעגל.

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
)

השוואת ערכי התצפית אחרי המיתון

לבסוף, שרטט את ערך התצפית של כל אובזרבבל. PEC + TREX משחזר את ערכי התצפית המדויקים בתוך גבולות השגיאה שלו, בעוד שהתוצאה ללא מיתון מוטה הרבה מתחתיה.

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()