Zero-noise extrapolation dengan gate folding dan ekstrapolasi kustom

Zero-noise extrapolation (ZNE) adalah teknik untuk memitigasi noise kuantum dalam perhitungan nilai ekspektasi observable. Caranya dengan menjalankan sirkuit kuantum pada skala noise yang makin besar, mengambil nilai ekspektasi untuk masing-masing skala noise tersebut, lalu mengekstrapolasi untuk memperkirakan nilai ekspektasi tanpa noise. Tahap amplifikasi noise maupun ekstrapolasi bisa diimplementasikan dengan berbagai cara. Panduan ini mengamplifikasi noise dengan "digital gate folding", yaitu gate dua-qubit diganti dengan urutan ekuivalen berupa gate dan inversnya. Misalnya, mengganti unitary $U$ dengan $U U^† U$ akan menghasilkan faktor amplifikasi noise 3. Kamu juga bisa mengimplementasikan amplifikasi dengan probabilistic error amplification (PEA), yang membutuhkan langkah tambahan berupa pembelajaran model noise dan dibahas di panduan probabilistic error amplification. Panduan ini mengekstrapolasi menggunakan beberapa strategi fitting bawaan (linear dan eksponensial), dan juga menunjukkan cara mengimplementasikan fungsi fitting kustom milikmu sendiri.

Perlu diingat, meskipun ZNE bisa efektif memitigasi noise kuantum dengan biaya sampling lebih rendah daripada metode lain, ZNE tidak punya batas untuk bias.

Alur kerjanya sebagai berikut:

1. Tentukan sirkuit dan observable
2. Hitung nilai ekspektasi eksak menggunakan simulator statevector Qiskit
3. Hitung nilai ekspektasi yang noisy dan belum dimitigasi
4. Siapkan program ZNE dengan digital gate folding dan sampel sirkuit dengan skala noise yang bervariasi
5. Ekstrapolasi ke batas zero-noise dengan fungsi fitting bawaan dan kustom
6. Bandingkan nilai ekspektasi hasil mitigasi

Tentukan sirkuit dan observable

Hitung magnetisasi rata-rata situs $\langle Z \rangle_{\text{avg}} = \frac{1}{N} \sum_{i=1}^{N} \langle Z_i \rangle$ dari model Ising medan transversal (TFIM) yang di-Trotter-kan pada delapan qubit dengan topologi cincin. Kamu akan menyapu beberapa nilai kekuatan medan transversal dan memplot magnetisasi sebagai fungsi kekuatan medan.

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

_needed = {"IPython": "IPython", "matplotlib": "matplotlib", "numpy": "numpy", "qiskit": "qiskit", "qiskit_aer": "qiskit-aer", "qiskit_ibm_runtime": "qiskit-ibm-runtime", "qiskit_mitigation": "qiskit-mitigation", "scipy": "scipy"}
_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")

import warnings

import numpy as np
from qiskit.circuit import Parameter, QuantumCircuit
from qiskit.quantum_info import SparsePauliOp

warnings.filterwarnings("ignore")

def trotter_circuit(num_qubits, theta, layers, num_steps):
    """Trotterized transverse-field Ising evolution: RX on every site, then a CZ layer per coupling layer."""
    circuit = QuantumCircuit(num_qubits)
    for _ in range(num_steps):
        circuit.rx(theta, range(num_qubits))
        for layer in layers:
            circuit.barrier()
            for edge in layer:
                circuit.cz(*edge)
                circuit.sdg(edge)
            circuit.barrier()
    return circuit

num_qubits = 8
trotter_steps = 15
ring_coupling = [(i, (i + 1) % num_qubits) for i in range(num_qubits)]
layers = [
    ring_coupling[0::2],
    ring_coupling[1::2],
]  # two disjoint layers of CZ gates
theta = Parameter("theta")
circuit = trotter_circuit(num_qubits, theta, layers, trotter_steps)

# Parameter values: theta in [0, 0.25 pi] - values for the transverse-field sweep
data_points = 6
parameter_values = np.linspace(0, 0.25 * np.pi, data_points).reshape(
    (data_points, 1)
)

# Define observables - Z on each site
observables = [
    SparsePauliOp("I" * (num_qubits - q - 1) + "Z" + "I" * q)
    for q in range(num_qubits)
]

print("Visualization of 1 Trotter step")
trotter_circuit(num_qubits, theta, layers, 1).draw("mpl", scale=0.7)

Visualization of 1 Trotter step

Hitung nilai ekspektasi eksak menggunakan simulator statevector Qiskit

Gunakan kelas Statevector Qiskit untuk menghitung nilai ekspektasi eksak $Z_i$ pada setiap situs untuk setiap nilai $\theta$.

import matplotlib.pyplot as plt
from qiskit.quantum_info import Statevector

exact_evs = np.array(
    [
        [
            Statevector(circuit.assign_parameters([value]))
            .expectation_value(obs)
            .real
            for value in parameter_values.ravel()
        ]
        for obs in observables
    ]
)  # shape (observables, parameter values)

x_vals = parameter_values.ravel() / np.pi
fig, ax = plt.subplots(figsize=(7, 4))
ax.plot(x_vals, exact_evs.mean(axis=0), "ko-", label="Exact")
ax.set_xlabel(r"$\theta/\pi$")
ax.set_ylabel(r"$\langle Z \rangle_{\text{avg}}$")
ax.legend()
plt.tight_layout()
plt.show()

Hitung nilai ekspektasi yang noisy dan belum dimitigasi

Sebuah QuantumProgram berisi sirkuit kuantum dan samplex dalam format yang cocok untuk primitif Executor. MitigationTask.prepare membangun QuantumProgram dari sirkuit, observable, dan nilai parameter, lalu membungkus sirkuit dalam box dengan opsi boxing yang diberikan. Kamu mengatur opsi boxing measure_annotations ke "change_basis", yang merotasi pengukuran ke basis yang dibutuhkan, dan menambahkan add_tags: "unique_box" agar simulator bisa memetakan model noise ke box-box pada sirkuit. Program ini di-Pauli-twirl secara default, yang pada hardware mengubah noise koheren menjadi noise Pauli yang diasumsikan oleh metode mitigasi. num_randomizations menentukan berapa banyak instans sirkuit hasil twirl yang disampel dan shots_per_randomization menentukan berapa banyak shot yang diterima tiap instans, sehingga hasil kali keduanya adalah anggaran shot per sirkuit.

Pengguna dengan akses QPU akan menjalankan program dengan Executor yang dibangun dari backend milik QiskitRuntimeService. Di sini, kamu akan memberikan AerSimulator ke Executor, yang menjalankan program secara lokal, dan memakai opsi simulator layer_noise_model-nya untuk menyuntikkan model noise Pauli-Lindblad acak pada setiap layer gate entangling, yang dihasilkan dengan random_pauli_lindblad_noise dari qiskit-mitigation dengan error per layered gate (EPLG) sebesar 0.003, tipikal untuk QPU Heron saat ini.

from qiskit.circuit import CircuitInstruction
from qiskit.quantum_info import PauliLindbladMap, pauli_basis
from qiskit_aer import AerSimulator
from qiskit_ibm_runtime import Executor
from qiskit_mitigation import MitigationTask, find_combined_unique_layers

def generate_random_pauli_lindblad_noise(
    layers,
    num_qubits: int,
    eplg: float = 3e-3,
    seed: int | None = None,
) -> list[tuple[CircuitInstruction, PauliLindbladMap]]:
    """Generate random Pauli-Lindblad noise over the full Pauli basis."""
    layer_noise_model = []
    for layer in layers:
        layer_instructions = layer.operation.params[0]
        contain_measurement = any(
            inst.name.startswith("measure") for inst in layer_instructions
        )
        if contain_measurement:  # do not add measurement noise
            layer_noise = PauliLindbladMap.identity(num_qubits)
        else:
            edges = [
                (inst.qubits[0]._index, inst.qubits[1]._index)
                for inst in layer_instructions
            ]  # extract the edges of the layer

            basis_paulis = [p for p in pauli_basis(2) if np.sum(p.x + p.z)]
            rng = np.random.default_rng(seed=seed)
            num_generators = len(basis_paulis)
            rates = rng.random(len(basis_paulis)) * eplg / num_generators

            layer_noise = PauliLindbladMap.from_sparse_list(
                [
                    (pauli.to_label(), edge, rate)
                    for pauli, rate in zip(basis_paulis, rates, strict=True)
                    for edge in edges
                ],
                num_qubits=num_qubits,
            )
        layer_noise_model.append((layer, layer_noise))
    return layer_noise_model

num_randomizations = 50
shots_per_randomization = 64

boxing_options = {
    "measure_annotations": "change_basis",
    "add_tags": "unique_box",
}

# find the unique layers to generate the layered noise for
box_layers = find_combined_unique_layers(
    [circuit], custom_boxing_options=boxing_options
)
layer_noise_model = generate_random_pauli_lindblad_noise(
    box_layers, num_qubits, seed=42
)

# a local noisy simulator in place of a QPU backend
noisy_executor = Executor(
    AerSimulator(),
    options={"simulator": {"layer_noise_model": layer_noise_model}},
)

task = MitigationTask()
program = task.prepare(
    circuit=circuit,
    observables=observables,
    parameters=parameter_values,
    num_randomizations=num_randomizations,
    shots_per_randomization=shots_per_randomization,
    custom_boxing_options=boxing_options,
)
program._semantic_role = "simulation"  # flag the program for the simulator
unmitigated_res = task.postprocess(noisy_executor.run(program).result())

fig, ax = plt.subplots(figsize=(7, 4))
ax.plot(x_vals, exact_evs.mean(axis=0), "k-", label="Exact")
ax.errorbar(
    x_vals,
    np.mean(unmitigated_res.data.evs, axis=0),
    np.mean(unmitigated_res.data.stds, axis=0),
    fmt="o-",
    label="Unmitigated",
)
ax.set_xlabel(r"$\theta/\pi$")
ax.set_ylabel(r"$\langle Z \rangle_{\text{avg}}$")
ax.legend()
plt.tight_layout()
plt.show()

Siapkan program ZNE dengan digital gate folding dan sampel sirkuit dengan skala noise yang bervariasi

Siapkan program ZNE yang berisi sirkuit dengan skala noise yang makin besar. Setiap faktor noise mewakili kelipatan gate entangling yang dimuat sirkuit (faktor noise 3.0 memuat gate entangling tiga kali lebih banyak daripada sirkuit aslinya). Sirkuit tiap faktor noise di-twirl dengan jumlah randomisasi dan shot yang sama seperti program yang tidak dimitigasi.

from qiskit_mitigation import GateFolding

task_gf = GateFolding()
program_zne = task_gf.prepare(
    circuit=circuit,
    observables=observables,
    parameters=parameter_values,
    num_randomizations=num_randomizations,
    shots_per_randomization=shots_per_randomization,
    custom_boxing_options=boxing_options,
    noise_factors=[1, 3, 5],  # ZNE option
)

Gambar Trotter step pertama dari template hasil twirl dengan faktor noise 1 dan faktor noise 3 untuk memvisualisasikan amplifikasinya. Untuk faktor noise 3.0, setiap layer sirkuit diganti dengan layer itu, inversnya, lalu layer itu lagi.

import IPython.display

def first_trotter_step(template, noise_factor):
    """The template up to the closing barrier of its first Trotter step (two layers, each folded ``noise_factor`` times)."""
    toy = template.copy_empty_like()
    for instruction in template.data:
        toy.append(instruction)
        if (
            instruction.operation.name == "barrier"
            and instruction.operation.label.startswith(
                f"R{2 * noise_factor - 1}@"
            )
        ):
            break
    return toy

IPython.display.display(
    first_trotter_step(program_zne.items[0].circuit, 1).draw(
        "mpl", scale=0.7, plot_barriers=False
    )
)
first_trotter_step(program_zne.items[1].circuit, 3).draw(
    "mpl", fold=-1, plot_barriers=False
)

Terakhir, jalankan program pada simulator yang noisy. Hasilnya menjadi masukan untuk ekstrapolasi di langkah berikutnya.

program_zne._semantic_role = (
    "simulation"  # flag the program for the simulator
)
job_zne_results = noisy_executor.run(program_zne).result()

Ekstrapolasi ke batas zero-noise dengan fungsi fitting bawaan dan kustom

Setelah punya data untuk tiap skala noise, kamu bisa mencocokkan sebuah fungsi ke data dan mengekstrapolasi kembali ke nilai zero-noise untuk memperkirakan nilai ekspektasi observable bebas noise. Ada banyak cara untuk fitting data, dan pilihan fungsi yang tepat bergantung pada sifat data yang disampel. GateFolding.postprocess menyediakan beberapa fungsi ekstrapolasi bawaan, termasuk linear, exponential, dan double_exponential. Metode postprocess menerima daftar metode ekstrapolasi, mengevaluasinya secara berurutan, dan mengembalikan nilai pertama dengan nilai yang valid secara fisik dan varians yang berhingga. Artinya urutan ekstrapolator yang kamu tentukan sangat berpengaruh pada ekstrapolator mana yang dipilih. Mengirim fallback sebagai ekstrapolator berarti kalau tidak ada jawaban valid yang ditemukan, nilai ekspektasi dari faktor noise terkecil akan dipilih dan dikembalikan. Kalau fallback tidak dikirim dan tidak ada ekstrapolator yang menghasilkan nilai valid secara fisik, hasil dengan std terkecil yang dikembalikan. Di bawah noise Pauli, sebuah observable meluruh secara eksponensial terhadap faktor noise, jadi dengan 30 layer (15 Trotter step, 2 layer per step) gate entangling, ekstrapolator eksponensial adalah pilihan yang wajar; kamu juga bisa mengekstrapolasi secara linear, yang umum dipilih ketika peluruhannya kecil.

Untuk memakai fungsi fitting kustom, kirim ekstrapolator "custom" bersama tuple custom_fit yang berisi fungsi fitting, daftar tebakan awal untuk tiap parameternya, dan batas-batasnya. Semua ini diteruskan ke scipy.optimize.curve_fit. Post-processing bisa diulang dengan ekstrapolator berbeda tanpa menjalankan ulang programnya. Di sini kamu memakai fungsi fitting akar kuadrat kustom. Untuk memvisualisasikan ekstrapolasinya, plot magnetisasi rata-rata situs yang dihitung pada tiap faktor noise untuk $\theta = 0.15\pi$, ketiga fit, dan rata-rata situs dari nilai ekspektasi hasil ekstrapolasi di titik zero-noise. Ketiga titik data tampak hampir linear, namun ekstrapolasi linear jatuh jelas di bawah nilai eksak. Seperti yang diprediksi, fit eksponensial paling efektif memperkirakan nilai ekspektasi tanpa noise.

from scipy.optimize import curve_fit

zne_mitigated_exp = task_gf.postprocess(
    job_zne_results, extrapolator=["exponential", "linear"]
)
zne_mitigated_lin = task_gf.postprocess(
    job_zne_results, extrapolator=["linear"]
)

def sqrt_fit(x, a, b):
    return a * np.sqrt(x) + b

custom_fit = (sqrt_fit, [1.0, 1.0], (-np.inf, np.inf))
zne_mitigated_custom = task_gf.postprocess(
    job_zne_results,
    extrapolator=["custom", "exponential", "linear"],
    custom_fit=custom_fit,
)

# The result holds the expectation value at each noise factor and the extrapolated value, but
# not the fit parameters, so we refit the same models here to draw the curves.
noise_factors = np.array([1, 3, 5])
angle = 3  # theta = 0.15 pi, the angle with the largest unmitigated error
amplifications = zne_mitigated_exp.data.evs_noise_factors[:, angle, :].mean(
    axis=0
)
x_fit = np.linspace(0, 5, 51)
fits = {
    "exponential": (
        lambda x, a, b: a * np.exp(b * x),
        zne_mitigated_exp,
        "black",
    ),
    "linear": (lambda x, a, b: a * x + b, zne_mitigated_lin, "blue"),
    "square-root": (sqrt_fit, zne_mitigated_custom, "green"),
}

plt.plot(
    noise_factors, amplifications, "o", color="red", label="amplifications"
)
for name, (fit, res, color) in fits.items():
    params, _ = curve_fit(fit, noise_factors, amplifications, p0=[1.0, -0.1])
    plt.plot(
        x_fit, fit(x_fit, *params), "--", color=color, label=f"{name} fit"
    )
    plt.plot(
        [0],
        res.data.evs[:, angle].mean(),
        "*",
        color=color,
        markersize=15,
        label=f"{name} extrapolation",
    )
plt.plot(
    [0],
    exact_evs[:, angle].mean(),
    "_",
    color="red",
    markersize=15,
    markeredgewidth=3,
    label="exact",
)
plt.xlabel("noise factor")
plt.ylabel(r"$\langle Z \rangle_{\text{avg}}$")
plt.legend()
plt.show()

Bandingkan nilai ekspektasi hasil mitigasi

Plot magnetisasi rata-rata situs dari tiap metode sepanjang sapuan sudut rotasi Ising. Ekstrapolasi eksponensial memulihkan magnetisasi eksak dalam batas error bar-nya, sedangkan ekstrapolasi linear menghilangkan sebagian besar bias tetapi konsisten meremehkan nilainya, karena peluruhan sebenarnya bersifat eksponensial. Fit akar kuadrat melenceng jauh ke atas.

results = {
    "Unmitigated": unmitigated_res,
    "ZNE exponential": zne_mitigated_exp,
    "ZNE linear": zne_mitigated_lin,
    "ZNE square-root": zne_mitigated_custom,
}

fig, ax = plt.subplots(figsize=(7, 4))
ax.plot(x_vals, exact_evs.mean(axis=0), "k-", label="Exact")
for label, res in results.items():
    ax.errorbar(
        x_vals,
        np.mean(res.data.evs, axis=0),
        np.mean(res.data.stds, axis=0),
        fmt="o-",
        label=label,
    )
ax.set_xlabel(r"$\theta/\pi$")
ax.set_ylabel(r"$\langle Z \rangle_{\text{avg}}$")
ax.legend()
ax.set_title(
    f"Error-Mitigated Magnetization [{num_qubits}-qubit, {trotter_steps} Trotter steps]"
)
plt.tight_layout()
plt.show()