שיפור ערכי ציפייה: ספיגת רעש מתפשט (PNA)

במדריך הזה נלמד איך למנף את הכלים האחרונים במערכת Qiskit כדי להפעיל תהליך עבודה מתאים לשגיאות, שניתן להתאמה מלאה. נציג את טכניקת PNA ונשתמש בה כדי להפחית שגיאות Gate. כמו כן, נשתמש ב-TREX להפחתת שגיאות קריאה ובבחירה-מוקדמת להפחתת שגיאות שלא נלכדות במודל הרעש הנלמד.

מתווה

• סקירה קצרה של PNA
• יצירת Circuit קוונטי מסוג Trotterized ושל Observable. Transpile שלו ל-Backend וכלילת מדידות לבחירה-מוקדמת.
• שימוש ב-samplomatic כדי לסובב (twirl) שכבות של שערי 2Q ומדידות. מציאת שכבות 2Q ייחודיות כדי לצמצם את עלות למידת הרעש.
• שימוש ב-NoiseLearnerV3 ללמידת מודל השגיאות המשפיע על שערי 2Q ומדידות.
• שימוש ב-qiskit-addon-pna ליצירת observable מפחית-רעש.
• שימוש ב-Executor primitive של qiskit-ibm-runtime ליצירת הדגימות הגולמיות מה-QPU, המשקפות כל shot לכל אקראות סיבוב וכל בסיס מדידה.
• שימוש ב-qiskit-addon-utils לעיבוד הנתונים לאחר מכן לערך ציפייה מופחת-שגיאות.

מהי ספיגת רעש מתפשט (PNA)?

טכניקה להפחתת שגיאות Gate על ידי הפצת ה-observable דרך ערוץ הרעש ההפוך המשפיע על שערי 2-qubit, ויצירת observable מפחית-רעש. שערי 2Q בניסוי שאנחנו רוצים להריץ יושפעו מרעש משמעותי. אם נלמד את מודל הרעש, נוכל להחיל את ההיפוך שלו ולבטל את הרעש. במקום לממש את ערוץ הרעש ההפוך על ידי דגימתו ב-QPU כמו ב-PEC, אפשר לממש אותו קלאסית ב-observable הנמדד באמצעות הפצת Pauli. התוצאה היא observable מורכב יותר שכאשר נמדד, יש לו את ההשפעה של הפחתת רעש ה-Gate הנלמד.

יצירת Circuit Trotter מוראי ו-observable

לניסוי הזה נחקור את דינמיקת הזמן של מודל Ising מובעט (kicked) עם 30 אתרים על שרשרת ספין 1D. ה-Hamiltonian הנבחר הוא:

$H = -J\sum\limits_{\langle i,j \rangle} Z_iZ_j + h\sum\limits_iX_i$,

כאשר $J>0$ מתאר את הצימוד של ספינים שכנים, $i<j$, ושדה הרוחב הגלובלי, $h$, מוגדר כ-$\frac{\pi}{8}$. ככל ש-$h$ רחוק יותר מזווית Clifford (כלומר $\theta=n\frac{\pi}{2}, n \in \mathbb{Z}$), כך קשה יותר להפיץ את מחוללי האנטי-רעש דרך ה-Circuit.

לבחירת ה-observable, נבחן את המגנטיזציה הממוצעת של אתר בודד, $\frac{1}{N} \sum_{i=1}^{N} \langle z_i \rangle$, כאשר $N$ הוא מספר האתרים.

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

import numpy as np
from qiskit import QuantumCircuit
from qiskit.quantum_info import Pauli, SparsePauliOp

num_qubits = 30
num_trotter_steps = 10
rx_angle = np.pi / 8

# Avg single-site magnetization
id_pauli = Pauli("I" * num_qubits)
observable = SparsePauliOp([id_pauli.dot(Pauli("Z"), [i]) for i in range(num_qubits)]) / num_qubits

# Implement Trotterized kicked-Ising model
circuit = QuantumCircuit(num_qubits)
for _step in range(num_trotter_steps):
    circuit.rx(rx_angle, range(num_qubits))
    for first_qubit in (1, 2):
        for idx in range(first_qubit, num_qubits, 2):
            # equivalent to Rzz(-pi/2):
            circuit.sdg([idx - 1, idx])
            circuit.cz(idx - 1, idx)
circuit.compose(circuit.inverse(), inplace=True)
circuit.measure_active()
circuit.draw("mpl", fold=-1)

בשלב הבא נבחר שרשרת של Qubits ב-ibm_kingston שמדווחות על שיעורי שגיאה נמוכים ונבצע Transpile של ה-Circuit ל-Backend.

from qiskit.transpiler import generate_preset_pass_manager
from qiskit_ibm_runtime import QiskitRuntimeService

backend_name = "ibm_kingston"
service = QiskitRuntimeService()
backend = service.backend(backend_name, use_fractional_gates=True)

# Use a chain of low-noise qubits
layout = [
    44,
    45,
    46,
    47,
    57,
    67,
    68,
    69,
    78,
    89,
    88,
    87,
    97,
    107,
    106,
    105,
    117,
    125,
    126,
    127,
    128,
    129,
    118,
    109,
    110,
    111,
    98,
    91,
    92,
    93,
]

pm = generate_preset_pass_manager(backend=backend, initial_layout=layout, optimization_level=0)
isa_circuit = pm.run(circuit)
isa_observable = observable.apply_layout(isa_circuit.layout)
isa_circuit.draw("mpl", fold=-1)

qiskit_runtime_service._discover_account:WARNING:2025-11-10 14:30:57,148: Loading account with the given token. A saved account will not be used.

סיבוב שכבות שערי 2-qubit ומדידות ומציאת שכבות ייחודיות

כאן אנחנו מוודאים שמנהל ה-pass מסמן את התיבות עם הערות Twirl ו-InjectNoise, שמאפשרות לנו ללמוד את הרעש שישפיע על ה-Circuit שלנו ולשייך את הרעש הזה לשכבת ה-Circuit המתאימה לו.

• enable_gates/enable_measure: True: מקיף (Box) את כל שכבות שערי 2Q ומדידות הקצה. שערים חד-ביטיים יירשמו כ-left-dressed בתוך התיבות.
• measure_annotations: all כולל הערות Twirl ו-ChangeBasis על תיבת המדידה.
• twirling_strategy: active: מסובב את כל ה-Qubits הפעילים בכל תיבה המכילה שערים מסבכים (entangling).
• inject_noise_targets: gates: הערות InjectNoise יתווספו לכל תיבות המסומנות ב-Twirl המכילות שערים מסבכים.
• inject_noise_strategy: uniform_modification: כל שכבות הרעש ישוקללו בצורה שווה.

from samplomatic.transpiler import generate_boxing_pass_manager

# Box up circuit with Twirl and InjectNoise annotations
pm = generate_boxing_pass_manager(
    enable_gates=True,
    enable_measures=True,
    measure_annotations="all",
    twirling_strategy="active",
    inject_noise_targets="gates",
    inject_noise_strategy="uniform_modification",
    remove_barriers=True,
)
boxed_circuit = pm.run(isa_circuit)

draw_circ = QuantumCircuit(boxed_circuit.num_qubits)
draw_circ.append(boxed_circuit.data[0], qargs=boxed_circuit.data[0].qubits)
draw_circ.append(boxed_circuit.data[1], qargs=boxed_circuit.data[1].qubits)
draw_circ.draw("mpl", fold=-1, scale=0.3, idle_wires=False)

יצירת ה-Circuit תבנית וה-samplex, הגדרת אופן הדגימה של ה-Circuit

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

import samplomatic
from qiskit.transpiler import PassManager
from qiskit_addon_utils.noise_management.post_selection.transpiler.passes import (
    AddPostSelectionMeasures,
    AddSpectatorMeasures,
)

# Build template circuit and samplex for later use with the "Executor"
template_circuit, samplex = samplomatic.build(boxed_circuit)

# Add post-selection instructions to the template circuit
post_selection_pm = PassManager(
    [
        AddSpectatorMeasures(backend.coupling_map),
        AddPostSelectionMeasures(x_pulse_type="rx"),
    ]
)
template_circuit = post_selection_pm.run(template_circuit)

draw_circ = template_circuit.copy_empty_like()
draw_circ.data = template_circuit.data[:324]
draw_circ.draw("mpl", fold=-1, scale=0.3, idle_wires=False)

למידת הרעש

לפני שמריצים את הניסויים, נלמד את מודל הרעש המשפיע על שערי ה-entangling ועל המדידות ב-Circuit. נחוץ מודל רעש מדויק כדי להפחית שגיאות ביעילות. למידת הרעש ממש לפני הרצת הניסויים נותנת את ההזדמנות הטובה ביותר שמודל הרעש יתאר נאמנה את הרעש האמיתי שמשפיע על השערים בזמן ההרצה.

לפני שלומדים את הרעש, עלינו למצוא את שכבות ה-2-qubit הייחודיות ב-Circuit שלנו, כדי שנוכל לצמצם את מספר ה-shots הדרושים ללמידת הרעש עבור כל ה-Circuit. אנחנו משתמשים ב-find_unique_box_instructions מ-samplomatic כדי לקבל את השכבות הייחודיות מה-Circuit המוקפן, כולל שכבת המדידה. אלה הן השכבות שמעבירים ללומד הרעש.

ברגע שידועות השכבות, אפשר ללמוד את הרעש. יש כמה פרמטרים שיש לשקול:

• num_randomizations: מספר ה-circuits האקראיים לשימוש לכל קונפיגורציה של Circuit ללמידה.
• shots_per_randomization: המספר הכולל של shots לשימוש לכל Circuit לימוד אקראי.
• layer_pair_depths: עומקי ה-Circuit (נמדדים במספר זוגות) לשימוש בניסויי הלמידה.
• post_selection: נשתמש בבחירה-מוקדמת מבוססת-קצה במהלך הלמידה עם שערי rx ליישום פעימות המדידה שלאחר מכן.

from qiskit_ibm_runtime.noise_learner_v3.noise_learner_v3 import NoiseLearnerV3
from qiskit_ibm_runtime.options import NoiseLearnerV3Options
from samplomatic.utils import find_unique_box_instructions

# Load noise learner data from a shared job
load_saved_nl_result = True

# Noise learning parameters
num_randomizations_nl = 64
shots_per_randomization_nl = 128
strategy = "edge"
enable_postsel = True
x_pulse_type = "rx"

# Find the unique instructions (layers) from boxed-up circuit
unique_2q_layers_and_meas = find_unique_box_instructions(
    boxed_circuit, normalize_annotations=None, undress_boxes=True
)

noise_learner_params = {
    "num_randomizations": num_randomizations_nl,
    "shots_per_randomization": shots_per_randomization_nl,
    "layer_pair_depths": [1, 2, 4, 8, 12, 16, 24, 32, 40, 48],
    "post_selection": {
        "enable": enable_postsel,
        "strategy": strategy,
        "x_pulse_type": x_pulse_type,
    },
    "experimental": {},
}
# set the options
noise_learner_options = NoiseLearnerV3Options(**noise_learner_params)

# run the noise learner job
noise_learner = NoiseLearnerV3(backend, noise_learner_options)
noise_learner_job = noise_learner.run(unique_2q_layers_and_meas)
noise_learner_result = noise_learner_job.result()

nl_metadata = noise_learner_params | {"layout": layout}

import matplotlib.pyplot as plt

hw_rates_1q = []
hw_rates_2q = []
for nlr in noise_learner_result[:2]:
    plm_list = nlr.to_pauli_lindblad_map().to_sparse_list()
    hw_rates_1q += [rate for (pstr, qubits, rate) in plm_list if len(pstr) == 1]
    hw_rates_2q += [rate for (pstr, qubits, rate) in plm_list if len(pstr) == 2]
hw_rates_1q = sorted(hw_rates_1q)
hw_rates_2q = sorted(hw_rates_2q)
median_1q = hw_rates_1q[len(hw_rates_1q) // 2]
median_2q = hw_rates_2q[len(hw_rates_2q) // 2]
fig, ax = plt.subplots(1, 1, figsize=(14, 5))
ax.scatter(
    (hw_rates_1q),
    [(i) / (len(hw_rates_1q) - 1) for i in range(len(hw_rates_1q))],
    color="red",
    label="1q rates",
)
ax.set_xscale("log")
ax.set_ylim(0, 1.1)
ax.vlines(median_1q, 0, 1, color="red")
ax.text(median_1q * 1.1, 0.1, f"{median_1q:.2e}")
ax.scatter(
    (hw_rates_2q),
    [(i) / (len(hw_rates_2q) - 1) for i in range(len(hw_rates_2q))],
    color="blue",
    label="2q rates",
)
ax.set_xscale("log")
ax.set_ylim(0, 1.1)
ax.vlines(median_2q, 0, 1, color="blue")
ax.text(median_2q * 1.1, 0.2, f"{median_2q:.2e}")
ax.set_title("Learned noise rates")
ax.set_xlabel("Noise rate")
ax.set_yticks([])
plt.legend()

<matplotlib.legend.Legend at 0x321dd63f0>

שיוך תיבות Circuit למודל הרעש שנלמד

כאן אנחנו יוצרים מיפוי בין מזהי ה-InjectNoise של כל תיבה לבין מודל הרעש שנלמד (PauliLindbladMap) המשפיע על שערי ה-Qubit הדו-ביטיים באותה תיבה.

from samplomatic.annotations import InjectNoise
from samplomatic.utils import get_annotation

# map inject noise refs to pauli lindblad maps
refs_to_noise_models = {}
for instruction, result in zip(unique_2q_layers_and_meas, noise_learner_result, strict=False):
    if inject_noise_annot := get_annotation(instruction.operation, InjectNoise):
        refs_to_noise_models[inject_noise_annot.ref] = result.to_pauli_lindblad_map()

הפצת ה-Observable דרך האנטי-רעש שנלמד לקבלת Observable מקזז-רעש

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

כאשר מריצים PNA, תמיד יש לספק Circuit רועש ו-Observable. אם ה-Circuit הרועש שלך הוא Circuit ממוסגר עם אנוטציות InjectNoise, יש לספק את המיפוי שיצרנו בשלב לעיל. ניתן גם להעביר Circuit לא-ממוסגר המכיל הוראות PauliLindbladError מ-qiskit-aer. במקרה זה, אין צורך לספק את refs_to_noise_models. בנוסף לקלטים הראשיים, המשתמשים כדאי שישקלו:

• max_err_terms: מספר האיברים לשמור בכל גנרטור אנטי-רעש כאשר הוא מופץ קדמית. מתן ערך גדול יותר בדרך כלל מגדיל את הדיוק, אך ההתנהגות הזו אינה מובטחת להיות מונוטונית.
• max_obs_terms: מספר האיברים לשמור ב-Observable מקזז-הרעש, $\tilde{O}$, כאשר הוא מופץ אחורית דרך האנטי-רעש המתפתח. ערכים גדולים יותר בדרך כלל מגדילים את הדיוק, אך לא מובטח שהם יעשו זאת בצורה מונוטונית.
• num_processes: מספר הליבות להקצות לתהליך. זכור, הגנרטורים מופצים קדמית ומוחלים על ה-Observable במקביל.
• search_step: שלב ההפצה-האחורית משתמש בשיטה חמדנית כדי לקרב שני אופרטורים בבסיס הפאולי. ניתן להאיץ שיטה זו על ידי הגדלת search_step. ראה את תיעוד pauli-prop למידע נוסף.
• num_to_measure: למרות שמשתנה זה אינו קלט ל-generate_noise_mitigating_observable, אנחנו משתמשים בו כדי לשלוט בכמה איברים מ-$\tilde{O}$ אנחנו רוצים למדוד בפועל. כאן נמדוד רק את 30 האיברים העליונים, שהם האיברים המקוריים ב-Observable שלנו. האיברים עוצבו מחדש כך שמדידתם מביאה לאפקט של הפחתת רעש השערים שנלמד. למרות שאנחנו מודדים רק 30 איברים מ-$\tilde{O}$, עדיין לעתים קרובות שימושי לאפשר לו לגדול, שכן זה מגדיל את הדיוק של גורמי ה-scaling של האיברים המובילים.

from qiskit_addon_pna import generate_noise_mitigating_observable

# PNA parameters
num_processes = 8
max_err_terms = 10_000
max_obs_terms = 10_000
num_to_measure = num_qubits

obs_tilde_isa = generate_noise_mitigating_observable(
    boxed_circuit,
    isa_observable,
    refs_to_noise_models,
    max_err_terms=max_err_terms,
    max_obs_terms=max_obs_terms,
    num_processes=num_processes,
    print_progress=True,
    search_step=8,
)
p_2_v = {p: v for v, p in enumerate(layout)}
obs_tilde_virtual = SparsePauliOp.from_sparse_list(
    [
        (pstr, [p_2_v[p] for p in p_qubits], coeff)
        for (pstr, p_qubits, coeff) in obs_tilde_isa.to_sparse_list()
    ],
    num_qubits=num_qubits,
)
obs_tilde_virtual = obs_tilde_virtual[np.argsort(np.abs(obs_tilde_virtual.coeffs))[::-1]][
    :num_to_measure
]

Finished! 13560 / 13560 generators propagated.

obs_tilde_isa = obs_tilde_isa[np.argsort(np.abs(obs_tilde_isa.coeffs))][::-1]
plt.xscale("log")
plt.yscale("log")
plt.title(r"$\tilde{O}$ coeff magnitudes")
plt.ylabel("Magnitude")
plt.xlabel("Pauli term index")
plt.plot(np.abs(obs_tilde_isa.coeffs), ".")

[<matplotlib.lines.Line2D at 0x16b69e840>]

המרת בסיסי המדידה לצורה קנונית

בשלב הבא, נמצא קבוצה מינימלית של בסיסים למדידה כך שנוכל לכסות כל איבר פאולי ב-Observable הנמדד (ניתן למדוד observables רבים בו-זמנית אם הם מתחלפים Qubit-wise). מאחר שאנחנו מודדים רק את האיברים ב-Observable המקורי שלנו, שהוא סכום כל הפאוליים הבודדי-Z, נדרש בסיס יחיד — בסיס ה-Z הכולל.

בנוסף למציאת קבוצת בסיסי מדידת פאולי, עלינו למפות איברים אלה לצורה הקנונית המצופה על ידי ה-Primitive Executor. למידע נוסף על סדר Qubit קנוני, כנס ל-תיעוד samplomatic.

from qiskit_addon_utils.exp_vals.measurement_bases import get_measurement_bases

meas_box = boxed_circuit.data[-1]
canonical_qubits = [
    idx for idx, qubit in enumerate(boxed_circuit.qubits) if qubit in meas_box.qubits
]
c_2_p = {c: p for c, p in enumerate(canonical_qubits)}  # canonical -> physical
p_2_v = {p: v for v, p in enumerate(layout)}  # physical -> virtual
c_2_v = {c: p_2_v[p] for c, p in c_2_p.items()}  # canonical -> virtual
meas_bases, bases_reverser = get_measurement_bases(obs_tilde_virtual)
meas_bases_canonical = [
    np.array([base[c_2_v[c]] for c in range(num_qubits)], dtype=np.uint8) for base in meas_bases
]

הגדרת אופן הדגימה ב-QuantumProgram

ה-QuantumProgram הוא המקום בו אנחנו מגדירים כיצד לדגום את הניסוי:

• template_circuit: ה-Circuit המכיל את כל השערים הנדרשים לביצוע כל האקראיות הרצויות (מאקראיות twirling, פרמטרים וכדומה).
• samplex: אובייקט המגדיר התפלגות הסתברות על פני כל אקראיות ה-Circuit האפשריות שממנה לדגום.
• samplex_arguments: קשירות הנדרשות להגדרה מלאה של ה-samplex
  • basis_changes: כאן אנחנו מציינים קבוצת בסיסים למדידה שתכסה את כל איברי הפאולי ב-Observable הנמדד.
  • noise_scales.ref: אנחנו מגדירים את scale של כל שכבת רעש ל-0.0 כדי למנוע הזרקת רעש נוסף לדגימות שלנו
  • pauli_lindblad_maps: נדרש אם מעבירים noise_scales. זה פשוט ממפה שכבות רעש למודל הרעש המשויך.
• shape: tuple של צורה להרחבת הצורה הסמויה המוגדרת על ידי samplex_arguments. צירים לא-טריביאליים שהוכנסו על ידי הרחבה זו מונים אקראיות.

from qiskit_ibm_runtime import QuantumProgram

# Control the # of shots during execution
shots_per_randomization_exec = 64
num_randomizations_exec = 6144

# Zero out the noise to prevent noise from being injected during execution.
# We only added InjectNoise annotations so PNA could associate the noise
# to layers in the circuit
samplex_inputs = {f"noise_scales.{ref}": 0.0 for ref in refs_to_noise_models}
samplex_inputs |= {"pauli_lindblad_maps": refs_to_noise_models}

# Specify the bases to measure
bases_broadcastable = np.expand_dims(np.array(meas_bases_canonical), axis=1)
samplex_inputs |= {"basis_changes": {"basis0": bases_broadcastable}}

# Convert samplex_inputs into a dict to pass to QuantumProgram
samplex_arguments = samplex.inputs().make_broadcastable().bind(**samplex_inputs)

# Instantiate the QuantumProgram with the specified parameters
program = QuantumProgram(shots=shots_per_randomization_exec)
program.append(
    circuit=template_circuit,
    samplex=samplex,
    samplex_arguments=samplex_arguments,
    shape=(num_randomizations_exec),
)

דגימת ה-Circuit באמצעות ה-Primitive Executor

עכשיו שהגדרנו את ה-QuantumProgram שלנו, ביצוע הניסוי הוא פשוט. אנחנו פשוט יוצרים מופע של אובייקט ה-Executor, מספקים לו את ה-Backend, ומריצים את התוכנית.

from qiskit_ibm_runtime import Executor

# Execute (sample) the circuit
executor = Executor(backend)
job_exec = executor.run(program)
exec_results = job_exec.result()

עיבוד-לאחר של הדגימות לחישוב ערך צפיה מופחת-שגיאות

כדי לחשב ערך צפיה מופחת-שגיאות, נבצע את הפעולות הבאות:

• חישוב גורמי ה-scaling של TREX בהתבסס על הרעש שנלמד המשפיע על המדידות
• יצירת מסכה לשמירת דגימות עם post-selection בלבד
• שימוש בפונקציה executor_expectation_values מ-qiskit-addon-utils לשילוב כל הנתונים לערך צפיה מופחת-שגיאות.

from qiskit_addon_utils.exp_vals.expectation_values import executor_expectation_values
from qiskit_addon_utils.noise_management import trex_factors
from qiskit_addon_utils.noise_management.post_selection import PostSelector

# Computing the TREX factors
measurement_noise_map = noise_learner_result[2].to_pauli_lindblad_map()
trex_rescale_factors = trex_factors(measurement_noise_map, bases_reverser)

# Post-select the results
post_selector = PostSelector.from_circuit(
    circuit=template_circuit, coupling_map=backend.coupling_map
)

# Compute the ps mask for filtering results
mask = post_selector.compute_mask(exec_results[0], strategy="edge")

# Compute expvals using post selected results
results = executor_expectation_values(
    exec_results[0]["meas"],
    bases_reverser,
    meas_basis_axis=0,
    avg_axis=1,
    measurement_flips=exec_results[0]["measurement_flips.meas"],
    pauli_signs=exec_results[0].get("pauli_signs", None),
    postselect_mask=mask,
    rescale_factors=trex_rescale_factors,
)

bases_reverser_unmit = {Pauli("Z" * num_qubits): [observable]}
args = [
    (bases_reverser_unmit, None, None),
    (bases_reverser, None, None),
    (bases_reverser, None, trex_rescale_factors),
    (bases_reverser, mask, None),
    (bases_reverser, mask, trex_rescale_factors),
]

evs = []
for reverser, postsel_mask, factors in args:
    # Compute expvals using post selected results
    res_ps = executor_expectation_values(
        exec_results[0]["meas"],
        reverser,
        meas_basis_axis=0,
        avg_axis=1,
        measurement_flips=exec_results[0]["measurement_flips.meas"],
        pauli_signs=exec_results[0].get("pauli_signs", None),
        postselect_mask=postsel_mask,
        rescale_factors=factors,
    )
    res_ps = np.array(res_ps)
    evs.append(res_ps[:, 0][0])

experiments = ["PNA", "PNA+TREX", "PNA+PS", "PNA+PS+TREX"]
colors = ["#d9d9d9", "#b0b0b0", "#7f7f7f", "#4c4c4c"]
plt.bar(experiments, evs[1:], color=colors)
plt.axhline(y=1, color="green", linestyle="--", linewidth=2, label="Ideal")
plt.axhline(y=evs[0], color="red", linestyle="--", linewidth=2, label="Unmitigated")
plt.ylabel("Expectation value", fontsize=14)

plt.title(r"30q Mirrored Ising, 10 Trotter steps, $\theta_{rx}=\frac{\pi}{8}$", fontsize=14)
plt.legend(loc="upper left", bbox_to_anchor=(1.05, 1), borderaxespad=0.0)
plt.xticks(rotation=45)
plt.tight_layout()
plt.show()