Saltar al contenido principal
Este proyecto está en fase beta. ¿Encontraste un error? Abre un issue en GitHub. ¿Tienes una pregunta o una idea? Pregunta en Discussions — ¡nos encantaría tu opinión!
¿Cómo es la calidad de la traducción en esta página?

Cancelación probabilística de errores con modelos de ruido lógico

Estimación de uso: 28 minutos en un procesador Heron r3 (NOTA: Esto es solo una estimación. Tu tiempo de ejecución puede variar.)

Hexagonal Ising lattice embedded on a heavy-hex qubit layout, with Ising sites on the degree-3 qubits and mediator qubits on the edges between them In this tutorial we will evaluate observables of a 22-site Ising model on a hexagonal lattice using a Heron QPU, which has a heavy-hexagonal qubit topology. Many demonstrations on Heron QPUs focus on systems defined on heavy-hexagonal lattices to match the system's connectivity; however, hexagonal models tend to be more interesting to study as they appear more commonly in nature and are more difficult to simulate classically due to their denser connectivity. Since the qubit topology of the QPU cannot directly support the connectivity of a hexagonal lattice, we must decide how to efficiently embed the hexagonal model onto the qubit topology. In this example we represent each site in the Ising model with a qubit lying on a vertex of the heavy-hex QPU lattice. We use the qubits on the edges (mediator qubits) to facilitate entanglement between the Ising qubits on the vertices. Additionally, the mediator qubits implement entanglement in such a way that they are always expected to return to the ground state ∣0⟩|0\rangle. If a mediator qubit measures ∣1⟩|1\rangle, it indicates that a logical error occurred during circuit execution; postselecting only samples with no detected errors on the mediator qubits yields a smaller distribution of logical samples that might have a higher fidelity than the raw, noisy distribution. Not only can we detect that an error occurred by measuring a mediator qubit, we can also determine exactly which logical errors that qubit is capable of detecting throughout the circuit. Removing the noise generators the mediator qubits' symmetry checks are capable of detecting from a noise model leaves a reduced noise model, which can be mitigated using techniques such as probabilistic error cancellation (PEC).

En este ejemplo de notebook, combinaremos la técnica de detección de errores descrita anteriormente con PEC para contrarrestar el ruido cuántico de forma más eficaz que cualquiera de las dos técnicas por separado. Usaremos 49 qubits para incorporar el modelo de Ising de 22 qubits y usaremos los 27 qubits mediadores adicionales para la detección de errores. Ejecutaremos PEC en el canal de ruido posseleccionado para mitigar el ruido que evade las comprobaciones de simetría. Además de combinar la detección de errores con PEC, usaremos técnicas como la mitigación de errores de lectura TREX y comprobaciones de errores no markovianos para combatir aún más el efecto del ruido cuántico.

El flujo de trabajo es el siguiente:

  1. Simular clásicamente los valores esperados de los observables para el modelo de Ising hexagonal de 22 qubits
  2. Implementar el modelo de Ising hexagonal detector de errores de 49 qubits con un circuito cuántico y transpilarlo al backend QPU
  3. Especificar las capas de entrelazamiento en el circuito con samplomatic. Aprenderemos y mitigaremos el ruido que afecta a estas capas de entrelazamiento.
  4. Agregar comprobaciones de errores no markovianos al circuito
  5. Aprender el ruido de puertas y de lectura que afecta al circuito
  6. Podar el modelo de ruido de los términos detectados por las comprobaciones de simetría
  7. Muestrear el circuito detector de errores. Posseleccionar únicamente las muestras que pasan todas las comprobaciones de simetría y de errores no markovianos, y usar la mitigación de lectura TREX para todos los cálculos de valores esperados
    • Muestrear el circuito detector de errores
    • Muestrear el circuito detector de errores con PEC, pero sin posseleccionar según las comprobaciones de simetría y mitigando todo el canal de ruido aprendido.
    • Muestrear el circuito detector de errores con PEC. Posseleccionar según las comprobaciones de simetría y mitigar únicamente el ruido que las comprobaciones no pueden detectar.
  8. Calcular los valores esperados y comparar las estrategias. Observar que QED+PEC cancela el ruido que afecta al experimento de forma más eficaz que cualquiera de los métodos por separado y produce valores esperados convergidos usando muchos menos shots de los que requeriría PEC solo.

Requisitos​

Antes de comenzar este tutorial, ejecuta el siguiente comando para obtener todos los paquetes necesarios:

%pip install networkx numpy "qiskit[visualization]" "qiskit-ibm-runtime[visualization]" samplomatic qiskit-mitigation matplotlib "qiskit-noise-learning @ git+https://github.com/Qiskit/qiskit-noise-learning.git@ieee-demo-2026" La celda colapsada de abajo define las funciones auxiliares de figuras usadas a lo largo de este notebook.

Haz clic para ver las funciones auxiliares de figuras
# Added by doQumentation — installs the packages this notebook needs if they are missing
import importlib.util

_needed = {"matplotlib": "matplotlib", "networkx": "networkx", "numpy": "numpy", "qiskit": "qiskit", "qiskit_ibm_runtime": "qiskit-ibm-runtime", "qiskit_mitigation": "qiskit-mitigation", "qiskit_noise_learning": "qiskit-noise-learning@git+https://github.com/Qiskit/qiskit-noise-learning.git@ieee-demo-2026", "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")
# Figure helpers for this notebook (collapsed): all plotting and styling lives here.
# The chip maps follow the style of Fig. 30 of arXiv:2607.25998.

from math import comb

import matplotlib.pyplot as plt
from matplotlib import cm
from matplotlib.colors import BoundaryNorm
from matplotlib.patches import Circle, Rectangle, Wedge

# one typography scheme for every figure in the notebook
plt.rcParams.update(
    {
        "font.size": 13,
        "axes.titlesize": 15,
        "axes.labelsize": 13,
        "xtick.labelsize": 11,
        "ytick.labelsize": 11,
        "legend.fontsize": 11,
    }
)

def x_labels(layout, n_data):
    """Per-observable tick labels carrying the physical qubit indices (site i = qubit i)."""
    return [f"$X_{{{layout[i]}}}$" for i in range(n_data)]

def plot_layout(backend, layout, n_data):
    """Chip cartoon of the embedding: data qubits green, check qubits orange."""
    from qiskit.visualization import plot_coupling_map

    return plot_coupling_map(
        num_qubits=backend.num_qubits,
        qubit_coordinates=None,
        coupling_map=list(backend.coupling_map.get_edges()),
        figsize=(9, 9),
        qubit_color=[
            "#4CAF50"
            if q in layout[:n_data]
            else "#FF9800"
            if q in layout[n_data:]
            else "#DDDDDD"
            for q in range(backend.num_qubits)
        ],
        qubit_size=220,
        line_width=2,
        font_size=90,
    )

def plot_exact(obs_exact, tick_labels, title):
    plt.figure(figsize=(12, 4))
    plt.plot(obs_exact, "o-")
    plt.title(title)
    plt.xticks(np.arange(len(obs_exact)), tick_labels)
    plt.xlabel("Observable")
    plt.ylabel(r"$\langle X \rangle$")
    plt.grid()

def draw_toy_circuit(
    generate_ed_ising, zz_coeff, x_coeff, include_checks=True
):
    """The boxing pipeline on a 3-plaquette-ring miniature (1 Trotter step) so the box
    structure is legible; every box carries the Twirl / InjectNoise annotations, and with
    ``include_checks`` the terminal xslow non-Markovian error check pattern is appended as in the
    production pipeline."""
    import networkx as nx
    from qiskit.circuit import ClassicalRegister
    from qiskit_mitigation.postselection import XSlowGate
    from samplomatic.transpiler import generate_boxing_pass_manager

    toy, _, _ = generate_ed_ising(nx.cycle_graph(3), 1, zz_coeff, x_coeff)
    toy.add_register(
        ClassicalRegister(3, "data"), ClassicalRegister(3, "check")
    )
    toy.barrier()
    toy.measure(range(3), range(3))
    toy.measure(range(3, 6), range(3, 6))
    toy_boxed = generate_boxing_pass_manager(
        enable_gates=True,
        enable_measures=True,
        inject_noise_targets="gates",
        inject_noise_strategy="individual_modification",
        inject_noise_site="after",
        twirling_strategy="active_circuit",
        measure_annotations="all",
    ).run(toy)
    if include_checks:
        toy_boxed.add_register(
            ClassicalRegister(3, "data_ps"), ClassicalRegister(3, "check_ps")
        )
        toy_boxed.barrier()
        for qb in range(6):
            toy_boxed.append(XSlowGate(), [qb])
        toy_boxed.measure(range(3), toy_boxed.cregs[2])
        toy_boxed.measure(range(3, 6), toy_boxed.cregs[3])
    return toy_boxed.draw("mpl", fold=-1, scale=0.6)

# --- chip-level noise maps ---------------------------------------------------------

BANDS = [1e-5, 2e-5, 3e-5, 4e-5, 6e-5, 1e-4, 2e-4, 3e-4, 4e-4, 6e-4, 1e-3]
CMAP = plt.get_cmap("YlOrBr")
NORM = BoundaryNorm(BANDS, CMAP.N, extend="both")

def _color(rate):
    return (
        "white"
        if rate < BANDS[0]
        else CMAP(NORM(min(rate, BANDS[-1] * 0.999)))
    )

def _layer_sparse(mit, weights):
    """Per-layer labelled terms from the saved run, box-local -> physical qubits."""
    phys = np.sort(mit["layout"])
    return [
        [
            (p, tuple(int(phys[q]) for q in qs if q >= 0), r)
            for p, qs, r in zip(
                mit[f"label_paulis_{i}"],
                mit[f"label_qubits_{i}"],
                w,
                strict=True,
            )
        ]
        for i, w in enumerate(weights)
    ]

def _aggregate(layers):
    """w1[qubit][P] and w2[(a,b)][PaPb]: rates summed over the 3 layers."""
    w1, w2 = {}, {}
    for terms in layers:
        for p, qs, r in terms:
            if len(qs) == 1:
                w1.setdefault(qs[0], dict.fromkeys("XYZ", 0.0))[p] += r
            else:
                (a, pa), (b, pb) = sorted(zip(qs, p, strict=True))
                w2.setdefault(
                    (a, b), {x + y: 0.0 for x in "XYZ" for y in "XYZ"}
                )[pa + pb] += r
    return w1, w2

def draw_noise_map(mit, backend, reduced=False):
    """Chip map of the learned model: X/Y/Z wheel per qubit, 3x3 two-qubit Pauli grid per
    coupler, log-banded colors, dashed outlines for unused hardware.

    With ``reduced=True``, keeps only the error terms the checks cannot see: detectability
    depends on circuit position, so each layer's 0/1 site scales are averaged over its uses.
    """
    from qiskit_ibm_runtime.visualization.embeddings import Embedding

    if reduced:
        scales = [
            mit["site_scales"][mit["site_layer"] == i].mean(axis=0)
            for i in range(3)
        ]
        weights = [mit[f"rates_{i}"] * scales[i] for i in range(3)]
        title = f"Reduced noise model ($\\gamma$ = {mit['gammas'][1]:.1f})"
    else:
        weights = [mit[f"rates_{i}"] for i in range(3)]
        title = f"Full noise model ($\\gamma$ = {mit['gammas'][0]:.0f})"
    w1, w2 = _aggregate(_layer_sparse(mit, weights))

    xy = np.array(
        [(c, -r) for r, c in Embedding.from_backend(backend).coordinates]
    )
    fig, ax = plt.subplots(figsize=(13, 7))
    for a, b in {tuple(sorted(e)) for e in backend.coupling_map.get_edges()}:
        if (
            (
                a,
                b,
            )
            in w2
        ):  # 3x3 Pauli grid laid along the bond (columns: qubit a, rows: qubit b)
            d = xy[b] - xy[a]
            u = d / np.hypot(*d)
            v = np.array([-u[1], u[0]])
            cl, cw = (np.hypot(*d) - 0.6) / 3, 0.17
            for i, Pa in enumerate("XYZ"):
                for j, Pb in enumerate("XYZ"):
                    ax.add_patch(
                        Rectangle(
                            xy[a] + u * (0.3 + i * cl) + v * ((j - 1.5) * cw),
                            cl,
                            cw,
                            angle=np.degrees(np.arctan2(u[1], u[0])),
                            facecolor=_color(w2[a, b][Pa + Pb]),
                            edgecolor="black",
                            lw=0.4,
                            zorder=2,
                        )
                    )
        else:
            ax.plot(
                *zip(xy[a], xy[b], strict=True),
                ls="--",
                lw=0.7,
                color="black",
                alpha=0.5,
                zorder=1,
            )
    for q, (xq, yq) in enumerate(xy):
        if q in w1:  # three-sector wheel: X top, Y lower left, Z lower right
            for P, t0 in (("X", 30), ("Y", 150), ("Z", 270)):
                ax.add_patch(
                    Wedge(
                        (xq, yq),
                        0.3,
                        t0,
                        t0 + 120,
                        facecolor=_color(w1[q][P]),
                        edgecolor="black",
                        lw=0.6,
                        zorder=3,
                    )
                )
            ax.annotate(
                str(q),
                (xq + 0.39, yq - 0.39),
                fontsize=9,
                color="gray",
                ha="left",
                va="top",
                zorder=4,
            )  # southeast, clear of the bond grids
        else:
            ax.add_patch(
                Circle(
                    (xq, yq),
                    0.24,
                    facecolor="none",
                    edgecolor="black",
                    ls="--",
                    lw=0.7,
                    alpha=0.5,
                    zorder=3,
                )
            )
    ux = xy[sorted(w1)]
    ax.set_xlim(ux[:, 0].min() - 2.2, ux[:, 0].max() + 2.2)
    ax.set_ylim(ux[:, 1].min() - 1.6, ux[:, 1].max() + 1.6)
    ax.set_aspect("equal")
    ax.axis("off")
    ax.set_title(title, fontsize=15)
    cb = fig.colorbar(
        cm.ScalarMappable(norm=NORM, cmap=CMAP),
        ax=ax,
        fraction=0.035,
        pad=0.02,
        extend="both",
        ticks=BANDS,
    )
    cb.set_label("coefficient (log bands; white < 1e-05)", fontsize=11)
    cb.ax.set_yticklabels(
        [f"{b:.0e}".replace("e-0", "e-") for b in BANDS], fontsize=10
    )
    _noise_map_legends(fig, ax)

def _noise_map_legends(fig, ax):
    """Weight-1 wheel and weight-2 grid keys, in a reserved band left of the lattice."""
    fig.subplots_adjust(left=0.17)
    axl = ax.inset_axes([-0.185, 0.70, 0.13, 0.24])
    for P, t0 in (("X", 30), ("Y", 150), ("Z", 270)):
        axl.add_patch(
            Wedge(
                (0.5, 0.45),
                0.38,
                t0,
                t0 + 120,
                facecolor="white",
                edgecolor="black",
                lw=0.8,
            )
        )
        axl.annotate(
            P,
            (
                0.5 + 0.2 * np.cos(np.radians(t0 + 60)),
                0.45 + 0.2 * np.sin(np.radians(t0 + 60)),
            ),
            ha="center",
            va="center",
            fontsize=9,
        )
    axl.set_title("weight-1", fontsize=10)
    axl.set_xlim(0, 1)
    axl.set_ylim(0, 1)
    axl.set_aspect("equal")
    axl.axis("off")
    axm = ax.inset_axes([-0.185, 0.32, 0.14, 0.30])
    for i, Pa in enumerate("XYZ"):
        for j, Pb in enumerate("XYZ"):
            axm.add_patch(
                Rectangle(
                    (i / 3, 1 - (j + 1) / 3),
                    1 / 3,
                    1 / 3,
                    facecolor="white",
                    edgecolor="black",
                    lw=0.6,
                )
            )
            axm.annotate(
                Pa + Pb,
                ((i + 0.5) / 3, 1 - (j + 0.5) / 3),
                ha="center",
                va="center",
                fontsize=7.5,
            )
        axm.annotate(Pa, ((i + 0.5) / 3, 1.08), ha="center", fontsize=8.5)
        axm.annotate(
            "XYZ"[i],
            (-0.13, 1 - (i + 0.5) / 3),
            ha="center",
            va="center",
            fontsize=8.5,
        )
    axm.annotate("qubit a", (0.5, 1.27), ha="center", fontsize=9)
    axm.annotate(
        "qubit b", (-0.33, 0.5), rotation=90, va="center", fontsize=9
    )
    axm.set_xlim(-0.35, 1.05)
    axm.set_ylim(-0.05, 1.35)
    axm.set_aspect("equal")
    axm.axis("off")

# --- run diagnostics ---------------------------------------------------------------

def plot_trex(trex_rescale, tick_labels):
    _fig, axt = plt.subplots(figsize=(12, 3))
    axt.stem(np.arange(len(trex_rescale)), (trex_rescale - 1) * 100)
    axt.set_xticks(np.arange(len(trex_rescale)), tick_labels)
    axt.set_xlabel("Observable")
    axt.set_ylabel("Readout correction (%)")
    axt.set_title("TREX rescale factors")
    plt.tight_layout()

def plot_postselection(mit):
    """Accepted shots per PEC randomization against a single binomial at the mean
    acceptance rate: agreement means acceptance is independent of the sampled circuit
    instance, the condition under which pooling accepted shots across randomizations
    is a consistent estimator."""
    _fig, ax = plt.subplots(figsize=(7.5, 3.8))
    counts = mit["acc_counts_post"]
    K, p = 64, counts.mean() / 64
    ks = np.arange(K + 1)
    pmf = np.array([comb(K, k) * p**k * (1 - p) ** (K - k) for k in ks])
    ax.hist(
        counts,
        bins=np.arange(-0.5, K + 1.5),
        density=True,
        alpha=0.6,
        color="#da1e28",
        label="measured",
    )
    ax.plot(ks, pmf, "k-", lw=1.5, label=f"Binomial(64, {p:.3f})")
    ax.set_xlim(-0.5, max(int(counts.max()) + 3, 20))
    ax.set_xlabel("Accepted shots per randomization")
    ax.set_ylabel("Probability")
    ax.legend()
    ax.grid(alpha=0.4)
    plt.tight_layout()

def plot_convergence(mit, obs_exact, n_data):
    """Running site-averaged estimate vs. randomizations for both PEC arms (the S5-consistent
    signed-ratio estimator, evaluated on growing prefixes of the sweep)."""

    def running(prefix):
        bits = np.squeeze(
            np.unpackbits(mit[f"data_{prefix}"], axis=-1)[..., :n_data]
            ^ mit[f"flips_{prefix}"]
        )
        mask = np.squeeze(mit[f"mask_{prefix}"])
        signs = 1 - 2 * (np.squeeze(mit[f"signs_{prefix}"]).sum(axis=-1) % 2)
        qv = (1 - 2 * bits.astype(int)) * mit["trex_rescale"]
        u = (
            (signs[:, None, None] * mask[..., None] * qv)
            .sum(axis=1)
            .mean(axis=1)
        )
        v = signs * mask.sum(axis=1)
        Rs = np.arange(500, len(u) + 1, 500)
        est, err = [], []
        for R in Rs:
            e = u[:R].sum() / v[:R].sum()
            est.append(e)
            err.append(
                np.sqrt(((u[:R] - e * v[:R]) ** 2).sum()) / abs(v[:R].sum())
            )
        return Rs, np.array(est), np.array(err)

    _fig, ax = plt.subplots(figsize=(12, 5))
    ideal_avg = float(np.mean(obs_exact))
    ax.axhline(ideal_avg, color="black", label="ideal")
    ax.fill_between(
        [-400, 24000],
        ideal_avg - 0.025,
        ideal_avg + 0.025,
        color="grey",
        alpha=0.22,
        label=r"$\pm 0.025$",
    )
    for prefix, label, color in (
        ("van", "vanilla PEC", "#8a3ffc"),
        ("post", "PEC + error detection", "#da1e28"),
    ):
        Rs, est, err = running(prefix)
        ax.errorbar(
            Rs,
            est,
            yerr=err,
            marker="o",
            linestyle="",
            markerfacecolor="none",
            color=color,
            alpha=0.85,
            capsize=3,
            label=label,
        )
    ax.set_xlim(-400, 24000)
    ax.set_ylim(ideal_avg - 0.08, ideal_avg + 0.08)
    ax.set_xlabel("# randomizations")
    ax.set_ylabel(r"Site-averaged $\langle X \rangle$")
    ax.legend(ncols=2)

def plot_final(
    obs_exact, baseline, ed, pec, post, gammas, tick_labels, title
):
    """Per-site <X> for every method, with an rms-deviation inset.
    baseline/ed/pec/post are (values, errors) pairs; gammas is (gamma, gamma_post)."""
    x = np.arange(len(obs_exact))
    _fig, ax = plt.subplots(figsize=(13, 5))
    h_ideal = ax.errorbar(
        x, obs_exact, fmt="-", capsize=4, label="ideal", color="black"
    )
    h_base = ax.errorbar(
        x,
        baseline[0],
        yerr=baseline[1],
        fmt=".--",
        capsize=4,
        label="baseline",
        color="#0f62fe",
    )
    h_ed = ax.errorbar(
        x,
        ed[0],
        yerr=ed[1],
        fmt="^",
        capsize=4,
        label="error detection",
        color="#009d9a",
        alpha=0.8,
    )
    h_pec = ax.errorbar(
        x,
        pec[0],
        yerr=pec[1],
        fmt="x",
        capsize=4,
        label=f"PEC ($\\gamma$={gammas[0]:.0f})",
        color="#8a3ffc",
        alpha=0.7,
    )
    h_post = ax.errorbar(
        x,
        post[0],
        yerr=post[1],
        fmt="d",
        capsize=4,
        label=f"PEC + error detection ($\\gamma$={gammas[1]:.0f})",
        color="#da1e28",
    )
    ax.set_xticks(x, tick_labels)
    ax.set_xlabel("Observable")
    ax.set_ylabel("Expectation value")
    # Set the y-axis range using the ideal, baseline, error detection, and PEC + error
    # detection curves, including error bars. Exclude PEC-only values when setting the
    # range so large fluctuations do not make the other curves difficult to distinguish.
    # PEC-only values may fall outside the visible range. Leave room above for the inset.
    series = [np.asarray(obs_exact)] + [
        np.asarray(v) + s * np.asarray(e)
        for v, e in (baseline, ed, post)
        for s in (-1, 1)
    ]
    lo, hi = min(a.min() for a in series), max(a.max() for a in series)
    span = max(hi - lo, 0.1)
    ax.set_ylim(lo - 0.1 * span, hi + 1.05 * span)
    # legend ordered to match the curves' vertical positions in the chart
    ax.legend(
        handles=[h_ideal, h_pec, h_post, h_ed, h_base],
        ncols=5,
        loc="lower center",
        bbox_to_anchor=(0.5, 1.02),
        frameon=False,
    )
    ax.set_title(title, pad=44)

    # inset: rows bottom-to-top so it reads top-to-bottom: PEC, PEC+ED, ED, baseline
    axi = ax.inset_axes([0.36, 0.68, 0.28, 0.26])
    methods = [
        ("baseline", baseline[0], "#0f62fe"),
        ("QED", ed[0], "#009d9a"),
        ("PEC+QED", post[0], "#da1e28"),
        ("PEC", pec[0], "#8a3ffc"),
    ]
    for k, (_nm, vals, color) in enumerate(methods):
        axi.barh(
            k,
            np.sqrt(np.mean((vals - np.array(obs_exact)) ** 2)),
            color=color,
            alpha=0.9,
        )
    axi.set_yticks(range(4), [m[0] for m in methods], fontsize=10)
    axi.set_title("RMS deviation from ideal", fontsize=11)
    axi.tick_params(labelsize=10)
    axi.patch.set_alpha(1.0)
    axi.set_zorder(5)

Simular clásicamente ⟨X⟩i\langle X \rangle_i para cada sitio en el modelo de Ising hexagonal de 22 qubits​

Primero, usamos la clase Statevector de Qiskit para simular los valores objetivo exactos de nuestro experimento. Aunque el problema de 22 qubits que estamos resolviendo es manejable clásicamente, el simulador de vector de estado de Qiskit no escalará a demostraciones mucho más grandes que esta, y para escalar más allá de unos 50 qubits, necesitaríamos usar técnicas de simulación no exactas, como la propagación de Pauli. Este experimento de 49 qubits nos permite investigar la combinación de la detección de errores con la mitigación de errores a escalas mayores, manteniendo al mismo tiempo el acceso a los valores esperados ideales.

import warnings

import networkx as nx
import numpy as np
from qiskit.circuit import QuantumCircuit
from qiskit.quantum_info import Pauli, Statevector

# Silence two harmless upstream warnings (a Samplomatic default-change notice and a
# numpy datetime timezone notice from the noise-learning circuit generator)
warnings.filterwarnings(
    "ignore", message="The default of the 'inject_noise_site'"
)
warnings.filterwarnings(
    "ignore", message="no explicit representation of timezones"
)

# Hexagonal lattice
data_graph = nx.convert_node_labels_to_integers(
    nx.hexagonal_lattice_graph(2, 3)
)

# Number of Trotter steps and rotation angles
depth = 4
zz_coeff = -np.pi / 4
x_coeff = 3 * np.pi / 8
n_data = data_graph.order()
n_checks = data_graph.size()
n_qubits = n_data + n_checks

qc_data = QuantumCircuit(n_data)
qc_data.h(range(n_data))
for _ in range(depth):
    for edge in data_graph.edges:
        qc_data.rzz(zz_coeff, *edge)
    qc_data.rx(x_coeff, range(n_data))
psi_exact = Statevector(qc_data)

# X on site i = qubit i
observables = ["I" * (n_data - 1 - i) + "X" + "I" * i for i in range(n_data)]
obs_exact = [psi_exact.expectation_value(Pauli(o)).real for o in observables]

plot_exact(
    obs_exact,
    x_labels(range(n_data), n_data),
    f"Statevector simulation, {n_data} qubit hex-Ising, {depth} Trotter steps",
)

Output of the previous code cell

Implementar el modelo de Ising hexagonal detector de errores de 49 qubits con un circuito cuántico y transpilarlo al backend QPU​

Este circuito simula 4 pasos de Trotter de un Hamiltoniano de Ising de campo transversal de 22 qubits en una red hexagonal. Los 22 qubits de datos están incorporados en un subgrafo heavy-hex de 49 qubits de la QPU Heron r3, ibm_boston. Los 27 qubits auxiliares (ancilla) se usan con dos fines: mediar el entrelazamiento entre los qubits de Ising y detectar errores lógicos durante la ejecución del circuito. Las capas de entrelazamiento del circuito se implementan de tal manera que los qubits mediadores regresan al estado fundamental ∣0⟩|0\rangle antes de sus mediciones finales. Uno o más qubits mediadores que midan ∣1⟩|1\rangle indican que la muestra se corrompió por un error lógico. Descartar estas muestras puede mejorar la fidelidad de la distribución muestreada.

En el gráfico de abajo, los qubits verdes representan los 22 qubits de Ising, y los qubits naranjas representan los 27 qubits mediadores.

from qiskit.circuit import ClassicalRegister
from qiskit.transpiler import generate_preset_pass_manager
from qiskit_ibm_runtime import QiskitRuntimeService

backend_name = "ibm_boston"

def generate_ed_ising(graph, depth, zz_coeff, x_coeff):
    """Build the mediated error-detecting Ising circuit for data-qubit `graph`; returns (circuit, CZ layers, hardware graph)."""
    hw_graph = nx.Graph()
    for i, (a, b) in enumerate(graph.edges()):
        hw_graph.add_edges_from(
            [(a, i + graph.order()), (b, i + graph.order())]
        )
    coloring = nx.coloring.greedy_color(
        nx.line_graph(hw_graph), strategy="DSATUR"
    )
    layers_coupling = [
        [e for e, c in coloring.items() if c == i]
        for i in set(coloring.values())
    ]
    circuit = QuantumCircuit(hw_graph.order())
    circuit.h(range(hw_graph.order()))
    for _ in range(depth):
        for angle, qubits in (
            (zz_coeff, range(graph.order(), hw_graph.order())),
            (x_coeff, range(graph.order())),
        ):
            circuit.barrier()
            for layer in layers_coupling:
                for edge in layer:
                    circuit.cz(*edge)
            circuit.barrier()
            circuit.rx(angle, qubits)
    circuit.barrier()
    circuit.h(range(graph.order(), hw_graph.order()))
    return circuit, layers_coupling, hw_graph

service = QiskitRuntimeService()
backend = service.backend(backend_name)

circuit, layers_coupling, hw_graph = generate_ed_ising(
    data_graph, depth, zz_coeff, x_coeff
)

circ_meas = circuit.copy()
data_reg, check_reg = (
    ClassicalRegister(n_data, "data"),
    ClassicalRegister(n_checks, "check"),
)
circ_meas.add_register(data_reg, check_reg)
circ_meas.barrier()
circ_meas.measure(range(n_data), data_reg)
circ_meas.measure(range(n_data, n_qubits), check_reg)

# Physical qubits hosting the 22 Ising sites on ibm_boston, one heavy-hex row per lattice chain;
# the mediator for each edge is the physical qubit sitting between its two data qubits
data_layout = [
    *[95, 93, 91, 89, 87],
    *[115, 113, 111, 109, 107, 105],
    *[135, 133, 131, 129, 127, 125],
    *[155, 153, 151, 149, 147],
]
coupling_graph = nx.Graph(list(backend.coupling_map.get_edges()))
layout = data_layout + [
    next(
        iter(
            nx.common_neighbors(
                coupling_graph, data_layout[a], data_layout[b]
            )
        )
    )
    for a, b in data_graph.edges()
]

circ_trans = generate_preset_pass_manager(
    backend=backend, optimization_level=1, initial_layout=layout
).run(circ_meas)

plot_layout(backend, layout, n_data)

Output of the previous code cell

Especificar las capas de entrelazamiento en el circuito con samplomatic.​

El paso generate_boxing_pass_manager de Samplomatic agrupa las capas de entrelazamiento del circuito en cajas anotadas para el twirling de Pauli y la inyección de ruido de PEC. El circuito plantilla resultante y el samplex (una distribución paramétrica sobre el circuito plantilla que especifica sus twirls e inyección de ruido) se usan para definir y ejecutar todos los experimentos de aprendizaje de ruido y muestreo en la QPU. Las cajas de medición también llevan una anotación ChangeBasis, de modo que las rotaciones que miden los qubits de datos en la base X se proporcionan al samplex como una entrada en el momento del muestreo, en lugar de estar integradas en el circuito.

A continuación, se aplica el mismo paso de boxing a una versión diminuta del circuito de Ising detector de errores para que podamos visualizar la estructura de un circuito en cajas.

from samplomatic.builders import build
from samplomatic.transpiler import generate_boxing_pass_manager
from samplomatic.utils import find_unique_box_instructions

boxed = generate_boxing_pass_manager(
    enable_gates=True,
    enable_measures=True,
    inject_noise_targets="gates",
    inject_noise_strategy="individual_modification",
    inject_noise_site="after",
    twirling_strategy="active_circuit",
    measure_annotations="all",
).run(circ_trans)
template, samplex = build(boxed)
unique_instructions = find_unique_box_instructions(
    boxed, normalize_annotations=None, undress_boxes=True
)

# Measure X on the data qubits and Z on the mediators
basis_key = next(
    s.name
    for s in samplex.inputs().get_specs()
    if s.name.startswith("basis_changes.")
)
meas_basis = np.array(
    [2 if q in set(layout[:n_data]) else 1 for q in sorted(layout)],
    dtype=np.uint8,
)

draw_toy_circuit(generate_ed_ising, zz_coeff, x_coeff, include_checks=False)

Output of the previous code cell

Agregar comprobaciones de errores no markovianos al circuito​

A continuación, agregamos comprobaciones de errores no markovianos al circuito para protegernos contra el ruido que no está modelado en nuestro protocolo de aprendizaje. Estas comprobaciones funcionan implementando una inversión de bit de pulso largo y verificando que el qubit se movió correctamente de un estado clásico a otro. Si ambos qubits de una arista fallan la comprobación, la muestra se descarta. Las comprobaciones de errores no markovianos se pueden usar al principio o al final del circuito (o en ambos); aquí las usamos solo al final. También colocamos comprobaciones en qubits "espectadores" no utilizados adyacentes al circuito, lo que proporciona más cobertura contra el ruido que estas comprobaciones están diseñadas para detectar.

Recuerda que, de forma similar a las comprobaciones de simetría, esto también es una forma de posselección, por lo que debemos asegurarnos de tener en cuenta la caída en la tasa de posselección al combinar estas técnicas. Concretamente, si P(non-Markovian error)=αP(\text{non-Markovian error}) = \alpha y P(symmetry error)=βP(\text{symmetry error}) = \beta, entonces necesitas aproximadamente 1(1−α)(1−β)\frac{1}{(1-\alpha)(1-\beta)} muestras ruidosas para recuperar una muestra lógica.

from qiskit.transpiler import PassManager
from qiskit_mitigation.postselection import PostSelector
from qiskit_mitigation.postselection.passes import (
    AddPostCircuitNonMarkovianErrorChecks,
    AddSpectatorPostCircuitNonMarkovianErrorChecks,
)

add_checks = PassManager(
    [
        AddPostCircuitNonMarkovianErrorChecks(x_pulse_type="xslow"),
        AddSpectatorPostCircuitNonMarkovianErrorChecks(
            backend.coupling_map, x_pulse_type="xslow"
        ),
    ]
)
template_checked = add_checks.run(template)
selector = PostSelector.from_circuit(template_checked, backend.coupling_map)

draw_toy_circuit(generate_ed_ising, zz_coeff, x_coeff, include_checks=True)

Output of the previous code cell

Aprender el ruido de puertas y de lectura para el circuito de Ising verificado de 49 qubits​

Para mitigar el ruido, debemos modelar cómo afecta el ruido a nuestras puertas de entrelazamiento. Aquí aprendemos un modelo de ruido de Pauli-Lindblad para cada una de las tres capas de entrelazamiento únicas del circuito usando qiskit-noise-learning. Más adelante, podaremos este modelo de ruido de los generadores de errores que son detectables por las comprobaciones de simetría y mitigaremos con PEC únicamente el canal de ruido restante. El mapa de abajo muestra las tres capas aprendidas en la QPU; cada qubit es una rueda de sus tasas X/Y/Z de peso 1, y cada acoplador es una cuadrícula de 3×3 de las tasas de Pauli de peso 2 en ese par. El valor de γ\gamma de 1010 implica un sobrecosto de muestreo de γ2=102=100\gamma^2=10^2=100.

Las correcciones de lectura (TREX) provienen directamente de las rutas SPAM del ajuste de aprendizaje de ruido. El gráfico de tallo debajo del mapa de ruido de la QPU muestra los factores de reescalado TREX por observable.

from qiskit.quantum_info import PauliLindbladMap, QubitSparsePauli
from qiskit_ibm_runtime import Executor, Session
from qiskit_ibm_runtime.quantum_program import QuantumProgram
from qiskit_mitigation.trex import TREX
from qiskit_noise_learning.analysis import (
    ComputeObservables,
    CurveFitObservables,
    FlipPostSelect,
    LeastSquaresSolve,
)
from qiskit_noise_learning.circuit_generator import ExecutorCircuitGenerator
from qiskit_noise_learning.experiment_builder import (
    BindFragmentDepths,
    CompleteSequences,
    EvenDepthVanillaPaths,
    Experiment,
    GenerateInstructionSequences,
    IdentifyRelations,
    MergeInstructionSequences,
    SPAMPaths,
    VanillaInstructionSequences,
)
from qiskit_noise_learning.gate_sets import QiskitGateSet
from qiskit_noise_learning.models import PauliLindbladModel
from qiskit_noise_learning.models.utils import split_pauli_lindblad_model
from samplomatic.annotations import InjectNoise
from samplomatic.utils import get_annotation

gate_set = QiskitGateSet(
    target=backend.target,
    qubit_subset=sorted(
        {
            boxed.find_bit(q).index
            for i in unique_instructions
            for q in i.qubits
        }
    ),
)
ref_to_qubits = {}
for inst in unique_instructions:
    if ann := get_annotation(inst.operation, InjectNoise):
        gate_set.add_box_as_gate(inst, name=ann.ref)
        ref_to_qubits[ann.ref] = sorted(
            boxed.find_bit(q).index for q in inst.qubits
        )

fidelity_model = PauliLindbladModel.k_local(
    gate_set, gate_k={**{r: 2 for r in ref_to_qubits}, "M": 1, "P": 1}
)
experiment = (
    EvenDepthVanillaPaths()
    + VanillaInstructionSequences()
    + IdentifyRelations()
    + SPAMPaths()
    + GenerateInstructionSequences()
    + MergeInstructionSequences()
    + CompleteSequences()
    + BindFragmentDepths([2, 4, 8, 12])
).run(Experiment(fidelity_model=fidelity_model, shots=384, randomizations=64))

circuit_generator = ExecutorCircuitGenerator(
    gate_set, pass_manager=add_checks
)
program_learn, data_mapper = circuit_generator.generate(experiment)

session = Session(backend)
executor_options = {"environment": {"job_tags": ["TUT_PECLNM"]}}
fit = circuit_generator.collect(
    Executor(session, options=executor_options).run(program_learn).result(),
    data_mapper,
)
fit = (
    FlipPostSelect()
    + ComputeObservables()
    + CurveFitObservables()
    + LeastSquaresSolve()
).run(fit)
print(
    "learning shots removed by checks:",
    f"{float(fit.raw_data.datatree['0']['data_mask'].mean()):.4f}",
)

# TREX factors from the fit's SPAM paths: the fit only identifies the product of
# state-prep and measurement error, so hand TREX the composition of the two maps
spam = fidelity_model.to_pauli_lindblad_maps(
    fit.model_data, include_spam=True
)
spam_map = spam["P"].compose(spam["M"])
z_terms = [
    QubitSparsePauli(("Z", [layout[i]]), num_qubits=backend.num_qubits)
    for i in range(n_data)
]
trex_rescale = np.array(
    [TREX.calculate_trex_factor(spam_map, z) for z in z_terms]
)

# learned maps come back in backend qubit indexing; samplex wants box-local order
plm = split_pauli_lindblad_model(fit.model).model
noise_maps = {}
for ref, m in plm.to_pauli_lindblad_maps(fit.model_data).items():
    box = ref_to_qubits[ref]
    noise_maps[ref] = PauliLindbladMap.from_sparse_list(
        [
            (p, tuple(box.index(q) for q in qs), r)
            for p, qs, r in m.to_sparse_list()
        ],
        num_qubits=len(box),
    )
# Full-PEC cost: each learned layer acts twice per Trotter step (compute + uncompute)
gamma = float(
    np.exp(2 * depth * sum(2 * sum(m.rates) for m in noise_maps.values()))
)

# Collect the learned model in the form the figure helpers expect
mit = dict(
    **{
        f"rates_{i}": np.asarray(m.rates)
        for i, m in enumerate(noise_maps.values())
    },
    **{
        f"label_paulis_{i}": np.array([p for p, _, _ in m.to_sparse_list()])
        for i, m in enumerate(noise_maps.values())
    },
    **{
        f"label_qubits_{i}": np.array(
            [
                list(qs) + [-1] * (2 - len(qs))
                for _, qs, _ in m.to_sparse_list()
            ]
        )
        for i, m in enumerate(noise_maps.values())
    },
    layout=np.array(layout),
    trex_rescale=trex_rescale,
    gammas=np.array([gamma, np.nan]),
)
draw_noise_map(mit, backend)
plot_trex(mit["trex_rescale"], x_labels(layout, n_data))
learning shots removed by checks: 0.2393

Output of the previous code cell

Output of the previous code cell

Podar el modelo de ruido de los términos detectables por las comprobaciones de simetría​

Cada comprobación de simetría es capaz de detectar algún subconjunto de generadores de errores en el modelo de ruido, lo que significa que no necesitas mitigar esos errores. Aquí creamos una máscara sobre los generadores de ruido detectables usando la función create_postselected_noise_mask de qiskit_mitigation. Esta máscara se usará para eliminar del modelo de ruido los términos detectables antes de realizar PEC. Realizar PEC con respecto al modelo de ruido podado puede producir valores esperados convergidos por una fracción del costo de muestreo requerido para realizar PEC en el modelo de ruido completo, como se ilustra en los mapas de modelos de ruido anteriores.

El mapa de abajo se dibuja en la misma escala de color que el modelo aprendido anterior, conservando únicamente los generadores de errores que las comprobaciones no pueden detectar. Vemos que se ha eliminado del modelo un gran número de generadores de errores, y el sobrecosto de muestreo se ha reducido significativamente. El sobrecosto de muestreo para ejecutar PEC en el modelo de ruido podado es γ2=2.52≈6.3\gamma^2=2.5^2\approx6.3, aproximadamente una reducción de 16 veces respecto al sobrecosto requerido para mitigar el canal de ruido completo.

from qiskit.quantum_info import Pauli
from qiskit_mitigation.noise import create_postselected_noise_mask

SHOTS_PER_RAND = 64  # shots per PEC randomization
N_RAND = 23_054  # PEC randomizations per experiment
MAX_PAIR_RATE = 0.02  # Max value of any coupler's summed two-qubit error rate

# The detectors are virtual Pauli-Z's representing the measurements on the mediator qubits
# We can find the set of detectable noise generators in the model for each detector by conjugating
# it backward through the circuit and calculating what generators it anti-commutes with.
detectors = [
    Pauli("I" * (boxed.num_qubits - 1 - q) + "Z" + "I" * q)
    for q in layout[n_data:]
]
# Detectable generators are handled by postselection; prune them from the model for PEC
local_scales, gamma2_post = create_postselected_noise_mask(
    boxed, noise_maps, detectors
)
gamma_post = gamma2_post**0.5
print(f"gamma PEC {gamma:.2f} | gamma QED+PEC {gamma_post:.2f}")

# Kill the job if any coupler's summed rate exceeds MAX_PAIR_RATE
pair_lam = {}
for m in noise_maps.values():
    for _, qs, r in m.to_sparse_list():
        if len(qs) == 2:
            pair_lam[tuple(sorted(qs))] = (
                pair_lam.get(tuple(sorted(qs)), 0.0) + r
            )
assert (
    max(pair_lam.values()) < MAX_PAIR_RATE
), f"kill: degraded coupler {max(pair_lam, key=pair_lam.get)}"

# Detectability depends on circuit position, so record each site's 0/1 scales by layer
site_ref = {
    a.modifier_ref: a.ref
    for inst in boxed.data
    if inst.operation.name == "box"
    and (a := get_annotation(inst.operation, InjectNoise))
    and a.ref
}
sites = sorted(local_scales, key=lambda s: int(s[1:]))
mit["site_scales"] = np.stack([local_scales[s] for s in sites])
mit["site_layer"] = np.array(
    [list(noise_maps).index(site_ref[s]) for s in sites]
)
mit["gammas"] = np.array([gamma, gamma_post])
draw_noise_map(mit, backend, reduced=True)
gamma PEC 10.08 | gamma QED+PEC 2.52

Output of the previous code cell

Muestrear el circuito detector de errores​

Aquí muestreamos el circuito de Ising detector de errores de 49 qubits. Habilitamos el twirling de Pauli y las comprobaciones de errores no markovianos, pero no realizamos muestreo con PEC. Usaremos estas muestras para calcular valores esperados de referencia, así como valores esperados solo con detección de errores.

# Baseline = the ED+PEC template itself at noise scale 0, i.e. twirling only
program_tw = QuantumProgram(shots=100)
program_tw.append_samplex_item(
    template_checked,
    samplex=samplex,
    shape=(1000, 1),
    samplex_arguments={
        "pauli_lindblad_maps": dict(noise_maps),
        basis_key: meas_basis,
        **{f"noise_scales.{k}": 0.0 for k in local_scales},
    },
)
job_tw = Executor(session, options=executor_options).run(program_tw)

Muestrear el circuito detector de errores con PEC​

Ahora muestreamos de nuevo e incluimos aleatorizaciones de PEC para mitigar el ruido que las comprobaciones de simetría no pueden detectar. Estas muestras se usarán para calcular valores esperados usando la detección de errores en combinación con PEC. A continuación, graficamos los shots aceptados por aleatorización de PEC frente a una binomial con la tasa de aceptación media. Vemos que la aceptación no está correlacionada con la instancia de circuito muestreada porque QED+PEC inyecta únicamente Paulis indetectables, que conmutan con las mediciones de todos los qubits mediadores.

def launch(session, template, samplex, samplex_args, n_rand, shots_per_rand):
    """Sample `n_rand` randomizations of `template` on `session` in <=100k-randomization jobs; returns the jobs."""
    jobs = []
    for start in range(0, n_rand, 100_000):
        program = QuantumProgram(shots=shots_per_rand)
        program.append_samplex_item(
            template,
            samplex=samplex,
            samplex_arguments=samplex_args,
            shape=(min(100_000, n_rand - start), 1),
        )
        jobs.append(Executor(session, options=executor_options).run(program))
    return jobs

# Reduced PEC
args_post = {
    "pauli_lindblad_maps": dict(noise_maps),
    basis_key: meas_basis,
    **{f"noise_scales.{k}": -1.0 for k in local_scales},
    **{f"local_scales.{k}": v for k, v in local_scales.items()},
}
# PEC on the full noise model
args_van = {
    "pauli_lindblad_maps": dict(noise_maps),
    basis_key: meas_basis,
    **{f"noise_scales.{k}": -1.0 for k in local_scales},
}

# Run sampling jobs
jobs = launch(
    session, template_checked, samplex, args_post, N_RAND, SHOTS_PER_RAND
)
jobs_van = launch(
    session, template_checked, samplex, args_van, N_RAND, SHOTS_PER_RAND
)
outs = [j.result()[0] for j in jobs]
outs_van = [j.result()[0] for j in jobs_van]
(out_tw,) = job_tw.result()
session.close()

Recopila las muestras y verifica que las aleatorizaciones de PEC conmuten con las mediciones en los qubits mediadores y no afecten a las estadísticas de posselección.

def bitflip_mask(out, selector):
    """Per-shot True/False for job result `out`: shot passes `selector`'s non-Markovian error checks."""
    regs = {
        k: np.asarray(out[k])
        for k in out
        if not k.startswith(("measurement_flips", "pauli_signs"))
    }
    return selector.compute_mask(regs, "edge", mode="post")

def symmetry_mask(out):
    """Per-shot True/False for job result `out`: every twirl-corrected symmetry check reads 0."""
    return ~(
        np.asarray(out["check"]) ^ np.asarray(out["measurement_flips.check"])
    ).any(axis=-1)

def keep_mask(out, selector):
    """Per-shot True/False for job result `out`: shot passes both check types."""
    return symmetry_mask(out) & bitflip_mask(out, selector)

def fracs(out_list, selector):
    """Fractions of shots across the job results in `out_list` passing [no, non-Markovian error, symmetry, both] checks."""
    bf = np.concatenate([bitflip_mask(o, selector) for o in out_list])
    sy = np.concatenate([symmetry_mask(o) for o in out_list])
    return [1.0, bf.mean(), sy.mean(), (bf & sy).mean()]

mask_post = np.concatenate([keep_mask(o, selector) for o in outs])

# Collect the sampled data alongside the learned model, in the form the figure helpers expect
mit.update(
    ps_fracs=np.array(
        [
            fracs([out_tw], selector),
            fracs(outs_van, selector),
            fracs(outs, selector),
        ]
    ),
    acc_counts_post=np.squeeze(mask_post).sum(axis=-1),
    data_tw=np.asarray(out_tw["data"]),
    flips_tw=np.asarray(out_tw["measurement_flips.data"]),
    mask_tw=keep_mask(out_tw, selector),
    data_post=np.packbits(
        np.concatenate([np.asarray(o["data"]) for o in outs]), axis=-1
    ),
    flips_post=np.concatenate(
        [np.asarray(o["measurement_flips.data"]) for o in outs]
    ),
    signs_post=np.concatenate([np.asarray(o["pauli_signs"]) for o in outs]),
    mask_post=mask_post,
    data_van=np.packbits(
        np.concatenate([np.asarray(o["data"]) for o in outs_van]), axis=-1
    ),
    flips_van=np.concatenate(
        [np.asarray(o["measurement_flips.data"]) for o in outs_van]
    ),
    signs_van=np.concatenate(
        [np.asarray(o["pauli_signs"]) for o in outs_van]
    ),
    mask_van=np.concatenate([bitflip_mask(o, selector) for o in outs_van]),
)

# Plot the PEC bias check
plot_postselection(mit)

Output of the previous code cell

Calcular valores esperados y comparar estrategias​

Finalmente, calculamos todos los valores esperados con la función auxiliar executor_expectation_values de qiskit-mitigation, que aplica por nosotros las inversiones de medición, las máscaras de posselección, los factores de reescalado TREX y los signos de cuasiprobabilidad.

Gráfico superior: Vemos que ambas variantes de PEC convergen, pero QED+PEC alcanza la banda de ±0.025\pm 0.025 con menos aleatorizaciones y barras de error más ajustadas. El sesgo residual en el cálculo de QED+PEC está por debajo de 0.01 y puede atribuirse a alguna combinación de discrepancia del modelo de ruido, deriva de los qubits a lo largo del tiempo, y efectos de fuentes de ruido no modeladas.

Gráfico inferior: La estimación acumulada del ⟨X⟩\langle X \rangle promediado por sitio a medida que se acumulan las aleatorizaciones, usando shots idénticos por aleatorización para ambas variantes de PEC. PEC + detección de errores se estabiliza dentro de la banda en unos pocos miles de aleatorizaciones; PEC solo, al pagar el sobrecosto de muestreo completo, requiere muchas más aleatorizaciones para converger.

from qiskit.quantum_info import SparsePauliOp
from qiskit_mitigation.utils import executor_expectation_values

gamma, gamma_post = mit["gammas"]
basis_map = {Pauli("X" * n_data): [SparsePauliOp(o) for o in observables]}
rescale = dict(zip(observables, mit["trex_rescale"], strict=True))

def evs(bits, basis_map, **kwargs):
    """Per-site (means, standard errors) from boolean shot data `bits` of shape (rands, 1, shots, n_data)."""
    out = executor_expectation_values(
        bits, basis_map, None, avg_axis=(0, 1), **kwargs
    )
    return np.array([m for m, _ in out]).ravel(), np.sqrt(
        [v for _, v in out]
    ).ravel()

def unpack(mit, prefix, n_data):
    """Restore `mit[f"data_\{prefix\}"]` from packed bytes to booleans of shape (rands, 1, shots, n_data)."""
    return np.unpackbits(mit[f"data_{prefix}"], axis=-1)[..., :n_data].astype(
        bool
    )

unmit_tw, unmit_tw_err = evs(
    mit["data_tw"], basis_map, measurement_flips=mit["flips_tw"]
)
ed_tw, ed_tw_err = evs(
    mit["data_tw"],
    basis_map,
    measurement_flips=mit["flips_tw"],
    postselect_mask=mit["mask_tw"],
    rescale_factors=rescale,
)
post, post_err = evs(
    unpack(mit, "post", n_data),
    basis_map,
    measurement_flips=mit["flips_post"],
    pauli_signs=mit["signs_post"],
    postselect_mask=mit["mask_post"],
    rescale_factors=rescale,
)
pec, pec_err = evs(
    unpack(mit, "van", n_data),
    basis_map,
    measurement_flips=mit["flips_van"],
    pauli_signs=mit["signs_van"],
    postselect_mask=mit["mask_van"],
    rescale_factors=rescale,
)

# effective ED+PEC overhead measured from the data: accepted shots per signed accepted shot
signs_post = 1 - 2 * (np.squeeze(mit["signs_post"]).sum(axis=-1) % 2)
gamma_eff = (
    mit["mask_post"].sum()
    / (signs_post * np.squeeze(mit["mask_post"]).sum(axis=-1)).sum()
)

print(
    f"gamma PEC {gamma:.2f} | gamma QED+PEC {gamma_post:.2f} (model) / {gamma_eff:.2f} (data) | "
    f"QED survival {mit['mask_tw'].mean():.3f} | QED+PEC survival {mit['mask_post'].mean():.3f} | "
    f"mean QED+PEC err {post_err.mean():.4f}"
)
plot_final(
    obs_exact,
    (unmit_tw, unmit_tw_err),
    (ed_tw, ed_tw_err),
    (pec, pec_err),
    (post, post_err),
    (gamma, gamma_post),
    x_labels(layout, n_data),
    "",
)
gamma PEC 10.08 | gamma QED+PEC 2.52 (model) / 2.51 (data) | QED survival 0.304 | QED+PEC survival 0.305 | mean QED+PEC err 0.0037

Output of the previous code cell

plot_convergence(mit, obs_exact, n_data)

Output of the previous code cell

Fuente: Documentación de IBM Quantum — actualizado el 13 oct 2026
Versión en inglés en doQumentation — actualizado el 13 oct 2026
¿Problema con el sitio o la traducción?¿Problema de contenido? Editar en IBM Quantum docs