Annulation probabiliste des erreurs à l'aide de modèles de bruit logiques
Estimation de la durée d'exécution : 28 minutes sur un processeur Heron r3 (REMARQUE : il s'agit uniquement d'une estimation. (Votre temps d'exécution peut varier.)
Dans ce tutoriel, nous allons évaluer les observables d'un modèle d'Ising à 22 sites sur un réseau hexagonal à l'aide d'un QPU Heron, doté d'une topologie de qubits « heavy-hex ». De nombreuses démonstrations réalisées sur les QPU Heron portent sur des systèmes définis sur des réseaux hexagonaux denses afin de s'adapter à la connectivité du système; cependant, les modèles hexagonaux sont généralement plus intéressants à étudier, car ils sont plus courants dans la nature et plus difficiles à simuler de manière classique en raison de leur connectivité plus dense. Étant donné que la topologie des qubits du QPU ne permet pas directement d’assurer la connectivité d’un réseau hexagonal, nous devons déterminer comment intégrer efficacement le modèle hexagonal à la topologie des qubits. Dans cet exemple, nous représentons chaque site du modèle d'Ising par un qubit situé à un sommet du réseau « heavy-hex » du QPU. Nous utilisons les qubits situés sur les arêtes ( qubits médiateurs ) pour faciliter l'intrication entre les qubits d'Ising situés aux sommets. De plus, les qubits médiateurs mettent en œuvre l'intrication de telle sorte qu'ils devraient toujours revenir à l'état fondamental . Si un qubit médiateur affiche la mesure , cela indique qu'une erreur logique s'est produite pendant l'exécution du circuit; le fait de ne retenir que les échantillons pour lesquels aucune erreur n'a été détectée sur les qubits médiateurs permet d'obtenir une distribution plus restreinte d' échantillons logiques, susceptible d'avoir une fidélité supérieure à celle de la distribution brute, entachée de bruit. Non seulement nous pouvons détecter qu'une erreur s'est produite en mesurant un qubit médiateur, mais nous pouvons également déterminer avec précision quelles erreurs logiques ce qubit est capable de détecter à n'importe quel endroit du circuit. En retirant d'un modèle de bruit les sources de bruit que les contrôles de symétrie des qubits médiateurs sont capables de détecter, on obtient un modèle de bruit réduit, qui peut être atténué à l'aide de techniques telles que l'annulation probabiliste des erreurs (PEC).
Dans cet exemple de cahier d'exercices, nous allons associer la technique de détection d'erreurs décrite ci-dessus au PEC afin de contrer le bruit quantique plus efficacement que ne le permettrait chacune de ces techniques prise isolément. Nous utiliserons 49 qubits pour implémenter le modèle d'Ising à 22 qubits et nous utiliserons les 27 qubits médiateurs supplémentaires pour la détection d'erreurs. Nous appliquerons la méthode PEC au canal de bruit présélectionné afin d'atténuer le bruit qui échappe aux contrôles de symétrie. Outre la combinaison de la détection d'erreurs avec le PEC, nous utiliserons des techniques telles que l'atténuation des erreurs de lecture TREX et les contrôles d'erreurs non markoviens afin de lutter davantage contre les effets du bruit quantique.
Le flux de travaux est le suivant :
- Simulation classique des valeurs attendues observables pour le modèle hex-Ising à 22 qubits
- Implémenter le modèle hex-Ising à 49 qubits avec détection d'erreurs à l'aide d'un circuit quantique, puis le transposer vers un backend QPU
- Indiquez les couches d'intrication du circuit à l'aide de
samplomatic. Nous allons étudier et atténuer le bruit qui affecte ces couches d’entrelacement. - Ajouter des contrôles d'erreurs non markoviennes au circuit
- Découvrez le bruit de grille et le bruit de lecture qui affectent le circuit
- Éliminer les termes du modèle de bruit détectés par les contrôles de symétrie
- Testez le circuit de détection d'erreurs. Ne sélectionner a posteriori que les échantillons qui satisfont à tous les contrôles de symétrie et d'erreurs non markoviennes, et utiliser la technique d'atténuation des erreurs de lecture TREX pour tous les calculs de valeurs attendues
- Exemple de circuit de détection d'erreurs
- Effectuez un échantillonnage du circuit de détection d'erreurs avec PEC, mais ne procédez pas à une post-sélection basée sur des contrôles de symétrie et atténuez l'intégralité du canal de bruit appris.
- Testez le circuit de détection d'erreurs à l'aide de la méthode PEC. Effectuer une sélection a posteriori sur la base de contrôles de symétrie et n'atténuer que le bruit qui n'est pas détectable par ces contrôles.
- Calculez les valeurs attendues et comparez les stratégies. On constate que la méthode QED+PEC élimine plus efficacement le bruit affectant l'expérience que chacune de ces deux méthodes prise isolément et permet d'obtenir des valeurs attendues convergentes en utilisant un nombre de mesures bien inférieur à celui requis par la méthode PEC seule.
Exigences
Avant de commencer ce tutoriel, exécutez la commande suivante pour installer tous les paquets nécessaires :
%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 cellule repliée ci-dessous définit les éléments d'aide à la création des figures utilisés tout au long de ce cahier.
# 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)
Simuler de manière classique l' e pour chaque site dans le modèle hex-Ising à 22 qubits
Tout d'abord, nous utilisons la classe Statevector Qiskit pour simuler les valeurs cibles exactes de notre expérience. Bien que le problème à 22 qubits que nous résolvons soit traitable de manière classique, le simulateur de vecteurs d’état de Qiskit ne permettra pas de traiter des démonstrations beaucoup plus volumineuses que celle-ci; pour dépasser environ 50 qubits, il faudrait recourir à des techniques de simulation inexactes, telles que la propagation de Pauli. Cette expérience à 49 qubits nous permet d'étudier la combinaison de la détection d'erreurs et de l'atténuation d'erreurs à plus grande échelle, tout en conservant l'accès aux valeurs d'espérance idéales.
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:
Implémenter le modèle hex-Ising à 49 qubits avec détection d'erreurs à l'aide d'un circuit quantique, puis le transposer vers un backend QPU
Ce circuit simule 4 étapes de Trotter d'un hamiltonien d'Ising à champ transversal de 22 qubits sur un réseau hexagonal. Les 22 qubits de données sont intégrés dans un sous-graphe « heavy-hex » de 49 qubits du processeur quantique (QPU) Heron r3ibm_boston. Les 27 qubits auxiliaires ont deux fonctions : assurer l'intrication entre les qubits d'Ising et détecter les erreurs logiques lors de l'exécution du circuit. Les couches d'intrication du circuit sont mises en œuvre de telle sorte que les qubits médiateurs reviennent à l'état fondamental avant les mesures finales. La détection d’ s par un ou plusieurs qubits médiateurs indique que l’échantillon a été altéré par une erreur logique. Le fait d'écarter ces échantillons peut améliorer la fidélité de la distribution échantillonnée.
Dans le graphique ci-dessous, les qubits verts représentent les 22 qubits d'Ising, tandis que les qubits orange représentent les 27 qubits médiateurs.
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:
Indiquez les couches d'intrication du circuit à l'aide de samplomatic.
Le « Samplomatic generate_boxing_pass_manager » regroupe les couches d'intrication du circuit dans des cases annotées destinées à la rotation de Pauli et à l'injection de bruit PEC. Le circuit modèle ainsi obtenu et le « samplex » (une distribution paramétrique sur le circuit modèle spécifiant ses « twirls » et l'injection de bruit) sont utilisés pour définir et mener toutes les expériences d'apprentissage du bruit et d'échantillonnage sur QPU. Les boîtes de mesure comportent également une ChangeBasis annotation; ainsi, les rotations permettant de mesurer les qubits de données dans la base X sont fournies au samplex en tant qu'entrée au moment de l'échantillonnage, plutôt que d'être intégrées au circuit.
Ci-dessous, cette même opération de « boxing » est appliquée à une version miniature du circuit d'Ising de détection d'erreurs, ce qui nous permet de visualiser la structure d'un circuit « boxé ».
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:
Ajouter des contrôles d'erreurs non markoviennes au circuit
Ensuite, nous ajoutons au circuit des contrôles d'erreurs non markoviens afin de le protéger contre le bruit qui n'est pas pris en compte dans notre protocole d'apprentissage. Ces vérifications consistent à appliquer une inversion de bits à impulsion longue et à s'assurer que le qubit est bien passé d'un état classique à un autre. Si les deux qubits d'une arête échouent au contrôle, l'échantillon est écarté. Les contrôles d'erreurs non markoviens peuvent être utilisés au début ou à la fin du circuit (ou aux deux endroits); ici, nous ne les utilisons qu'à la fin. Nous plaçons également des qubits « spectateurs » inutilisés à proximité du circuit, ce qui renforce la protection contre le bruit que ces qubits sont censés détecter.
N'oubliez pas que, tout comme les contrôles de symétrie, il s'agit là aussi d'une forme de post-sélection; nous devons donc veiller à tenir compte de la baisse du taux de post-sélection lorsque nous combinons ces techniques. Concrètement, si et , il faut alors environ échantillons bruités pour récupérer un échantillon logique.
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:
Découvrez le bruit de porte et le bruit de lecture du circuit Ising vérifié à 49 qubits
Pour atténuer le bruit, nous devons modéliser la manière dont celui-ci affecte nos portes d'intrication. Nous allons ici apprendre un modèle de bruit de Pauli-Lindblad pour chacune des trois couches d'intrication distinctes du circuit à l'aide de qiskit-noise-learning. Par la suite, nous éliminerons de ce modèle de bruit les générateurs d'erreurs détectables par les contrôles de symétrie et n'atténuerons que le canal de bruit restant à l'aide de la PEC. La carte ci-dessous représente les trois couches apprises sur le QPU; chaque qubit correspond à une roue indiquant ses taux d' weight-1 s X/Y/Z, et chaque coupleur à une grille 3×3 représentant les taux de Pauli weight-2 pour cette paire. La valeur « » de « » implique une surcharge d'échantillonnage de « ».
Les corrections de lecture (TREX) proviennent directement des trajectoires SPAM issues de l'ajustement par apprentissage du bruit. Le diagramme en tiges situé sous la carte de bruit du QPU présente les facteurs de mise à l'échelle TREX par 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)
fit = circuit_generator.collect(
Executor(session).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))Output:
learning shots removed by checks: 0.2393
Éliminer du modèle de bruit les termes détectables par les contrôles de symétrie
Chaque contrôle de symétrie permet de détecter un sous-ensemble de sources d'erreurs dans le modèle de bruit, ce qui signifie que vous n'avez pas besoin de corriger ces erreurs. Ici, nous créons un masque recouvrant les sources de bruit détectables à l'aide de la fonction create_postselected_noise_mask issue de qiskit_mitigation. Ce masque servira à éliminer les termes détectables du modèle de bruit avant d'effectuer le PEC. L'application de la méthode PEC au modèle de bruit allégé peut permettre d'obtenir des valeurs d'espérance convergentes pour une fraction du coût d'échantillonnage nécessaire à l'application de la méthode PEC au modèle de bruit complet, comme l'illustrent les cartes du modèle de bruit ci-dessus.
La carte ci-dessous est représentée selon la même échelle de couleurs que le modèle entraîné ci-dessus, en ne conservant que les sources d'erreur que les contrôles ne permettent pas de détecter. On constate qu'un grand nombre de sources d'erreurs ont été supprimées du modèle et que la charge liée à l'échantillonnage a considérablement diminué. La surcharge d'échantillonnage liée à l'exécution de la PEC sur le modèle de bruit élagué est de l' , ce qui représente une réduction d'environ l' 16x e par rapport à la surcharge nécessaire pour atténuer le canal de bruit complet.
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)Output:
gamma PEC 10.08 | gamma QED+PEC 2.52
Exemple de circuit de détection d'erreurs
Nous présentons ici un exemple de circuit d'Ising à 49 qubits permettant la détection d'erreurs. Nous prenons en charge le « Pauli twirling » et les contrôles d'erreurs non markoviens, mais nous n'effectuons pas d'échantillonnage PEC. Nous utiliserons ces échantillons pour calculer les valeurs attendues de référence ainsi que les valeurs attendues relatives à la seule détection d'erreurs.
# 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).run(program_tw)Tester le circuit de détection d'erreurs avec PEC
Nous procédons à un nouvel échantillonnage et intégrons des randomisations PEC afin d'atténuer le bruit que les contrôles de symétrie ne permettent pas de détecter. Ces échantillons serviront à calculer les valeurs attendues en combinant la détection d'erreurs et la méthode PEC. Nous représentons ci-dessous le nombre de prises acceptées par randomisation PEC en fonction d'une distribution binomiale correspondant au taux d'acceptation moyen. On constate que l'acceptation n'est pas corrélée à l'instance de circuit échantillonnée, car QED+PEC n'injecte que des particules de Pauli indétectables, qui commutent avec les mesures de chaque qubit médiateur.
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).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()Prélevez les échantillons et vérifiez que les randomisations PEC commutent avec les mesures effectuées sur les qubits médiateurs et n'affectent pas les statistiques de post-sélection.
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:
Calculer les valeurs attendues et comparer les stratégies
Enfin, nous calculons toutes les valeurs d'espérance à l'aide de la fonction d'aide executor_expectation_values de qiskit-mitigation, qui se charge pour nous d'appliquer les inversions de mesure, les masques de post-sélection, les facteurs de mise à l'échelle TREX et les signes de quasi-probabilité.
Graphique du haut : on constate que les deux variantes du modèle PEC convergent, mais que la méthode QED+PEC atteint la bande d’ s avec moins de randomisations et des barres d’erreur plus étroites. Le biais résiduel dans le calcul QED+PEC est inférieur à 0.01 et peut être attribué à une combinaison de divergences entre les modèles de bruit, de dérive des qubits au fil du temps et d'effets provenant de sources de bruit non modélisées.
Graphique du bas : estimation glissante de l' e moyenne par site à mesure que les randomisations s'accumulent, en utilisant le même nombre de tirs par randomisation pour les deux variantes du PEC. La détection d'erreurs PEC + se stabilise dans la bande en quelques milliers d'itérations aléatoires; la méthode PEC seule, qui supporte l'intégralité de la charge d'échantillonnage, nécessite beaucoup plus d'itérations aléatoires pour 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),
"",
)Output:
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
plot_convergence(mit, obs_exact, n_data)Output: