Simulace 2D Isingova modelu s naklopeným polem pomocí funkce QESEM

poznámka

Qiskit Functions jsou experimentální funkce dostupné pouze uživatelům IBM Quantum® Premium Plan, Flex Plan a On-Prem (přes IBM Quantum Platform API) Plan. Jsou ve stavu náhledu a mohou se měnit.

Odhadovaná doba použití: 20 minut na procesoru Heron r2. (POZNÁMKA: Toto je pouze odhad. Skutečná doba běhu se může lišit.)

Pozadí

Tento tutoriál ukazuje, jak používat QESEM, Qiskit Function od Qedmy, k simulaci dynamiky kanonického kvantového spinového modelu – 2D Isingova modelu s naklopeným polem (TFI) s non-Cliffordovými úhly:

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

kde $\langle i,j \rangle$ označuje nejbližší sousedy na mřížce. Simulace časového vývoje mnohatělových kvantových systémů je pro klasické počítače výpočetně náročný úkol. Kvantové počítače jsou naopak přirozeně navrženy k efektivnímu provádění tohoto úkolu. TFI model se zejména stal populárním benchmarkem na kvantovém hardwaru díky svému bohatému fyzikálnímu chování a implementaci přátelské k hardwaru.

Místo simulace dynamiky v kontinuálním čase používáme úzce příbuzný model „kicked Ising". Dynamiku lze vyjádřit přesně jako periodický kvantový Circuit, kde každý krok vývoje se skládá ze tří vrstev zlomkových dvoQubitových Gate $R_{ZZ} (\alpha_{ZZ})$, proložených vrstvami jednoQubitových Gate $R_X (\alpha_X)$ a $R_Z (\alpha_Z)$.

Použijeme obecné úhly, které jsou náročné jak pro klasickou simulaci, tak pro chybovou mitigaci. Konkrétně jsme zvolili $\alpha_{ZZ} = 1.0$, $\alpha_X = 0.53$ a $\alpha_Z = 0.1$, čímž model umísťujeme daleko od jakéhokoli integrabilního bodu.

V tomto tutoriálu provedeme následující:

• Odhadneme očekávanou dobu běhu QPU pro plnou chybovou mitigaci pomocí analytických a empirických funkcí odhadování času QESEM.
• Sestavíme a nasimulujeme obvod 2D Isingova modelu s naklopeným polem s využitím rozložení Qubitů inspirovaných hardwarem a vrstev Gate.
• Vizualizujeme konektivitu Qubitů zařízení a vybrané podgrafy pro náš experiment.
• Ukážeme použití zpětné propagace operátorů (OBP) ke snížení hloubky obvodu. Tato technika ořezává operace z konce obvodu za cenu více měření operátorů.
• Provedeme nestrannou chybovou mitigaci (EM) pro více observabilů současně pomocí QESEM a porovnáme ideální, zašuměné a mitigované výsledky.
• Analyzujeme a vykreslíme dopad chybové mitigace na magnetizaci při různých hloubkách obvodů.

Poznámka: OBP obecně vrátí sadu pozorovatelných veličin, které spolu nemusí komutovat. QESEM automaticky optimalizuje měřicí báze, když cílové pozorovatelné veličiny obsahují nekomutující členy. Generuje kandidátní sady měřicích bází pomocí několika heuristických algoritmů a vybere sadu, která minimalizuje počet různých bází. To znamená, že QESEM seskupuje kompatibilní pozorovatelné veličiny do společných bází, aby snížil celkový počet potřebných konfigurací měření a zlepšil efektivitu.

O QESEM

QESEM je spolehlivý, vysoce přesný software založený na charakterizaci, implementující efektivní, nestrannou kvaziprobabilistickou chybovou mitigaci. Je navržen pro mitigaci chyb v obecných kvantových Circuit a je aplikačně nezávislý. Byl ověřen na různých hardwarových platformách, včetně experimentů v utility měřítku na zařízeních IBM® Eagle a Heron. Fáze pracovního postupu QESEM jsou následující:

1. Charakterizace zařízení – mapuje věrnosti Gate a identifikuje koherentní chyby, poskytuje kalibrační data v reálném čase. Tato fáze zajišťuje, že mitigace využívá operace s nejvyšší dostupnou věrností.

2. Transpilace zohledňující šum – generuje a vyhodnocuje alternativní mapování Qubitů, sady operací a měřicí báze, vybírá variantu minimalizující odhadovanou dobu běhu QPU s volitelnou paralelizací pro urychlení sběru dat.

3. Potlačení chyb – předefinuje nativní Gate, aplikuje Pauliho twirling a optimalizuje řízení na úrovni pulzů (na podporovaných platformách) pro zlepšení věrnosti.

4. Charakterizace obvodu – sestaví přizpůsobený lokální chybový model a přizpůsobí ho měřením QPU pro kvantifikaci reziduálního šumu.

5. Chybová mitigace – konstruuje kvaziprobabilistické dekompozice více typů a vzorkuje z nich v adaptivním procesu, který minimalizuje dobu QPU pro mitigaci a citlivost na výkyvy hardwaru, dosahuje vysoké přesnosti při velkých objemech obvodu.

Další informace o QESEM a experimentu v utility měřítku tohoto modelu na 103-Qubitovém podgrafu s vysokou konektivitou nativní geometrie heavy-hex zařízení ibm_marrakesh najdeš v Reliable high-accuracy error mitigation for utility-scale quantum circuits.

Požadavky

Před spuštěním notebooku nainstaluj následující Python balíčky:

• Qiskit SDK v2.0.0 nebo novější (pip install qiskit)

• Qiskit Runtime v0.40.0 nebo novější (pip install qiskit-ibm-runtime)

• Qiskit Functions Catalog v0.8.0 nebo novější ( pip install qiskit-ibm-catalog )

• Qiskit addon pro zpětnou propagaci operátorů v0.3.0 nebo novější ( pip install qiskit-addon-obp )

• Qiskit Utils addon v0.1.1 nebo novější ( pip install qiskit-addon-utils )

• Qiskit Aer simulátor v0.17.1 nebo novější ( pip install qiskit-aer )

• Matplotlib v3.10.3 nebo novější ( pip install matplotlib )

Nastavení

Nejprve importuj příslušné knihovny:

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

_needed = {"matplotlib": "matplotlib", "numpy": "numpy", "qiskit": "qiskit", "qiskit_addon_obp": "qiskit-addon-obp", "qiskit_addon_utils": "qiskit-addon-utils", "qiskit_aer": "qiskit-aer", "qiskit_ibm_catalog": "qiskit-ibm-catalog", "qiskit_ibm_runtime": "qiskit-ibm-runtime"}
_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")

%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.visualization import (
    plot_gate_map,
)

Dále se ověř pomocí svého API klíče z řídicího panelu IBM Quantum Platform. Pak vyber Qiskit Function takto. (Poznámka: z bezpečnostních důvodů je nejlepší uložit si přihlašovací údaje k účtu do místního prostředí, pokud jsi na důvěryhodném počítači, abys nemusel/a zadávat API klíč při každém ověření.)

# 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")

Krok 1: Mapování klasických vstupů na kvantový problém

Začneme definováním funkce, která vytvoří Trotterův Circuit:

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

Dále vytvoříme funkci pro výpočet ideálních očekávaných hodnot pomocí AerSimulator.

Poznámka: pro velké Circuit (30 a více Qubitů) doporučujeme použít předpočítané hodnoty ze simulací belief-propagation (BP) PEPS. Tento kód zahrnuje předpočítané hodnoty pro 35 Qubitů jako příklad, založené na přístupu BP pro vývoj tenzorové sítě PEPS představeném v tomto článku (který označujeme jako PEPS-BP), s využitím Python balíčku pro tenzorové sítě 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

Použijeme mapování vrstvy $R_{ZZ}$ na hardwarové bázi převzaté ze zařízení Heron, ze kterého vystřihneme vrstvy podle počtu Qubitů, které chceme simulovat. Definujeme podgrafy pro 10, 21, 28 a 35 Qubitů, které zachovávají 2D strukturu (neváhej změnit na svůj oblíbený podgraf):

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

Nyní vizualizujeme rozmístění Qubitů na zařízení Heron pro vybraný podgraf:

catalog = QiskitFunctionsCatalog(
    channel=channel,
    token=token,
    instance=instance,
)
backend = catalog.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()

Všimni si, že konektivita zvoleného rozmístění Qubitů nemusí být nutně lineární a může pokrývat velké oblasti zařízení Heron v závislosti na zvoleném počtu Qubitů.

Nyní vygenerujeme Trotterův Circuit a pozorovatelnou veličinu průměrné magnetizace pro zvolený počet Qubitů a parametry:

# 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: "
    f"{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])

Krok 2: Optimalizace problému pro spuštění na kvantovém hardwaru

Odhad času QPU s OBP a bez OBP

Uživatelé obvykle chtějí vědět, kolik času QPU jejich experiment vyžaduje. Jde však o problém, který je pro klasické počítače obtížně řešitelný.

QESEM nabízí dva režimy odhadu času, které tě informují o proveditelnosti tvých experimentů:

1. Analytický odhad času – poskytuje velmi hrubý odhad a nevyžaduje žádný čas QPU. Lze ho využít k testování, zda by transpilační průchod potenciálně snížil čas QPU.
2. Empirický odhad času (ukázaný zde) – poskytuje docela dobrý odhad a spotřebuje několik minut času QPU.

V obou případech QESEM vrátí odhad času potřebného k dosažení požadované přesnosti pro všechny pozorovatelné veličiny.

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={
        # "empirical" - gets actual time estimates without running full mitigation
        "estimate_time_only": "empirical",
        "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

Nyní použijeme zpětnou propagaci operátorů (OBP). (Více podrobností o doplňku OBP pro Qiskit najdeš v dokumentaci OBP.) Vytvoříme funkci, která generuje řezy obvodu pro zpětnou propagaci:

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, "
            f"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 "
            f"{metadata.accumulated_error(0):.3e}"
        )
        print("-----------------")
    return bp_circuit_vec, bp_observable_vec

Zavoláme tuto funkci:

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:

Vidíme, že zpětná propagace zkrátila Circuit o dvě vrstvy. Nyní, když máme zkrácený Circuit a rozšířené pozorovatelné veličiny, provedeme odhad času pro Circuit po zpětné propagaci:

# 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

Vidíme, že OBP snižuje časové náklady na mitigaci obvodu.

Krok 3: Spuštění pomocí Qiskit primitiv

Spuštění na skutečném backendu

Nyní spustíme celý experiment na několika Trotterových krocích. Počet qubitů, požadovanou přesnost a maximální čas QPU lze upravit podle dostupných zdrojů QPU. Všimni si, že omezení maximálního času QPU ovlivní výslednou přesnost, jak uvidíš v závěrečném grafu níže.

Analyzujeme čtyři Circuit s 5, 7 a 9 Trotterovými kroky při přesnosti 0,05 a porovnáme jejich ideální, zašuměné a chybově zmírněné očekávané hodnoty:

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)

Opět provádíme OBP na každém Circuit, abychom zkrátili dobu běhu:

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

Nyní spustíme dávku úplných úloh QESEM. Pro lepší kontrolu nad rozpočtem QPU omezíme maximální dobu běhu QPU pro každý z bodů.

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)

Zde zkontrolujeme stav každé úlohy:

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

DONE
DONE
DONE
DONE

Krok 4: Následné zpracování a vrácení výsledku v požadovaném klasickém formátu

Jakmile jsou všechny úlohy dokončeny, můžeme porovnat jejich zašuměné a mitigované očekávané hodnoty.

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

Nakonec můžeme vykreslit magnetizaci v závislosti na počtu kroků. To shrnuje přínos použití funkce QESEM Qiskit pro nezkreslené zmírňování chyb na zašuměných kvantových zařízeních.

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

poznámka

Devátý krok má velkou statistickou chybovou úsečku, protože jsme omezili čas QPU na 5 minut. Pokud tento krok spustíš po dobu 15 minut (jak naznačuje empirický odhad času), získáš menší chybovou úsečku. Zmírněná hodnota se tak přiblíží ideální hodnotě.