Mga multi-product formula para mabawasan ang Trotter error

Tinatayang paggamit: Apat na minuto sa isang Heron r2 processor (TANDAAN: Tantiya lamang ito. Maaaring mag-iba ang iyong runtime.)

Mga learning outcome

Kung paano binabawasan ng multi-product formulas (MPFs) ang Trotter error sa Hamiltonian simulation sa pamamagitan ng pagsasama-sama ng expectation values mula sa maraming shallow circuits
Kailan kapaki-pakinabang ang MPFs kaysa sa standard product formulas at kailan ito hindi ang tamang tool
Kung paano kompyutin ang static at dynamic na MPF coefficients gamit ang qiskit_addon_mpf package
Kung paano isagawa ang isang MPF workflow end-to-end sa IBM Quantum® hardware, kabilang ang transpilation, error mitigation, at post-processing

Prerequisites

Compilation methods for Hamiltonian simulation circuits — nagpapakilala sa Trotter (product formula) circuits sa Qiskit.
Product formulas sa Qiskit, partikular ang mga synthesis class na SuzukiTrotter at LieTrotter.
Qiskit primitives and the Estimator interface.

Pangkalahatang-ideya

Ano ang mga multi-product formula?

Kapag ginagaya ang mga quantum system sa isang quantum computer, isang sentral na gawain ay ang pag-approximate ng time-evolution operator $e^{-iHt}$ para sa isang Hamiltonian $H$ . Ang standard approach ay gumagamit ng product formulas (PFs), na kilala rin bilang Trotter-Suzuki decompositions. Dinedekompos nito ang $H = \sum_{a=1}^d F_a$ sa mga terms kung saan ang kani-kaniyang unitaries na $e^{-iF_a t}$ ay episyenteng maipapatupad, at pagkatapos ay ina-approximate ang buong evolution bilang isang ordered product ng mga simpleng unitary na ito.

Ang first-order product formula (Lie-Trotter) ay:

S_1(t) := \prod_{a=1}^d e^{-i F_a t},

na nagdudulot ng quadratic error: $S_1(t) = e^{-iHt} + \mathcal{O}(t^2)$ . Ang mga higher-order symmetric formula na $S_{2\chi}(t)$ , kung saan ang $\chi$ ay nagla-label sa order ng symmetric product formula (tingnan ang Ref. [1]), ay kumokombergehe nang mas mabilis bilang $e^{-iHt} + \mathcal{O}(t^{2\chi+1})$ , ngunit sa halaga ng mas malalim na circuits bawat hakbang.

Upang mabawasan ang error sa isang fixed na order $\chi$ , karaniwang hinahati ang kabuuang evolution time $t$ sa $k$ na mas maliliit na Trotter steps. Bawat hakbang ay ina-approximate ang $e^{-iHt/k}$ ng isang product formula at ang mga hakbang ay pinagsasama-sama:

e^{-iHt} \approx \left[S_{2\chi}(t/k)\right]^k.

Para sa isang $2\chi$ -order symmetric formula, ang residual Trotter error ay nagsi-scale bilang $\mathcal{O}\!\left(t^{2\chi+1} / k^{2\chi}\right)$ . Kaya ang pagtaas ng $k$ ay mabilis na sumusupress sa Trotter error — ngunit ito rin ay linear na nagpapalalim sa circuit, at sa noisy hardware nangangahulugan ito ng mas maraming naiipong gate noise. Ang tension na ito sa pagitan ng Trotter error (pumapabor sa mas malaking $k$ ) at hardware noise (pumapabor sa mas maliit na $k$ ) ang eksaktong disenyo ng multi-product formulas upang lutasin. Tandaan na ang MPFs ay tungkol sa pagsasama-sama ng mga resulta mula sa magkakaibang pagpili ng $k$ sa isang fixed order $\chi$ — hindi nito binabago ang order ng underlying product formula.

Ang Multi-product formulas (MPFs) [1] ay bumubuo ng isang weighted linear combination ng expectation values na nakuha mula sa ilang mas mabababaw na Trotter circuits, bawat isa ay gumagamit ng ibang bilang ng Trotter steps na $k_1, k_2, \ldots, k_r$ (isang set ng $r$ na step counts):

\langle A \rangle_{\text{MPF}}(t) = \sum_{j=1}^r x_j \, \langle A \rangle_{k_j}(t),

kung saan ang $\langle A \rangle_{k_j}(t)$ ay ang expectation value ng isang observable $A$ sa oras $t$ na na-estimate mula sa isang Trotter circuit na may $k_j$ steps, at ang mga coefficient na $\{x_j\}_{j=1}^r$ ay pinipili upang mag-cancel ang mga leading Trotter-error terms sa combination. Babalikan natin ang expression na ito sa Step 4, kung saan susuriin natin ito nang explicit upang pagsamahin ang ating mga resultang Trotter. Ang pangunahing praktikal na punto ay ang pinakamalalim na circuit sa MPF ay nangangailangan lamang ng $k_{\max}$ na steps, na mas maliit kaysa sa iisang $k$ na kakailanganin upang maabot direkta ang parehong effective Trotter error. Ang mas mabababaw na circuits ang gumagawa sa MPF approach na mas angkop para sa noisy hardware.

Paano tinutukoy ang mga coefficient?

May dalawang pamilya ng MPF coefficients:

Ang Static coefficients ay independiyente sa Hamiltonian, sa initial state, at sa evolution time. Nakukuha ang mga ito sa pamamagitan ng paglutas ng isang linear system na $Ax = b$ na nag-eenforce ng cancellation ng leading Trotter error terms. Para sa isang set ng Trotter steps na $\{k_j\}_{j=1}^r$ na ginamit kasama ng isang $2\chi$ -order symmetric product formula, ang pag-expand ng Trotter error sa inverse powers ng $k_j$ ay humahantong sa constraint equations ng anyo:

\sum_{j=1}^r x_j = 1, \quad \sum_{j=1}^r \frac{x_j}{k_j^{\eta_n}} = 0 \quad (n = 0, \ldots, r-2),

kung saan ang mga integer exponent na $\{\eta_n\}$ ay ang mga order ng magkakasunod na Trotter-error terms para sa napiling product formula. Para sa isang symmetric $2\chi$ -order PF, ang leading error sa $\left[S_{2\chi}(t/k)\right]^k$ ay nagsi-scale bilang $1/k^{2\chi}$ , na may kasunod na mga correction sa $1/k^{2\chi+2}, 1/k^{2\chi+4}, \ldots$ — kaya ang mga exponent ay $\eta_n = 2\chi + 2n$ . Para sa non-symmetric PFs, kapwa ang odd at even powers ay nag-aambag at $\eta_n = 2\chi + n$ . Tingnan ang Ref. [1] para sa buong derivation. Ang unang equation sa system sa itaas ay tinitiyak ang unbiasedness (ginagaya ng MPF ang eksaktong expectation value sa $k_j \to \infty$ limit), at ang natitirang $r-1$ na equations ay magkakasunod na kinakansela ang unang $r-1$ na Trotter-error terms. Kapag ang resultang $L_1$ -norm na $\|x\|_1$ ay masyadong malaki (na nagpapalakas sa sampling noise), maaari kang lumutas sa halip ng isang approximate optimization na nagli-limit sa $\|x\|_1$ habang minamaliit ang $\|Ax - b\|$ .

Ang Dynamic coefficients [2], [3] ay nakadepende rin sa Hamiltonian, sa initial state, at sa evolution time na $t$ . Minimize nila ang Frobenius-norm distance sa pagitan ng totoong time-evolved state at ng MPF approximation:

\|\rho(t) - \mu^D(t)\|_F^2 = 1 + \sum_{i,j} M_{ij}(t)\, x_i(t)\, x_j(t) - 2\sum_i L_i(t)\, x_i(t),

kung saan ang $M_{ij}(t) = \mathrm{Tr}[\rho_{k_i}(t)\,\rho_{k_j}(t)]$ ay ang Gram matrix ng overlaps sa pagitan ng mga Trotter-evolved states para sa iba't ibang step counts na $k_i, k_j$ , at ang $L_i(t) = \mathrm{Tr}[\rho(t)\,\rho_{k_i}(t)]$ ay sumusukat ng overlap sa (approximate) eksaktong state. Sa tutorial na ito, ang mga dami na ito ay kinokompyut nang episyente gamit ang tensor-network methods, partikular ang mga TeNPy-based backend sa qiskit_addon_mpf.

Kailan gagamitin ang mga MPF

Ang MPFs ay pinaka-kapaki-pakinabang kapag:

Ang circuit depth ang bottleneck. Kung ang hardware noise ang naglilimita sa depth na maaari mong patakbuhin, gamitin ang MPFs upang makamit ang mas mataas na effective Trotter accuracy mula sa mas mabababaw na circuits.
Kailangan mo ng tumpak na expectation values, hindi buong state preparation. Ang MPFs ay gumagana sa antas ng expectation values — pinagsasama nila ang mga classical number, hindi quantum states. Kaya't ideal ang mga ito para sa observable estimation kapag gumagamit ng Estimator primitive.
Pinagsasama mo ang katamtamang bilang ng Trotter-step counts. Karaniwan, ang pagsasama ng $r = 3$ – $5$ na magkakaibang step counts na $k_j$ ay sapat na upang kanselahin ang ilang leading Trotter-error terms habang pinapanatiling manageable ang $\|x\|_1$ .

Kapag hindi nakakatulong ang MPFs

Napakaikling oras ng ebolusyon. Kapag maliit na ang $t$ na sapat na ang isang mababang-order na Trotter formula ay tama na, ang overhead ng pagpapatakbo ng maraming circuit ay hindi na kailangan.
Mga gawain sa paghahanda ng estado. Ang mga MPF ay gumagawa ng naitamang expectation value, hindi naitamang quantum state. Kung kailangan mo ang aktwal na time-evolved na estado (halimbawa, bilang input sa isa pang quantum subroutine), hindi angkop ang MPFs.
Mga bilang ng Trotter-step na lumalabag sa convergence regime. Ang static-coefficient derivation ay nagpapalawak ng bawat indibidwal na $\left[S_{2\chi}(t/k_j)\right]^{k_j}$ bilang isang series sa $t/k_j$ ; ang expansion na ito ay mahusay lang nagko-converge kapag $t/k_{\min} \lesssim 1$ . Kung masyadong maliit ang $k_{\min}$ na napili para sa ibinigay na $t$ , ang pinakamababaw na circuit ay malayo sa perturbative regime, ang mga higher-order error terms na hindi na-cancel ng MPF ay nagiging malaki, at ang cancellation ay maaaring mangailangan ng malalaking coefficients. Ang $L_1$ -norm na $\|x\|_1$ ay ang praktikal na diagnostic: kapag $\|x\|_1 \gg 1$ , ang sampling overhead na $\propto \|x\|_1^2$ ay maaaring mas malaki kaysa sa pagbawas ng Trotter-error. Tingnan ang gabay sa pagpili ng Trotter steps para sa detalye.

Ano ang sinasaklaw ng tutorial na ito

Ang tutorial na ito ay dumadaan sa isang end-to-end MPF workflow sa dalawang yugto. Una, isang maliit na halimbawa gamit ang simulator (10-qubit Heisenberg chain) na nagpapakita kung paano itakda ang problema, kalkulahin ang static at dynamic na MPF coefficients, at ikumpara ang resultang expectation values laban sa exact diagonalization. Pagkatapos, isang malaking halimbawa gamit ang hardware (50-qubit XXZ chain) na nagpapakita kung paano mag-transpile, magpatakbo sa IBM Quantum hardware gamit ang error mitigation, at i-post-process ang mga resulta gamit ang MPF coefficients. Sa buong tutorial, ginagamit natin ang qiskit_addon_mpf package kasama ng standard na Qiskit tools.

Mga Kinakailangan

Bago simulan ang tutorial na ito, tiyakin na naka-install ang mga sumusunod:

Qiskit SDK v2.0 o mas bago, na may suporta sa visualization
Qiskit Runtime v0.22 o mas bago (pip install qiskit-ibm-runtime)
Qiskit Aer simulator (pip install qiskit-aer)
MPF Qiskit addon na may TeNPy backend (pip install "qiskit-addon-mpf[tenpy]")
Qiskit addon utilities (pip install qiskit-addon-utils)
SciPy (pip install scipy)

Setup

Sa ibaba ay tinipon natin ang lahat ng package imports na ginamit sa buong tutorial na ito sa isang cell. Nagtatakda rin tayo ng CollectAndCollapse transpiler pass na pinagsasama-sama ang magkatabing rxx at ryy rotations sa iisang XXPlusYYGate. Ang pass na ito ay inilalapat parehong sa panahon ng pagbuo ng circuit sa Hakbang 1 (upang mapanatiling mababa ang bilang ng gate) at hindi tuwiran kapag kinuha natin ang layered structure para sa dynamic MPF sa Hakbang 4 (inaasahan ng TeNPy ang two-qubit gates, hindi pares ng unfused rotations).

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

_needed = {"matplotlib": "matplotlib", "numpy": "numpy", "qiskit": "qiskit", "qiskit_addon_mpf": "qiskit-addon-mpf", "qiskit_addon_utils": "qiskit-addon-utils", "qiskit_aer": "qiskit-aer", "qiskit_ibm_runtime": "qiskit-ibm-runtime", "scipy": "scipy", "tenpy": "physics-tenpy"}
_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
import matplotlib.pyplot as plt
from functools import partial
from copy import deepcopy

from qiskit import QuantumCircuit
from qiskit.quantum_info import Pauli, SparsePauliOp, Statevector
from qiskit.synthesis import SuzukiTrotter
from qiskit.transpiler import CouplingMap, PassManager
from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager
from qiskit.circuit.library import XXPlusYYGate
from qiskit.transpiler.passes.optimization.collect_and_collapse import (
    CollectAndCollapse,
    collect_using_filter_function,
    collapse_to_operation,
)

from qiskit_aer import AerSimulator
from qiskit_ibm_runtime import EstimatorV2 as Estimator, QiskitRuntimeService

from qiskit_addon_utils.problem_generators import (
    generate_xyz_hamiltonian,
    generate_time_evolution_circuit,
)
from qiskit_addon_utils.slicing import slice_by_depth
from qiskit_addon_mpf.static import setup_static_lse
from qiskit_addon_mpf.dynamic import setup_dynamic_lse
from qiskit_addon_mpf.costs import (
    setup_exact_problem,
    setup_sum_of_squares_problem,
    setup_frobenius_problem,
)
from qiskit_addon_mpf.backends.tenpy_layers import (
    LayerModel,
    LayerwiseEvolver,
)
from qiskit_addon_mpf.backends.tenpy_tebd import MPOState, MPS_neel_state

from scipy.linalg import expm

# Suppress TeNPy's `unit_cell_width` future-API warning. The default
# (`unit_cell_width=len(sites)`) is correct for Chain lattices, which is what
# `CouplingMap.from_line(...)` produces here, so the warning is informational.
warnings.filterwarnings(
    "ignore",
    message=r".*unit_cell_width.*",
    category=UserWarning,
)

# --- Helper: collect XX + YY rotations into a single gate ---
def filter_function(node):
    return node.op.name in {"rxx", "ryy"}

collect_function = partial(
    collect_using_filter_function,
    filter_function=filter_function,
    split_blocks=True,
    min_block_size=1,
)

def collapse_to_xx_plus_yy(block):
    param = 0.0
    for node in block.data:
        param += node.operation.params[0]
    return XXPlusYYGate(param)

collapse_function = partial(
    collapse_to_operation,
    collapse_function=collapse_to_xx_plus_yy,
)

pm = PassManager()
pm.append(CollectAndCollapse(collect_function, collapse_function))

Maliit na halimbawa gamit ang simulator

Hakbang 1: I-map ang mga classical input sa isang quantum na problema

Nagsisimula tayo sa isang 10-qubit Heisenberg model sa isang linya, gamit ang Néel state na $\vert 0101\ldots01 \rangle$ bilang paunang estado. Ang Hamiltonian ay:

\hat{\mathcal{H}}_{\text{Heis}} = J \sum_{i=1}^{L-1} \left(X_i X_{i+1} + Y_i Y_{i+1} + Z_i Z_{i+1}\right),

kung saan ang $J$ ay ang nearest-neighbor coupling strength. Sinusukat natin ang ZZ correlator na $Z_{L/2-1} Z_{L/2}$ sa isang pares ng qubits sa gitna ng chain, at gumagamit ng Trotter steps na $k_j = [1, 2, 4]$ na may second-order product formula.

L = 10

# Generate coupling map and Hamiltonian
coupling_map = CouplingMap.from_line(L, bidirectional=False)

hamiltonian = generate_xyz_hamiltonian(
    coupling_map,
    coupling_constants=(1.0, 1.0, 1.0),
    ext_magnetic_field=(0.0, 0.0, 0.0),
)
print(hamiltonian)

SparsePauliOp(['IIIIIIIXXI', 'IIIIIIIYYI', 'IIIIIIIZZI', 'IIIIIXXIII', 'IIIIIYYIII', 'IIIIIZZIII', 'IIIXXIIIII', 'IIIYYIIIII', 'IIIZZIIIII', 'IXXIIIIIII', 'IYYIIIIIII', 'IZZIIIIIII', 'IIIIIIIIXX', 'IIIIIIIIYY', 'IIIIIIIIZZ', 'IIIIIIXXII', 'IIIIIIYYII', 'IIIIIIZZII', 'IIIIXXIIII', 'IIIIYYIIII', 'IIIIZZIIII', 'IIXXIIIIII', 'IIYYIIIIII', 'IIZZIIIIII', 'XXIIIIIIII', 'YYIIIIIIII', 'ZZIIIIIIII'],
              coeffs=[1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j,
 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j,
 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j, 1.+0.j])

# Observable: ZZ on the middle pair of qubits
observable = SparsePauliOp.from_sparse_list(
    [("ZZ", (L // 2 - 1, L // 2), 1.0)], num_qubits=L
)
print(observable)

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

# MPF parameters
mpf_trotter_steps = [1, 2, 4]
order = 2
symmetric = False

trotter_times = np.arange(0.5, 1.55, 0.1)
exact_evolution_times = np.arange(trotter_times[0], 1.55, 0.05)

Bumuo ng mga Trotter circuit

Gumagawa tayo ng mga circuit na nagpapatupad ng tinatayang Trotter time-evolutions para sa bawat time point at bawat Trotter step count. Ang CollectAndCollapse pass na itinakda sa seksyon ng Setup ay pinagsasama-sama ang XX at YY rotations sa iisang XX+YY gates, upang maghanda para sa mas mahusay na tensor-network simulation sa susunod.

# Initial Neel state preparation
initial_state_circ = QuantumCircuit(L)
initial_state_circ.x([i for i in range(L) if i % 2 != 0])

all_circs = []
for total_time in trotter_times:
    mpf_trotter_circs = [
        generate_time_evolution_circuit(
            hamiltonian,
            time=total_time,
            synthesis=SuzukiTrotter(reps=num_steps, order=order),
        )
        for num_steps in mpf_trotter_steps
    ]

    mpf_trotter_circs = pm.run(
        mpf_trotter_circs
    )  # Collect XX and YY into XX + YY

    mpf_circuits = [
        initial_state_circ.compose(circuit) for circuit in mpf_trotter_circs
    ]
    all_circs.append(mpf_circuits)

mpf_circuits[-1].draw("mpl", fold=-1)

Hakbang 2: I-optimize ang problema para sa pagpapatupad sa quantum hardware

Para sa maliit na halimbawa, tina-target natin ang Aer simulator. Dalawang transformation ang nangyayari bago handa nang isagawa ang mga circuit:

Pagkolekta ng gate sa antas ng Hamiltonian-simulation. Sa Setup cell binuo natin ang isang CollectAndCollapse pass na pinagsasama-sama ang magkatabing rxx at ryy rotations sa iisang XXPlusYYGate. Nailapat na natin ang pass na ito noong binuo natin ang mga Trotter circuit sa Hakbang 1 (ang pm.run(...) call). Binabawasan nito ang bilang ng two-qubit gate at gumagawa ng structure na mas angkop para sa tensor-network simulation para sa dynamic-coefficient computation sa susunod.
Pagbaba sa ISA ng simulator. Sa ibaba, pinapatakbo natin ang Qiskit preset pass manager sa optimization_level=3 upang ibaba ang bawat Trotter circuit sa instruction-set architecture (ISA) ng simulator.

aer_sim = AerSimulator()
pm_sim = generate_preset_pass_manager(backend=aer_sim, optimization_level=3)

isa_circs_all_times = [
    pm_sim.run([deepcopy(c) for c in mpf_circuits])
    for mpf_circuits in all_circs
]

Hakbang 3: Isagawa gamit ang Qiskit primitives

Para sa maliit na halimbawa, pinapatakbo natin ang mga ISA-lowered Trotter circuits sa pamamagitan ng EstimatorV2 primitive na suportado ng Aer. Sa paggawa nito, nakakakuha tayo ng walang ingay na reference value para sa bawat $(k_j, t)$ pares — ito ang mga $\langle A \rangle_{k_j}(t)$ values na pagsasamahin ng MPF sa Hakbang 4. Sinusuri natin ang iba't ibang oras ng ebolusyon upang mas maaari nating i-plot ang buong time-series curve ng bawat indibidwal na product formula at ng MPF.

estimator = Estimator(mode=aer_sim)

mpf_expvals_all_times, mpf_stds_all_times = [], []
for isa_circuits in isa_circs_all_times:
    result = estimator.run(
        [(circuit, observable) for circuit in isa_circuits], precision=0.005
    ).result()
    mpf_expvals_all_times.append([res.data.evs for res in result])
    mpf_stds_all_times.append([res.data.stds for res in result])

Hakbang 4: I-post-process at ibalik ang resulta sa nais na klasikal na format

Ang Hakbang 4 ay kung saan aktwal na binubuo ang MPF. Kahit na kinakalkula dito ang mga coefficients na $x_j$ (at, para sa dynamic variant, maaaring maging masinsinan ang kalkulasyong ito), sa konsepto ito ay isang klasikal na recipe para pagsamahin ang mga quantum measurements mula sa Hakbang 3 sa isang naitamang expectation value — kaya't itinuturing natin ang buong coefficient at combination workflow bilang post-processing.

Upang tasahin kung gaano kahusay sinusubaybayan ng MPF ang tunay na dynamics, kalkulahin muna natin ang eksaktong time-evolved expectation values sa pamamagitan ng direktang pag-exponentiate ng Hamiltonian. Ito ay magagawa lamang dahil $L = 10$ ; sa malaking halimbawa gamit ang hardware sa ibaba, kailangan nating umasa sa tensor-network estimates sa halip.

exact_expvals = []
for t in exact_evolution_times:
    exp_H = expm(-1j * t * hamiltonian.to_matrix())
    initial_state = Statevector(initial_state_circ).data
    time_evolved_state = exp_H @ initial_state

    exact_obs = (
        time_evolved_state.conj()
        @ observable.to_matrix()
        @ time_evolved_state
    ).real
    exact_expvals.append(exact_obs)

Mga static na MPF coefficient

Ang static MPFs ay gumagamit ng coefficients na $x_j$ na independiyente sa oras ng ebolusyon, ang Hamiltonian, at ang paunang estado. Itinatakda natin ang linear system na $Ax = b$ na inilarawan sa Background at nireresolba para sa coefficients. Ang matrix na $A$ ay tinutukoy ng Trotter step counts na $k_j$ , ang order $\chi$ ng product formula, at kung ang formula ay symmetric (na kumokontrol sa exponents na $\eta_n$ ).

Para sa ating maliit na halimbawa, ginagamit natin ang $k_j = [1, 2, 4]$ na may non-symmetric order- $2\chi=2$ Suzuki-Trotter formula (kaya $\chi=1$ at $\eta_n = 2 + n$ , ibinibigay ang $\eta_0 = 2,\, \eta_1 = 3$ ). Ang system ay nagiging:

A = \begin{bmatrix} 1 & 1 & 1\\ 1 & \frac{1}{2^2} & \frac{1}{4^2} \\ 1 & \frac{1}{2^3} & \frac{1}{4^3} \\ \end{bmatrix}, \quad b = \begin{bmatrix} 1 \\ 0 \\ 0 \end{bmatrix}.

Ang unang row ay nagpapatupad ng unbiasedness ( $\sum_j x_j = 1$ ); ang ikalawa at ikatlong row ay kina-cancel ang unang-orden na $1/k^2$ at next-order na $1/k^3$ Trotter-error terms, ayon sa pagkakasunod-sunod.

Itakda ang LSE

Ginagamit natin ang setup_static_lse mula sa qiskit_addon_mpf.static upang buuin ang matrix na $A$ at right-hand-side vector na $b$ na inilarawan sa itaas. Ang matrix na $A$ ay nakadepende hindi lang sa $k_j$ kundi pati na rin sa ating pinili na product formula — partikular na sa order nitong $\chi$ at kung ito ay symmetric. Ang symmetric flag ay kumokontrol sa exponent pattern na $\eta_n$ (ang symmetric formulas ay gumagawa lamang ng even-power Trotter-error terms; tingnan ang Ref. [1]). Tandaan na, tulad ng ipinapakita sa Ref. [2], ang pagtatakda ng symmetric=True ay hindi mahigpit na kinakailangan kahit symmetric ang underlying PF — ang non-symmetric LSE ay nananatiling valid (nagpapatupad ito ng ekstrang di-kailangang constraints).

Para sa ating halimbawa, naitakda na natin ang order = 2 at symmetric = False sa Hakbang 1.

lse = setup_static_lse(mpf_trotter_steps, order=order, symmetric=symmetric)

Suriin ang binuong matrix na $A$ at vector na $b$ upang kumpirmahin na tumutugma ang mga ito sa system na isinulat sa itaas.

lse.A

array([[1.      , 1.      , 1.      ],
       [1.      , 0.25    , 0.0625  ],
       [1.      , 0.125   , 0.015625]])

lse.b

array([1., 0., 0.])

Taglay ang LSE, nireresolba natin ang static coefficients na $x_j$ sa pamamagitan ng lse.solve() (ito ang direktang solusyong $x = A^{-1}b$ ).

mpf_coeffs = lse.solve()
print(
    f"The static coefficients associated with the ansatze are: {mpf_coeffs}"
)

The static coefficients associated with the ansatze are: [ 0.04761905 -0.57142857  1.52380952]

I-optimize ang $x$ gamit ang isang eksaktong modelo

Bilang alternatibo sa pagkalkula ng $x = A^{-1}b$ , maaari mo ring gamitin ang setup_exact_model upang itayo ang isang cvxpy.Problem instance na gumagamit ng LSE bilang mga constraint at ang optimal na solusyon nito ay magbubunga ng $x$ .

model_exact, coeffs_exact = setup_exact_problem(lse)
model_exact.solve()
print(coeffs_exact.value)

[ 0.04761905 -0.57142857  1.52380952]

print(
    "L1 norm of the exact coefficients:",
    np.linalg.norm(coeffs_exact.value, ord=1),
)

L1 norm of the exact coefficients: 2.1428571428556378

I-optimize ang $x$ gamit ang isang tinantyang modelo

Maaaring mangyari na ang $L_1$ norm para sa piniling hanay ng mga halaga ng $k_j$ ay masyadong mataas. Kung ganoon ang kaso at hindi ka makapipili ng ibang hanay ng mga halaga ng $k_j$ , maaari kang gumamit ng tinantyang solusyon na naghihigpit sa $L_1$ -norm sa isang piniling threshold habang minimisado ang $\|Ax - b\|$ . Tingnan ang gabay sa Paano gamitin ang approximate model.

model_approx, coeffs_approx = setup_sum_of_squares_problem(
    lse, max_l1_norm=1.5
)
model_approx.solve()
print(coeffs_approx.value)
print(
    "L1 norm of the approximate coefficients:",
    np.linalg.norm(coeffs_approx.value, ord=1),
)

[-1.10294118e-03 -2.48897059e-01  1.25000000e+00]
L1 norm of the approximate coefficients: 1.5

Mga dynamic na MPF coefficient

Ang static MPF ay kina-cancel ang Trotter-error terms sa paraang agnostiko sa Hamiltonian at estado, kaya hindi nito palaging naiibigay ang pinakamaliit na posibleng approximation error para sa isang ibinigay na Hamiltonian at paunang estado. Sa halip, ang dynamic MPF (Refs. [2], [3]) ay naghahanap ng time-dependent coefficients na $x_i(t)$ na minimisado ang Frobenius-norm distance na $\|\rho(t) - \mu^D(t)\|_F^2$ sa bawat oras na $t$ . Tulad ng ipinapakita sa Background, ito ay nangangailangan ng overlap matrix na $M_{ij}(t)$ sa pagitan ng mga Trotter-evolved states at ang overlap na $L_i(t)$ sa eksaktong estado — kapwa tinatantiya natin gamit ang tensor-network (TeNPy) backends sa qiskit_addon_mpf.

Upang itakda ang dynamic LSE, kailangan natin ng tatlong sangkap:

Isang approximate evolver factory na papatakbuhin ng addon para sa bawat $k_j$ upang bumuo ng $\rho_{k_j}(t)$ bilang MPS/MPO. Binubuo natin ito mula sa layered structure ng order- $2$ Trotter circuit (isang layer bawat slice_by_depth), na na-wrap bilang LayerwiseEvolver na may TeNPy truncation parameters.
Isang exact evolver factory na gumagawa ng mataas-ang-accuracy na reference na $\rho(t)$ . Ginagamit natin ang isang small-time-step fourth-order Suzuki-Trotter circuit (dt=0.1, order=4) bilang proxy para sa eksaktong ebolusyon.
Isang identity factory at isang initial-state MPS na nagsisimula ng TeNPy simulation.

Ang cell sa ibaba ay bumubuo ng approximate evolver factory.

# Create approximate time-evolution circuits
single_2nd_order_circ = generate_time_evolution_circuit(
    hamiltonian, time=1.0, synthesis=SuzukiTrotter(reps=1, order=order)
)
single_2nd_order_circ = pm.run(single_2nd_order_circ)  # collect XX and YY

# Find layers in the circuit
layers = slice_by_depth(single_2nd_order_circ, max_slice_depth=1)

# Create tensor network models
models = [
    LayerModel.from_quantum_circuit(layer, conserve="Sz") for layer in layers
]

# Create the time-evolution object
approx_factory = partial(
    LayerwiseEvolver,
    layers=models,
    options={
        "preserve_norm": False,
        "trunc_params": {
            "chi_max": 64,
            "svd_min": 1e-8,
            "trunc_cut": None,
        },
        "max_delta_t": 2,
    },
)

babala

Ang mga option ng LayerwiseEvolver na tumutukoy sa mga detalye ng tensor network simulation ay dapat piliin nang maingat para maiwasan ang pagtatakda ng ill-defined na optimization problem.

Inaproximate natin ang eksaktong time-evolved state gamit ang fourth-order Suzuki-Trotter formula gamit ang maliit na time step na dt=0.1. Ang TeNPy truncation parameters ay maaaring makaapekto sa accuracy, kaya mahalagang tuklasin ang saklaw ng mga halaga.

single_4th_order_circ = generate_time_evolution_circuit(
    hamiltonian, time=1.0, synthesis=SuzukiTrotter(reps=1, order=4)
)
single_4th_order_circ = pm.run(single_4th_order_circ)
exact_model_layers = [
    LayerModel.from_quantum_circuit(layer, conserve="Sz")
    for layer in slice_by_depth(single_4th_order_circ, max_slice_depth=1)
]

exact_factory = partial(
    LayerwiseEvolver,
    layers=exact_model_layers,
    dt=0.1,
    options={
        "preserve_norm": False,
        "trunc_params": {
            "chi_max": 64,
            "svd_min": 1e-8,
            "trunc_cut": None,
        },
        "max_delta_t": 2,
    },
)

Sa wakas, tinatakda natin ang isang identity_factory na gumagawa ng paunang MPO state at inihahanda ang Néel initial state bilang MPS na tumutugma sa lattice na ginamit ng layered Trotter model.

def identity_factory():
    return MPOState.initialize_from_lattice(models[0].lat, conserve=True)

mps_initial_state = MPS_neel_state(models[0].lat)

Taglay ang mga factory, kinakalkula natin ngayon ang dynamic coefficients sa bawat oras ng ebolusyon. Para sa bawat $t$ , ang setup_dynamic_lse ay bumubuo ng kaugnay na overlap matrices sa pamamagitan ng TeNPy, at ang setup_frobenius_problem ay nagbabalik ng cvxpy.Problem na minimisado ang Frobenius-norm cost. Ang solver ay nagbabalik ng coefficients na $x_j(t)$ na naangkop para sa oras na iyon; tinitipon natin ang mga ito sa mpf_dynamic_coeffs_list. Kung nabigo ang solver para sa isang ibinigay na $t$ , bumabalik tayo sa zero coefficients upang magpatuloy ang loop.

mpf_dynamic_coeffs_list = []
for t in trotter_times:
    print(f"Computing dynamic coefficients for time={t}")
    lse = setup_dynamic_lse(
        mpf_trotter_steps,
        t,
        identity_factory,
        exact_factory,
        approx_factory,
        mps_initial_state,
    )
    problem, coeffs = setup_frobenius_problem(lse)
    try:
        problem.solve()
        mpf_dynamic_coeffs_list.append(coeffs.value)
    except Exception as error:
        mpf_dynamic_coeffs_list.append(np.zeros(len(mpf_trotter_steps)))
        print(error, "Calculation Failed for time", t)
    print("")

Computing dynamic coefficients for time=0.5

Computing dynamic coefficients for time=0.6

Computing dynamic coefficients for time=0.7

Computing dynamic coefficients for time=0.7999999999999999

Computing dynamic coefficients for time=0.8999999999999999

Computing dynamic coefficients for time=0.9999999999999999

Computing dynamic coefficients for time=1.0999999999999999

Computing dynamic coefficients for time=1.1999999999999997

Computing dynamic coefficients for time=1.2999999999999998

Computing dynamic coefficients for time=1.4

Computing dynamic coefficients for time=1.4999999999999998

Pagsamahin ang Trotter expectation values sa MPF coefficients

Ngayon ay ini-evaluate natin ang $\langle A \rangle_{\text{MPF}}(t) = \sum_j x_j \, \langle A \rangle_{k_j}(t)$ para sa bawat hanay ng coefficients (static-exact, static-approximate, at dynamic), ipinapalaganap ang per-circuit standard errors, at ini-plot ang resultang time series laban sa exact-diagonalization curve.

sym = {1: "^", 2: "s", 4: "p"}
# Get expectation values at all times for each Trotter step
for k, step in enumerate(mpf_trotter_steps):
    trotter_curve, trotter_curve_error = [], []
    for trotter_expvals, trotter_stds in zip(
        mpf_expvals_all_times, mpf_stds_all_times
    ):
        trotter_curve.append(trotter_expvals[k])
        trotter_curve_error.append(trotter_stds[k])

    plt.errorbar(
        trotter_times,
        trotter_curve,
        yerr=trotter_curve_error,
        alpha=0.5,
        markersize=4,
        marker=sym[step],
        color="grey",
        label=f"{mpf_trotter_steps[k]} Trotter steps",
    )

# Get expectation values at all times for the static MPF with exact coeffs
exact_mpf_curve, exact_mpf_curve_error = [], []
for trotter_expvals, trotter_stds in zip(
    mpf_expvals_all_times, mpf_stds_all_times
):
    mpf_std = np.sqrt(
        sum(
            [
                (coeff**2) * (std**2)
                for coeff, std in zip(coeffs_exact.value, trotter_stds)
            ]
        )
    )
    exact_mpf_curve_error.append(mpf_std)
    exact_mpf_curve.append(trotter_expvals @ coeffs_exact.value)

plt.errorbar(
    trotter_times,
    exact_mpf_curve,
    yerr=exact_mpf_curve_error,
    markersize=4,
    marker="o",
    label="Static MPF - Exact",
    color="purple",
)

# Get expectation values at all times for the static MPF with approximate coeffs
approx_mpf_curve, approx_mpf_curve_error = [], []
for trotter_expvals, trotter_stds in zip(
    mpf_expvals_all_times, mpf_stds_all_times
):
    mpf_std = np.sqrt(
        sum(
            [
                (coeff**2) * (std**2)
                for coeff, std in zip(coeffs_approx.value, trotter_stds)
            ]
        )
    )
    approx_mpf_curve_error.append(mpf_std)
    approx_mpf_curve.append(trotter_expvals @ coeffs_approx.value)

plt.errorbar(
    trotter_times,
    approx_mpf_curve,
    yerr=approx_mpf_curve_error,
    markersize=4,
    marker="o",
    label="Static MPF - Approx",
    color="orange",
)

# Get expectation values at all times for the dynamic MPF
dynamic_mpf_curve, dynamic_mpf_curve_error = [], []
for trotter_expvals, trotter_stds, dynamic_coeffs in zip(
    mpf_expvals_all_times, mpf_stds_all_times, mpf_dynamic_coeffs_list
):
    mpf_std = np.sqrt(
        sum(
            [
                (coeff**2) * (std**2)
                for coeff, std in zip(dynamic_coeffs, trotter_stds)
            ]
        )
    )
    dynamic_mpf_curve_error.append(mpf_std)
    dynamic_mpf_curve.append(trotter_expvals @ dynamic_coeffs)

plt.errorbar(
    trotter_times,
    dynamic_mpf_curve,
    yerr=dynamic_mpf_curve_error,
    markersize=4,
    marker="o",
    label="Dynamic MPF",
    color="pink",
)

# Exact expectation values
plt.plot(
    exact_evolution_times,
    exact_expvals,
    color="red",
    linestyle="--",
    label="Exact time-evolution",
)

plt.title(f"$\\langle Z_{{{L//2-1}}} Z_{{{L//2}}} \\rangle$ vs time")
plt.xlabel("Time")
plt.ylabel("Expectation Value")
plt.legend(loc="upper center", bbox_to_anchor=(0.5, -0.2), ncol=2)
plt.grid(alpha=0.1)
plt.tight_layout()
plt.show()

Output of the previous code cell

Ipinapakita ng plot sa itaas ang interplay sa pagitan ng Trotter error at sampling error.

Trotter error. Ang mga indibidwal na product formulas (kulay-abong markers) ay lumihis sa eksaktong curve nang lumihis nang lumihis habang tumataas ang oras. Ang $k=1$ circuit ay may pinakamalaking paglihis at ito ang pinakamababaw, ngunit nasa regime na rin kung saan $t/k \gtrsim 1$ , kaya malaki ang unang $1/k^{2}$ error term. Ang mga kombinasyon ng MPF (mga markers na may kulay) ay kina-cancel ang ilan sa mga leading Trotter-error terms na ito, kaya't mas malapit silang sumusubaybay sa eksaktong curve kaysa sa alinmang solong circuit na $k_j$ . Ang natitirang agwat ay sumasalamin sa higher-order Trotter terms na hindi kina-cancel ng MPF: ang isang order- $2$ , $r=3$ static MPF ay pinapatay lamang ang unang dalawang error orders, at sa malaking $t/k_{\min}$ ay unti-unting nangingibabaw ang di-nakansel na tail — kaya't hindi ginagarantiya ng MPF na mananatiling tumpak ang napakababaw na circuits sa kahit anong oras.
Sampling error. Ang mas malalapad na error bars sa MPF curves ay direktang bunga ng linear combination: ang pagpapalaganap ng independiyenteng per-circuit standard errors na $\sigma_{k_j}$ ay nagbibigay ng kabuuang variance na $\sigma_{\text{MPF}}^2 = \sum_j x_j^2 \, \sigma_{k_j}^2$ . Kaya, kung mas malaki ang $\|x\|_2$ (at sa praktika ang $\|x\|_1$ , na siyang kinokontrol natin), mas maraming shots ang kailangan upang maabot ang isang target na uncertainty. Ito ang trade-off sa likod ng opsyong approximate-solver sa Background: hinihigpitan natin ang $\|x\|_1$ upang mapanatiling manageable ang overhead na ito. Mahalagang tandaan, hindi tulad ng Trotter error, ang sampling error ay bumababa kasama ng $1/\sqrt{N_{\text{shots}}}$ , kaya't palagi itong maaaring bawasan sa pamamagitan ng paggamit ng mas maraming shots.

Sa malaking halimbawa gamit ang hardware sa ibaba, ang hardware noise ay pumapasok bilang karagdagang error source sa bawat $\langle A \rangle_{k_j}$ , na katulad ding pinapalakas ng MPF coefficients. Makikita natin kung paano nag-i-interact ang error mitigation sa MPFs sa seksyong iyon.

Malaking halimbawa gamit ang hardware

Sa seksyong ito, pinalalaki natin ang problema lampas sa kung ano ang posibleng simulahin nang eksakto. Ginagaya natin ang ilan sa mga resulta na ipinakita sa Ref. [3], gamit ang isang 50-qubit XXZ chain sa oras na $t = 3$ . Sinusundan natin ang parehong four-step workflow ng maliit na halimbawa, ngayon ay tina-target ang tunay na quantum hardware na may error mitigation. Tulad ng sa template, bawat hakbang ay minarkahan inline sa code, at maaaring sumaklaw ang isang solong hakbang sa maraming cells kapag karapat-dapat suriin ang mga intermediate outputs nito. Ang mapping ay sumasalamin sa maliit na halimbawa: magtakda ng Hamiltonian, pumili ng Trotter parameters, kalkulahin ang MPF coefficients (static at dynamic), at bumuo ng circuits. Ang mga pangunahing pagkakaiba ay:

Isang XXZ Hamiltonian sa 50 sites na may random na couplings na kinuha mula sa $\mathcal{U}(0.5, 1.5)$ (Ref. [3]).
Isang symmetric na second-order Trotter formula na may $k_j = [3, 4, 6]$ (kaya $\chi=1$ , symmetric=True).
Isang solong nakapirming oras ng ebolusyon na $t = 3$ . Sa $k_{\min}=3$ ito ay nagbibigay ng $t/k_{\min}=1$ , na nagpapanatili sa mababaw na sangkap sa loob ng Trotter-convergence regime kung saan valid ang leading-error model na inaasahan ng MPF.
Isang karagdagang single-circuit comparison run na may $k = 10$ Trotter steps, na ginamit bilang baseline. Pinili natin ang $k = 10$ dahil ang two-qubit depth nito sa hardware ay mas malalim kaysa sa pinakamalalim na sangkap ng MPF ( $k_{\max}=6$ ) kasama ang overhead ng pagpapatakbo ng maraming MPF circuits — sapat na malalim upang maging noise-limited, na siyang regime kung saan inaasahang mas magiging mahusay ang kombinasyon ng MPF kaysa sa single-circuit baseline. Ito ay isang paghahambing na "single deep circuit" laban sa kombinasyon ng MPF, hindi circuit na tina-target ang epektibong Trotter error ng MPF (na mangangailangan ng mas maraming steps).

Tandaan na kahit nasa Hakbang 1 pa tayo dito (mapping at pagbuo ng circuit), kinakalkula na rin natin nang paunang ang dynamic coefficients kasabay ng static ones sa cell na ito. Ang dynamic coefficients ay nakadepende sa $H$ at $t$ ngunit hindi sa mga quantum measurements, kaya maaari silang kalkulahin anumang oras bago ang Hakbang 4. Ginagawa natin ito ngayon upang mapanatili ang lahat ng MPF-specific na setup sa isang lugar.

# -------------------------Step 1-------------------------
L = 50
coupling_map = CouplingMap.from_line(L, bidirectional=False)

# XXZ Hamiltonian with random couplings (Ref. [3])
np.random.seed(0)
even_edges = list(coupling_map.get_edges())[::2]
odd_edges = list(coupling_map.get_edges())[1::2]

Js = np.random.uniform(0.5, 1.5, size=L)
hamiltonian = SparsePauliOp(Pauli("I" * L))
for i, edge in enumerate(even_edges + odd_edges):
    hamiltonian += SparsePauliOp.from_sparse_list(
        [
            ("XX", (edge), 2 * Js[i]),
            ("YY", (edge), 2 * Js[i]),
            ("ZZ", (edge), 4 * Js[i]),
        ],
        num_qubits=L,
    )

observable = SparsePauliOp.from_sparse_list(
    [("ZZ", (L // 2 - 1, L // 2), 1.0)], num_qubits=L
)

total_time = 3
mpf_trotter_steps = [3, 4, 6]
order = 2
symmetric = True

# Static coefficients
lse = setup_static_lse(mpf_trotter_steps, order=order, symmetric=symmetric)
mpf_coeffs = lse.solve()
print(f"Static coefficients: {mpf_coeffs}")
print(f"L1 norm: {np.linalg.norm(mpf_coeffs, ord=1)}")

model_approx, coeffs_approx = setup_sum_of_squares_problem(
    lse, max_l1_norm=2.0
)
model_approx.solve()
print(f"Approximate coefficients: {coeffs_approx.value}")
print(f"L1 norm (approx): {np.linalg.norm(coeffs_approx.value, ord=1)}")

# -------------------------Dynamic coefficients-------------------------
single_2nd_order_circ = generate_time_evolution_circuit(
    hamiltonian, time=1.0, synthesis=SuzukiTrotter(reps=1, order=order)
)
single_2nd_order_circ = pm.run(single_2nd_order_circ)

layers = slice_by_depth(single_2nd_order_circ, max_slice_depth=1)
models = [
    LayerModel.from_quantum_circuit(layer, conserve="Sz") for layer in layers
]

approx_factory = partial(
    LayerwiseEvolver,
    layers=models,
    options={
        "preserve_norm": False,
        "trunc_params": {"chi_max": 64, "svd_min": 1e-8, "trunc_cut": None},
        "max_delta_t": 4,
    },
)

single_4th_order_circ = generate_time_evolution_circuit(
    hamiltonian, time=1.0, synthesis=SuzukiTrotter(reps=1, order=4)
)
single_4th_order_circ = pm.run(single_4th_order_circ)
exact_model_layers = [
    LayerModel.from_quantum_circuit(layer, conserve="Sz")
    for layer in slice_by_depth(single_4th_order_circ, max_slice_depth=1)
]

exact_factory = partial(
    LayerwiseEvolver,
    layers=exact_model_layers,
    dt=0.1,
    options={
        "preserve_norm": False,
        "trunc_params": {"chi_max": 64, "svd_min": 1e-8, "trunc_cut": None},
        "max_delta_t": 3,
    },
)

def identity_factory():
    return MPOState.initialize_from_lattice(models[0].lat, conserve=True)

mps_initial_state = MPS_neel_state(models[0].lat)

print(f"Computing dynamic coefficients for time={total_time}")
lse_dyn = setup_dynamic_lse(
    mpf_trotter_steps,
    total_time,
    identity_factory,
    exact_factory,
    approx_factory,
    mps_initial_state,
)
problem, coeffs_dyn = setup_frobenius_problem(lse_dyn)
try:
    problem.solve()
    mpf_dynamic_coeffs = coeffs_dyn.value
except Exception as error:
    mpf_dynamic_coeffs = np.zeros(len(mpf_trotter_steps))
    print(error, "Calculation Failed")

# -------------------------Step 1 (cont): Build circuits-------------------------
mpf_circuits = []
for k in mpf_trotter_steps:
    circuit = QuantumCircuit(L)
    circuit.x([i for i in range(L) if i % 2])
    trotter_circ = generate_time_evolution_circuit(
        hamiltonian,
        synthesis=SuzukiTrotter(reps=k, order=order),
        time=total_time,
    )
    circuit.compose(trotter_circ, qubits=range(L), inplace=True)
    mpf_circuits.append(circuit)

# Baseline "single deep circuit" comparison run with k=10 Trotter steps.
# Its two-qubit depth is deeper than the deepest MPF constituent (k_max=6) plus
# the overhead of running multiple circuits, pushing it into the noise-limited
# regime where MPF is expected to outperform. It does NOT target the MPF's effective
# Trotter error (which would require many more steps).
comp_circuit = QuantumCircuit(L)
comp_circuit.x([i for i in range(L) if i % 2])
trotter_circ = generate_time_evolution_circuit(
    hamiltonian,
    synthesis=SuzukiTrotter(reps=10, order=order),
    time=total_time,
)
comp_circuit.compose(trotter_circ, qubits=range(L), inplace=True)
mpf_circuits.append(comp_circuit)

Static coefficients: [ 0.42857143 -1.82857143  2.4       ]
L1 norm: 4.65714285714286
Approximate coefficients: [-0.4942491   0.40206845  1.09218065]
L1 norm (approx): 1.9884981979026675
Computing dynamic coefficients for time=3

Ngayon ay ino-optimize natin ang mga circuit para sa piniling backend. Ginagamit natin ang Qiskit preset pass manager sa optimization_level=3, na awtomatikong pumipili ng magandang hanay ng physical qubits at ni-ru-route ang bawat circuit patungo sa device topology.

# -------------------------Step 2-------------------------
service = QiskitRuntimeService()
# backend = service.least_busy(operational=True, simulator=False, min_num_qubits=L)
backend = service.backend("ibm_fez")
print(backend)

transpiler = generate_preset_pass_manager(
    optimization_level=3, backend=backend
)
transpiled_circuits = [transpiler.run(circ) for circ in mpf_circuits]

isa_observables = [
    observable.apply_layout(circ.layout) for circ in transpiled_circuits
]

<IBMBackend('ibm_fez')>

Ang pagpapatakbo ng mas malalim na circuits sa tunay na hardware ay nangangailangan ng agresibong error mitigation. Pinapagana natin ang dynamical decoupling, gate at measurement twirling, measurement error mitigation, at zero-noise extrapolation (ZNE). Tandaan na ang mga ZNE noise factors na ginamit dito (1, 1.2, 1.4) ay mas maliit kaysa sa isang shallow-circuit scenario, dahil ang mas malalim na sangkap ng MPF ay malapit na sa noise threshold at ang malalaking noise amplifications ay maaaring itulak sila lampas sa punto kung saan reliable pa ang ZNE extrapolation.

Isinusumite natin ang lahat ng apat na circuits (tatlong MPF constituents sa $k_j = [3, 4, 6]$ kasama ang $k = 10$ baseline) sa isang solong Estimator job.

# -------------------------Step 3-------------------------
estimator = Estimator(mode=backend)
estimator.options.default_shots = 30000

# Error suppression/mitigation
estimator.options.dynamical_decoupling.enable = True
estimator.options.twirling.enable_gates = True
estimator.options.twirling.enable_measure = True
estimator.options.twirling.num_randomizations = "auto"
estimator.options.twirling.strategy = "active-accum"
estimator.options.resilience.measure_mitigation = True
estimator.options.experimental.execution_path = "gen3-turbo"

estimator.options.resilience.zne_mitigation = True
estimator.options.resilience.zne.noise_factors = (1, 1.2, 1.4)
estimator.options.resilience.zne.extrapolator = "linear"

estimator.options.environment.job_tags = ["TUT_MPF"]

job_50 = estimator.run(
    [
        (circ, observable)
        for circ, observable in zip(transpiled_circuits, isa_observables)
    ]
)

Kinukuha natin ang per-circuit expectation values at standard deviations mula sa resulta ng job, at pagkatapos ay pinagsasama natin ang mga ito sa bawat hanay ng MPF coefficients nang eksaktong tulad ng sa maliit na halimbawa: $\langle A \rangle_{\text{MPF}} = \sum_j x_j \, \langle A \rangle_{k_j}$ , na may pinalaganap na variance na $\sigma^2 = \sum_j x_j^2 \sigma_{k_j}^2$ .

# -------------------------Step 4-------------------------
result = job_50.result()
evs = [res.data.evs for res in result]
std = [res.data.stds for res in result]

print(evs)
print(std)

[array(-0.07916195), array(-0.04479681), array(-0.2560756), array(-0.06045848)]
[array(0.04605538), array(0.10056336), array(0.14426151), array(0.04059092)]

exact_mpf_std = np.sqrt(
    sum([(coeff**2) * (std**2) for coeff, std in zip(mpf_coeffs, std[:3])])
)
print(
    "Exact static MPF expectation value: ",
    evs[:3] @ mpf_coeffs,
    "+-",
    exact_mpf_std,
)
approx_mpf_std = np.sqrt(
    sum(
        [
            (coeff**2) * (std**2)
            for coeff, std in zip(coeffs_approx.value, std[:3])
        ]
    )
)
print(
    "Approximate static MPF expectation value: ",
    evs[:3] @ coeffs_approx.value,
    "+-",
    approx_mpf_std,
)
dynamic_mpf_std = np.sqrt(
    sum(
        [
            (coeff**2) * (std**2)
            for coeff, std in zip(mpf_dynamic_coeffs, std[:3])
        ]
    )
)
print(
    "Dynamic MPF expectation value: ",
    evs[:3] @ mpf_dynamic_coeffs,
    "+-",
    dynamic_mpf_std,
)

Exact static MPF expectation value:  -0.5665938395816946 +- 0.3925273058119915
Approximate static MPF expectation value:  -0.25856647611537903 +- 0.164249927266166
Dynamic MPF expectation value:  -0.12667812062949296 +- 0.06059471006973169

sym = {3: "^", 4: "s", 6: "p"}
for k, step in enumerate(mpf_trotter_steps):
    plt.errorbar(
        k,
        evs[k],
        yerr=std[k],
        alpha=0.5,
        markersize=4,
        marker=sym[step],
        color="grey",
        label=f"{mpf_trotter_steps[k]} Trotter steps",
    )

plt.errorbar(
    3,
    evs[-1],
    yerr=std[-1],
    alpha=0.5,
    markersize=8,
    marker="x",
    color="blue",
    label="10 Trotter steps",
)

plt.errorbar(
    4,
    evs[:3] @ mpf_coeffs,
    yerr=exact_mpf_std,
    markersize=4,
    marker="o",
    color="purple",
    label="Static MPF",
)

plt.errorbar(
    5,
    evs[:3] @ coeffs_approx.value,
    yerr=approx_mpf_std,
    markersize=4,
    marker="o",
    color="orange",
    label="Approximate static MPF",
)

plt.errorbar(
    6,
    evs[:3] @ mpf_dynamic_coeffs,
    yerr=dynamic_mpf_std,
    markersize=4,
    marker="o",
    color="pink",
    label="Dynamic MPF",
)

exact_obs = -0.24384471447172074  # Calculated via Tensor Network calculation
plt.axhline(
    y=exact_obs, linestyle="--", color="red", label="Exact time-evolution"
)

plt.title(
    f"$\\langle Z_{{{L//2-1}}} Z_{{{L//2}}} \\rangle$ at time {total_time} for the different methods"
)
plt.xlabel("Method")
plt.ylabel("Expectation Value")
plt.legend(loc="upper center", bbox_to_anchor=(0.5, -0.2), ncol=2)
plt.grid(alpha=0.1)
plt.tight_layout()
plt.show()

Output of the previous code cell

Ilang obserbasyon tungkol sa mga resulta ng hardware sa itaas:

Ang paglalim ay hindi libre sa hardware. Direktang ipinapakita ito ng single-circuit baselines: ang $k = 6$ circuit ay halos eksakto ( $-0.256$ laban sa reference na $-0.244$ ), gayunpaman ang mas malalim na $k = 10$ baseline ay mas masahol ( $-0.061$ , off by $\sim 0.18$ ), hindi mas mahusay. Kapag maliit na ang Trotter error, ang pagdagdag ng steps ay pangunahing nagpapalalim lamang ng circuit at nagpapatipon ng mas maraming gate noise at decoherence. Ito mismo ang regime na para dito binuo ang MPFs: maabot ang katumpakan ng isang malalim na circuit gamit lamang ang mababaw na sangkap.
Isang maliit-ang-norm na MPF ay natatalo ang malalim na single circuit. Ang approximate-static MPF (hinihigpitan sa $\|x\|_1 \approx 2$ ) ay napunta sa $-0.259$ , sa loob ng $\sim 0.015$ mula sa reference at mas malapit kaysa sa $k = 10$ baseline. Ang dynamic MPF ( $-0.127$ ) ay komportable ring natatalo ang baseline na iyon. Pareho silang pinagsasama lamang ang mababaw na $k_j = [3, 4, 6]$ circuits, ngunit nabawi ang isang sagot na hindi kayang makuha ng malalim na single circuit.
Ang norm ng coefficient ay mas mahalaga kaysa sa mathematical optimality. Ang exact-static MPF ay may $\|x\|_1 = 4.66$ at ito ang pinakamasahol na estimator sa lahat ( $-0.567$ , off by more than $0.3$ ): ang malaking coefficient norm ay pinapalakas ang natitirang gate noise, decoherence, at ZNE error sa bawat $\langle A \rangle_{k_j}$ nang halos parehong salik, na natatabunan ang pagkansel ng Trotter-error na binibili nito. Ang paghihigpit sa norm (ang approximate-static solver, $\|x\|_1 \approx 2$ ) ay inaalis ang pagkatabunan na ito at nagbibigay ng pinakamahusay na estimate — kahit hindi na eksaktong kina-cancel ng coefficients nito ang leading Trotter error.
Ang mga indibidwal na mababaw na circuit ay maaari pa ring maging competitive. Ang nag-iisang $k = 6$ na sangkap ( $-0.256$ ) ay halos eksakto na rin dito — sa run na ito, ito ay bahagyang mas malapit pa kaysa sa approximate-static MPF. Ang problema ay hindi mo alam nang maaga kung aling solong $k$ ang nasa sweet spot ng "converged pero hindi pa noise-limited," at ang mukhang-ligtas na pagpipiliang basta pinalalim lang ( $k = 10$ ) upang garantiyahan ang Trotter convergence ay eksaktong ang nabibigo. Ang MPF ay nagbibigay ng makatwirang kombinasyon ng mababaw na circuits na hindi nangangailangan ng paghula sa tamang lalim.

Ang praktikal na aral ay sa hardware, ang MPFs ay dapat ipares sa malakas na error mitigation sa bawat indibidwal na $\langle A \rangle_{k_j}$ , ang $L_1$ -norm ng coefficient ay dapat panatilihing katamtaman (gumamit ng approximate solver, o ng dynamic MPF), at ang Trotter steps na $k_j$ ay dapat piliin upang $t/k_{\min} \lesssim 1$ — dito ang $k_{\min} = 3$ sa $t = 3$ ay nagbibigay ng $t/k_{\min} = 1$ , na pinananatili ang mga sangkap sa loob ng convergent regime kung saan valid ang leading-error model na inaasahan ng static MPF. Sa mga pagpipiliang iyon, ang small-norm MPFs dito ay tumutugma sa isang converged single circuit habang ang naive na "basta pinalalim lang" baseline ay hindi, na nabawi ang bentaha sa depth-versus-accuracy na ipinakita sa Ref. [3]. Tandaan din na ang mga indibidwal na run ay maingay — sa ibang submission ng parehong job (o ibang backend), maaaring magbago ang eksaktong pagkakasunod-sunod; ang matatag na trends ay ang small- $\|x\|_1$ MPFs ay gumaganda, ang large- $\|x\|_1$ exact-static MPF ay pinapalakas ng hardware noise, at ang over-deep na single circuit ay noise-limited.

Mga susunod na hakbang

Mga rekomendasyon

Kung nasiyahan ka sa gawaing ito, maaaring interesado ka sa mga sumusunod na materyal:

Paano piliin ang Trotter steps para sa MPF — praktikal na gabay sa pagpili ng mga halaga ng $k_j$ upang maiwasan ang mga instability
Paano gamitin ang approximate model — pag-tune ng $L_1$ -norm constraint at solver options para sa approximate static MPF
qiskit-addon-mpf API reference — buong dokumentasyon para sa static, dynamic, at backend modules

Mga Sanggunian

[1] Vazquez, A. C., Egger, D. J., Ochsner, D., & Woerner, S. Well-conditioned multi-product formulas for hardware-friendly Hamiltonian simulation. Quantum, 7, 1067 (2023)

[2] Zhuk, S., Robertson, N. F., & Bravyi, S. Trotter error bounds and dynamic multi-product formulas for Hamiltonian simulation. Physical Review Research, 6(3), 033309 (2024)

[3] Robertson, N. F., et al. Tensor network enhanced dynamic multiproduct formulas. arXiv:2407.17405 (2024)

Mga learning outcome​

Prerequisites​

Pangkalahatang-ideya​

Ano ang mga multi-product formula?​

Paano tinutukoy ang mga coefficient?​

Kailan gagamitin ang mga MPF​

Kapag hindi nakakatulong ang MPFs​

Ano ang sinasaklaw ng tutorial na ito​

Mga Kinakailangan​

Setup​

Maliit na halimbawa gamit ang simulator​

Hakbang 1: I-map ang mga classical input sa isang quantum na problema​

Bumuo ng mga Trotter circuit​

Hakbang 2: I-optimize ang problema para sa pagpapatupad sa quantum hardware​

Hakbang 3: Isagawa gamit ang Qiskit primitives​

Hakbang 4: I-post-process at ibalik ang resulta sa nais na klasikal na format​

Mga static na MPF coefficient​

Itakda ang LSE​

I-optimize ang xxx gamit ang isang eksaktong modelo​

I-optimize ang xxx gamit ang isang tinantyang modelo​

Mga dynamic na MPF coefficient​

Pagsamahin ang Trotter expectation values sa MPF coefficients​

Malaking halimbawa gamit ang hardware​

Mga susunod na hakbang​

Mga Sanggunian​