Przejdź do głównej treści

Symulacja rozpraszania neutronów za pomocą workflow Serverless AQC + dynamika Trottera

Szacowany czas użycia: 18 minut na procesorze Heron r3 (UWAGA: To tylko szacunek. Twój czas wykonania może się różnić.)

Efekty uczenia się​

  • Jak nieelastyczne widmo rozpraszania neutronów mapuje się na dynamiczny czynnik struktury S(q,ω)S(q, \omega) jednowymiarowego magnesu kwantowego.

  • Jak przygotować stan podstawowy KCuF3_3 (izotropowy Heisenberg) za pomocą grupy renormalizacji macierzy gęstości (DMRG) i maksymalizacji wierności matrix product state (MPS).

  • Jak uruchomić ewolucję czasową Trottera, kompresję circuit przybliżonej kompilacji kwantowej (AQC) i wykonanie z mitygacją błędów jako pojedyncze wywołanie funkcji.

  • Jak przetworzyć szereg czasowy ⟨σz⟩(t)\langle \sigma_z \rangle(t) per site na S(q,ω)S(q, \omega) i zidentyfikować kontinuum dwuspinonowe.

Wymagania wstępne​

Tło​

Nieelastyczne rozpraszanie neutronów mierzy dynamiczny czynnik struktury S(q,ω)S(q, \omega), przekształcenie Fouriera w przestrzeni i czasie funkcji korelacji spin-spin, więc odtworzenie S(q,ω)S(q, \omega) na podstawie mikroskopowego modelu spinowego stanowi bezpośredni, falsyfikowalny test symulacji kwantowej. Ten samouczek analizuje KCuF3_3, antyferromagnetyczny łańcuch Heisenberga o spinie 12\frac{1}{2}, którego wzbudzenia to nie pojedyncze przerzuty spinu, lecz pary sfrakcjonowanych spinonów: zamiast ostrej dyspersji magnonu, S(q,ω)S(q, \omega) pokazuje szerokie kontinuum dwuspinonowe, ograniczone od dołu przez π2∣sin⁡q∣\tfrac{\pi}{2}|\sin q|, a od góry przez π∣sin⁡(q/2)∣\pi|\sin(q/2)|. To są przerywane krzywe na poniższych wykresach. Pełna fizyka oraz porównanie z zmierzonymi danymi neutronowymi są omówione w oryginalnym samouczku oraz w pracy Lee i in., arXiv:2603.15608.

Workflow kwantowy odzwierciedla eksperyment rozpraszania:

  1. Przygotuj stan podstawowy łańcucha ∣ψ0⟩|\psi_0\rangle.

  2. Wprowadź lokalne zaburzenie w centralnym site w postaci rotacji π/2\pi/2 ZZ, naśladując transfer pędu i energii od neutronu.

  3. Ewoluuj czasowo pod hamiltonianem Heisenberga, e−iHte^{-iHt}, za pomocą wzoru produktowego Trottera.

  4. Zmierz magnetyzację per site ⟨σzj⟩(t)\langle \sigma_z^j \rangle(t). Jako funkcja site jj i czasu tt, jest to dokładnie retardowana funkcja Greena GR(j,jc,t)G^R(j, j_c, t), więc przed przekształceniem Fouriera w kroku 5 nie jest potrzebna żadna konwersja.

  5. Przekształć Fourierem GRG^R na S(q,ω)S(q, \omega).

Problemy mogą pojawić się w kroku 3, gdy dokładne circuits Trottera dla długich ewolucji stają się zbyt głębokie dla sprzętu. AQC z sieciami tensorowymi rozwiązuje to poprzez skompresowanie bloku kroków Trottera w stały, płytki sparametryzowany ansatz, którego wierność stanu w stosunku do dokładnej ewolucji jest maksymalizowana klasycznie za pomocą symulatora MPS (arXiv:2301.08609). AQC Dynamics Template pakuje cały ten kwantowy rdzeń (synteza Trottera, kompresja AQC i wykonanie z mitygacją) w jedno wywołanie:

PRZED (ten notebook)FUNKCJA (aqc-dynamics-function)PO (ten notebook)
Stan podstawowy z DMRG plus maksymalizacja wierności MPS, z wbudowanym kopnięciem neutronowym w tym samym circuitSynteza Trottera → kompresja AQC → wykonanie na statevector, fake lub runtime, zwracające ⟨σzj⟩(t)\langle \sigma_z^j \rangle(t) per siteS(q,ω)S(q, \omega), dynamiczny czynnik struktury

Praca specyficzna dla eksperymentu pozostaje tutaj, w notebooku: przygotowanie stanu podstawowego (PRZED) i post-processing S(q,ω)S(q, \omega) (PO). Dwa najbardziej kwantowo intensywne kroki, kompresja i wykonanie, uruchamiane są wewnątrz funkcji.

Ten samouczek jest uzupełnieniem Symulacji rozpraszania neutronów w materiałach kwantowych za pomocą circuits kwantowych, który buduje ten sam eksperyment w linii: ten sam model KCuF3_3, przygotowanie stanu podstawowego, kopnięcie neutronowe i post-processing, z syntezą Trottera, kompresją AQC i wykonaniem z mitygacją błędów zapisanymi krok po kroku. Przeczytaj tamten samouczek, aby dowiedzieć się, jak działa kompresja AQC. Przeczytaj ten, aby uruchomić ten sam eksperyment poprzez wdrożony szablon funkcji: kwantowy rdzeń staje się pojedynczym wywołaniem funkcji, a wielogodzinna kompresja AQC uruchamiana jest wewnątrz workera Serverless zamiast na Twoim komputerze, więc podczas jej wykonywania nie potrzebujesz systemu HPC ani otwartego kernela. To samo wywołanie napędza również inne jednowymiarowe eksperymenty dynamiki.

Wymagania​

Zanim rozpoczniesz ten samouczek, upewnij się, że masz następujące elementy:

  • Funkcja wdrożona na Twoim koncie Qiskit Serverless. Najpierw uruchom towarzyszący szablon funkcji: Wdrażanie i uruchamianie szablonu funkcji dynamiki AQC + Trotter. Ten przewodnik przeprowadzi Cię przez pobranie plików źródłowych i przesłanie funkcji na Twoje konto. Ten samouczek jedynie wywołuje wdrożoną funkcję.

  • Zapisane poświadczenia IBM Quantum® dla QiskitServerless (zobacz szablon funkcji). Oba przykłady w tym samouczku wywołują wdrożoną funkcję, więc oba ich potrzebują.

  • Qiskit SDK w wersji 2.0 lub nowszej (pip install qiskit).

  • Klient Qiskit IBM Catalog (pip install qiskit-ibm-catalog).

  • NumPy, SciPy i Matplotlib (pip install numpy scipy matplotlib). SciPy 1.14 lub nowszy jest wymagany dla optymalizatora COBYQA używanego w przygotowaniu stanu podstawowego.

  • Stos sieci tensorowych AQC, ponieważ przygotowanie stanu podstawowego w kroku 1 uruchamia się lokalnie w tym notebooku: pip install 'qiskit-addon-aqc-tensor[quimb-jax]==0.3.1'.

Pierwsze wywołanie nowo wdrożonej funkcji czeka, podczas gdy worker Serverless instaluje swoje zależności, więc spodziewaj się dodatkowego opóźnienia przy tym uruchomieniu.

Konfiguracja​

Zaimportuj biblioteki i zdefiniuj pomocnicze funkcje specyficzne dla eksperymentu, używane później: build_gs_ansatz (wariacyjny ansatz hamiltonianu, czyli HVA, do przygotowania stanu podstawowego), prepare_ground_state (DMRG plus maksymalizacja wierności MPS) oraz get_spectrum, plot_green i plot_spectrum (post-processing S(q,ω)S(q, \omega)). Zostały one zaadaptowane z oryginalnego samouczka o rozpraszaniu neutronów.

# Added by doQumentation — required packages for this notebook
!pip install -q matplotlib numpy qiskit qiskit-addon-aqc-tensor qiskit-ibm-catalog quimb scipy
from functools import partial

import matplotlib.pyplot as plt
import numpy as np
import scipy.optimize

import quimb.tensor as qtn
from qiskit import QuantumCircuit
from qiskit.quantum_info import SparsePauliOp
from qiskit_addon_aqc_tensor.simulation import tensornetwork_from_circuit
from qiskit_addon_aqc_tensor.simulation.quimb import QuimbSimulator
from qiskit_ibm_catalog import QiskitServerless
#  Dynamical structure factor via discrete Fourier transform

def get_spectrum(n, Gjjc, dt, time_steps, q_steps, w_steps):
    """Compute the dynamical structure factor from the retarded Green's function.

    Uses the center-site approximation and a discrete Fourier transform.
    """
    green = Gjjc / 4  # sigma -> S=1/2
    omega_max = np.pi / dt
    qpoints = np.arange(0, 2 * np.pi, 2 * np.pi / q_steps)
    omegas = np.arange(0, omega_max, omega_max / w_steps)
    green_map = np.zeros((omegas.shape[0], qpoints.shape[0]))
    center = n // 2 - 1
    for iw, w in enumerate(omegas):
        exponent = np.exp(1j * w * dt * np.arange(1, time_steps + 1))
        S_w = np.dot(green.T, exponent) * dt
        for iq, q in enumerate(qpoints):
            q_matrix = np.exp(-1j * q * np.arange(-center, center + 2, 1))
            green_map[iw, iq] = np.imag(np.dot(S_w, q_matrix))
    return green_map

#  Plotting helpers

def plot_spectrum(
    dsf,
    dt,
    q_steps,
    w_steps,
    lower_bound=False,
    upper_bound=False,
    title=None,
):
    """Heat-map of the dynamical structure factor."""
    omega_max = np.pi / dt
    qpoints = np.arange(0, 2 * np.pi, 2 * np.pi / q_steps)
    omegas = np.arange(0, omega_max, omega_max / w_steps)
    x, y = np.meshgrid(qpoints, omegas)
    fig, ax = plt.subplots(figsize=(8, 5))
    c = ax.pcolormesh(x, y, dsf / np.max(dsf), cmap="viridis", shading="auto")
    fig.colorbar(c, ax=ax, label="Normalized intensity")
    if lower_bound:
        ax.plot(
            qpoints,
            np.pi * np.abs(np.sin(qpoints)) / 2,
            "--",
            color="white",
            lw=1.5,
            label="Lower bound",
        )
    if upper_bound:
        ax.plot(
            qpoints,
            np.pi * np.abs(np.sin(qpoints / 2)),
            "--",
            color="red",
            lw=1.5,
            label="Upper bound",
        )
    ax.set_ylim(0, 3.6)
    ax.set_xlim(0, 2 * np.pi - 2 * np.pi / q_steps)
    ax.set_xlabel(r"$q$", fontsize=16)
    ax.set_ylabel(r"$\tilde{\omega} = \omega / J$", fontsize=16)
    ax.set_xticks([0, np.pi / 2, np.pi, 3 * np.pi / 2, 2 * np.pi])
    ax.set_xticklabels(["0", r"$\pi/2$", r"$\pi$", r"$3\pi/2$", r"$2\pi$"])
    if lower_bound or upper_bound:
        ax.legend(loc="upper right", fontsize=11)
    if title:
        ax.set_title(title, fontsize=14)
    plt.tight_layout()
    plt.show()

def plot_green(n, Gjjc, time_steps, dt, title=None):
    """Heat-map of the retarded Green's function in real space and time."""
    fig, ax = plt.subplots(figsize=(8, 6))
    t_axis = np.arange(1, time_steps + 1) * dt
    site_axis = np.arange(n)
    x, y = np.meshgrid(t_axis, site_axis)
    c = ax.pcolormesh(
        x,
        y,
        np.real(Gjjc).T,
        cmap="RdBu",
        vmax=0.5,
        vmin=-0.5,
        shading="auto",
    )
    fig.colorbar(c, ax=ax, label=r"Re $G^R(j, j_c, t)$")
    ax.set_xlabel(r"Time  ($t / J^{-1}$)", fontsize=16)
    ax.set_ylabel("Site index $j$", fontsize=16)
    if title:
        ax.set_title(title, fontsize=14)
    plt.tight_layout()
    plt.show()

#  Variational ground-state ansatz (HVA)

def _apply_xxz_pair_gate(qc, q0, q1, theta):
    """Apply the parameterized XXZ-type two-qubit gate used in the HVA."""
    qc.cx(q0, q1)
    qc.rz(theta, q1)
    qc.h(q0)
    qc.rz(theta + np.pi / 2, q0)
    qc.cx(q0, q1)
    qc.rz(-theta, q1)
    qc.h(q1)
    qc.cx(q1, q0)
    qc.rz(np.pi / 2, q1)
    qc.rz(-np.pi / 2, q0)
    qc.h(q1)
    qc.h(q0)

def build_gs_ansatz(n, params, layers):
    """Build the Hamiltonian variational ansatz (HVA) circuit for
    ground-state preparation of the 1D Heisenberg model.

    Starts from a product of singlet pairs and applies alternating
    odd/even layers of parameterized XXZ gates. For layer r,
    params[2 * r] is the odd-layer (inter-pair) angle and
    params[2 * r + 1] is the even-layer (intra-pair) angle.
    """
    qc = QuantumCircuit(n)
    # Initial singlet product state
    for i in range(n // 2):
        qc.x(2 * i)
        qc.x(2 * i + 1)
        qc.h(2 * i + 1)
        qc.cx(2 * i + 1, 2 * i)
    # Variational layers
    for r in range(layers):
        for i in range(1, (n + 1) // 2):  # odd layer
            _apply_xxz_pair_gate(qc, 2 * i - 1, 2 * i, params[2 * r])
        for i in range(n // 2):  # even layer
            _apply_xxz_pair_gate(qc, 2 * i, 2 * i + 1, params[2 * r + 1])
    return qc

def prepare_ground_state(n, gs_layers=5, max_bond=128, cutoff=1e-8):
    """Prepare the KCuF3 (isotropic Heisenberg) ground state as a QuantumCircuit.

    Runs DMRG (quimb MPO + DMRG2) to get the chain's ground state, then optimizes
    the HVA angles to maximize the MPS overlap |<psi_ansatz|psi_DMRG>|^2. No exact
    diagonalization, so it scales to larger n.
    """
    J = Jz = 1.0
    builder = qtn.SpinHam1D(S=1 / 2)
    builder += J * 0.5, "+", "-"
    builder += J * 0.5, "-", "+"
    builder += Jz, "Z", "Z"
    H_mpo = builder.build_mpo(L=n)
    dmrg = qtn.DMRG2(H_mpo)
    dmrg.solve(tol=1e-8, verbosity=0)

    gs_sim = QuimbSimulator(
        quimb_circuit_factory=partial(
            qtn.CircuitMPS, gate_opts=dict(cutoff=cutoff, max_bond=max_bond)
        ),
        autodiff_backend="jax",
    )

    def gs_infidelity(params):
        psi = tensornetwork_from_circuit(
            build_gs_ansatz(n, params, gs_layers), gs_sim
        ).psi
        return 1 - abs(psi.H @ dmrg.state) ** 2

    # Seed and optimizer match the original tutorial. Each layer starts at
    # [0, pi/2]: an odd-layer angle of 0 makes the inter-pair gate the identity,
    # and an even-layer angle of pi/2 makes the intra-pair gate a SWAP (since
    # 0.5 * (XX + YY + ZZ) = SWAP - I/2). That puts the seed at the singlet-pair
    # product limit, which is already a decent approximation to the Heisenberg
    # ground state, so the optimizer only has to refine it. The small jitter
    # (fixed RNG seed, so runs are reproducible) breaks the exact symmetry
    # between layers; COBYQA then runs for up to 100 iterations.
    rng = np.random.default_rng(12345)
    x0 = np.tile([0.0, np.pi / 2], gs_layers) + rng.normal(
        scale=0.1, size=2 * gs_layers
    )
    result_gs = scipy.optimize.minimize(
        gs_infidelity, x0, method="COBYQA", options={"maxiter": 100}
    )
    print(f"DMRG ground-state energy: {dmrg.energy:.6f}")
    print(f"GS fidelity: {1 - result_gs.fun:.4f}")
    return build_gs_ansatz(n, result_gs.x, gs_layers)

print("Setup complete - helpers defined.")
Setup complete - helpers defined.

Wczytaj szablon funkcji​

Połącz się z Qiskit Serverless i wczytaj wdrożoną funkcję aqc-dynamics-function. Oba przykłady w tym samouczku wywołują ten sam uchwyt fn, więc funkcja jest wczytywana tylko raz, tutaj.

# Credentials are read from the account saved once via QiskitServerless.save_account(...)
serverless = QiskitServerless()
fn = serverless.load("aqc-dynamics-function")

Przykład na małą skalę na symulatorze​

Najpierw uruchamiamy pełny workflow na małym 10-site łańcuchu, używając dokładnego backendu statevector. Waliduje to potok PRZED → FUNKCJA → PO przed poświęceniem jakiegokolwiek czasu QPU.

Krok 1: Zmapuj klasyczne dane wejściowe na problem kwantowy​

Zbuduj hamiltonian KCuF3_3 jako SparsePauliOp (izotropowy Heisenberg: XX+YY+ZZXX + YY + ZZ przy sprzężeniu 14\tfrac14 na każdym wiązaniu najbliższych sąsiadów; ciągi są operatorami Pauliego, więc 14\tfrac14 daje sprzężenie spinu-12\frac{1}{2}). Przygotuj stan podstawowy za pomocą DMRG plus maksymalizacji wierności MPS, następnie wbuduj kopnięcie neutronowe: rotację π/2\pi/2 ZZ w centralnym site. Przygotowany circuit to to, co przekazujemy funkcji jako initial_state. Pozostawiamy observables na wartości domyślnej (Z per site), co jest dokładnie odczytem ⟨σzj⟩(t)\langle \sigma_z^j \rangle(t) potrzebnym w workflow neutronowym.

n = 10
dt = 0.6  # physical time per Trotter step (also the omega-axis unit in POST)
time_steps = 10
center = n // 2 - 1

# MPS-simulator settings, shared by the ground-state prep here and the AQC
# compression inside the function (matches the original tutorial).
mps_max_bond = 32
mps_cutoff = 1e-8

# 1D isotropic Heisenberg (KCuF3) Hamiltonian on n qubits
H = SparsePauliOp.from_sparse_list(
    [(p, [i, i + 1], 0.25) for i in range(n - 1) for p in ("XX", "YY", "ZZ")],
    num_qubits=n,
)

# Ground state (DMRG + fidelity max) + neutron kick baked into the same circuit
gs_circuit = prepare_ground_state(
    n, gs_layers=3, max_bond=mps_max_bond, cutoff=mps_cutoff
)
gs_circuit.rz(
    np.pi / 2, center
)  # exp(-i (pi/2)/2 Z_center): the neutron perturbation
print(
    f"Prepared {n}-qubit ground state with the neutron kick at site {center}."
)
DMRG ground-state energy: -4.258035
GS fidelity: 0.9841
Prepared 10-qubit ground state with the neutron kick at site 4.

Kroki 2 i 3: Kompresja i wykonanie za pomocą szablonu funkcji​

W ręcznie pisanym workflow są to dwa oddzielne etapy: optymalizacja circuits pod sprzęt (Krok 2) i ich wykonanie (Krok 3). Szablon funkcji łączy oba w jedno wywołanie. Wykonuje syntezę Trottera, kompresję AQC i transpilację pod sprzęt, a następnie uruchamia circuits (tutaj na dokładnym symulatorze, później z wbudowaną mitygacją błędów na sprzęcie). Dwa parametry dostrajające to aqc_segments (plan kompresji) i aqc_options (ustawienia MPS i optymalizatora). Każdy segment {"n_steps": k, "ansatz_steps": m} kompresuje k kolejnych kroków Trottera w ansatz zbudowany z celu Trottera o m krokach, a wszelkie kroki poza sum(n_steps) uruchamiane są jako zwykły Trotter. Wczesne kroki o niskim splątaniu kompresują się dobrze do płytkiego (ansatz_steps=1) ansatzu, więc tutaj kompresujemy pierwsze trzy kroki w jednowarstwowy ansatz, a kolejne dwa w głębszy dwuwarstwowy ansatz; pozostałych pięć z 10 kroków Trottera uruchamiane jest jako zwykły Trotter. Dla aqc_options odwzorowujemy oryginalny samouczek: wymiar wiązania MPS max_bond=32, cutoff=1e-8 i optymalizator L-BFGS-B ograniczony do 100 iteracji.

Wywołaj funkcję załadowaną w Setup. backend="statevector" uruchamia dokładną ścieżkę referencyjną: bez czasu QPU, z obwodami działającymi na dokładnym symulatorze wektora stanu wewnątrz workera serverless (nadal potrzebne jest zapisane konto Qiskit Serverless, aby ją wywołać). initial_state przenosi przygotowany stan podstawowy (w tym kopnięcie); observables jest pominięte, więc funkcja mierzy domyślne ZZ na miejsce.

job = fn.run(
    t_steps=time_steps,
    aqc_segments=[
        {
            "n_steps": 3,
            "ansatz_steps": 1,
        },  # early steps -> shallow 1-layer ansatz
        {
            "n_steps": 2,
            "ansatz_steps": 2,
        },  # later steps -> deeper 2-layer ansatz
    ],
    aqc_options={
        "max_bond": mps_max_bond,  # MPS bond dimension for AQC compression
        "cutoff": mps_cutoff,
        "optimizer_settings": {
            "method": "L-BFGS-B",
            "jac": True,
            "options": {"maxiter": 100},
        },
    },
    dt=dt,
    hamiltonian=H,
    initial_state=gs_circuit,  # prepared ground state including the neutron kick
    # observables omitted -> default per-site Z (the neutron sigma_z readout)
    backend="statevector",
)
print(job.status())  # rerun this cell until status says DONE
DONE
# The per-site <sigma_z>(t) the function returns is the retarded Green's function
# G(j, j_c, t). The workflow samples t = 1..time_steps, so drop the t = 0 row (the
# prepared+kicked state before any evolution) before post-processing.
result = job.result()
print(
    "AQC fidelities:",
    {k: round(v, 4) for k, v in result["metadata"]["aqc_fidelities"].items()},
)

ev = np.array(result["expectation_values"])
Gjjc = ev[1:]  # shape (time_steps, n)
print("Green's function shape:", Gjjc.shape)
AQC fidelities: {'1': 1.0, '2': 0.9999, '3': 0.9992, '4': 0.9998, '5': 0.9995}
Green's function shape: (10, 10)

Krok 4: Post-processing i zwrot wyniku w żądanym formacie klasycznym​

Przekształć Fourierem funkcję Greena na S(q,ω)S(q, \omega), symetryzuj lustrzanie i przytnij wartości ujemne: standardowy neutronowy post-processing. Odbicie lustrzane jest dokładne, ponieważ S(q,ω)=S(−q,ω)S(q, \omega) = S(-q, \omega) dla tego modelu, a wartości ujemne, które pozostają, są artefaktami przekształcenia Fouriera skończonego, dyskretnie próbkowanego szeregu czasowego, więc są przycinane do zera. W tym małym dokładnym uruchomieniu kontinuum dwuspinonowe jest rozpoznane tylko z grubsza, ale mechanizm jest identyczny jak w następującym uruchomieniu sprzętowym.

q_res, w_res = 100, 100
spectrum = get_spectrum(n, Gjjc, dt, time_steps, q_res, w_res)
spectrum = -(spectrum + spectrum[:, ::-1]) / 2  # mirror symmetry
spectrum = np.clip(spectrum, a_min=0, a_max=None)  # clip negatives

plot_green(
    n,
    Gjjc,
    time_steps,
    dt,
    title=f"Retarded Green's function - {n} qubits (AQC, statevector)",
)
plot_spectrum(
    spectrum,
    dt,
    q_res,
    w_res,
    lower_bound=True,
    upper_bound=True,
    title=f"Dynamical structure factor - {n} qubits (AQC, statevector)",
)

Output of the previous code cell

Output of the previous code cell

Przykład na dużą skalę na sprzęcie​

Ten sam workflow skaluje się bez zmiany jakiegokolwiek kodu naukowego: 30-site łańcuch, dwukrotna głębokość Trottera (20 kroków), plan kompresji, który różnicuje głębokość ansatzu (głębszy ansatz dla późniejszych, bardziej splątanych kroków), oraz wykonanie na procesorze IBM Quantum z wbudowaną mitygacją błędów funkcji (dynamiczne rozsprzęganie, Pauli twirling i twirled readout error extinction (TREX)). Przechodzimy przez te same cztery kroki co w przykładzie symulatora, ponownie wykorzystując uchwyt fn z Konfiguracji.

Mała skalaDuża skala
Qubity1030
Kroki Trottera1020
Kroki skompresowane AQC (1-warstwa + 2-warstwy)3 + 2 = 56 + 4 = 10
Warstwy ansatzu stanu podstawowego35
Maksymalny wymiar wiązania MPS32128
BackendstatevectorQPU z DD, Pauli twirling i TREX

Krok 1: Zmapuj klasyczne dane wejściowe na problem kwantowy​

Zbuduj ten sam SparsePauliOp Heisenberga KCuF3_3 i przygotuj stan podstawowy, teraz z głębszym ansatzem gs_layers=5 dla dłuższego łańcucha, następnie wbuduj kopnięcie neutronowe π/2\pi/2 ZZ w centralnym site. Jest to identyczne z mapowaniem na małą skalę, ale przy n=30n = 30.

Spodziewaj się niższej wierności stanu podstawowego niż w przebiegu 10-site: około 0,82 tutaj wobec 0,98 dla mniejszego łańcucha, ponieważ pięć warstw HVA nie może w pełni uchwycić stanu podstawowego 30-site. To spodziewane, a nie porażka, a oryginalny samouczek akceptuje w przybliżeniu 0,65 przy 50 site z tego samego powodu. Zwiększenie gs_layers lub limitu iteracji COBYQA poprawia to kosztem dodatkowej mocy obliczeniowej.

n = 30
dt = 0.6
time_steps = 20
center = n // 2 - 1

# Same MPS settings as the original large-scale run: a larger bond for the
# longer, more-entangled chain (shared by GS prep and AQC compression).
mps_max_bond = 128
mps_cutoff = 1e-8

# Same KCuF3 Hamiltonian and ground-state prep, on a larger chain
H = SparsePauliOp.from_sparse_list(
    [(p, [i, i + 1], 0.25) for i in range(n - 1) for p in ("XX", "YY", "ZZ")],
    num_qubits=n,
)
gs_circuit = prepare_ground_state(
    n, gs_layers=5, max_bond=mps_max_bond, cutoff=mps_cutoff
)
gs_circuit.rz(np.pi / 2, center)  # neutron kick at the center site
print(
    f"Prepared {n}-qubit ground state with the neutron kick at site {center}."
)
DMRG ground-state energy: -13.111355
GS fidelity: 0.8201
Prepared 30-qubit ground state with the neutron kick at site 14.

Kroki 2 i 3: Kompresja i wykonanie za pomocą szablonu funkcji​

To samo pojedyncze wywołanie co w przykładzie symulatora, teraz z backend_name wskazującym na procesor IBM Quantum, więc funkcja transpiluje i wykonuje tam. Plan kompresji różnicuje głębokość ansatzu: pierwsze sześć (niskosplątanych) kroków Trottera kompresuje się do płytkiego jednowarstwowego ansatzu, kolejne cztery do głębszego dwuwarstwowego ansatzu, a pozostałe 10 z 20 kroków uruchamianych jest jako zwykły Trotter. aqc_options podnosi wymiar wiązania MPS do max_bond=128 dla dłuższego, bardziej splątanego łańcucha (zgodnie z oryginałem), zachowując ten sam optymalizator L-BFGS-B ograniczony do 100 iteracji. estimator_options włącza wbudowaną mitygację błędów: dynamiczne rozsprzęganie (XY4), twirling bramek i mitygację odczytu TREX. Domyślne wartości funkcji już odpowiadają oryginalnemu samouczkowi we wszystkich tych aspektach oprócz budżetu uczenia TREX (measure_noise_learning). Cały blok jest mimo to wypisany, ponieważ estimator_options dostarczone przez wywołującego zastępuje domyślne wartości funkcji w całości zamiast się z nimi łączyć, więc pominięcie klucza spowodowałoby powrót do domyślnej wartości IBM Quantum Compute, a nie funkcji.

# Steps 2 + 3: the function compresses (varied ansatz) and executes on hardware.
job = fn.run(
    t_steps=time_steps,
    aqc_segments=[
        {
            "n_steps": 6,
            "ansatz_steps": 1,
        },  # early steps -> shallow 1-layer ansatz
        {
            "n_steps": 4,
            "ansatz_steps": 2,
        },  # later steps -> deeper 2-layer ansatz
    ],
    aqc_options={
        "max_bond": mps_max_bond,  # 128 for the longer chain
        "cutoff": mps_cutoff,
        "optimizer_settings": {
            "method": "L-BFGS-B",
            "jac": True,
            "options": {"maxiter": 100},
        },
    },
    dt=dt,
    hamiltonian=H,
    initial_state=gs_circuit,
    backend_name="ibm_pittsburgh",
    # Mitigation settings from the original tutorial. Only the two
    # measure_noise_learning values differ from the function's defaults; the rest
    # restates them, because a caller-supplied estimator_options dict replaces the
    # function's defaults wholesale rather than merging into them.
    estimator_options={
        "environment": {"job_tags": ["TUT-SNS"]},
        "dynamical_decoupling": {"enable": True, "sequence_type": "XY4"},
        "twirling": {
            "enable_gates": True,
            "num_randomizations": 1000,
            "shots_per_randomization": 128,
        },
        "resilience": {
            "measure_mitigation": True,
            "measure_noise_learning": {
                "num_randomizations": 32,
                "shots_per_randomization": 100,
            },
        },
    },
)
print("job ID (save this to reconnect later):", job.job_id)
job ID (save this to reconnect later): 43ed8d07-6d7d-4f33-b70a-7f31b765b310
Ponowne łączenie się z długo działającym zadaniem

Przebieg na dużą skalę nie jest szybki, a większość czasu jest klasyczna, a nie na QPU. Kompresja AQC uruchamiana jest wewnątrz funkcji, zanim cokolwiek dotrze do QPU: przy 30 site z max_bond=128 zajęło to prawie cztery godziny w naszym przebiegu, wobec około 18 minut czasu QPU podanego w Szacowanym czasie użycia na początku tego samouczka. Czas oczekiwania w kolejce dochodzi do obu tych wartości. Nie musisz utrzymywać otwartego tego notebooka ani kernela podczas jego działania.

Skopiuj identyfikator zadania wypisany przez poprzednią komórkę i zapisz go. Kolejne trzy komórki pozwalają Ci później do niego wrócić:

  1. Ponowne połączenie, potrzebne tylko w nowej sesji kernela: ponownie uruchom komórki Konfiguracji, aby odtworzyć serverless, a następnie odbuduj uchwyt job z zapisanego identyfikatora. Pomiń tę komórkę, jeśli nadal jesteś w sesji, w której złożyłeś zadanie, ponieważ uchwyt jest już aktywny.

  2. Sprawdzenie statusu: uruchamiaj ponownie, aż zgłosi DONE.

  3. Pobranie wyniku: uruchom tylko wtedy, gdy status to DONE.

Poniższa komórka do ponownego łączenia zawiera symbol zastępczy. Zastąp go swoim własnym job_id:

# Reconnect to a previously submitted job by its ID. Only needed in a NEW kernel
# session; if you are still in the session where you submitted, the `job` handle
# from the preceding cell is already live, so skip this cell. Replace the ID that follows with your own.
job = serverless.get_job_by_id("<your job ID>")
# Check where the job is. Re-run this until it reports DONE before fetching the
# result in the following cell: QUEUED -> INITIALIZING -> RUNNING: OPTIMIZING_FOR_HARDWARE ->
# RUNNING: WAITING_FOR_QPU -> RUNNING: EXECUTING_QPU -> RUNNING: POST_PROCESSING
# -> DONE.
print(job.status())
DONE
# Run this only once the preceding status cell reports DONE. result() blocks until
# the job finishes, so calling it earlier just waits (possibly for hours).
result = job.result()
print(
    "AQC fidelities:",
    {k: round(v, 4) for k, v in result["metadata"]["aqc_fidelities"].items()},
)

ev = np.array(result["expectation_values"])
Gjjc = ev[1:]  # drop the t = 0 row -> shape (time_steps, n)
AQC fidelities: {'1': 1.0, '2': 0.9994, '3': 0.9944, '4': 0.9853, '5': 0.9747, '6': 0.959, '7': 0.9495, '8': 0.9542, '9': 0.9533, '10': 0.9451}

Krok 4: Post-processing i zwrot wyniku w żądanym formacie klasycznym​

Identyczny post-processing jak w przebiegu symulatora: przekształć Fourierem funkcję Greena na S(q,ω)S(q, \omega), symetryzuj lustrzanie i przytnij wartości ujemne. Przy dłuższym łańcuchu i ewolucji kontinuum dwuspinonowe jest znacznie lepiej rozpoznane. Powinno wypełniać pasmo między przerywanymi granicami, najjaśniejsze w pobliżu q=πq = \pi.

n = result["metadata"]["n"]
q_res, w_res = 100, 100
spectrum = get_spectrum(n, Gjjc, dt, time_steps, q_res, w_res)
spectrum = -(spectrum + spectrum[:, ::-1]) / 2  # mirror symmetry
spectrum = np.clip(spectrum, a_min=0, a_max=None)  # clip negatives

plot_green(
    n,
    Gjjc,
    time_steps,
    dt,
    title=f"Retarded Green's function - {n} qubits (AQC, hardware)",
)
plot_spectrum(
    spectrum,
    dt,
    q_res,
    w_res,
    lower_bound=True,
    upper_bound=True,
    title=f"Dynamical structure factor - {n} qubits (AQC, hardware)",
)

Output of the previous code cell

Output of the previous code cell

Dodatek​

Poprzedni przykład sprzętowy uruchamia pojedynczą długość łańcucha. Trzy widma, które następują, pochodzą z wcześniejszych przebiegów sprzętowych tego samego workflow na ibm_pittsburgh przy 10, 20 i 30 site, ze wszystkimi pozostałymi danymi wejściowymi utrzymanymi na stałym poziomie: 20 kroków Trottera przy dt = 0.6, plan kompresji sześciu jednowarstwowych i czterech dwuwarstwowych kroków skompresowanych AQC oraz max_bond = 128. Są to zarejestrowane wyniki, a nie wyjście z poprzednich komórek.

Te same ustawienia są używane przy wszystkich trzech rozmiarach, więc widma są bezpośrednio porównywalne. Dostrajanie ich dla każdej długości łańcucha, na przykład z większą liczbą warstw ansatzu stanu podstawowego lub większym max_bond, może dać lepsze wyniki niż jakiekolwiek pokazane tutaj.

Dynamical structure factor at 10 sites, a single sharp bright peak at q = pi near the lower bound

Dynamical structure factor at 20 sites, spectral weight filling the band between the two dashed two-spinon bounds

Dynamical structure factor at 30 sites, the continuum resolved more finely with fainter contrast and some weight outside the bounds

Kolejne kroki​

Rekomendacje
Źródło: IBM Quantum docs — zaktualizowano 9 wrz 2026
Wersja angielska na doQumentation — zaktualizowano 9 wrz 2026
Problem ze stroną lub tłumaczeniem?Problem z treścią? Edytuj w IBM Quantum docs