זיהוי שגיאות בעומס נמוך באמצעות קודי מרחב-זמן

במדריך הזה נממש תהליך עבודה לדגימה ממעגל Clifford ונשתמש בחבילה qiskit-paulice כדי לשפר את הנאמנות (fidelity) של ההתפלגות הנמדדת על ידי הוספת בדיקות Pauli במרחב-זמן לאורך המעגל. בדיקות Pauli במרחב-זמן דומות ל-בדיקות Pauli קוהרנטיות (CPC) בכך שמעגל Clifford "מטען" (payload) שזור עם כמה qubits עזר כדי "לבדוק" אינווריאנטים מסוימים במעגל. מדידת ה-ancillas מפיקה סינדרום שמשמש לקבוע אם זוהתה שגיאה לוגית במהלך הרצת המעגל. השארת דגימות שבהן לא זוהתה שגיאה בלבד (postselection) יכולה לשפר את הנאמנות של ההתפלגויות הנדגמות, במחיר של ירידה בשיעורי ה-postselection. ההבדל המרכזי בין CPC לבדיקות מרחב-זמן הוא שמדידת סינדרומים של CPC דורשת מדידת אופרטורים בעלי משקל גבוה ומקומיים בזמן. בטופולוגיות qubits עם קישוריות מוגבלת, כמו heavy hex, מימוש הבדיקות האלה דורש הרבה שערי SWAP, ולעיתים קרובות המעגל נעשה עמוק מכדי להריץ אותו בפועל. מימוש בדיקות Pauli כקודי מרחב-זמן מאפשר לפזר אותן לאורך מעגל המטען, ולקבל קידוד qubits יעיל מבחינת חומרה שנשאר אפקטיבי בזיהוי שגיאות לוגיות במהלך הרצת המעגל. התפקיד העיקרי של החבילה qiskit-paulice הוא למצוא ולהכניס בדיקות Pauli תקפות במרחב-זמן למקומות במעגל שממקסמים את יעילות זיהוי השגיאות וממזערים את עומס ה-qubits.

במדריך הזה נממש תהליך עבודה לשיפור הנאמנות של התפלגות שנדגמה ממעגל קוונטי באמצעות בדיקות מרחב-זמן. מכיוון שנשתמש במעגל Clifford, אפשר להעריך את הנאמנות בקלות מממוצע ערכי התוחלת של היציבות (stabilizers) של המעגל. נבצע את הצעדים הבאים:

1. יצירת מעגל: מעגל Clifford אקראי רדוד בצורת brickwork חד-ממדי, ששימש להמחשת הבדיקות.

2. זיהוי qubits כיעדים לבדיקות: מיפוי המעגל ל-qubits פיזיים ב-QPU מסוג Heron בשם ibm_boston וזיהוי qubits עם ancillas שכנים.

3. טרנספילציה של המעגל ל-Backend: בחירת פריסת qubits ותרגום שערי המעגל לקבוצת השערים הטבעית של ibm_boston.

4. יצירת מודל רעש מנתוני ה-Backend: יצירת מודל רעש גס על סמך נתוני benchmark של ה-Backend. נתוני הרעש האלה ישמשו לדירוג היעילות של כל בדיקה נוספת.

5. הצגה חזותית של בדיקות Pauli: הוספת כמה בדיקות Pauli למעגל הרדוד, כדי שנוכל לראות איך הבדיקות נוספות.

6. קביעת היציבות (stabilizers) למדידה: הגדלת עומק המעגל כדי להגביר את ההשפעה של רעש השערים. דגימה אקראית של תת-קבוצה של יציבות מהמעגל העמוק יותר.

7. מציאת בדיקות טובות לכל המצב: מציאת קבוצה אחת של בדיקות ליציבות של המצב שהמעגל מכין. בדיקות כאלה תקפות לא משנה מה נמדד אחר כך, ולכן אותו מעגל עם בדיקות משרת כל יציבות שנדגמה.

8. דגימה מהמעגל הרועש: מדידת המעגל החשוף והמעגל עם הבדיקות בבסיס של כל יציבות שנדגמה באמצעות samplomatic, ודגימה שלהם עם הפרימיטיב Executor על סימולטור רועש שמשתמש באותו רעש depolarizing ששימש למציאת הבדיקות.

9. Postselection של דגימות ממעגל עם בדיקות ללא שגיאות שזוהו: הסרת דגימות שנכשלו בלפחות בדיקה אחת. חישוב ערכי תוחלת עבור המעגל הרועש החשוף והמעגל שמכיל בדיקות.

10. השוואת נאמנות של המצבים הרועשים המקוריים ושל אלה אחרי postselection: הצגה שהשארת דגימות ללא שגיאות שזוהו בלבד משפרת את הנאמנות המוערכת של המצב הרועש.

1. יצירת המעגל

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

# 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_paulice": "qiskit-paulice", "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 numpy as np
from qiskit import QuantumCircuit

def random_clifford_circuit(
    num_qubits: int, depth: int, rng: np.random.Generator
) -> QuantumCircuit:
    """Brickwork random Clifford on `num_qubits`, with `depth` CZ layers."""
    qc = QuantumCircuit(num_qubits)
    qc.h(range(num_qubits))
    for d in range(depth):
        for i in range(d % 2, num_qubits - 1, 2):
            qc.cz(i, i + 1)
        for q in range(num_qubits):
            if rng.integers(0, 2):
                qc.sx(q)
            if rng.integers(0, 2):
                qc.s(q)
            if rng.integers(0, 2):
                qc.sx(q)
    return qc

num_qubits = 12
depth = 4
seed = 1764
rng = np.random.default_rng(seed)
np.random.seed(seed)
circuit = random_clifford_circuit(num_qubits, depth, rng)
circuit.draw("mpl", fold=-1, scale=0.6)

2. זיהוי qubits כיעדים לבדיקות

בשלב הבא נבחר פריסת qubits חד-ממדית על ה-QPU מסוג Heron r3 בשם ibm_boston. נבחר פריסה שמתפתלת באלכסון בגרף הצימוד וממקסמת את מספר qubits העזר הזמינים. כדי למצוא ancillas זמינים נשתמש בפונקציה get_check_qubits, שמחזירה רשימות של זוגות qubit יעד/עזר. בדיקת Pauli על target_qubits[i] תמומש באמצעות ancilla_qubits[i].

בגרף הצימוד שלהלן, ה-qubits הירוקים הם qubits המטען, וה-qubits הכתומים הם ה-ancillas שבהם נשתמש כדי לממש את הבדיקות. ל-qubits 68, 69, 89, 91, 111, 113 ו-133 יש ancillas סמוכים, ולכן נשתמש ב-qubits האלה (שמתאימים ל-qubits וירטואליים 0, 1, 3, 5, 7, 9 ו-11) כ-qubits יעד למימוש בדיקות. שים לב ש-qubit 134 סמוך ל-qubit 133 במטען אבל לא מסומן כ-ancilla. הסיבה לכך היא ש-133 כבר נבחר (באופן שרירותי) לממש בדיקה עם ancilla 132.

from qiskit.visualization import plot_coupling_map
from qiskit_ibm_runtime import QiskitRuntimeService
from qiskit_paulice.layout import get_check_qubits

service = QiskitRuntimeService()
backend = service.backend("ibm_boston")

# Choose a layout and pair each available target qubit with a neighboring ancilla
layout = [68, 69, 78, 89, 90, 91, 98, 111, 112, 113, 119, 133]
target_qubits, ancilla_qubits = get_check_qubits(backend, layout)

print(f"Target qubits:  {target_qubits}")
print(f"Ancilla qubits: {ancilla_qubits}")
plot_coupling_map(
    num_qubits=backend.num_qubits,
    qubit_coordinates=getattr(
        backend.configuration(), "qubit_coordinates", None
    ),
    coupling_map=backend.configuration().coupling_map,
    figsize=(12, 12),
    qubit_color=[
        "#4CAF50"
        if i in set(layout)
        else "#FF9800"
        if i in set(ancilla_qubits)
        else "#DDDDDD"
        for i in backend.coupling_map.graph.node_indices()
    ],
    qubit_size=220,
    line_width=2,
    font_size=90,
)

Target qubits:  [68, 69, 89, 91, 111, 113, 133]
Ancilla qubits: [67, 70, 88, 92, 110, 114, 132]

3. טרנספילציה של המעגל ל-Backend

מכיוון שה-Backend כבר נבחר, אפשר להמשיך ולבצע טרנספילציה של המעגל הווירטואלי שלנו למעגל ISA. כאן אנחנו רק קובעים את פריסת ה-qubits ומתרגמים את השערים לקבוצת השערים הטבעית של ibm_boston.

from qiskit.transpiler import generate_preset_pass_manager

pm = generate_preset_pass_manager(
    optimization_level=0, backend=backend, initial_layout=layout
)
circuit_isa = pm.run(circuit)
circuit_isa.draw("mpl", fold=-1, scale=0.6)

4. יצירת מודל רעש מנתוני ה-Backend

עכשיו, כשאנחנו יודעים איך המעגל ייפרס על ה-QPU, אפשר לדמות איך רעש השערים והקריאה ב-Backend ישפיע על הרצת המעגל. מודל הרעש משמש לקביעת המקומות האופטימליים במעגל להצבת בדיקות Pauli. אנחנו רוצים להציב בדיקות במקומות שבהם הן לוכדות את כמות השגיאה המרבית. מודל רעש מדויק יותר יכול לשפר את יכולות הזיהוי של הבדיקות, אבל בדרך כלל אין צורך ללמוד מודל רעש על ידי דגימת ה-QPU. כאן אנחנו מסיקים מודל רעש depolarizing אחיד לשערים באמצעות נתוני benchmark של ה-Backend מ-qiskit-ibm-runtime בלבד, ומשמיטים את שגיאת הקריאה שהוא מדווח גם כן: הבדיקות שנמצא בשלב 7 קוראות רק את ה-ancilla שלהן, ולכן שגיאות קריאה של המטען בלתי נראות להן והיו רק מורידות את הנאמנות אחרי postselection בשלב 10. למידע נוסף על הגדרת רעש שערים וקריאה, ראה את התיעוד של qiskit_paulice.NoiseModel. למי שרוצה להשתמש במודל רעש Pauli-Lindblad נלמד, ראה NoiseModel.from_pauli_lindblad_maps.

from dataclasses import replace

from qiskit_paulice.noise_models import NoiseModel

# Model the gate noise as a uniform depolarizing channel based on backend benchmark data
noise_model = replace(
    NoiseModel.from_backend(backend, layout, uniform_gate_noise=True),
    readout_noise=None,
)
print(noise_model)

NoiseModel(gate_noise=0.0015808382794757693, readout_noise=None, idling_noise=None)

5. הצגה חזותית של בדיקות Pauli

בוא נוסיף כמה בדיקות למעגל ונראה איך הן נראות. אפשר להשתמש בפונקציה qiskit_paulice.add_pauli_checks כדי למצוא קבוצה טובה של בדיקות תקפות. אנחנו מספקים כקלט את מעגל ה-Clifford של המטען, יחד עם רשימת qubits של המטען שעליהם יש להוסיף בדיקות (qubits יעד) ואת מודל הרעש שיצרנו בשלב הקודם. בדוגמה הזו נוסיף בדיקות לכל 7 ה-qubits היעד שיש להם ancillas זמינים. הבדיקות מתווספות ל-qubits לפי הסדר שבו הם מופיעים ב-target_qubits, וזה מקל להבין את הפריסה הסופית של המעגל עם הבדיקות: final_layout = layout + ancilla_qubits. אם משתמש מחליט שהוא מעדיף להריץ אחד ממעגלי הפלט עם פחות בדיקות (i), הפריסה הסופית תהיה final_layout = layout + ancilla_qubits[:i+1].

ההצגה החזותית של המעגל מראה שהבדיקות מומשו נכון באמצעות זוגות היעד/עזר שהוגדרו.

פרטים על qiskit_paulice.add_pauli_checks:

• הקלט הראשון צריך להיות מעגל Clifford. כאן אין בו מדידות ואנחנו מעבירים stabilizers="all", ולכן הבדיקות נמצאות ליציבות של המצב שהמעגל מכין: בדיקה תקפה אם מכפלת Pauli שלה במרחב-זמן מודדת יציבות כזו, הסינדרום שלה הוא סיבית ה-ancilla שלה בלבד, והיא נשארת תקפה לא משנה באיזה בסיס המטען נמדד אחר כך. זה מה שהערכת נאמנות צריכה, כי כל יציבות שנדגמה נמדדת בבסיס אחר. (השמטת stabilizers מוצאת במקום זאת בדיקות למדידות הסופיות של מעגל.)
• הארגומנט המיקומי השני הוא רשימה של אינדקסים של qubits יעד. כל qubit יעד ישתף שערי הסתבכות עם qubit עזר כדי ליצור בדיקה. הבדיקות יימצאו ויתווספו למעגל לפי הסדר שבו הן נתונות בארגומנט הזה.
• הארגומנט המיקומי האחרון הוא מודל רעש שמייצג את הרעש הצפוי להשפיע על המעגל בזמן ההרצה. המודל הזה משמש להערכת היעילות של בדיקה נתונה.
• ancilla_qubits: אם למעגל הקלט יש פריסה, ה-kwarg הזה מורה ל-add_pauli_checks איך להקצות qubits פיזיים ל-ancillas שנוספו. ancilla_qubits צריך להיות באותו אורך כמו target_qubits (הארגומנט המיקומי השני), ו-ancilla_qubits[i] ישמש למימוש בדיקות על target_qubits[i].
• cost: פונקציית העלות ששימשה להערכת בדיקות. יכולה להיות gamma או LER (שיעור שגיאה לוגית), שכולן דרכים לתאר כמה מערוץ הרעש נשאר לא מזוהה על ידי הבדיקות. עם stabilizers="all" שגיאה נחשבת לוגית אם היא משנה את המצב המוכן, וזה בדיוק מה שהנאמנות מודדת.
• method: השיטה ששימשה למציאת בדיקות תקפות. אנחנו משתמשים ב-windowed שדוגמת באקראי הרבה תת-קבוצות של חלונות מרחביים מתוך התומך (החוטים) של qubit יעד. ערכים תקפים אחרים הם genetic ו-windowed_genetic.
• הפלט של add_pauli_checks מכיל מעגלים עם מספר גדל והולך של בדיקות, ממעגל ללא בדיקות ועד מעגל עם בדיקות על כל qubit יעד (בהנחה שנמצאו בדיקות תקפות).

ראה סעיפים II-IV של המידע המשלים למידע מפורט על מציאת בדיקות Pauli טובות במרחב-זמן

import matplotlib.pyplot as plt
from qiskit_paulice import add_pauli_checks

# Find checks on circuit
checked = add_pauli_checks(
    circuit_isa,
    target_qubits,
    noise_model,
    ancilla_qubits=ancilla_qubits,
    cost="gamma",
    method="windowed",
    seed=seed,
    stabilizers="all",
)

print(f"Physical layout of payload+ancillas: {layout + ancilla_qubits}")
print("Checked circuit:")
checked[-1].circuit.draw("mpl", fold=-1, idle_wires=False)

Physical layout of payload+ancillas: [68, 69, 78, 89, 90, 91, 98, 111, 112, 113, 119, 133, 67, 70, 88, 92, 110, 114, 132]
Checked circuit:

6. קביעת היציבות (stabilizers) למדידה

כדי לבדוק את ביצועי זיהוי השגיאות, נעריך את הנאמנות של מצב היציבות, $|\psi\rangle = U|0\rangle^{\otimes n}$, שהמעגל מכין באופן אידיאלי, מול המצב הרועש, $\rho$, שהחומרה מוציאה בפועל. ההטלה על מצב היציבות הטהור, $|\psi\rangle$, מוגדרת כממוצע אחיד על $2^n$ האיברים של חבורת היציבות שלו $\mathcal{S}$:

$|\psi\rangle\langle\psi| = \frac{1}{2^n}\sum_{G \in \mathcal{S}} G$

אם נציב זאת במשוואת הנאמנות, נראה שהנאמנות של המצב הרועש, $\rho$, שווה לערך התוחלת הממוצע של כל היציבות, $G \in \mathcal{S}$, ביחס ל-$\rho$:

$F = \mathrm{Tr}(\rho|\psi\rangle\langle\psi|) = \mathrm{Tr}({\frac{1}{2^n}\sum_{G \in \mathcal{S}} \rho G}) = \frac{1}{2^n} \sum_{G \in \mathcal{S}} \mathrm{Tr}(\rho G) = \frac{1}{2^n} \sum_{G \in \mathcal{S}} \langle G \rangle_\rho$

בבעיות גדולות יותר, מדידה (או אפילו ספירה) של כל $2^n$ היציבות אינה אפשרית, ולכן הנאמנות מוערכת לרוב על ידי מדידת מדגם אקראי של יציבות. דגימה של $M$ יציבות $G_1, \ldots, G_M$ באופן אחיד ואקראי מתוך $\mathcal{S}$ נותנת אומד חסר הטיה לנאמנות, $\hat F_M$,

$\hat F_M = \frac{1}{M} \sum_{i=1}^M \langle G_i \rangle_\rho$ בניסוי הזה נעמיק את המעגל משמעותית כדי להגביר את ההשפעה של רעש השערים, ונוסיף בדיקות לכל qubit שיש לו qubit עזר שכן. יש דרכים מתוחכמות יותר לדגום יציבות להערכת נאמנות, אבל בדוגמה הזו נדגום $100$ יציבות בהסתברות אחידה.

from qiskit.quantum_info import Clifford, Pauli, PauliList

num_stabilizers = 100
num_shots = 1_000
num_checks = len(target_qubits)
depth = 72

circuit = random_clifford_circuit(num_qubits, depth, rng)
# Create all 2^N stabilizers
stabilizers = PauliList([Pauli("I" * num_qubits)])
for generator in (
    Pauli(label) for label in Clifford(circuit).to_labels(mode="S")
):
    stabilizers = stabilizers + stabilizers.compose(generator)

# Choose some random stabilizers from the set
keep = np.where(stabilizers.x.any(axis=1) | stabilizers.z.any(axis=1))[0]
chosen = np.random.default_rng(seed).choice(
    keep, size=min(num_stabilizers, len(keep)), replace=False
)
stabilizers = stabilizers[chosen]
depth = circuit.depth(lambda x: x.operation.num_qubits == 2)
print(
    f"Sampled {len(stabilizers)} stabilizers of {circuit.num_qubits}-qubit circuit with two-qubit depth of {depth}: {{{stabilizers[0]}, {stabilizers[1]}, ...}}"
)

Sampled 100 stabilizers of 12-qubit circuit with two-qubit depth of 72: {-YIIYZZZIZIYY, -XIXYIIZXYIZX, ...}

7. מציאת בדיקות טובות לכל המצב

נמצא קבוצה אחת של בדיקות למעגל העמוק בדיוק כמו בשלב 5. מכיוון שהבדיקות קוראות רק את ה-ancillas שלהן, אותו מעגל עם בדיקות משרת כל יציבות שנדגמה. המעגלים המוחזרים מודדים רק את ה-ancillas של הבדיקות, ולכן נוסיף את מדידות המטען בשלב הבא.

import time

# Transpile the deep payload without measurements and find one set of checks for the whole state
payload_isa = pm.run(circuit)
t0 = time.time()
checked = add_pauli_checks(
    payload_isa,
    target_qubits,
    noise_model,
    ancilla_qubits=ancilla_qubits,
    cost="gamma",
    method="windowed",
    seed=seed,
    stabilizers="all",
)
print(
    f"Added {len(checked[-1].check_qubits)} checks in {time.time() - t0:.0f}s."
)
print(
    f"2-qubit depth increased from {payload_isa.depth(lambda x: len(x.qubits) == 2)} to "
    f"{checked[-1].circuit.depth(lambda x: len(x.qubits) == 2)} when adding {num_checks} checks."
)

Added 7 checks in 3s.
2-qubit depth increased from 72 to 85 when adding 7 checks.

8. דגימה מהמעגל הרועש

אנחנו דוגמים כל וריאנט, מהמטען החשוף ועד המעגל עם כל 7 הבדיקות, עבור כל יציבות שנדגמה. הווריאנט עם אפס בדיקות הוא המטען החשוף עם אותן מדידות, ולכן הוא משמש כהתייחסות הרועשת. במקום לצרף את סיבוב הבסיס של כל יציבות לעותק של כל וריאנט, אנחנו נותנים ל-samplomatic לבחור את בסיס המדידה בזמן הדגימה. לכל וריאנט אנחנו:

• מודדים את המטען ומקבצים את המדידות הסופיות לתוך box באמצעות ה-boxing pass manager של samplomatic. ה-box נושא הערת ChangeBasis, ולכן build מחזיר מעגל תבנית ששערי שינוי הבסיס שלו מפורמטרים, יחד עם samplex שמחשב את הפרמטרים עבור בסיס מבוקש.

• מוסיפים את התבנית ואת ה-samplex ל-QuantumProgram. ה-samplex מקבל בסיס מדידה אחד לכל qubit ב-box, בסדר פיזי עולה, מקודד כ-$I=0$, $Z=1$, $X=2$, $Y=3$. העברת מערך של בסיסים עם ציר מוביל באורך 100 גורמת לפריט יחיד למדוד את כל 100 היציבות: המטען בבסיס של כל יציבות וה-ancillas ב-$Z$. qubit מחוץ לתומך של יציבות נמדד גם הוא ב-$Z$, והסיבית שלו מתעלמת בעת חישוב ערך התוחלת של היציבות.

הפרימיטיב Executor מ-qiskit-ibm-runtime מריץ את כל התוכנית בעבודה אחת. כשנותנים לו AerSimulator במקום QPU, הוא רץ מקומית על הסימולטור הזה, ולכן אנחנו נותנים לו סימולטור עם אותו רעש depolarizing בשערים ששימש אותנו למציאת הבדיקות, ומשתמשים בשיטת היציבות כי כל המעגלים הם Clifford. Qiskit Aer נחוץ רק לבניית הסימולטור הזה: המצב המקומי לוקח את ה-seed מאפשרויות ה-Executor אבל מתעלם מה-noise_model שלהן, ולכן הרעש חייב להיקבע על הסימולטור עצמו. העברת Backend של QPU במקום זאת תריץ את אותה תוכנית על חומרה. אחרי שהסימולציות יסתיימו, נחשב את שיעור ה-postselection ואת ערכי התוחלת כדי לקבל תחושה עד כמה זיהוי השגיאות עבד, לפני שנציג את הנתונים חזותית.

from qiskit import ClassicalRegister
from qiskit_aer import AerSimulator
from qiskit_aer.noise import NoiseModel as AerNoiseModel
from qiskit_aer.noise import depolarizing_error
from qiskit_ibm_runtime import Executor
from qiskit_ibm_runtime.quantum_program import QuantumProgram
from samplomatic.builders import build
from samplomatic.transpiler import generate_boxing_pass_manager

# Group the terminal measurements into a box whose measurement basis is chosen at sampling time
boxing_pm = generate_boxing_pass_manager(
    enable_gates=False, measure_annotations="change_basis"
)

# bases[i, q]: basis in which physical qubit q is measured for stabilizer i (I=0, Z=1, X=2, Y=3)
bases = np.ones((len(stabilizers), backend.num_qubits), dtype=np.uint8)
bases[:, layout] = stabilizers.z + 2 * stabilizers.x

program = QuantumProgram(shots=num_shots)
for variant in checked:
    circ = variant.circuit.copy()
    circ.add_register(meas := ClassicalRegister(num_qubits, "meas"))
    circ.measure(layout, meas)
    template, samplex = build(boxing_pm.run(circ))
    (basis_input,) = samplex.inputs().get_specs()
    measured_qubits = sorted(layout + list(variant.check_qubits))
    program.append_samplex_item(
        template,
        samplex=samplex,
        samplex_arguments={basis_input.name: bases[:, measured_qubits]},
        shape=(len(stabilizers),),
    )

# Simulate the program with the same gate noise we used to find the checks
aer_nm = AerNoiseModel()
aer_nm.add_all_qubit_quantum_error(
    depolarizing_error(noise_model.gate_noise, 2), ["cz"]
)
noisy_sim = AerSimulator(method="stabilizer", noise_model=aer_nm)

# outs[k][register]: bits of the variant with k checks, shape (stabilizers, shots, bits)
t0 = time.time()
outs = (
    Executor(noisy_sim, options={"simulator": {"seed_simulator": seed}})
    .run(program)
    .result()
)
print(
    f"Sampled {len(checked)} circuits for {len(stabilizers)} stabilizers in {time.time() - t0:.0f}s."
)

Sampled 8 circuits for 100 stabilizers in 77s.

9. Postselection של דגימות ממעגל עם בדיקות

בחבילה הזו בדיקה ממומשת באמצעות שערי הסתבכות בין qubit עזר אחד ל-qubit יעד אחד. כל ancilla מתחיל ב-$|0\rangle$, ולכן $Z_\text{anc}$ מייצב את מצב הקלט שלו. הפצה קדימה של $Z_\text{anc}$ מתחילת ה-ancilla דרך כל המעגל עם הבדיקות מניבה $Z_\text{anc}$ כפול יציבות של המצב המוכן, ולכן הבדיקה עוברת אם סיבית ה-ancilla היא $0$. דגימה נשמרת אם כל הבדיקות עוברות.

המעגלים עם הבדיקות מודדים את ה-ancillas אל האוגר checks_c, ולכן דגימה נשמרת אם כל סיבית באוגר הזה היא 0. באופן כללי, השיטה get_postselection_method של מופעי CheckedCircuit שמוחזרים מ-add_pauli_checks ממפה מחרוזת סיביות של shot לסינדרום של כל הבדיקות בבת אחת; עבור בדיקות שנמצאו ליציבות של המצב המוכן, כל סינדרום הוא פשוט סיבית ה-ancilla שלו, ולכן אנחנו קוראים את האוגר ישירות.

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

def keep_mask(out):
    """Per-shot True/False: every ancilla read 0, so every Pauli check passed."""
    if "checks_c" not in out:  # the bare payload has no checks
        return np.ones(out["meas"].shape[:-1], dtype=bool)
    return ~out["checks_c"].any(axis=-1)

# masks[k][i]: kept shots of stabilizer i measured on the variant with k checks
masks = [keep_mask(out) for out in outs]
psrs_arr = np.array([mask.mean(axis=-1) for mask in masks]).T
ks = np.arange(psrs_arr.shape[1])

fig, ax = plt.subplots(figsize=(8, 4))
ax.plot(ks, psrs_arr.T, color="darkorange", alpha=0.15, linewidth=1)
ax.plot(
    ks,
    np.median(psrs_arr, axis=0),
    color="black",
    linewidth=1,
    linestyle="--",
    label="median",
)
ax.set(
    xlabel="Checks committed",
    ylabel="Post-selection rate",
    ylim=(0, 1.05),
    title=f"Per-stabilizer post-selection rate ({len(psrs_arr)} stabilizers)",
)
ax.legend()
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.show()

10. השוואת נאמנות של המצבים הרועשים המקוריים ושל אלה אחרי postselection

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

def expectation_values(out, mask):
    """EV of each stabilizer from the shots of ``out`` kept by ``mask``.

    Bit q of ``meas`` is virtual qubit q, so a stabilizer's eigenvalue in a shot is the parity of
    the bits in its support, times its sign.
    """
    support = stabilizers.x | stabilizers.z
    parity = (out["meas"] & support[:, None, :]).sum(axis=-1) % 2
    signs = np.where(stabilizers.phase == 2, -1.0, 1.0)
    eigenvalues = np.where(parity, -1.0, 1.0) * signs[:, None]
    kept = mask.sum(axis=-1)
    return np.where(
        kept > 0,
        (eigenvalues * mask).sum(axis=-1) / np.maximum(kept, 1),
        np.nan,
    )

evs_noisy = expectation_values(outs[0], masks[0])
evs_post = expectation_values(outs[-1], masks[-1])
print(
    f"ideal fidelity: 1.0\n"
    f"noisy fidelity: {np.nanmean(evs_noisy):.4f}\n"
    f"postselected fidelity: {np.nanmean(evs_post):.4f}\n"
    f"mean post-sel rate: {masks[-1].mean():.3f}\n"
)

evs_ideal = np.ones(len(stabilizers))
idx = np.arange(len(stabilizers))

def cum_mean_sem(ys):
    valid = ~np.isnan(ys)
    s = np.cumsum(np.where(valid, ys, 0.0))
    s2 = np.cumsum(np.where(valid, ys**2, 0.0))
    n = np.maximum(np.cumsum(valid).astype(float), 1)
    mean = s / n
    sem = np.sqrt(np.maximum(s2 / n - mean**2, 0) / n)
    return np.where(np.cumsum(valid) > 0, mean, np.nan), sem

def strip(ax, ys, color, label):
    m, s = np.nanmean(ys), np.nanstd(ys)
    ax.axhspan(
        m - s, m + s, color=color, alpha=0.15, label=f"{label} mean ± std"
    )
    ax.axhline(m, color=color, linewidth=1, linestyle="--")

fig, ax = plt.subplots(figsize=(8, 4))
ax.axhline(np.nanmean(evs_ideal), color="black", linewidth=1.5, label="ideal")
strip(ax, evs_noisy, "red", "noisy")
strip(ax, evs_post, "green", "post-sel")
ax.scatter(
    idx, evs_noisy, color="red", s=22, alpha=0.7, zorder=3, label="noisy <G>"
)
ax.scatter(
    idx,
    evs_post,
    color="green",
    s=22,
    alpha=0.7,
    zorder=3,
    label="post-sel <G>",
)
ax.set(
    xlabel="stabilizer index",
    ylabel=r"$\langle G \rangle$",
    ylim=(-0.1, 1.1),
    title="Per-stabilizer EVs",
)
ax.legend(loc="lower left", ncol=2)
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.show()

M = np.arange(1, len(stabilizers) + 1)
fig, ax = plt.subplots(figsize=(8, 4))
for ys, color, label in [
    (evs_ideal, "black", "ideal"),
    (evs_noisy, "red", "noisy"),
    (evs_post, "green", "post-sel"),
]:
    cm, sem = cum_mean_sem(ys)
    ax.plot(M, cm, color=color, linewidth=1.5, label=label)
    ax.fill_between(M, cm - sem, cm + sem, color=color, alpha=0.15)
ax.set(
    xlabel="number of stabilizers averaged",
    ylabel="running fidelity estimate",
    title="Fidelity convergence vs. number of stabilizers",
)
ax.legend()
ax.grid(True, alpha=0.3)
plt.tight_layout()
plt.show()

ideal fidelity: 1.0
noisy fidelity: 0.5582
postselected fidelity: 0.9747
mean post-sel rate: 0.535

מכיוון שקבוצה אחת של בדיקות משרתת כל יציבות, יש עקומת gamma אחת להסתכל עליה: ערך פונקציית העלות אחרי שכל בדיקה נקבעה. בבעיות אחרות על קבוצה שונה של qubits פיזיים, ועם מודל רעש לא מושלם (זכור, מודל הרעש במדריך הזה דימה בדיוק את הרעש שהשפיע על הדגימות), לסדר שבו הבדיקות מתווספות ל-qubits היעד יכולה להיות השפעה בולטת יותר על האופן שבו העקומה מתכנסת.

עקומת gamma זו שמתכנסת ל-1.0 מרמזת שבדיקות נוספות לא יוכלו לזהות עוד הרבה שגיאות, כי יש מעט מאוד שגיאה מדומה שהבדיקות לא מכסות. זה מחזק את התצפית הקודמת שלנו שההתכנסות בשיעור ה-postselection נובעת מכך שאין במעגל שגיאה לא מכוסה.

costs = [variant.cost for variant in checked]

fig, ax = plt.subplots(figsize=(8, 4))
ax.plot(range(len(costs)), costs, color="steelblue", marker="o")
ax.set(
    xlabel="Checks committed",
    ylabel="Gamma",
    yscale="log",
    title="Gamma of the noise left uncovered by the checks",
)
ax.grid(True, alpha=0.3, which="both")
plt.tight_layout()
plt.show()