量子回路を用いて量子材料中の中性子散乱をシミュレートする
使用量の見積もり: Heron r2プロセッサで13分(注: これはあくまで見積もりです。実際の実行時間は異なる場合があります。)
どのチュートリアルを使うべきか?
このチュートリアルを使用すると、古典計算をノートパソコン上でローカルに実行し、ハードウェア実行をIBM Quantum® QPU上で行いながら、実装をステップごとに学ぶことができます。大規模な例では大量のメモリが必要であり、近似量子コンパイル(AQC)は一般的なノートパソコンのリソースで数時間かかることがあります。圧縮をQiskit Serverlessの計算リソースとメモリリソースにオフロードするには、対応するServerlessチュートリアルを使用してください。
学習成果
-
非弾性中性子散乱(INS)スペクトルが、量子スピン模型の動的構造因子(DSF)とどのように結びついているか。
-
基底状態を準備し、局所的な摂動を加え、量子回路上でトロッター時間発展を実行する方法。
-
qiskit-addon-aqc-tensorを使った近似量子コンパイル(AQC)により、ハードウェア実行のために深いトロッター回路を圧縮する方法。 -
キュービットの期待値から遅延グリーン関数(RGF)を抽出し、それをフーリエ変換してDSFを得る方法。
前提条件
背景
このチュートリアルでは、Lee他, arXiv:2603.15608の結果を再現します。
中性子散乱は、バッテリー電極を含む持続可能性に関連する材料の特性評価に役立ちます。このチュートリアルでは、量子回路を用いて磁気励起をシミュレートし、理論モデルを散乱測定に結び付ける方法を探ります。その量子磁性体の例では、材料挙動を研究する手法を発展させ、材料発見のためのより広範な研究ツールキットに貢献します。
非弾性中性子散乱と動的構造因子
非弾性中性子散乱(INS)は、量子材料中の磁気励起を調べる最も強力な実験プローブの1つです。熱中性子または冷中性子のビームが結晶に入射すると、個々の中性子は運動量とエネルギーの両方を磁気サブシステムと交換します。測定される散乱強度は、動的構造因子(DSF)に比例します。
これはスピン自由度の完全な時空相関を符号化します。
KCuF: 典型的なラッティンジャー液体磁性体
フッ化カリウム銅(KCuF)は準1次元反強磁性体であり、スピンのCuイオンの鎖が最近接ハイゼンベルグ交換を介して相互作用し、鎖間結合はのわずかです。INSデータが得られるでは、スペクトルはトモナガ・ラッティンジャー液体に特徴的な分数化されたスピノン励起に支配されます。鎖内ダイナミクスは、等方点()における1次元スピン XXZハミルトニアンによってよく記述されるため、
KCuFは量子シミュレーションの理想的なベンチマークとなります。ハミルトニアンは量子プロセッサ上で実装できるほど単純である一方、基底状態は強くエンタングルしており、励起スペクトルは広い2スピノン連続体を特徴とします。
注: このチュートリアルではエネルギー単位としてを設定し、規格化を採用しています。これは、論文の回路実装(補足図S3)で使われている局所ハミルトニアンに対応します。論文の完全なハミルトニアン(式3)には全体に追加の因子2が掛かっているため、論文のはここで使うの2倍です。
シミュレートして測定する内容
計算する物理量は遅延グリーン関数(RGF)であり、時間依存のスピン-スピン相関関数として定義されます。
ここでは基準サイト(鎖の中心)であり、はハイゼンベルグ描像のスピン演算子です。このチュートリアルでは成分()に焦点を当てます。量子コンピュータ上では、基底状態を準備し、に局所的な摂動を加え、摂動された状態を時間発展させ、各時間ステップですべてのサイトについて単一キュービットの期待値を測定することでRGFにアクセスします。重要な洞察は、各が基底状態磁化からの差を与えるということです。等方的なハイゼンベルグ反強磁性体はサイトあたりの正味の磁化がゼロ()であるため、測定された生の値は明示的な差し引きなしに直接を与えます。
すべてのサイトと時間ステップにわたってを収集することで、2次元のデータセットを構築し、それを空間と時間の両方でフーリエ変換して動的構造因子を得ます。DSFはINS実験で直接測定される量であり、各運動量とエネルギーにどのような磁気励起が存在するかを教えてくれます。等方的なハイゼンベルグ鎖の場合、厳密な励起スペクトルは2スピノン連続体であり、散乱強度の広い帯状の形状は量子シミュレーションの厳格なエンドツーエンドのベンチマークとして機能します。基底状態の準備、摂動、トロッター時間発展、測定プロトコルのすべてを一度に検証します。
量子シミュレーション・ワークフロー
このワークフローはINSイベントの物理を反映しています。(1) 個のキュービット上で多体基底状態を準備し、(2) 中性子のスピン移動を模倣するために鎖の中心に局所的なスピンフリップ摂動を加え、(3) 2次トロッター化を使っての下で離散的な時間ステップで発展させ、(4) 各ステップですべてのキュービット上でを測定してRGFを得ます。次に、2次元離散フーリエ変換によってDSF が得られます。
各時間ステップで測定する観測量は、各キュービット上のです。Qiskitでは、これはSparsePauliOp演算子のリストとして表現されます。各サイトに対応する、キュービットの恒等文字列に埋め込まれた単一キュービットのです。これらの観測量は、問題マッピングのステップ(ステップ1)で問題インスタンスごとに1回構築され、そのスケールのすべての回路で再利用されます。
近似量子コンパイル(AQC)
深いTrotter回路は**近似量子コンパイル(AQC)**によって圧縮できます。これは、最初の数層のTrotter層を、パラメータが古典的に最適化されて元の深い回路とのMPSレベルの忠実度を最大化するより短いパラメータ化されたアンザッツに置き換えます。残りのTrotterステップは正確に付加され、2キュービット・ゲートが大幅に少ない「混合」AQC+Trotter回路が生成されます。
MPSシミュレーション
1次元系の場合、行列積状態(MPS)法は、基底状態の準備(密度行列繰り込み群、すなわちDMRGを用いる)と回路レベルの時間発展の両方を効率的にシミュレートできます。ボンド次元を制御することで、精度と計算コストのトレードオフを調整できます。このチュートリアルでは、qiskit-addon-aqc-tensor内でMPSシミュレーションを使用して、ハードウェア実行のために深いトロッター回路を圧縮する高忠実度のAQCアンザッツを計算します。
要件
このチュートリアルを始める前に、以下がインストールされていることを確認してください。
-
可視化サポート付きのQiskit SDK
-
Qiskit Runtime(
pip install qiskit-ibm-runtime) -
quimbとJAXのエクストラ付きのqiskit-addon-aqc-tensor(pip install 'qiskit-addon-aqc-tensor[quimb-jax]')
セットアップ
# 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_aqc_tensor": "qiskit-addon-aqc-tensor", "qiskit_ibm_runtime": "qiskit-ibm-runtime", "qiskit_quimb": "qiskit-quimb", "quimb": "quimb", "scipy": "scipy"}
_missing = [pip for module, pip in _needed.items()
if importlib.util.find_spec(module) is None]
# One at a time, so a package that fails to install does not block the others
for _pip in _missing:
%pip install -q {_pip}
if not _missing:
print("\u2713 All required packages are installed")
import timeit
import warnings
from collections.abc import Iterator, Sequence
from functools import partial
import matplotlib.pyplot as plt
import numpy as np
import quimb.tensor as qtn
import scipy.optimize
from numpy.typing import NDArray
from qiskit import QuantumCircuit
from qiskit.circuit import CircuitInstruction, ParameterVector, Qubit
from qiskit.circuit.library import PauliEvolutionGate
from qiskit.primitives import StatevectorEstimator
from qiskit.quantum_info import SparsePauliOp
from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager
from qiskit_addon_aqc_tensor import generate_ansatz_from_circuit
from qiskit_addon_aqc_tensor.objective import MaximizeStateFidelity
from qiskit_addon_aqc_tensor.simulation import tensornetwork_from_circuit
from qiskit_addon_aqc_tensor.simulation.quimb import QuimbSimulator
from qiskit_ibm_runtime import EstimatorV2 as Estimator
from qiskit_ibm_runtime import QiskitRuntimeService
from qiskit_quimb import quimb_circuit
from scipy.sparse import SparseEfficiencyWarning
# scipy.sparse.linalg.expm (reached via PauliEvolutionGate.to_matrix) prefers
# CSC and warns when handed the CSR matrix from SparsePauliOp.to_matrix; the
# result is unaffected, so silence the cosmetic warning.
warnings.filterwarnings("ignore", category=SparseEfficiencyWarning)
def xxz_hamiltonian_mpo(
n_qubits: int, interaction: float = 1.0, anisotropy: float = 1.0
) -> qtn.MatrixProductOperator:
"""1D XXZ Hamiltonian as a quimb MPO.
Builds the Hamiltonian using ``qtn.SpinHam1D``.
Args:
n_qubits: Number of sites.
interaction: Overall interaction strength (J in the paper).
anisotropy: XY/Z anisotropy (ε in the paper). ``ε = 1`` is the
isotropic Heisenberg point; ``ε = 0`` is the Ising limit.
Returns:
The Hamiltonian as a matrix product operator.
"""
builder = qtn.SpinHam1D(S=1 / 2)
builder += interaction * anisotropy * 0.5, "+", "-"
builder += interaction * anisotropy * 0.5, "-", "+"
builder += interaction, "Z", "Z"
return builder.build_mpo(L=n_qubits)
def build_ground_state_ansatz(n_qubits: int, n_layers: int) -> QuantumCircuit:
"""Hamiltonian variational ansatz (HVA) circuit for ground-state preparation.
Starts from a product of singlet pairs and applies alternating
odd/even layers of parameterized XXZ pair-evolution gates.
The returned circuit is parameterized: it carries a ``ParameterVector``
named ``"theta"`` of length ``2 * n_layers`` whose values must be
assigned (e.g. via ``circuit.assign_parameters``) before simulation.
Element ``2 * r`` is the odd-layer angle and element ``2 * r + 1`` is the
even-layer angle of layer ``r``.
Args:
n_qubits: Number of qubits (must be even).
n_layers: Number of HVA layers.
Returns:
The parameterized HVA preparation circuit.
"""
theta = ParameterVector("theta", 2 * n_layers)
circuit = QuantumCircuit(n_qubits)
# Initial singlet product state
for i in range(n_qubits // 2):
circuit.x(2 * i)
circuit.x(2 * i + 1)
circuit.h(2 * i + 1)
circuit.cx(2 * i + 1, 2 * i)
# Variational layers: each pair gate is exp(-i theta (XX+YY+ZZ)/2)
pair_ham = SparsePauliOp(
["XX", "YY", "ZZ"], coeffs=[0.5, 0.5, 0.5]
) # H_pair (HVA form)
for r in range(n_layers):
for i in range(1, (n_qubits + 1) // 2): # odd layer
circuit.append(
PauliEvolutionGate(pair_ham, time=theta[2 * r]),
[2 * i - 1, 2 * i],
)
for i in range(n_qubits // 2): # even layer
circuit.append(
PauliEvolutionGate(pair_ham, time=theta[2 * r + 1]),
[2 * i, 2 * i + 1],
)
return circuit
def optimize_ground_state_ansatz(
ansatz: QuantumCircuit,
x0: NDArray[np.floating],
target_mps: qtn.MatrixProductState,
*,
max_bond: int | None = None,
cutoff: float = 1e-10,
method: str = "COBYQA",
options: dict | None = None,
) -> scipy.optimize.OptimizeResult:
"""Optimize HVA parameters by maximizing fidelity with a target MPS.
The HVA circuit is simulated as a matrix product state with the given
bond-dimension truncation, and the parameters are optimized to maximize
the state fidelity ``|<psi_HVA | target_mps>|**2`` with the DMRG ground
state. Both states are normalized, so the minimized objective is the
infidelity ``1 - |<psi_HVA | target_mps>|**2``.
Args:
ansatz: Parameterized HVA circuit from ``build_ground_state_ansatz``.
The length of ``x0`` must equal ``ansatz.num_parameters``.
x0: Initial parameters.
target_mps: Target MPS (DMRG ground state) to maximize fidelity with.
max_bond: Maximum MPS bond dimension during gate application.
cutoff: Singular-value cutoff during gate application.
method: ``scipy.optimize.minimize`` method.
options: Options dict forwarded to ``scipy.optimize.minimize``.
Returns:
The Scipy OptimizeResult.
"""
def infidelity(params: NDArray[np.floating]) -> float:
circuit = ansatz.assign_parameters(params)
circuit_mps = quimb_circuit(
circuit.decompose(["PauliEvolution"]),
quimb_circuit_class=qtn.CircuitMPS,
max_bond=max_bond,
cutoff=cutoff,
)
return 1 - abs(circuit_mps.psi.H @ target_mps) ** 2
return scipy.optimize.minimize(
infidelity, np.asarray(x0), method=method, options=options
)
def trotter_evolution(
qubits: Sequence[Qubit],
interaction: float,
anisotropy: float,
time_step: float,
n_steps: int,
) -> Iterator[CircuitInstruction]:
"""Second-order Trotter steps of the XXZ pair Hamiltonian.
While the paper used a hand-optimized circuit for the Trotter steps, we use
PauliEvolutionGate here for simplicity and generality. The final two-qubit gate
count and gate depth are equivalent when transpiled with ``optimization_level=3``.
Args:
qubits: Qubits to act on (length ``n_qubits``).
interaction: Overall interaction strength (J in the paper).
anisotropy: XY/Z anisotropy (ε in the paper). ``ε = 1`` is the
isotropic Heisenberg point; ``ε = 0`` is the Ising limit.
time_step: Per-step Trotter time.
n_steps: Number of Trotter steps.
Yields:
``CircuitInstruction``s implementing the Trotter steps.
"""
if n_steps == 0:
return
n_qubits = len(qubits)
pair_ham = SparsePauliOp(
["XX", "YY", "ZZ"],
coeffs=[
0.25 * interaction * anisotropy,
0.25 * interaction * anisotropy,
0.25 * interaction,
],
)
half_evo = PauliEvolutionGate(pair_ham, time=time_step / 2)
full_evo = PauliEvolutionGate(pair_ham, time=time_step)
for i in range(n_qubits // 2): # half even layer
yield CircuitInstruction(half_evo, (qubits[2 * i], qubits[2 * i + 1]))
for i in range(n_qubits // 2 - 1): # full odd layer
yield CircuitInstruction(
full_evo, (qubits[2 * i + 1], qubits[2 * i + 2])
)
for _ in range(n_steps - 1): # interior steps
for i in range(n_qubits // 2):
yield CircuitInstruction(
full_evo, (qubits[2 * i], qubits[2 * i + 1])
)
for i in range(n_qubits // 2 - 1):
yield CircuitInstruction(
full_evo, (qubits[2 * i + 1], qubits[2 * i + 2])
)
for i in range(n_qubits // 2): # half even layer
yield CircuitInstruction(half_evo, (qubits[2 * i], qubits[2 * i + 1]))
def get_dsf(
n_qubits: int,
rgf_mat: NDArray[np.floating],
time_step: float,
n_steps: int,
n_points_momentum: int,
n_points_frequency: int,
) -> NDArray[np.floating]:
"""Compute the dynamical structure factor from the retarded Green's function.
Uses the center-site approximation and a discrete Fourier transform.
The result is symmetrized about the momentum axis and clipped to
non-negative values, ready for plotting.
Args:
n_qubits: Number of qubits (sites).
rgf_mat: RGF matrix of shape ``(n_steps, n_qubits)``.
time_step: Trotter time-step size.
n_steps: Number of time steps.
n_points_momentum: Number of momentum points.
n_points_frequency: Number of frequency points.
Returns:
DSF array of shape ``(n_points_frequency, n_points_momentum)``,
symmetrized about the momentum axis and clipped to non-negative values.
"""
max_frequency = np.pi / time_step
momentum_range = np.linspace(0, 2 * np.pi, n_points_momentum)
frequency_range = np.linspace(0, max_frequency, n_points_frequency)
result = np.zeros((frequency_range.shape[0], momentum_range.shape[0]))
center = n_qubits // 2 - 1
for iw, w in enumerate(frequency_range):
exponent = np.exp(1j * w * time_step * np.arange(1, n_steps + 1))
# S = sigma/2, so two spin operators contribute a factor of (1/2)(1/2) = 1/4.
rgf_omega = (
np.dot(rgf_mat.T, exponent) * time_step / 4
) # S(omega): time Fourier slice of the Green's function
for iq, q in enumerate(momentum_range):
momentum_phases = np.exp(
-1j * q * np.arange(-center, center + 2, 1)
)
result[iw, iq] = np.imag(np.dot(rgf_omega, momentum_phases))
result = -(result + result[:, ::-1]) / 2
result = np.clip(result, a_min=0, a_max=None)
return result
def plot_dsf(
dsf: NDArray[np.floating],
time_step: float,
n_points_momentum: int,
n_points_frequency: int,
title: str | None = None,
) -> None:
"""Heat-map of the dynamical structure factor.
Args:
dsf: DSF array of shape ``(n_points_frequency, n_points_momentum)``.
time_step: Trotter time-step size.
n_points_momentum: Number of momentum points.
n_points_frequency: Number of frequency points.
title: Optional plot title.
"""
max_frequency = np.pi / time_step
momentum_range = np.linspace(0, 2 * np.pi, n_points_momentum)
frequency_range = np.linspace(0, max_frequency, n_points_frequency)
x, y = np.meshgrid(momentum_range, frequency_range)
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")
ax.set_ylim(0, 3.6)
ax.set_xlim(0, 2 * np.pi)
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 title:
ax.set_title(title, fontsize=14)
plt.tight_layout()
plt.show()
def plot_rgf(
n_qubits: int,
rgf_mat: NDArray[np.floating],
time_step: float,
n_steps: int,
title: str | None = None,
) -> None:
"""Heat-map of the retarded Green's function in real space and time.
Args:
n_qubits: Number of qubits (sites).
rgf_mat: RGF matrix of shape ``(n_steps, n_qubits)``.
time_step: Trotter time-step size.
n_steps: Number of time steps.
title: Optional plot title.
"""
fig, ax = plt.subplots(figsize=(8, 6))
qubit_axis = np.arange(n_qubits)
t_axis = np.arange(1, n_steps + 1) * time_step
x, y = np.meshgrid(qubit_axis, t_axis)
c = ax.pcolormesh(
x,
y,
np.real(rgf_mat),
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("Qubit", fontsize=16)
ax.xaxis.set_major_locator(
plt.matplotlib.ticker.MaxNLocator(integer=True)
)
ax.set_ylabel(r"Time", fontsize=16)
if title:
ax.set_title(title, fontsize=14)
plt.tight_layout()
plt.show()
def uniform_2q_depth(circuit: QuantumCircuit) -> int:
"""Two-qubit gate depth in a standardized basis."""
pass_manager = generate_preset_pass_manager(
optimization_level=0, basis_gates=["cz", "id", "rz", "sx", "x"]
)
return pass_manager.run(circuit).depth(
lambda inst: inst.operation.num_qubits == 2
)
小規模シミュレータの例
まず、10キュービットでの完全なワークフローを示します。MPSシミュレーションを使ってHVA基底状態アンザッツを最適化し、時間発展にはQiskitステートベクトル・シミュレータを使用します。DMRGは基準となる基底状態エネルギーと基準MPSを与えます。この小規模の例により、スケールアップする前にすべてのステップを検証できます。
ステップ1: 古典入力を量子問題にマッピングする
まず、物理モデルを定義し、量子回路を構築することから始めます。
ハミルトニアン。 KCuFは等方点(、、エネルギー単位としてを設定)における1次元XXZハミルトニアンでモデル化されます。
基底状態。 基底状態準備回路として、ハミルトニアン変分アンザッツ(HVA)回路build_ground_state_ansatzを使用します。HVAパラメータは、DMRG基底状態MPSとの状態忠実度を最大化することで古典的に最適化されます。ここで、HVA状態は回路をquimbのCircuitMPSで行列積状態としてシミュレートすることで評価されます。最適化にはscipy.optimize.minimizeを使用します。最適化されたアンザッツのエネルギーも参考のため計算されます。
トロッター・ゲート。 各最近接相互作用項()は、PauliEvolutionGate(H_pair, time=time_step)で構築されます。Qiskitはトランスパイル中にこれを最適な3CNOT分解に合成します。
摂動。 中心キュービット上のゲートがを実装し、散乱された中性子によって生じる局所的なスピンフリップを模倣します。
観測量。 各キュービット・サイトについて観測量を構築します。これらのSparsePauliOpオブジェクトは、ステップ3でEstimatorプリミティブに渡され、各時間ステップでを抽出します。
# -- Physical parameters --
n_qubits = 10
interaction = 1.0 # J
anisotropy = 1.0 # ε (isotropic point)
time_step = 0.6
n_steps = 10
mps_max_bond = 32
mps_cutoff = 1e-8
center = n_qubits // 2 - 1
# -- Hamiltonian MPO --
ham_mpo = xxz_hamiltonian_mpo(n_qubits, interaction, anisotropy)
# -- Reference ground-state energy via DMRG --
dmrg = qtn.DMRG2(ham_mpo)
dmrg.solve(tol=1e-8)
print(f"Ground-state energy (DMRG): {dmrg.energy:.6f}")
# -- Build ground state ansatz circuit --
gs_n_layers = 3
gs_ansatz = build_ground_state_ansatz(n_qubits, gs_n_layers)
# -- Optimize ground state ansatz parameters to maximize fidelity with the DMRG MPS --
# Initialize odd-layer angles near 0 (where the inter-pair gate is the
# identity) and even-layer angles near pi/2 (where the intra-pair gate
# is a SWAP, since 0.5*(XX + YY + ZZ) = SWAP - I/2).
rng = np.random.default_rng(12345)
x0 = np.tile([0, np.pi / 2], gs_n_layers) + rng.normal(
scale=0.1, size=2 * gs_n_layers
)
print("Optimizing ground state ansatz...")
t0 = timeit.default_timer()
result = optimize_ground_state_ansatz(
gs_ansatz,
x0,
dmrg.state,
max_bond=mps_max_bond,
cutoff=mps_cutoff,
options=dict(maxiter=100),
)
t1 = timeit.default_timer()
print(f"Finished optimizing ground state ansatz in {t1 - t0} seconds.")
print(f"Ground state ansatz fidelity: {1 - result.fun:.6f}")
gs_circuit = gs_ansatz.assign_parameters(result.x)
gs_circuit_mps = quimb_circuit(
gs_circuit.decompose(["PauliEvolution"]),
quimb_circuit_class=qtn.CircuitMPS,
max_bond=mps_max_bond,
cutoff=mps_cutoff,
)
gs_ansatz_energy = qtn.expec_TN_1D(
gs_circuit_mps.psi.H, ham_mpo, gs_circuit_mps.psi
)
print(f"Ground state ansatz energy: {gs_ansatz_energy:.6f}")
# -- Build circuits for each time step --
perturbed = gs_circuit.copy()
perturbed.rz(np.pi / 2, center)
circuits = []
for t in range(1, n_steps + 1):
circuit = perturbed.copy()
for instr in trotter_evolution(
circuit.qubits, interaction, anisotropy, time_step, t
):
circuit.append(instr)
circuits.append(circuit)
print(
f"Built {len(circuits)} circuits, deepest 2q depth (uniform basis) = "
f"{uniform_2q_depth(circuits[-1])}"
)
# -- Observables: Z on each qubit site --
observables = [
SparsePauliOp.from_sparse_list([("Z", [i], 1)], num_qubits=n_qubits)
for i in range(n_qubits)
]
print(f"Defined {len(observables)} Z observables.")
Ground-state energy (DMRG): -4.258035
Optimizing ground state ansatz...
Finished optimizing ground state ansatz in 8.850687621976249 seconds.
Ground state ansatz fidelity: 0.984277
Ground state ansatz energy: -4.232565
Built 10 circuits, deepest 2q depth (uniform basis) = 163
Defined 10 Z observables.
ステップ2: 量子ハードウェア実行のために問題を最適化する
実際のハードウェアでは、上記のトロッター回路は深すぎます。**近似量子コンパイル(AQC)**は、最初の層のトロッター層(基底状態回路を含む)を、元の深い回路とのMPSレベルの忠実度を最大化するように最適化されたより短いパラメータ化されたアンザッツに置き換えることでこれに対処します。残りのトロッター・ステップは正確に付加され、より浅い「AQC+トロッター」回路が生成されます。
表現力と回路の深さのバランスを取るため、2つのアンザッツが使用されます。単層アンザッツ(単一のトロッター・ステップから生成)は最も初期の時間ステップを圧縮し、より深い2層アンザッツ(2つのトロッター・ステップから生成)は、より高い忠実度が必要な次のいくつかのステップを圧縮します。
AQCワークフローには4つのサブステップがあります。
-
目標回路を構築する — ステップ1からの最初の個の回路が直接AQCのターゲットとして機能します。
-
目標MPSを計算する —
quimb.tensor.CircuitMPSを使用して、各目標回路を行列積状態としてシミュレートします。 -
アンザッツを生成して最適化する —
generate_ansatz_from_circuitが単層と2層のパラメータ化されたアンザッツを作成します。パラメータは、を最小化するために、JAXで高速化された勾配を使ったL-BFGS-Bによって最適化されます。最初のステップは単層アンザッツを使用し、次のステップは2層アンザッツを使用します。各ステージ内では、各ステップは前のステップの最適化されたパラメータからウォームスタートし、ステージの境界でパラメータはステージのデフォルトにリセットされます。 -
混合回路を組み立てる — AQCチェックポイントを超える時間ステップについては、
trotter_evolutionを使用して、最適化された2層のAQC回路に正確なトロッター層を追加します。
# Number of time steps to compress into an AQC ansatz with one layer
aqc_n_steps_1 = 3
# Number of time steps to compress into an AQC ansatz with two layers
aqc_n_steps_2 = 2
aqc_n_steps_total = aqc_n_steps_1 + aqc_n_steps_2
# ── Step 2a: Target circuits (first aqc_n_steps_total circuits from Step 1) ──
target_circuits = {
k: circuits[k - 1] for k in range(1, aqc_n_steps_total + 1)
}
# ── Step 2b: Compute target MPS ──
# Decompose PauliEvolutionGate to RXX/RYY/RZZ before passing to the AQC MPS
# backend, which does not understand PauliEvolutionGate natively.
aqc_sim = QuimbSimulator(
quimb_circuit_factory=partial(
qtn.CircuitMPS,
gate_opts=dict(max_bond=mps_max_bond, cutoff=mps_cutoff),
),
autodiff_backend="jax",
)
print("\nStep 2b — target MPS:")
target_mps = {}
for k in range(1, aqc_n_steps_total + 1):
target_mps[k] = tensornetwork_from_circuit(
target_circuits[k].decompose(["PauliEvolution"]), aqc_sim
)
print(f" k={k}: max bond = {target_mps[k].psi.max_bond()}")
# ── Step 2c: Generate ansätze (one-layer and two-layer) and optimize parameters ──
ansatz_1, initial_params_1 = generate_ansatz_from_circuit(
target_circuits[1].decompose(["PauliEvolution"]),
qubits_initially_zero=True,
)
initial_params_1 = np.array(initial_params_1)
ansatz_2, initial_params_2 = generate_ansatz_from_circuit(
target_circuits[2].decompose(["PauliEvolution"]),
qubits_initially_zero=True,
)
initial_params_2 = np.array(initial_params_2)
print(
f"\nStep 2c — one-layer ansatz: {ansatz_1.num_parameters} parameters, "
f"2q depth (uniform basis) = {uniform_2q_depth(ansatz_1)}"
)
print(
f"Step 2c — two-layer ansatz: {ansatz_2.num_parameters} parameters, "
f"2q depth (uniform basis) = {uniform_2q_depth(ansatz_2)}"
)
aqc_circuits = {}
aqc_params = {}
for k in range(1, aqc_n_steps_total + 1):
if k <= aqc_n_steps_1:
ansatz, base_params = ansatz_1, initial_params_1
else:
ansatz, base_params = ansatz_2, initial_params_2
# Warm-start from the previous step only within the same stage
same_stage = (k - 1 >= 1) and (
(k - 1 <= aqc_n_steps_1) == (k <= aqc_n_steps_1)
)
x0 = aqc_params[k - 1] if same_stage else base_params
obj = MaximizeStateFidelity(target_mps[k], ansatz, aqc_sim)
t0 = timeit.default_timer()
result = scipy.optimize.minimize(
obj.loss_function,
x0,
method="L-BFGS-B",
jac=True,
options=dict(maxiter=100),
)
elapsed = timeit.default_timer() - t0
aqc_params[k] = result.x
aqc_circuits[k] = ansatz.assign_parameters(result.x)
print(
f" k={k}: fidelity = {1 - result.fun:.4f}, "
f"2q depth (uniform basis) = "
f"{uniform_2q_depth(aqc_circuits[k])}, "
f"{elapsed:.1f}s"
)
# ── Step 2d: Assemble full circuit set (AQC + Trotter) ──
all_circuits = []
for k in range(1, aqc_n_steps_total + 1):
all_circuits.append(aqc_circuits[k])
base = aqc_circuits[aqc_n_steps_total]
for k in range(1, n_steps - aqc_n_steps_total + 1):
circuit = base.copy()
for instr in trotter_evolution(
circuit.qubits, interaction, anisotropy, time_step, k
):
circuit.append(instr)
all_circuits.append(circuit)
full_depths = [
uniform_2q_depth(circuits[k - 1]) for k in range(1, n_steps + 1)
]
aqc_depths = [uniform_2q_depth(circuit) for circuit in all_circuits]
full_2q = [
circuits[k - 1].decompose(["PauliEvolution"]).num_nonlocal_gates()
for k in range(1, n_steps + 1)
]
aqc_2q = [
circuit.decompose(["PauliEvolution"]).num_nonlocal_gates()
for circuit in all_circuits
]
print(f"\nStep 2d — assembled {len(all_circuits)} circuits")
print(
f" At step {n_steps} (uniform basis): "
f"Trotter 2q depth = {full_depths[-1]}, AQC+Trotter 2q depth = {aqc_depths[-1]}; "
f"2q gates: Trotter = {full_2q[-1]}, AQC+Trotter = {aqc_2q[-1]}"
)
steps_axis = np.arange(1, n_steps + 1)
fig, ax = plt.subplots(figsize=(7, 4))
ax.plot(steps_axis, full_depths, "-o", color="black", label="Full Trotter")
ax.plot(
steps_axis, aqc_depths, "-o", color="cadetblue", label="AQC + Trotter"
)
ax.set_xlabel("Trotter steps", fontsize=13)
ax.xaxis.set_major_locator(plt.matplotlib.ticker.MaxNLocator(integer=True))
ax.set_ylabel("2q gate depth (uniform basis)", fontsize=13)
ax.set_title(f"AQC circuit-depth reduction ({n_qubits} qubits)", fontsize=13)
ax.legend(fontsize=11)
plt.tight_layout()
plt.show()
Step 2b — target MPS:
k=1: max bond = 22
k=2: max bond = 22
k=3: max bond = 26
k=4: max bond = 27
k=5: max bond = 30
Step 2c — one-layer ansatz: 389 parameters, 2q depth (uniform basis) = 27
Step 2c — two-layer ansatz: 470 parameters, 2q depth (uniform basis) = 33
k=1: fidelity = 1.0000, 2q depth (uniform basis) = 27, 11.8s
k=2: fidelity = 0.9980, 2q depth (uniform basis) = 27, 13.2s
k=3: fidelity = 0.9890, 2q depth (uniform basis) = 27, 13.7s
k=4: fidelity = 0.9983, 2q depth (uniform basis) = 33, 16.4s
k=5: fidelity = 0.9957, 2q depth (uniform basis) = 33, 16.4s
Step 2d — assembled 10 circuits
At step 10 (uniform basis): Trotter 2q depth = 163, AQC+Trotter 2q depth = 99; 2q gates: Trotter = 371, AQC+Trotter = 200
ステップ3: Qiskitプリミティブを使って実行する
StatevectorEstimatorプリミティブを使用して、各AQCコンパイル済み回路をシミュレートします。
estimator = StatevectorEstimator()
pubs = [(circuit, observables) for circuit in all_circuits]
job = estimator.run(pubs)
result = job.result()
ステップ4: 望ましい古典形式で結果を後処理して返す
ここで、各時間ステップにおける各キュービットの期待値を抽出します。これらの値が遅延グリーン関数行列を形成します。次に、RGFをフーリエ変換して動的構造因子を得て、RGFとDSFの両方をプロットします。get_dsfは内部で鏡映対称性を適用し、負の値を切り取ります。
rgf_mat = np.stack([pub_result.data.evs for pub_result in result])
# -- Compute DSF --
n_points_momentum, n_points_frequency = 100, 100
spectrum = get_dsf(
n_qubits,
rgf_mat,
time_step,
n_steps,
n_points_momentum,
n_points_frequency,
)
# -- Plot retarded Green's function --
plot_rgf(
n_qubits,
rgf_mat,
time_step,
n_steps,
title=f"Retarded Green's function — {n_qubits} qubits (simulation)",
)
# -- Plot DSF --
plot_dsf(
spectrum,
time_step,
n_points_momentum,
n_points_frequency,
title=f"Dynamical structure factor — {n_qubits} qubits (simulation)",
)
大規模ハードウェア実行
ここで50キュービットへとスケールアップします。この規模では、最適化は小規模の例よりも低い忠実度に達します。基底状態アンザッツの忠実度はおよそ0.65まで低下し、AQCの忠実度は後半のチェックポイントでおよそ0.7まで低下します。これは想定内であり、追加の古典コストを払って基底状態アンザッツ層の数(gs_n_layers)や最適化の反復回数(maxiter)を増やすことでこれらの忠実度を改善できます。また、2層のAQC忠実度が単層のものよりも低くなることにも注意してください。これは劣化ではありません。後半の時間ステップはよりエンタングルメントを生成するため単純に圧縮が難しく、そのため、より表現力の高い2層アンザッツをそれらのステップに使用しています。
AQCの最適化には数時間の古典計算時間がかかることもあります(ここで示す実行ではおよそ6時間で、その大部分はステップ2cの2層チェックポイントに費やされています)。実行時間を短縮するには、このノートブックを高性能計算(HPC)システムなど、より強力な古典ハードウェア上で実行することを検討してください。あるいは、より小さな問題インスタンス(例えばキュービット数や時間ステップ数を減らす)にスケールダウンすることもできますが、その場合、結果はここで示したものとは異なります。
以下のコードは、小規模の例と同じ4ステップの構造に従います。HVA基底状態パラメータは、再びMPSシミュレーションを使って最適化されます。QPU上では、エラー抑制・緩和のために動的デカップリング(DD)、パウリ・トワイリング、トワイリング読み出しエラー消去(TREX)を有効にします。次の表は、大規模実験が小規模実験とどのように異なるかをまとめたものです。
| 小規模 | 大規模 | |
|---|---|---|
| キュービット数 | 10 | 50 |
| 時間ステップ数 | 10 | 20 |
| AQCチェックポイント(1層+2層) | 3 + 2 = 5 | 6 + 4 = 10 |
| 基底状態アンザッツ層数 | 3 | 5 |
| MPS最大ボンド次元 | 32 | 128 |
| Estimator | StatevectorEstimator | DD、パウリ・トワイリング、TREXを使用するQPU |
XLA slow-compile messages
AQC最適化中(以下のステップ2c)に、次のようなstderrメッセージが表示されることがあります。
[Compiling module jit_MakeArrayFn for CPU] Very slow compile? ...
The operation took 2m14s
これらは、qiskit-addon-aqc-tensorのJAX自動微分を支えるコンパイラであるXLAからの無害な診断メッセージです。50キュービットでMPSボンド次元が128の場合、勾配関数が初めてトレースされる際にXLAがコンパイルにおよそ2、3分を費やします。コンパイルは成功し、最適化の結果には影響しません。
# ── Parameters ──────────────────────────────────────────────────────────────
n_qubits = 50 # 10 → 50
interaction = 1.0 # J
anisotropy = 1.0 # ε (isotropic point)
time_step = 0.6
n_steps = 20 # 10 → 20
aqc_n_steps_1 = 6 # 3 → 6
aqc_n_steps_2 = 4 # 2 → 4
aqc_n_steps_total = aqc_n_steps_1 + aqc_n_steps_2
gs_n_layers = 5 # 3 → 5
mps_max_bond = 128 # 32 -> 128
mps_cutoff = 1e-8
center = n_qubits // 2 - 1
# ── Step 1: Map ──────────────────────────────────────────────────────────────
ham_mpo = xxz_hamiltonian_mpo(n_qubits, interaction, anisotropy)
dmrg = qtn.DMRG2(ham_mpo)
dmrg.solve(tol=1e-8)
print(f"Ground-state energy (DMRG): {dmrg.energy:.6f}")
gs_ansatz = build_ground_state_ansatz(n_qubits, gs_n_layers)
rng = np.random.default_rng(12345)
x0 = np.tile([0.0, np.pi / 2], gs_n_layers) + rng.normal(
scale=0.1, size=2 * gs_n_layers
)
print("Optimizing ground state ansatz...")
t0 = timeit.default_timer()
result = optimize_ground_state_ansatz(
gs_ansatz,
x0,
dmrg.state,
max_bond=mps_max_bond,
cutoff=mps_cutoff,
options=dict(maxiter=100),
)
t1 = timeit.default_timer()
print(f"Finished optimizing ground state ansatz in {t1 - t0} seconds.")
print(f"Ground state ansatz fidelity: {1 - result.fun:.6f}")
gs_circuit = gs_ansatz.assign_parameters(result.x)
gs_circuit_mps = quimb_circuit(
gs_circuit.decompose(["PauliEvolution"]),
quimb_circuit_class=qtn.CircuitMPS,
max_bond=mps_max_bond,
cutoff=mps_cutoff,
)
gs_ansatz_energy = qtn.expec_TN_1D(
gs_circuit_mps.psi.H, ham_mpo, gs_circuit_mps.psi
)
print(f"Ground state ansatz energy: {gs_ansatz_energy:.6f}")
perturbed = gs_circuit.copy()
perturbed.rz(np.pi / 2, center)
circuits = []
for t in range(1, n_steps + 1):
circuit = perturbed.copy()
for instr in trotter_evolution(
circuit.qubits, interaction, anisotropy, time_step, t
):
circuit.append(instr)
circuits.append(circuit)
print(
f"Built {len(circuits)} circuits, deepest 2q depth (uniform basis) = "
f"{uniform_2q_depth(circuits[-1])}"
)
observables = [
SparsePauliOp.from_sparse_list([("Z", [i], 1)], num_qubits=n_qubits)
for i in range(n_qubits)
]
print(f"Defined {len(observables)} Z observables.")
# ── Step 2: AQC ──────────────────────────────────────────────────────────────
target_circuits = {
k: circuits[k - 1] for k in range(1, aqc_n_steps_total + 1)
}
aqc_sim = QuimbSimulator(
quimb_circuit_factory=partial(
qtn.CircuitMPS,
gate_opts=dict(max_bond=mps_max_bond, cutoff=mps_cutoff),
),
autodiff_backend="jax",
)
print("\nStep 2b — target MPS:")
target_mps = {}
for k in range(1, aqc_n_steps_total + 1):
target_mps[k] = tensornetwork_from_circuit(
target_circuits[k].decompose(["PauliEvolution"]), aqc_sim
)
print(f" k={k}: max bond = {target_mps[k].psi.max_bond()}")
ansatz_1, initial_params_1 = generate_ansatz_from_circuit(
target_circuits[1].decompose(["PauliEvolution"]),
qubits_initially_zero=True,
)
initial_params_1 = np.array(initial_params_1)
ansatz_2, initial_params_2 = generate_ansatz_from_circuit(
target_circuits[2].decompose(["PauliEvolution"]),
qubits_initially_zero=True,
)
initial_params_2 = np.array(initial_params_2)
print(
f"\nStep 2c — one-layer ansatz: {ansatz_1.num_parameters} parameters, "
f"2q depth (uniform basis) = {uniform_2q_depth(ansatz_1)}"
)
print(
f"Step 2c — two-layer ansatz: {ansatz_2.num_parameters} parameters, "
f"2q depth (uniform basis) = {uniform_2q_depth(ansatz_2)}"
)
aqc_circuits = {}
aqc_params = {}
for k in range(1, aqc_n_steps_total + 1):
if k <= aqc_n_steps_1:
ansatz, base_params = ansatz_1, initial_params_1
else:
ansatz, base_params = ansatz_2, initial_params_2
# Warm-start from the previous step only within the same stage
same_stage = (k - 1 >= 1) and (
(k - 1 <= aqc_n_steps_1) == (k <= aqc_n_steps_1)
)
x0 = aqc_params[k - 1] if same_stage else base_params
obj = MaximizeStateFidelity(target_mps[k], ansatz, aqc_sim)
t0 = timeit.default_timer()
result = scipy.optimize.minimize(
obj.loss_function,
x0,
method="L-BFGS-B",
jac=True,
options=dict(maxiter=100),
)
elapsed = timeit.default_timer() - t0
aqc_params[k] = result.x
aqc_circuits[k] = ansatz.assign_parameters(result.x)
print(
f" k={k}: fidelity = {1 - result.fun:.4f}, "
f"2q depth (uniform basis) = "
f"{uniform_2q_depth(aqc_circuits[k])}, "
f"{elapsed:.1f}s"
)
all_circuits = []
for k in range(1, aqc_n_steps_total + 1):
all_circuits.append(aqc_circuits[k])
base = aqc_circuits[aqc_n_steps_total]
for k in range(1, n_steps - aqc_n_steps_total + 1):
circuit = base.copy()
for instr in trotter_evolution(
circuit.qubits, interaction, anisotropy, time_step, k
):
circuit.append(instr)
all_circuits.append(circuit)
full_depths = [
uniform_2q_depth(circuits[k - 1]) for k in range(1, n_steps + 1)
]
aqc_depths = [uniform_2q_depth(circuit) for circuit in all_circuits]
full_2q = [
circuits[k - 1].decompose(["PauliEvolution"]).num_nonlocal_gates()
for k in range(1, n_steps + 1)
]
aqc_2q = [
circuit.decompose(["PauliEvolution"]).num_nonlocal_gates()
for circuit in all_circuits
]
print(f"\nStep 2d — assembled {len(all_circuits)} circuits")
print(
f" At step {n_steps} (uniform basis): "
f"Trotter 2q depth = {full_depths[-1]}, AQC+Trotter 2q depth = {aqc_depths[-1]}; "
f"2q gates: Trotter = {full_2q[-1]}, AQC+Trotter = {aqc_2q[-1]}"
)
steps_axis = np.arange(1, n_steps + 1)
fig, ax = plt.subplots(figsize=(7, 4))
ax.plot(steps_axis, full_depths, "-o", color="black", label="Full Trotter")
ax.plot(
steps_axis, aqc_depths, "-o", color="cadetblue", label="AQC + Trotter"
)
ax.set_xlabel("Trotter steps", fontsize=13)
ax.xaxis.set_major_locator(plt.matplotlib.ticker.MaxNLocator(integer=True))
ax.set_ylabel("2q gate depth (uniform basis)", fontsize=13)
ax.set_title(f"AQC circuit-depth reduction ({n_qubits} qubits)", fontsize=13)
ax.legend(fontsize=11)
plt.tight_layout()
plt.show()
# ── Step 3: Execute on IBM Quantum hardware ───────────────────────────────────
# (replaces StatevectorEstimator)
service = QiskitRuntimeService()
backend = service.least_busy(
min_num_qubits=n_qubits,
operational=True,
simulator=False,
filters=lambda x: x.configuration().processor_type["family"] == "Heron",
)
print(f"Backend: {backend.name}")
pm = generate_preset_pass_manager(optimization_level=3, backend=backend)
isa_circuits = pm.run(all_circuits, num_processes=1)
isa_2q_depths = [
isa_circuit.depth(lambda inst: inst.operation.num_qubits == 2)
for isa_circuit in isa_circuits
]
print(
f"Transpiled 2q depth (deepest, ISA on {backend.name}): "
f"{max(isa_2q_depths)} "
f"(2q gates: {max(isa_circuit.num_nonlocal_gates() for isa_circuit in isa_circuits)})"
)
estimator = Estimator(backend)
estimator.options.environment.job_tags = ["TUT_SNS"]
estimator.options.dynamical_decoupling.enable = True
estimator.options.dynamical_decoupling.sequence_type = "XY4"
estimator.options.twirling.enable_gates = True
estimator.options.twirling.num_randomizations = 1000
estimator.options.twirling.shots_per_randomization = 128
estimator.options.resilience.measure_mitigation = True
estimator.options.resilience.measure_noise_learning.num_randomizations = 32
estimator.options.resilience.measure_noise_learning.shots_per_randomization = 100
pubs = [
(
isa_circuit,
[obs.apply_layout(isa_circuit.layout) for obs in observables],
)
for isa_circuit in isa_circuits
]
job = estimator.run(pubs)
print(f"Job ID: {job.job_id()}")
result = job.result()
# ── Step 4: Post-process ──────────────────────────────────────────────────────
rgf_mat = np.stack([pub_result.data.evs for pub_result in result])
n_points_momentum, n_points_frequency = 100, 100
spectrum = get_dsf(
n_qubits,
rgf_mat,
time_step,
n_steps,
n_points_momentum,
n_points_frequency,
)
plot_rgf(
n_qubits,
rgf_mat,
time_step,
n_steps,
title=f"Retarded Green's function — {n_qubits} qubits (QPU)",
)
plot_dsf(
spectrum,
time_step,
n_points_momentum,
n_points_frequency,
title=rf"KCuF$_3$ DSF — {n_qubits} qubits (QPU)"
"\n(AQC + DD + Pauli twirling + TREX)",
)
Ground-state energy (DMRG): -21.972109
Optimizing ground state ansatz...
Finished optimizing ground state ansatz in 132.2076231740648 seconds.
Ground state ansatz fidelity: 0.645956
Ground state ansatz energy: -21.616744
Built 20 circuits, deepest 2q depth (uniform basis) = 307
Defined 50 Z observables.
Step 2b — target MPS:
k=1: max bond = 44
k=2: max bond = 46
k=3: max bond = 53
k=4: max bond = 62
k=5: max bond = 75
k=6: max bond = 96
k=7: max bond = 118
k=8: max bond = 128
k=9: max bond = 128
k=10: max bond = 128
Step 2c — one-layer ansatz: 2971 parameters, 2q depth (uniform basis) = 39
Step 2c — two-layer ansatz: 3412 parameters, 2q depth (uniform basis) = 45
k=1: fidelity = 1.0000, 2q depth (uniform basis) = 39, 215.6s
k=2: fidelity = 0.9676, 2q depth (uniform basis) = 39, 269.8s
k=3: fidelity = 0.9014, 2q depth (uniform basis) = 39, 293.5s
k=4: fidelity = 0.8755, 2q depth (uniform basis) = 39, 279.6s
k=5: fidelity = 0.8450, 2q depth (uniform basis) = 39, 266.4s
k=6: fidelity = 0.8146, 2q depth (uniform basis) = 39, 371.0s
E0626 02:59:07.229070 1508496 slow_operation_alarm.cc:73]
********************************
[Compiling module jit_MakeArrayFn for CPU] Very slow compile? If you want to file a bug, run with envvar XLA_FLAGS=--xla_dump_to=/tmp/foo and attach the results.
********************************
E0626 02:59:21.713796 1508477 slow_operation_alarm.cc:140] The operation took 2m14.484932074s
********************************
[Compiling module jit_MakeArrayFn for CPU] Very slow compile? If you want to file a bug, run with envvar XLA_FLAGS=--xla_dump_to=/tmp/foo and attach the results.
********************************
k=7: fidelity = 0.7312, 2q depth (uniform basis) = 45, 10865.3s
k=8: fidelity = 0.7720, 2q depth (uniform basis) = 45, 1450.3s
k=9: fidelity = 0.7316, 2q depth (uniform basis) = 45, 3388.6s
E0626 07:21:05.148724 1508496 slow_operation_alarm.cc:73]
********************************
[Compiling module jit_MakeArrayFn for CPU] Very slow compile? If you want to file a bug, run with envvar XLA_FLAGS=--xla_dump_to=/tmp/foo and attach the results.
********************************
E0626 07:21:24.286095 1508477 slow_operation_alarm.cc:140] The operation took 2m19.137566498s
********************************
[Compiling module jit_MakeArrayFn for CPU] Very slow compile? If you want to file a bug, run with envvar XLA_FLAGS=--xla_dump_to=/tmp/foo and attach the results.
********************************
k=10: fidelity = 0.6726, 2q depth (uniform basis) = 45, 3396.1s
Step 2d — assembled 20 circuits
At step 20 (uniform basis): Trotter 2q depth = 307, AQC+Trotter 2q depth = 171; 2q gates: Trotter = 3775, AQC+Trotter = 1913
Backend: ibm_fez
Transpiled 2q depth (deepest, ISA on ibm_fez): 105 (2q gates: 2574)
Job ID: d8v39vhropqc738biotg
ハードウェアの結果は2スピノン連続体の主要な特徴を再現しています。散乱強度は低エネルギーの反強磁性波数ベクトル付近に集中しており、下側は正弦波状のスピノン分散によって境界付けられ、単一の鋭いモードではなくその上に広いスペクトル重みの連続体が広がっています。これは、KCuF上の非弾性中性子散乱によって測定されたものと同じ構造であり、50キュービット規模における基底状態準備、摂動、AQC圧縮されたトロッター発展、およびエラー緩和された測定の全ワークフローを検証しています。
次のステップ
この内容に興味を持った場合は、以下の資料も参考になるでしょう。
Recommendations
-
AQC+トロッター・ダイナミクスのServerlessワークフローで中性子散乱をシミュレートする — デプロイされたQiskit Serverless関数テンプレートを通じてこの同じ実験を実行する関連チュートリアル
-
Lee他、「中性子散乱実験による量子シミュレーションのベンチマーク」(arXiv:2603.15608) — このチュートリアルの基になった参考論文
-
エラー緩和・抑制手法 — ハードウェア実験で使用されるDD、パウリ・トワイリング、TREX
-
時間発展回路のための近似量子コンパイル — AQC-Tensorに関するチュートリアル