Kvantový šum a zmírnění chyb

poznámka

Toshinari Itoko (28. června 2024)

Stáhni si pdf původní přednášky. Pozor, některé úryvky kódu mohou být zastaralé, protože se jedná o statické obrázky.

Přibližný čas QPU pro spuštění tohoto experimentu je 1 m 40 s.

1. Úvod

V této lekci prozkoumáme šum a způsoby, jak jej lze na kvantových počítačích zmírnit. Začneme tím, že se podíváme na účinky šumu pomocí simulátoru, který dokáže šum simulovat několika způsoby, včetně použití profilů šumu ze skutečných kvantových počítačů. Poté přejdeme ke skutečným kvantovým počítačům, v nichž je šum neodmyslitelný. Podíváme se na účinky zmírnění chyb, včetně kombinací metod jako extrapolace k nulovému šumu (ZNE) a gate-twirling.

Začneme tím, že načteme několik balíčků.

# Added by doQumentation — required packages for this notebook
!pip install -q matplotlib qiskit qiskit-aer qiskit-ibm-runtime

# !pip install qiskit qiskit_aer qiskit_ibm_runtime
# !pip install jupyter
# !pip install matplotlib pylatexenc

import qiskit

qiskit.__version__

'2.0.2'

import qiskit_aer

qiskit_aer.__version__

'0.17.1'

import qiskit_ibm_runtime

qiskit_ibm_runtime.__version__

'0.40.1'

2. Zašuměná simulace bez zmírnění chyb

Qiskit Aer je klasický simulátor pro kvantové výpočty. Dokáže simulovat nejen ideální, ale i zašuměné provádění kvantových obvodů. Tento notebook ukazuje, jak spustit zašuměnou simulaci pomocí Qiskit Aer:

1. Sestav model šumu
2. Sestav zašuměný Sampler (simulátor) s modelem šumu
3. Spusť kvantový obvod na zašuměném Sampleru

noise_model = NoiseModel()
...
noisy_sampler = Sampler(options={"backend_options": {"noise_model": noise_model}})
job = noisy_sampler.run([circuit])

2.1 Sestavení testovacího obvodu

Uvažujeme hračkové 1-qubitové obvody, které jednoduše opakují hradla X d-krát (d=0 ... 100) a měří pozorovatelnou veličinu Z.

from qiskit.circuit import QuantumCircuit

MAX_DEPTH = 100
circuits = []
for d in range(MAX_DEPTH + 1):
    circ = QuantumCircuit(1)
    for _ in range(d):
        circ.x(0)
        circ.barrier(0)
    circ.measure_all()
    circuits.append(circ)

display(circuits[3].draw(output="mpl"))

from qiskit.quantum_info import SparsePauliOp

obs = SparsePauliOp.from_list([("Z", 1.0)])
obs

SparsePauliOp(['Z'],
              coeffs=[1.+0.j])

2.2 Sestavení modelu šumu

Abychom mohli provádět zašuměnou simulaci, potřebujeme zadat NoiseModel. V této části ukážeme, jak NoiseModel sestavit. Nejprve musíme definovat kvantové (nebo měřicí) chyby, které do modelu šumu přidáme.

from qiskit_aer.noise.errors import (
    coherent_unitary_error,
    amplitude_damping_error,
    ReadoutError,
)
from qiskit.circuit.library import RXGate

# Coherent (unitary) error: Over X-rotation error
# https://qiskit.github.io/qiskit-aer/stubs/qiskit_aer.noise.coherent_unitary_error.html#qiskit_aer.noise.coherent_unitary_error
OVER_ROTATION_ANGLE = 0.05
coherent_error = coherent_unitary_error(RXGate(OVER_ROTATION_ANGLE).to_matrix())

# Incoherent error: Amplitude dumping error
# https://qiskit.github.io/qiskit-aer/stubs/qiskit_aer.noise.amplitude_damping_error.html#qiskit_aer.noise.amplitude_damping_error
AMPLITUDE_DAMPING_PARAM = 0.02  # in [0, 1] (0: no error)
incoherent_error = amplitude_damping_error(AMPLITUDE_DAMPING_PARAM)

# Readout (measurement) error: Readout error
# https://qiskit.github.io/qiskit-aer/stubs/qiskit_aer.noise.ReadoutError.html#qiskit_aer.noise.ReadoutError
PREP0_MEAS1 = 0.03  # P(1|0): Probability of preparing 0 and measuring 1
PREP1_MEAS0 = 0.08  # P(0|1): Probability of preparing 1 and measuring 0
readout_error = ReadoutError(
    [[1 - PREP0_MEAS1, PREP0_MEAS1], [PREP1_MEAS0, 1 - PREP1_MEAS0]]
)

from qiskit_aer.noise import NoiseModel

noise_model = NoiseModel()
noise_model.add_quantum_error(coherent_error.compose(incoherent_error), "x", (0,))
noise_model.add_readout_error(readout_error, (0,))

2.3 Sestavení zašuměného sampleru s modelem šumu

from qiskit_aer.primitives import SamplerV2 as Sampler

noisy_sampler = Sampler(options={"backend_options": {"noise_model": noise_model}})

2.4 Spuštění kvantových obvodů na zašuměném sampleru

job = noisy_sampler.run(circuits, shots=400)

result = job.result()

result[0].data.meas.get_counts()

{'0': 389, '1': 11}

2.5 Vykreslení výsledků

import matplotlib.pyplot as plt

plt.title("Noisy simulation")
ds = list(range(MAX_DEPTH + 1))
plt.plot(
    ds,
    [result[d].data.meas.expectation_values(["Z"]) for d in ds],
    color="gray",
    linestyle="-",
)
plt.scatter(ds, [result[d].data.meas.expectation_values(["Z"]) for d in ds], marker="o")
plt.hlines(0, xmin=0, xmax=MAX_DEPTH, colors="black")
plt.ylim(-1, 1)
plt.xlabel("Circuit depth")
plt.ylabel("Measured <Z>")
plt.show()

2.6 Ideální simulace

ideal_sampler = Sampler()
job_ideal = ideal_sampler.run(circuits)
result_ideal = job_ideal.result()
plt.title("Ideal simulation")
ds = list(range(MAX_DEPTH + 1))
plt.plot(
    ds,
    [result_ideal[d].data.meas.expectation_values(["Z"]) for d in ds],
    color="gray",
    linestyle="-",
)
plt.scatter(
    ds, [result_ideal[d].data.meas.expectation_values(["Z"]) for d in ds], marker="o"
)
plt.hlines(0, xmin=0, xmax=MAX_DEPTH, colors="black")
plt.xlabel("Circuit depth")
plt.ylabel("Measured <Z>")
plt.show()

2.7 Cvičení

Úpravou níže uvedeného kódu:

• Zkus 25násobný počet vzorků (= 10_000 shots) a ověř, že získáš hladší graf
• Změň parametry šumu (OVER_ROTATION_ANGLE, AMPLITUDE_DAMPING_PARAM, PREP0_MEAS1 nebo PREP1_MEAS0) a sleduj, jak se graf mění

OVER_ROTATION_ANGLE = 0.05
coherent_error = coherent_unitary_error(RXGate(OVER_ROTATION_ANGLE).to_matrix())
AMPLITUDE_DAMPING_PARAM = 0.02  # in [0, 1] (0: no error)
incoherent_error = amplitude_damping_error(AMPLITUDE_DAMPING_PARAM)
PREP0_MEAS1 = 0.1  # P(1|0): Probability of preparing 0 and measuring 1
PREP1_MEAS0 = 0.05  # P(0|1): Probability of preparing 1 and measuring 0
readout_error = ReadoutError(
    [[1 - PREP0_MEAS1, PREP0_MEAS1], [PREP1_MEAS0, 1 - PREP1_MEAS0]]
)
noise_model = NoiseModel()
noise_model.add_quantum_error(coherent_error.compose(incoherent_error), "x", (0,))
noise_model.add_readout_error(readout_error, (0,))
options = {
    "backend_options": {"noise_model": noise_model},
}
noisy_sampler = Sampler(options=options)
job = noisy_sampler.run(circuits, shots=400)
result = job.result()
plt.title("Noisy simulation")
ds = list(range(MAX_DEPTH + 1))
plt.plot(
    ds,
    [result[d].data.meas.expectation_values(["Z"]) for d in ds],
    marker="o",
    linestyle="-",
)
plt.hlines(0, xmin=0, xmax=MAX_DEPTH, colors="black")
plt.ylim(-1, 1)
plt.xlabel("Depth")
plt.ylabel("Measured <Z>")
plt.show()

2.8 Realističtější zašuměná simulace

from qiskit_aer import AerSimulator
from qiskit_ibm_runtime import SamplerV2 as Sampler, QiskitRuntimeService

service = QiskitRuntimeService()
real_backend = service.least_busy(
    operational=True, simulator=False, min_num_qubits=127
)  # Eagle

<IBMBackend('ibm_strasbourg')>

aer = AerSimulator.from_backend(real_backend)
noisy_sampler = Sampler(mode=aer)
job = noisy_sampler.run(circuits)
result = job.result()
plt.title("Noisy simulation with noise model from real backend")
ds = list(range(MAX_DEPTH + 1))
plt.plot(
    ds,
    [result[d].data.meas.expectation_values(["Z"]) for d in ds],
    marker="o",
    linestyle="-",
)
plt.hlines(0, xmin=0, xmax=MAX_DEPTH, colors="black")
plt.ylim(-1, 1)
plt.xlabel("Depth")
plt.ylabel("Measured <Z>")
plt.show()

3. Reálný kvantový výpočet se zmírněním chyb

V této části si ukážeme, jak pomocí Qiskit Estimatoru získat výsledky (střední hodnoty) se zmírněním chyb. Uvažujeme 6qubitové Trotterizované obvody pro simulaci časového vývoje jednorozměrného Isingova modelu a podíváme se, jak chyba škáluje v závislosti na počtu časových kroků.

backend = service.least_busy(
    operational=True, simulator=False, min_num_qubits=127
)  # Eagle
backend

<IBMBackend('ibm_strasbourg')>

NUM_QUBITS = 6
NUM_TIME_STEPS = list(range(8))
RX_ANGLE = 0.1
RZZ_ANGLE = 0.1

3.1 Sestavení obvodů

# Build circuits with different number of time steps
circuits = []
for n_steps in NUM_TIME_STEPS:
    circ = QuantumCircuit(NUM_QUBITS)
    for i in range(n_steps):
        # rx layer
        for q in range(NUM_QUBITS):
            circ.rx(RX_ANGLE, q)
        # 1st rzz layer
        for q in range(1, NUM_QUBITS - 1, 2):
            circ.rzz(RZZ_ANGLE, q, q + 1)
        # 2nd rzz layer
        for q in range(0, NUM_QUBITS - 1, 2):
            circ.rzz(RZZ_ANGLE, q, q + 1)
    circ.barrier()  # need not to optimize the circuit
    # Uncompute stage
    for i in range(n_steps):
        for q in range(0, NUM_QUBITS - 1, 2):
            circ.rzz(-RZZ_ANGLE, q, q + 1)
        for q in range(1, NUM_QUBITS - 1, 2):
            circ.rzz(-RZZ_ANGLE, q, q + 1)
        for q in range(NUM_QUBITS):
            circ.rx(-RX_ANGLE, q)
    circuits.append(circ)

Abychom znali ideální výstup předem, používáme obvody typu compute-uncompute, které se skládají z první fáze, kde je aplikován původní obvod $U$, a druhé fáze, kde je obrácen $U^\dagger$. Všimni si, že ideální výsledek takových obvodů bude triviálně vstupní stav $|000000\rangle$, který má triviální střední hodnoty pro libovolné Pauliho pozorovatelné, například $\langle IIIIIZ \rangle = 1$.

# Print the circuit with 2 time steps
circuits[2].draw(output="mpl")

Poznámka: Jak je ukázáno výše, obvod s $k$ časovými kroky bude mít $4k$ vrstev dvouqubitových Gate.

obs = SparsePauliOp.from_sparse_list([("Z", [0], 1.0)], num_qubits=NUM_QUBITS)
obs

SparsePauliOp(['IIIIIZ'],
              coeffs=[1.+0.j])

3.2 Transpilace obvodů

Obvody transpilujeme pro Backend s optimalizací (optimization_level=1).

from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager

pm = generate_preset_pass_manager(optimization_level=1, backend=backend)
isa_circuits = pm.run(circuits)
display(isa_circuits[2].draw("mpl", idle_wires=False, fold=-1))

3.3 Spuštění pomocí Estimatoru (s různými úrovněmi odolnosti)

Nastavení úrovně odolnosti (estimator.options.resilience_level) je nejjednodušší způsob, jak při použití Qiskit Estimatoru aplikovat zmírnění chyb. Estimator podporuje následující úrovně odolnosti (ke dni 28. 6. 2024). Další podrobnosti najdeš v příručce Configure error mitigation.

from qiskit_ibm_runtime import Batch
from qiskit_ibm_runtime import EstimatorV2 as Estimator

jobs = []
job_ids = []
with Batch(backend=backend):
    for resilience_level in [0, 1, 2]:
        estimator = Estimator()
        estimator.options.resilience_level = resilience_level
        job = estimator.run(
            [(circ, obs.apply_layout(circ.layout)) for circ in isa_circuits]
        )
        job_ids.append(job.job_id())
        print(f"Job ID (rl={resilience_level}): {job.job_id()}")
        jobs.append(job)

Job ID (rl=0): d146vcnmya70008emprg
Job ID (rl=1): d146vdnqf56g0081sva0
Job ID (rl=2): d146ven5z6q00087c61g

# check job status
for job in jobs:
    print(job.status())

DONE
DONE
DONE

# REPLACE WITH YOUR OWN JOB IDS
jobs = [service.job(job_id) for job_id in job_ids]

# Get results
results = [job.result() for job in jobs]

3.4 Vykreslení výsledků

plt.title("Error mitigation with different resilience levels")
labels = ["0 (No mitigation)", "1 (TREX)", "2 (ZNE + Gate twirling)"]
steps = NUM_TIME_STEPS
for result, label in zip(results, labels):
    plt.errorbar(
        x=steps,
        y=[result[s].data.evs for s in steps],
        yerr=[result[s].data.stds for s in steps],
        marker="o",
        linestyle="-",
        capsize=4,
        label=label,
    )
plt.hlines(
    1.0, min(steps), max(steps), linestyle="dashed", label="Ideal", colors="black"
)
plt.xlabel("Time steps")
plt.ylabel("Mitigated <IIIIIZ>")
plt.legend()
plt.show()

4. (Volitelné) Přizpůsobení možností zmírnění chyb

Aplikaci technik zmírnění chyb můžeme přizpůsobit pomocí možností, jak je ukázáno níže.

# TREX
estimator.options.twirling.enable_measure = True
estimator.options.twirling.num_randomizations = "auto"
estimator.options.twirling.shots_per_randomization = "auto"

# Gate twirling
estimator.options.twirling.enable_gates = True
# ZNE
estimator.options.resilience.zne_mitigation = True
estimator.options.resilience.zne.noise_factors = [1, 3, 5]
estimator.options.resilience.zne.extrapolator = ("exponential", "linear")

# Dynamical decoupling
estimator.options.dynamical_decoupling.enable = True  # Default: False
estimator.options.dynamical_decoupling.sequence_type = "XX"

# Other options
estimator.options.default_shots = 10_000

Podrobnosti o možnostech zmírnění chyb najdeš v následujících průvodcích a v referenci API.

• Konfigurace zmírnění chyb
• Úvod do možností
• EstimatorOptions
• SamplerOptions