סימולציה של מודל Ising דו-ממדי עם שדה מוטה באמצעות פונקציית QESEM

הערה

Qiskit Functions are an experimental feature available only to IBM Quantum® Premium Plan, Flex Plan, and On-Prem (via IBM Quantum Platform API) Plan users. They are in preview release status and subject to change.

הערכת שימוש: 20 דקות על מעבד Heron r2. (הערה: זוהי הערכה בלבד. זמן הריצה שלך עשוי להשתנות.)

רקע

מדריך זה מראה כיצד להשתמש ב-QESEM, פונקציית Qiskit של Qedma, כדי לדמות את הדינמיקה של מודל ספין קוונטי קנוני, מודל Ising דו-ממדי עם שדה מוטה (TFI) עם זוויות שאינן Clifford:

$$H = J \sum_{\langle i,j \rangle} Z_i Z_j + g_x \sum_i X_i + g_z \sum_i Z_i ,$$

כאשר $\langle i,j \rangle$ מציין שכנים קרובים על סריג. סימולציה של התפתחות הזמן של מערכות קוונטיות רבות-גופיים היא משימה קשה מבחינה חישובית עבור מחשבים קלאסיים. מחשבים קוונטיים, לעומת זאת, מתוכננים באופן טבעי לבצע משימה זו ביעילות. מודל TFI, בפרט, הפך למדד ביצועים פופולרי על חומרה קוונטית בשל ההתנהגות הפיזיקלית העשירה שלו והיישום הידידותי לחומרה.

במקום לדמות דינמיקה בזמן רציף, אנו מאמצים את מודל Ising המופעל הקשור קרוב. הדינמיקה יכולה להתבטא בדיוק כמעגל קוונטי תקופתי, כאשר כל שלב התפתחות מורכב משלוש שכבות של שערים של שני קיוביטים שבריריים $R_{ZZ} (\alpha_{ZZ})$, משולבות עם שכבות של שערים של קיוביט בודד $R_X (\alpha_X)$ ו-$R_Z (\alpha_Z)$.

נשתמש בזוויות גנריות שמהוות אתגר גם לסימולציה קלאסית וגם להפחתת שגיאות. באופן ספציפי, בחרנו $\alpha_{ZZ} = 1.0$, $\alpha_X = 0.53$ ו-$\alpha_Z = 0.1$, ממקמים את המודל רחוק מכל נקודה ניתנת לאינטגרציה.

במדריך זה נבצע את הפעולות הבאות:

• הערכת זמן הריצה הצפוי של QPU להפחתת שגיאות מלאה באמצעות תכונות ההערכה האנליטיות והאמפיריות של QESEM.
• בניית וסימולציה של מעגל מודל Ising דו-ממדי עם שדה מוטה באמצעות פריסות קיוביט ושכבות שער בהשראת חומרה.
• הצגה חזותית של קישוריות קיוביטים של המכשיר ותתי-גרפים נבחרים עבור הניסוי שלך.
• הדגמת השימוש ב-הפצה לאחור של אופרטורים (OBP) להפחתת עומק המעגל. טכניקה זו קוצצת פעולות מסוף המעגל במחיר של יותר מדידות אופרטור.
• ביצוע הפחתת שגיאות (EM) לא מוטה עבור מספר ניתנים למדידה בו-זמנית באמצעות QESEM, תוך השוואה בין תוצאות אידיאליות, רועשות ומופחתות שגיאות.
• ניתוח ושרטוט של השפעת הפחתת השגיאות על המגנטיזציה בעומקי מעגל שונים.

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

אודות QESEM

QESEM הוא תוכנה אמינה, בדיוק גבוה, מבוססת-אפיון המיישמת הפחתת שגיאות קוואזי-הסתברותית יעילה ולא מוטה. היא מיועדת להפחית שגיאות במעגלים קוונטיים גנריים והיא אגנוסטית ליישום. היא אומתה על פני פלטפורמות חומרה מגוונות, כולל ניסויים בקנה מידה שימושיות על התקני Eagle ו-Heron של IBM®. שלבי זרימת העבודה של QESEM הם כדלקמן:

1. אפיון התקן - ממפה נאמנויות שער ומזהה שגיאות קוהרנטיות, מספק נתוני כיול בזמן אמת. שלב זה מבטיח שההפחתה ממנפת את הפעולות בעלות הנאמנות הגבוהה ביותר הזמינות.
2. טרנספיליציה מודעת-רעש - מייצר ומעריך מיפויי קיוביט חלופיים, קבוצות פעולות ובסיסי מדידה, בוחר את הגרסה שממזערת את זמן הריצה המוערך של QPU, עם הקבלה אופציונלית לאיסוף נתונים מואץ.
3. דיכוי שגיאות - מגדיר מחדש שערים מקוריים, מיישם סיבוב Pauli ומייעל בקרה ברמת הדופק (בפלטפורמות נתמכות) לשיפור הנאמנות.
4. אפיון מעגל - בונה מודל שגיאה מקומי מותאם ומתאים אותו למדידות QPU כדי לכמת רעש שיורי.
5. הפחתת שגיאות - בונה פירוקים קוואזי-הסתברותיים רב-סוגיים, ודוגם מהם בתהליך אדפטיבי שממזער את זמן QPU של ההפחתה ואת הרגישות לתנודות חומרה, ומשיג דיוקים גבוהים בנפחי מעגל גדולים.

למידע נוסף על QESEM וניסוי בקנה מידה שימושיות של מודל זה על תת-גרף של 103 קיוביטים, בעל קישוריות גבוהה, של הגיאומטריה המקורית heavy-hex של ibm_marrakesh, ראה Reliable high-accuracy error mitigation for utility-scale quantum circuits.

דרישות

התקן את חבילות Python הבאות לפני הרצת ה-notebook:

• Qiskit SDK v2.0.0 או מאוחר יותר (pip install qiskit)
• Qiskit Runtime v0.40.0 או מאוחר יותר (pip install qiskit-ibm-runtime)
• Qiskit Functions Catalog v0.8.0 או מאוחר יותר ( pip install qiskit-ibm-catalog )
• Operator Backpropagation Qiskit addon v0.3.0 או מאוחר יותר ( pip install qiskit-addon-obp )
• Qiskit Utils addon v0.1.1 או מאוחר יותר ( pip install qiskit-addon-utils )
• Qiskit Aer simulator v0.17.1 או מאוחר יותר ( pip install qiskit-aer )
• Matplotlib v3.10.3 או מאוחר יותר ( pip install matplotlib )

הגדרה

ראשית, ייבא את הספריות הרלוונטיות:

# Added by doQumentation — required packages for this notebook
!pip install -q matplotlib numpy qiskit qiskit-addon-obp qiskit-addon-utils qiskit-aer qiskit-ibm-catalog qiskit-ibm-runtime

%matplotlib inline

from typing import Sequence

import matplotlib.pyplot as plt
import numpy as np

import qiskit
from qiskit_ibm_runtime import EstimatorV2 as Estimator
from qiskit_ibm_catalog import QiskitFunctionsCatalog
from qiskit_aer import AerSimulator
from qiskit_addon_utils.slicing import combine_slices, slice_by_gate_types
from qiskit_addon_obp import backpropagate
from qiskit_addon_obp.utils.simplify import OperatorBudget
from qiskit_ibm_runtime import QiskitRuntimeService
from qiskit.visualization import (
    plot_gate_map,
)

לאחר מכן, אמת את הזהות באמצעות מפתח ה-API שלך מלוח המחוונים של IBM Quantum Platform. לאחר מכן, בחר את פונקציית Qiskit באופן הבא. (שים לב שלמען האבטחה, עדיף לשמור את אישורי החשבון שלך בסביבה המקומית שלך, אם אתה על מכונה מהימנה, כדי שלא תצטרך להזין את מפתח ה-API שלך בכל פעם שאתה מאמת זהות.)

# Paste here your instance and token strings

instance = "YOUR_INSTANCE"
token = "YOUR_TOKEN"
channel = "ibm_quantum_platform"

catalog = QiskitFunctionsCatalog(
    channel=channel, token=token, instance=instance
)
qesem_function = catalog.load("qedma/qesem")

שלב 1: מיפוי קלטים קלאסיים לבעיה קוונטית

אנו מתחילים על ידי הגדרת פונקציה שיוצרת את מעגל Trotter:

def trotter_circuit_from_layers(
    steps: int,
    theta_x: float,
    theta_z: float,
    theta_zz: float,
    layers: Sequence[Sequence[tuple[int, int]]],
    init_state: str | None = None,
) -> qiskit.QuantumCircuit:
    """
    Generates an ising trotter circuit
    :param steps: trotter steps
    :param theta_x: RX angle
    :param theta_z: RZ angle
    :param theta_zz: RZZ angle
    :param layers: list of layers (can be list of layers in device)
    :param init_state: Initial state to prepare. If None, will not prepare any state. If "+", will
     add Hadamard gates to all qubits.
    :return: QuantumCircuit
    """
    qubits = sorted({i for layer in layers for edge in layer for i in edge})
    circ = qiskit.QuantumCircuit(max(qubits) + 1)

    if init_state == "+":
        print("init_state = +")
        for q in qubits:
            circ.h(q)

    for _ in range(steps):
        for q in qubits:
            circ.rx(theta_x, q)
            circ.rz(theta_z, q)

        for layer in layers:
            for edge in layer:
                circ.rzz(theta_zz, *edge)
        circ.barrier(qubits)

    return circ

לאחר מכן אנו יוצרים פונקציה לחישוב ערכי תוחלת אידיאליים באמצעות AerSimulator.

שים לב שעבור מעגלים גדולים (30 קיוביטים או יותר) אנו ממליצים להשתמש בערכים מחושבים מראש מסימולציות PEPS של הפצה-אמונה (BP). קוד זה כולל ערכים מחושבים מראש עבור 35 קיוביטים כדוגמה, בהתבסס על גישת BP לאבולוציה של רשת טנזור PEPS המוצגת במאמר זה (שאליה אנו מתייחסים כ-PEPS-BP), תוך שימוש בחבילת Python לרשתות-טנזור quimb.

def calculate_ideal_evs(circ, obs, num_qubits, step):
    # Predefined results for large circuits - calculated using bppeps for 3, 5, 7, 9 trotter steps
    predefined_35 = [
        0.79537,
        0.78653,
        0.79699,
    ]

    if num_qubits == 35:
        print(
            "Using precalculated ideal values for large circuits calculated with belief propagation PEPS. Currently only for 35 qubits."
        )
        return predefined_35[step]

    else:
        simulator = AerSimulator()

        # Use Estimator primitive to get expectation value
        estimator = Estimator(simulator)
        sim_result = estimator.run([(circ, [obs])], precision=0.0001).result()

        # Extracting the result
        ideal_values = sim_result[0].data.evs[0]
        return ideal_values

אנו משתמשים במיפוי שכבות $R_{ZZ}$ מבוסס-חומרה שנלקח מהתקן Heron, שממנו אנו חותכים את השכבות בהתאם למספר הקיוביטים שברצוננו לדמות. אנו מגדירים תתי-גרפים עבור 10, 21, 28 ו-35 קיוביטים השומרים על מבנה דו-ממדי (אל תהסס לשנות לתת-גרף האהוב עליך):

LAYERS_HERON_R2 = [  # the full set of hardware layers for Heron r2
    [
        (2, 3),
        (6, 7),
        (10, 11),
        (14, 15),
        (20, 21),
        (16, 23),
        (24, 25),
        (17, 27),
        (28, 29),
        (18, 31),
        (32, 33),
        (19, 35),
        (36, 41),
        (42, 43),
        (37, 45),
        (46, 47),
        (38, 49),
        (50, 51),
        (39, 53),
        (60, 61),
        (56, 63),
        (64, 65),
        (57, 67),
        (68, 69),
        (58, 71),
        (72, 73),
        (59, 75),
        (76, 81),
        (82, 83),
        (77, 85),
        (86, 87),
        (78, 89),
        (90, 91),
        (79, 93),
        (94, 95),
        (100, 101),
        (96, 103),
        (104, 105),
        (97, 107),
        (108, 109),
        (98, 111),
        (112, 113),
        (99, 115),
        (116, 121),
        (122, 123),
        (117, 125),
        (126, 127),
        (118, 129),
        (130, 131),
        (119, 133),
        (134, 135),
        (140, 141),
        (136, 143),
        (144, 145),
        (137, 147),
        (148, 149),
        (138, 151),
        (152, 153),
        (139, 155),
    ],
    [
        (1, 2),
        (3, 4),
        (5, 6),
        (7, 8),
        (9, 10),
        (11, 12),
        (13, 14),
        (21, 22),
        (23, 24),
        (25, 26),
        (27, 28),
        (29, 30),
        (31, 32),
        (33, 34),
        (40, 41),
        (43, 44),
        (45, 46),
        (47, 48),
        (49, 50),
        (51, 52),
        (53, 54),
        (55, 59),
        (61, 62),
        (63, 64),
        (65, 66),
        (67, 68),
        (69, 70),
        (71, 72),
        (73, 74),
        (80, 81),
        (83, 84),
        (85, 86),
        (87, 88),
        (89, 90),
        (91, 92),
        (93, 94),
        (95, 99),
        (101, 102),
        (103, 104),
        (105, 106),
        (107, 108),
        (109, 110),
        (111, 112),
        (113, 114),
        (120, 121),
        (123, 124),
        (125, 126),
        (127, 128),
        (129, 130),
        (131, 132),
        (133, 134),
        (135, 139),
        (141, 142),
        (143, 144),
        (145, 146),
        (147, 148),
        (149, 150),
        (151, 152),
        (153, 154),
    ],
    [
        (3, 16),
        (7, 17),
        (11, 18),
        (22, 23),
        (26, 27),
        (30, 31),
        (34, 35),
        (21, 36),
        (25, 37),
        (29, 38),
        (33, 39),
        (41, 42),
        (44, 45),
        (48, 49),
        (52, 53),
        (43, 56),
        (47, 57),
        (51, 58),
        (62, 63),
        (66, 67),
        (70, 71),
        (74, 75),
        (61, 76),
        (65, 77),
        (69, 78),
        (73, 79),
        (81, 82),
        (84, 85),
        (88, 89),
        (92, 93),
        (83, 96),
        (87, 97),
        (91, 98),
        (102, 103),
        (106, 107),
        (110, 111),
        (114, 115),
        (101, 116),
        (105, 117),
        (109, 118),
        (113, 119),
        (121, 122),
        (124, 125),
        (128, 129),
        (132, 133),
        (123, 136),
        (127, 137),
        (131, 138),
        (142, 143),
        (146, 147),
        (150, 151),
        (154, 155),
        (0, 1),
        (4, 5),
        (8, 9),
        (12, 13),
        (54, 55),
        (15, 19),
    ],
]

subgraphs = {  # the subgraphs for the different qubit counts such that it's 2D
    10: list(range(22, 29)) + [16, 17, 37],
    21: list(range(3, 12)) + list(range(23, 32)) + [16, 17, 18],
    28: list(range(3, 12))
    + list(range(23, 32))
    + list(range(45, 50))
    + [16, 17, 18, 37, 38],
    35: list(range(3, 12))
    + list(range(21, 32))
    + list(range(41, 50))
    + [16, 17, 18, 36, 37, 38],
    42: list(range(3, 12))
    + list(range(21, 32))
    + list(range(41, 50))
    + list(range(63, 68))
    + [16, 17, 18, 36, 37, 38, 56, 57],
}

n_qubits = 35  # 21, 28, 35, 42

layers = [
    [
        edge
        for edge in layer
        if edge[0] in subgraphs[n_qubits] and edge[1] in subgraphs[n_qubits]
    ]
    for layer in LAYERS_HERON_R2
]

print(layers)

[[(6, 7), (10, 11), (16, 23), (24, 25), (17, 27), (28, 29), (18, 31), (36, 41), (42, 43), (37, 45), (46, 47), (38, 49)], [(3, 4), (5, 6), (7, 8), (9, 10), (21, 22), (23, 24), (25, 26), (27, 28), (29, 30), (43, 44), (45, 46), (47, 48)], [(3, 16), (7, 17), (11, 18), (22, 23), (26, 27), (30, 31), (21, 36), (25, 37), (29, 38), (41, 42), (44, 45), (48, 49), (4, 5), (8, 9)]]

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

service = QiskitRuntimeService(
    channel=channel,
    token=token,
    instance=instance,
)
backend = service.backend("ibm_fez")  # or any available device

selected_qubits = subgraphs[n_qubits]
num_qubits = backend.configuration().num_qubits
qubit_color = [
    "#ff7f0e" if i in selected_qubits else "#d3d3d3"
    for i in range(num_qubits)
]

plot_gate_map(
    backend=backend,
    figsize=(15, 10),
    qubit_color=qubit_color,
)
plt.show()

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

כעת אנו יוצרים את מעגל Trotter ואת התצפית ממוצע המגנטיזציה עבור מספר הקיוביטים והפרמטרים שנבחרו:

# Chosen parameters:
theta_x = 0.53
theta_z = 0.1
theta_zz = 1.0
steps = 9

circ = trotter_circuit_from_layers(steps, theta_x, theta_z, theta_zz, layers)
print(
    f"Circuit 2q layers: {circ.depth(filter_function=lambda instr: len(instr.qubits) == 2)}"
)
print("\nCircuit structure:")

circ.draw("mpl", scale=0.8, fold=-1, idle_wires=False)
plt.show()

observable = qiskit.quantum_info.SparsePauliOp.from_sparse_list(
    [("Z", [q], 1 / n_qubits) for q in subgraphs[n_qubits]],
    np.max(subgraphs[n_qubits]) + 1,
)  # Average magnetization observable

print(observable)
obs_list = [observable]

Circuit 2q layers: 27

Circuit structure:

SparsePauliOp(['IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIZIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIZIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIZIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIZIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIZIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIZIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIZIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIZIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIZIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'ZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIZIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII', 'IIIIIIIIIIIZIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIIII'],
              coeffs=[0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j,
 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j,
 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j,
 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j,
 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j,
 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j,
 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j,
 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j,
 0.02857143+0.j, 0.02857143+0.j, 0.02857143+0.j])

שלב 2: ייעול הבעיה לביצוע על חומרה קוונטית

הערכת זמן QPU עם ובלי OBP

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

QESEM מציעה שני מצבים של הערכת זמן כדי ליידע את המשתמשים על היתכנות הניסויים שלהם:

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

בשני המקרים, QESEM מפיקה את הערכת הזמן להגעה לדיוק הנדרש עבור כל הניתנים למדידה.

run_on_real_hardware = True

precision = 0.05
if run_on_real_hardware:
    backend_name = "ibm_fez"
else:
    backend_name = "fake_fez"

# Start a job for empirical time estimation
estimation_job_wo_obp = qesem_function.run(
    pubs=[(circ, obs_list)],
    instance=instance,
    backend_name=backend_name,  # E.g. "ibm_brisbane"
    options={
        "estimate_time_only": "empirical",  # "empirical" - gets actual time estimates without running full mitigation
        "max_execution_time": 120,  # Limits the QPU time, specified in seconds.
        "default_precision": precision,
    },
)

print(estimation_job_wo_obp.job_id)
print(estimation_job_wo_obp.status())

17d3828e-9fdb-482e-8e9b-392f3eefe313
DONE

# Get the result object (blocking method). Use job.status() in a loop for non-blocking.
# This takes 1-3 minutes
result = estimation_job_wo_obp.result()

print(
    f"Empirical time estimation (sec): {result[0].metadata['time_estimation_sec']}"
)

Empirical time estimation (sec): 1200

כעת נשתמש בהפצה לאחור של אופרטורים (OBP). (ראה את מדריך OBP לפרטים נוספים על תוסף Qiskit של OBP.) נצור פונקציה שמייצרת את פרוסות המעגל להפצה לאחור:

def run_backpropagation(circ_vec, observable, steps_vec, max_qwc_groups=8):
    """
    Runs backpropagation for a list of circuits and observables.
    Returns lists of backpropagated circuits and observables.
    """
    op_budget = OperatorBudget(max_qwc_groups=max_qwc_groups)
    bp_circuit_vec = []
    bp_observable_vec = []

    for i, circ in enumerate(circ_vec):
        slices = slice_by_gate_types(circ)
        bp_observable, remaining_slices, metadata = backpropagate(
            observable,
            slices,
            operator_budget=op_budget,
        )
        bp_circuit = combine_slices(remaining_slices, include_barriers=True)
        bp_circuit_vec.append(bp_circuit)
        bp_observable_vec.append(bp_observable)
        print(f"n.o. steps: {steps_vec[i]}")
        print(f"Backpropagated {metadata.num_backpropagated_slices} slices.")
        print(
            f"New observable has {len(bp_observable.paulis)} terms, which can be combined into "
            f"{len(bp_observable.group_commuting(qubit_wise=True))} groups.\n"
            f"After truncation, the error in our observable is bounded by {metadata.accumulated_error(0):.3e}"
        )
        print("-----------------")
    return bp_circuit_vec, bp_observable_vec

אנו קוראים לפונקציה:

bp_circ_vec, bp_obs_vec = run_backpropagation([circ], observable, [steps])

n.o. steps: 9
Backpropagated 11 slices.
New observable has 363 terms, which can be combined into 4 groups.
After truncation, the error in our observable is bounded by 0.000e+00
-----------------

print("The remaining circuit after backpropagation looks as follows:")
bp_circ_vec[-1].draw("mpl", scale=0.8, fold=-1, idle_wires=False)
None

The remaining circuit after backpropagation looks as follows:

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

# Start a job for empirical time estimation
estimation_job_obp = qesem_function.run(
    pubs=[(bp_circ_vec[-1], [bp_obs_vec[-1]])],
    instance=instance,
    backend_name=backend_name,
    options={
        "estimate_time_only": "empirical",
        "max_execution_time": 120,
        "default_precision": precision,
    },
)

print(estimation_job_obp.job_id)
print(estimation_job_obp.status())

8bae699d-a16b-4d39-bbd9-d123fbcce55d
DONE

result_obp = estimation_job_obp.result()

print(
    f"Empirical time estimation (sec): {result_obp[0].metadata['time_estimation_sec']}"
)

Empirical time estimation (sec): 900

אנו רואים ש-OBP מפחיתה את עלות הזמן להפחתה של המעגל.

שלב 3: ביצוע באמצעות פרימיטיבים של Qiskit

הרצה עם backend אמיתי

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

אנו מנתחים ארבעה מעגלים עם 5, 7 ו-9 שלבי Trotter בדיוק של 0.05, תוך השוואה בין ערכי התוחלת האידיאליים, הרועשים והמופחתי-שגיאות שלהם:

steps_vec = [5, 7, 9]

circ_vec = []
for steps in steps_vec:
    circ = trotter_circuit_from_layers(
        steps, theta_x, theta_z, theta_zz, layers
    )
    circ_vec.append(circ)

שוב, אנו מבצעים OBP על כל מעגל כדי להפחית את זמן הריצה:

bp_circ_vec_35, bp_obs_vec_35 = run_backpropagation(
    circ_vec, observable, steps_vec
)

n.o. steps: 5
Backpropagated 11 slices.
New observable has 363 terms, which can be combined into 4 groups.
After truncation, the error in our observable is bounded by 0.000e+00
-----------------
n.o. steps: 7
Backpropagated 11 slices.
New observable has 363 terms, which can be combined into 4 groups.
After truncation, the error in our observable is bounded by 0.000e+00
-----------------
n.o. steps: 9
Backpropagated 11 slices.
New observable has 363 terms, which can be combined into 4 groups.
After truncation, the error in our observable is bounded by 0.000e+00
-----------------

כעת אנו מריצים אצווה של עבודות QESEM מלאות. אנו מגבילים את זמן הריצה המקסימלי של QPU עבור כל אחת מהנקודות לשליטה טובה יותר בתקציב QPU.

run_on_real_hardware = True

precision = 0.05
if run_on_real_hardware:
    backend_name = "ibm_marrakesh"
else:
    backend_name = "fake_fez"

# Running full jobs for:
pubs_list = [
    [(bp_circ_vec_35[i], bp_obs_vec_35[i])] for i in range(len(bp_obs_vec_35))
]

# Initiating multiple jobs for different lengths
job_list = []
for pubs in pubs_list:
    job_obp = qesem_function.run(
        pubs=pubs,
        instance=instance,
        backend_name=backend_name,  # E.g. "ibm_brisbane"
        options={
            "max_execution_time": 300,  # Limits the QPU time, specified in seconds.
            "default_precision": 0.05,
        },
    )
    job_list.append(job_obp)

כאן אנו בודקים את הסטטוס של כל עבודה:

for job in job_list:
    print(job.status())

DONE
DONE
DONE
DONE

שלב 4: עיבוד-לאחר והחזרת תוצאה בפורמט קלאסי רצוי

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

ideal_values = []
noisy_values = []
error_mitigated_values = []
error_mitigated_stds = []

for i in range(len(job_list)):
    job = job_list[i]
    result = job.result()  # Blocking - takes 3-5 minutes
    noisy_results = result[0].metadata["noisy_results"]

    ideal_val = calculate_ideal_evs(circ_vec[i], observable, n_qubits, i)
    print("---------------------------------")
    print(f"Ideal: {ideal_val}")
    print(f"Noisy: {noisy_results.evs}")
    print(f"QESEM: {result[0].data.evs} \u00b1 {result[0].data.stds}")

    ideal_values.append(ideal_val)
    noisy_values.append(noisy_results.evs)
    error_mitigated_values.append(result[0].data.evs)
    error_mitigated_stds.append(result[0].data.stds)

Using precalculated ideal values for large circuits calculated with belief propagation PEPS. Currently only for 35 qubits.
---------------------------------
Ideal: 0.79537
Noisy: 0.7039237951821501
QESEM: 0.7828018244130982 ± 0.013257266977728376
Using precalculated ideal values for large circuits calculated with belief propagation PEPS. Currently only for 35 qubits.
---------------------------------
Ideal: 0.78653
Noisy: 0.6478583812958806
QESEM: 0.7875259197423828 ± 0.02703045139248604
Using precalculated ideal values for large circuits calculated with belief propagation PEPS. Currently only for 35 qubits.
---------------------------------
Ideal: 0.79699
Noisy: 0.6171787879868142
QESEM: 0.6918791909168913 ± 0.0740873782039517

לבסוף, אנו יכולים לשרטט את המגנטיזציה מול מספר השלבים. זה מסכם את היתרון של שימוש בפונקציית Qiskit של QESEM להפחתת שגיאות ללא הטיה על מכשירים קוונטיים רועשים.

plt.plot(steps_vec, ideal_values, "--", label="ideal")
plt.scatter(steps_vec, noisy_values, label="noisy")
plt.errorbar(
    steps_vec,
    error_mitigated_values,
    yerr=error_mitigated_stds,
    fmt="o",
    capsize=5,
    label="QESEM mitigation",
)
plt.legend()
plt.xlabel("n.o. steps")
plt.ylabel("Magnetization")

Text(0, 0.5, 'Magnetization')

הערה

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