{ "cells": [ { "cell_type": "markdown", "id": "title-cell", "metadata": {}, "source": [ "---\n", "title: \"ノイズの多い量子プロセッサにおける、堅牢かつコヒーレントな非アベル型ハドロンダイナミクスの観測\"\n", "description: \"IBM 量子ハードウェア上で、LSHフレームワークを用いてSU(2)格子ゲージ理論のハドロンダイナミクスをシミュレーションする。\"\n", "---\n", "\n", "{/* cspell:ignore Kogut Susskind expvals Pstep Nstep vmax vmin imshow fontsize cbar Ilčić */}\n", "\n", "\n", "\n", "# ノイズの多い量子プロセッサにおける、堅牢かつコヒーレントな非アベル型ハドロンダイナミクスの観測\n", "\n", "*推定実行時間：Heronプロセッサ（ibm\\_boston または同等のもの）で 6 分（注：これはあくまで推定値です。 （実行時間は状況によって異なる場合があります。）*\n", "\n" ] }, { "cell_type": "markdown", "id": "learning-outcomes", "metadata": {}, "source": [ "\n", "\n", "## 学習成果\n", "\n", "このチュートリアルを完了すると、以下のことを習得できます：\n", "\n", "* 非アーベル格子ゲージ理論（具体的にはSU(2)）を、効率的な量子シミュレーションのためにループ・ストリング・ハドロン（LSH）フレームワークを用いてどのように再定式化できるか\n", "* 近似SU(2)ゲージ理論のハミルトニアンに対するトロッター化時間発展回路の構築方法と、それらを量子ビットに写像する方法\n", "* IBM Quantum® ハードウェア上で、読み出し誤差の低減機能を備えたQiskit Estimatorプリミティブを使用して、これらの回路を実行する方法\n", "\n", "\n", "\n", "## 前提条件\n", "\n", "以下のトピックについて、あらかじめ理解しておくことをお勧めします：\n", "\n", "* [量子回路と量子ゲートの基礎](/learning/courses/basics-of-quantum-information)\n", "* [Qiskit Estimator プリミティブの概要](/docs/guides/get-started-with-estimator)\n", "* 量子場理論の概念について基本的な知識があること（あれば望ましいが必須ではない。背景のセクションで要点を解説している）\n", "\n" ] }, { "cell_type": "markdown", "id": "background", "metadata": {}, "source": [ "\n", "\n", "## 背景\n", "\n", "\n", "\n", "### モチベーション\n", "\n", "量子色力学（QCD）は、強い力を記述するSU(3)ゲージ理論であり、クォークをハドロンに束縛し、閉じ込めや弦の破れを支配している。 古典格子QCDの手法は静的性質の解析には優れているが、符号問題のため、リアルタイムのダイナミクスをシミュレーションすることはできない。 量子コンピュータは、ゲージ場の自由度をクビットに直接エンコードすることで、この障壁を乗り越える道筋を提供する。\n", "\n", "このチュートリアルでは、そのようなシミュレーションの手法を紹介します。 IBM Quantum ハードウェアを用いて、（1+1）次元のSU(2)格子ゲージ理論におけるハドロン伝播をリアルタイムでシミュレートします。これは、最も単純な非アーベルゲージ理論であり、完全なQCDへの足がかりとなるものです。\n", "\n", "\n", "\n", "### コグート・サスキンド・ハミルトニアン\n", "\n", "この理論は、 1D の空間格子上で、サイト上にスタッガード配置のフェルミオン（物質）が、リンク上にSU(2)ゲージ場が存在するという設定に基づいて定式化されている。 無次元形に換算すると、ハミルトニアンは次のようになる：\n", "\n", "$W = H_E^{\\text{(KS)}} + \\mu H_M + x H_I^{\\text{(KS)}},$\n", "\n", "ここで、 $H_E$ はクロモ電気場エネルギー、 $H_M$ はスタッガード質量項、 $H_I$ は物質・ゲージ相互作用（ホッピング）項、 $\\mu = 2\\frac{m}{g}\\sqrt{x}$ はフェルミオン質量を表し、 $x = \\frac{1}{g^2 a^2}$ は相互作用の強さを表す。 この理論の連続極限は、 $N \\to \\infty$ および $x \\to \\infty$ にある。\n", "\n", "\n", "\n", "### ループ・ストリング・ハドロン（LSH）フレームワーク\n", "\n", "重要な課題の一つは、各リンク上のゲージ場のヒルベルト空間が無限次元であるという点である。 **ループ・ストリング・ハドロン（LSH）** フレームワークは、この問題に対処するため、理論をゲージ不変の変数――フラックスのループ、分離した電荷を結ぶストリング、およびハドロン（サイトにおけるゲージシングレットフェルミオン対）――を用いて再定式化している。 LSH基底では、ガウスの法則は構成上自動的に満たされるため、すべての基底状態は物理的な状態である。 各格子点は、ループ数、入射ストリング、および出射ストリングを表す3つの量子数 $(n_l, n_i, n_o)$ によって特徴づけられる。ここで、 $n_i, n_o \\in \\{0,1\\}$ はフェルミオン的であり、 $n_l \\geq 0$ はボソン的である。 これらから、局所フェルミオン数は、サイト数が偶数の場合は $n_f(r) = n_i(r) + n_o(r)$、奇数の場合は $n_f(r) = 2 - [n_i(r) + n_o(r)]$ と定義される。\n", "\n", "\n", "\n", "### 完全ハミルトニアンから量子回路へ：3つの重要な近似\n", "\n", "この量子回路は、SU(2)ハミルトニアン全体を厳密にシミュレートする**ものではない**。 その代わりに、この手法では、 **弱結合領域** （ $x \\gg 1$ ）において有効な、制御された一連の近似を実装しています。何が近似され、何が近似されていないかを理解することが不可欠です：\n", "\n", "**近似 1 — $H_I$ における弱結合極限：** 完全相互作用ハミルトニアン $H_I^{\\text{(LSH)}}$ （式 [\\[1\\]](#references) の式(16)には、 $1/\\sqrt{n_l+1}$ のような項を通じて、ボソン量子数 $n_l$ に依存する係数が含まれている。弱結合領域（ $x \\gg 1$ ）では、ダイナミクスは電気項 $H_E$ によって支配され、これにより $n_l$ が大きい状態が優先される。 $n_l \\gg 1$ の場合、比 $n_l/(n_l+1) \\to 1$ となり、これらすべての係数は 1 に簡約される。 これにより、相互作用ハミルトニアンは、純粋に局所的な最近接ホッピングに還元される：\n", "\n", "$H_I^{\\text{approx}} = -\\sum_r \\left[\\sigma^-(r)\\sigma^+(r+1) + \\sigma^+(r)\\sigma^-(r+1)\\right],$\n", "\n", "これは $n_l$ とは独立しており、フェルミオン系の $(n_i, n_o)$ 量子ビットにのみ作用する。\n", "\n", "**近似 2 — $H_E$ における全球平均フラックス：** 電気エネルギーは、各リンクにおける $n_l$ に依存する。 弱結合の真空状態では、 $n_l$ は大きく、ほぼ一様である。 サイトごとの $n_l$ の値を、単一のグローバル平均値 $\\bar{n}_l$ に置き換え、 $H_E$ を、各サイトにおけるフェルミオン配置に比例する対角位相とする：\n", "\n", "$H_E^{\\text{approx}} = N h_E^0 + \\sum_{\\{r'\\}} \\left(\\frac{\\bar{n}_l}{2} + \\frac{3}{4}\\right)$\n", "\n", "ここで、 $\\{r'\\}$ は、フェルミオン配置 $(n_i=0, n_o=1)$ におけるサイトについて和をとったものであり、 $h_E^0$ は無視できるグローバル位相である。\n", "\n", "**近似 3 — トロッター化：** 持続時間 $\\delta_\\tau$ の 1 ステップに対する時間発展演算子は、次のように分解される：\n", "\n", "$e^{-i\\delta_\\tau W} \\approx e^{-i\\tilde{m} H_M} \\, e^{-i\\delta_\\tau H_E^{\\text{approx}}} \\, e^{-ic H_I^{\\text{approx}}}$\n", "\n", "ここで、 $c = \\delta_\\tau x$、 $\\tilde{m} = \\delta_\\tau \\mu$、および $\\theta = -\\delta_\\tau(\\bar{n}_l/2 + 3/4)$ である。この1次トロッター分解では誤差が生じるが、 $\\delta_\\tau \\to 0$ となるにつれてこの誤差は消える。ここでは、 $\\delta_\\tau = 0.0015$ を全域で固定する。\n", "\n", "これら3つの近似**の結果**、サイトあたり2つのフェルミオン量子ビット $(n_i, n_o)$ のみが動的であり、ボソン量子ビット $n_l$ の自由度は有効パラメータに吸収されている。 これにより、 $N$ 個の格子点に対して $2N$ 個の量子ビットからなるコンパクトな回路が得られる。ここで、各トロッターステップの2量子ビットゲートの深さは一定（1ステップあたり13）である。\n", "\n", "\n", "\n", "### このチュートリアルでシミュレートする内容\n", "\n", "このチュートリアルでは、 **ハドロンの伝播を**シミュレートします。強結合真空（積状態）から始め、格子の中心にメゾンを配置し、時間を進めます。 差分測定プロトコル――中心のメゾンを含む場合と含まない場合で回路を動作させ、その差を算出する――により、ハードウェアノイズや境界効果からコヒーレントなハドロン信号を分離することができる。 その結果、閉じ込められたメソン呼吸モードに特徴的な、フェルミオン密度の振動による光錐パターンが得られる。\n", "\n" ] }, { "cell_type": "markdown", "id": "requirements", "metadata": {}, "source": [ "\n", "\n", "## 要件\n", "\n", "このチュートリアルを始める前に、以下のものをインストールしてください：\n", "\n", "* Qiskit SDK v2.0 またはそれ以降のバージョンで、 [可視化](/docs/api/qiskit/visualization)機能をサポートしたもの\n", "* Qiskit Runtime v0.22 またはそれ以降 (`pip install qiskit-ibm-runtime`)\n", "* Pauli Propagation パッケージ (`pip install pauli-prop`)\n", "* NumPy (`pip install numpy`)\n", "* Matplotlib (`pip install matplotlib`)\n", "\n" ] }, { "cell_type": "markdown", "id": "setup-header", "metadata": {}, "source": [ "\n", "\n", "## セットアップ\n", "\n", "まず、必要なライブラリをインポートし、LSH時間発展のための量子回路を構築するヘルパー関数を定義します。 回路構築には、主に3つの機能があります：\n", "\n", "1. **`pair_hamiltonian_circuit`**: 隣接するサイト間の近似相互作用ハミルトニアンに対して、2量子ビットのユニタリー演算 $U_I$ を実装する。 ゲートの分解は次のとおりです： $\\text{CNOT} \\to H \\to R_z(-c) \\to \\text{CNOT} \\to R_z(c) \\to \\text{CNOT} \\to H \\to \\text{CNOT}$。\n", "\n", "2. **`electric_hamiltonian_circuit`**: 各サイトにおける近似電界エネルギーに対して、2量子ビットのユニタリー演算 $U_E$ を実装する。 ゲートの分解は次のとおりです： $X \\to R_z(\\theta/2) \\to \\text{CNOT} \\to R_z(-\\theta/2) \\to \\text{CNOT} \\to R_z(\\theta/2) \\to X$。\n", "\n", "3. **`construct_circuit`**: インタラクション項、電気項、質量項をSWAPゲートを用いて重ね合わせ、量子ビット間の接続を制御することで、完全な「トロッター化」回路を構築する。\n", "\n" ] }, { "cell_type": "code", "execution_count": null, "id": "setup-imports", "metadata": {}, "outputs": [], "source": [ "# Import libraries\n", "\n", "import numpy as np\n", "import matplotlib.pyplot as plt\n", "from matplotlib.colors import TwoSlopeNorm\n", "from qiskit.circuit import QuantumCircuit\n", "from qiskit.quantum_info import SparsePauliOp\n", "from typing import Optional\n", "\n", "import warnings\n", "\n", "warnings.filterwarnings(\"ignore\")" ] }, { "cell_type": "code", "execution_count": null, "id": "setup-functions", "metadata": {}, "outputs": [], "source": [ "def pair_hamiltonian_circuit(c: float) -> QuantumCircuit:\n", " \"\"\"Two-qubit unitary for the approximate interaction Hamiltonian H_I.\n", "\n", " Implements exp(-i * c * H_I^approx) for one pair of neighboring sites,\n", " where c = delta_tau * x.\n", " \"\"\"\n", " qc_temp = QuantumCircuit(2)\n", " qc_temp.cx(1, 0)\n", " qc_temp.h(1)\n", " qc_temp.rz(-c, 1)\n", " qc_temp.cx(0, 1)\n", " qc_temp.rz(c, 1)\n", " qc_temp.cx(0, 1)\n", " qc_temp.h(1)\n", " qc_temp.cx(1, 0)\n", " return qc_temp\n", "\n", "\n", "def electric_hamiltonian_circuit(theta: float) -> QuantumCircuit:\n", " \"\"\"Two-qubit unitary for the approximate electric field Hamiltonian H_E.\n", "\n", " Implements exp(-i * theta * H_E^approx) for one lattice site,\n", " where theta = -delta_tau * (n_bar_l / 2 + 3/4).\n", " \"\"\"\n", " qc_temp = QuantumCircuit(2)\n", " qc_temp.x(0)\n", " qc_temp.rz(theta / 2, 0)\n", " qc_temp.cx(0, 1)\n", " qc_temp.rz(-theta / 2, 1)\n", " qc_temp.cx(0, 1)\n", " qc_temp.rz(theta / 2, 1)\n", " qc_temp.x(0)\n", " return qc_temp\n", "\n", "\n", "def construct_circuit(\n", " num_lattice_point: int,\n", " num_trotter_steps: int,\n", " c: float,\n", " theta: float,\n", " m: float,\n", " theory: Optional[int] = 2,\n", " barriers: Optional[bool] = False,\n", " measurement: Optional[bool] = False,\n", " add_init_state: Optional[bool] = True,\n", " inverse_mid: Optional[bool] = False,\n", ") -> QuantumCircuit:\n", " \"\"\"Construct the full Trotterized time-evolution circuit.\n", "\n", " Builds a circuit implementing n Trotter steps of the approximate SU(2)\n", " LSH Hamiltonian evolution. The qubit layout uses a zigzag ordering:\n", " n_i(0), n_i(1), n_o(0), n_o(1), n_i(2), n_i(3), n_o(2), n_o(3), ...\n", " which minimizes the number of SWAP layers needed.\n", "\n", " Args:\n", " num_lattice_point: Number of lattice sites\n", " (num_qubits = 2 * num_lattice_point).\n", " num_trotter_steps: Number of Trotter steps.\n", " c: Interaction parameter (delta_tau * x).\n", " theta: Electric field phase parameter.\n", " m: Mass parameter (m_tilde = delta_tau * mu).\n", " theory: 1 for single chain, 2 for SU(2). Default 2.\n", " barriers: Insert barriers between Trotter layers for\n", " visualization.\n", " measurement: Append measurements at the end.\n", " add_init_state: Prepare the half-filled (strong-coupling vacuum)\n", " initial state.\n", " inverse_mid: Swap the central sites\n", " (for differential measurement protocol).\n", " \"\"\"\n", " num_qubits = theory * num_lattice_point\n", " qc = QuantumCircuit(num_qubits)\n", "\n", " if num_trotter_steps <= 0:\n", " return qc\n", "\n", " # --- Initial state preparation ---\n", " if add_init_state:\n", " i = 1\n", " while i < num_lattice_point:\n", " for j in range(theory):\n", " qc.x(i + j * num_lattice_point)\n", " i = i + 2\n", " if inverse_mid:\n", " mid_lattice_qubits = [num_qubits // 2 - 1, num_qubits // 2]\n", " qc.x(mid_lattice_qubits)\n", " else:\n", " i = 1\n", " while i < num_qubits - 1:\n", " qc.swap(i, i + 1)\n", " i = i + 4\n", "\n", " # --- Trotter steps ---\n", " for step in range(num_trotter_steps):\n", " if barriers:\n", " qc.barrier()\n", "\n", " # First SWAP layer (skipped at step 0 — absorbed into initial state mapping)\n", " if step > 0:\n", " i = 1\n", " while i < num_qubits - 1:\n", " qc.swap(i, i + 1)\n", " i = i + 4\n", "\n", " # First layer of pair interactions\n", " j = 0\n", " while j < num_qubits - 2:\n", " circ = pair_hamiltonian_circuit(c)\n", " qc.compose(circ, [j, j + 1], inplace=True)\n", " j = j + 2\n", " if num_lattice_point % 2 == 0:\n", " circ = pair_hamiltonian_circuit(c)\n", " qc.compose(circ, [j, j + 1], inplace=True)\n", "\n", " # Second SWAP layer\n", " i = 1\n", " while i < num_qubits - 1:\n", " qc.swap(i, i + 1)\n", " i = i + theory\n", "\n", " # Second layer of pair interactions\n", " j = 2\n", " while j < num_qubits - 3:\n", " circ = pair_hamiltonian_circuit(c)\n", " qc.compose(circ, [j, j + 1], inplace=True)\n", " j = j + 2\n", " if num_lattice_point % 2 != 0:\n", " circ = pair_hamiltonian_circuit(c)\n", " qc.compose(circ, [j, j + 1], inplace=True)\n", "\n", " # Third SWAP layer\n", " i = 3\n", " while i < num_qubits - 1:\n", " qc.swap(i, i + 1)\n", " i = i + 2 * theory\n", "\n", " # Electric field term\n", " if theta != 0:\n", " e_circ = electric_hamiltonian_circuit(theta)\n", " for j in range(num_lattice_point):\n", " qc.compose(e_circ, [2 * j, 2 * j + 1], inplace=True)\n", "\n", " # Mass term: Rz(-m_tilde) for even sites, Rz(m_tilde) for odd sites\n", " for q in range(num_qubits):\n", " if q % 2 == 0:\n", " qc.rz(-1 * m, q)\n", " else:\n", " qc.rz(m, q)\n", "\n", " if measurement:\n", " qc.measure_all()\n", "\n", " return qc" ] }, { "cell_type": "code", "execution_count": 3, "id": "setup-postprocess", "metadata": {}, "outputs": [], "source": [ "def get_probabilities(expval: float):\n", " \"\"\"Convert a Z-expectation value to site occupation probability.\n", "\n", " Since = p(0) - p(1), the occupation probability is p(1) = (1 - ) / 2.\n", " \"\"\"\n", " p1 = round((1 - expval) / 2, 3)\n", " return p1\n", "\n", "\n", "def get_number(expval_data, num_lattice_point):\n", " \"\"\"Convert raw Z-expectation values to staggered fermion number n_f at each site.\n", "\n", " n_f(r) = n_i(r) + n_o(r) for even r\n", " n_f(r) = 2 - [n_i(r) + n_o(r)] for odd r\n", "\n", " The two qubits per site encode (n_i, n_o), and occupation probabilities\n", " give us and .\n", " \"\"\"\n", " N = []\n", " for expvals in expval_data:\n", " Pstep = [get_probabilities(expval) for expval in expvals]\n", " Nstep = []\n", " for k in range(num_lattice_point):\n", " val = Pstep[2 * k] + Pstep[2 * k + 1]\n", " a = 2 * (k % 2) + (1 - 2 * (k % 2)) * val\n", " Nstep.append(float(a))\n", " N.append(Nstep)\n", " return N\n", "\n", "\n", "def calculate_difference(N, N_mid, num_lattice_point):\n", " \"\"\"Differential measurement protocol: |n_f(meson) - n_f(vacuum)|.\n", "\n", " Subtracting the vacuum (SCV) evolution from the meson evolution\n", " isolates the coherent hadron signal from symmetric noise and boundary effects.\n", " \"\"\"\n", " N_diff = []\n", " for i in range(len(N)):\n", " Nstep_diff = []\n", " for j in range(num_lattice_point):\n", " Nstep_diff.append(abs(N[i][j] - N_mid[i][j]))\n", " N_diff.append(Nstep_diff)\n", " return N_diff" ] }, { "cell_type": "markdown", "id": "sim-header", "metadata": {}, "source": [ "\n", "\n", "## 小規模シミュレータの例\n", "\n", "まず、6サイト格子（12キュービット）を用いて小規模なワークフローを実証し、ハードウェア上で実行する前に、回路の構築を確認し、物理的観測量を理解できるようにします。\n", "\n" ] }, { "cell_type": "markdown", "id": "step1-header", "metadata": {}, "source": [ "\n", "\n", "### ステップ 1：古典的な入力を量子問題に写像する\n", "\n", "本論文（ $x = 100$、 $m/g = 1$ ）で検討された弱結合領域に一致する物理パラメータを定義する。導出された回路パラメータは以下の通りである：\n", "\n", "* $c = \\delta_\\tau \\cdot x = 0.15$ (相互作用パラメータ)\n", "* $\\theta = -\\delta_\\tau (\\bar{n}_l/2 + 3/4) = 0.01$ （電界の位相）\n", "* $\\tilde{m} = \\delta_\\tau \\cdot \\mu = 0.03$ （質量パラメータ）\n", "\n", "トロッターのステップ数ごとに、 **2つの回路を**構築する。1つは中心でメゾンを初期化する回路（`inverse_mid=True`）であり、もう1つは強結合真空を準備する回路（`inverse_mid=False`）である。微分測定プロトコルでは、真空の進化を差し引くことで、ハドロン信号を分離する。\n", "\n" ] }, { "cell_type": "code", "execution_count": 4, "id": "step1-params", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Lattice sites: 6, Qubits: 12\n", "Parameters: c=0.15, theta=0.01, m_tilde=0.03\n" ] } ], "source": [ "# Physical / circuit parameters\n", "num_lattice_point = 6 # 6 lattice sites -> 12 qubits for SU(2)\n", "num_qubits = 2 * num_lattice_point\n", "c = 0.15 # delta_tau * x\n", "theta = 0.01 # electric field phase\n", "m = 0.03 # m_tilde = delta_tau * mu\n", "trotter_steps = range(1, 11) # 10 Trotter steps\n", "\n", "print(f\"Lattice sites: {num_lattice_point}, Qubits: {num_qubits}\")\n", "print(f\"Parameters: c={c}, theta={theta}, m_tilde={m}\")" ] }, { "cell_type": "code", "execution_count": 5, "id": "step1-circuits", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Circuit for 1 Trotter step: 12 qubits, depth 26\n" ] }, { "data": { "text/plain": [ "\"Output" ] }, "execution_count": 5, "metadata": {}, "output_type": "execute_result" } ], "source": [ "# Build circuits: meson initial state and vacuum (SCV) initial state\n", "circuits_mid = [\n", " construct_circuit(\n", " num_lattice_point,\n", " d,\n", " c,\n", " theta,\n", " m,\n", " barriers=False,\n", " measurement=False,\n", " add_init_state=True,\n", " inverse_mid=True,\n", " )\n", " for d in trotter_steps\n", "]\n", "\n", "circuits = [\n", " construct_circuit(\n", " num_lattice_point,\n", " d,\n", " c,\n", " theta,\n", " m,\n", " barriers=False,\n", " measurement=False,\n", " add_init_state=True,\n", " inverse_mid=False,\n", " )\n", " for d in trotter_steps\n", "]\n", "\n", "# Visualize a single Trotter step\n", "print(\n", " f\"Circuit for 1 Trotter step: {circuits[0].num_qubits} qubits, depth {circuits[0].depth()}\"\n", ")\n", "circuits[0].draw(\"mpl\", fold=-1)" ] }, { "cell_type": "markdown", "id": "step2-header", "metadata": {}, "source": [ "\n", "\n", "### ステップ 2：量子ハードウェアでの実行に向けて問題を最適化する\n", "\n", "観測量として、各量子ビットに対する単一量子ビットの $Z$ 測定を定義する。 $\\langle Z \\rangle$ から、占有確率を抽出し、各格子点 $r$ におけるスタッガードフェルミオン数 $n_f(r)$ を算出することができます。\n", "\n" ] }, { "cell_type": "code", "execution_count": 6, "id": "step2-observables", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Number of observables: 12\n" ] } ], "source": [ "# Z observable on each qubit\n", "observables = [\n", " SparsePauliOp(\"I\" * i + \"Z\" + \"I\" * (num_qubits - i - 1))\n", " for i in range(num_qubits)\n", "]\n", "\n", "print(f\"Number of observables: {len(observables)}\")" ] }, { "cell_type": "markdown", "id": "step3-header", "metadata": {}, "source": [ "\n", "\n", "### ステップ 3: `Qiskit primitives` を使用して実行する\n", "\n", "小規模でのノイズのない正確なシミュレーションには、 を使用 `StatevectorEstimator` してください。\n", "\n" ] }, { "cell_type": "code", "execution_count": 7, "id": "step3-simulate", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "Computed expectation values for 10 Trotter steps\n" ] } ], "source": [ "from qiskit.primitives import StatevectorEstimator\n", "\n", "estimator = StatevectorEstimator()\n", "\n", "# Run meson circuits\n", "pubs_mid = [(circuit, observables) for circuit in circuits_mid]\n", "result_mid = estimator.run(pubs_mid).result()\n", "\n", "# Run vacuum (SCV) circuits\n", "pubs = [(circuit, observables) for circuit in circuits]\n", "result = estimator.run(pubs).result()\n", "\n", "# Extract expectation values\n", "raw_expvals_mid = [\n", " result_mid[i].data.evs[::-1] for i in range(len(circuits_mid))\n", "]\n", "raw_expvals = [result[i].data.evs[::-1] for i in range(len(circuits))]\n", "\n", "print(f\"Computed expectation values for {len(raw_expvals)} Trotter steps\")" ] }, { "cell_type": "markdown", "id": "step4-header", "metadata": {}, "source": [ "\n", "\n", "### ステップ4：後処理を行い、希望する従来の形式で結果を返す\n", "\n", "期待値をスタッガードフェルミオン数 $n_f(r, t)$ に変換し、微分測定プロトコル（メソン $-$ 真空）を適用して、ハドロン伝播ヒートマップを生成する。 これは、参考論文の図3の構成を再現したものです。x軸には格子点 $r$、y軸にはトロッターステップ（時間） $t$、そして色スケールとして $n_f(r,t)$ が用いられています。\n", "\n" ] }, { "cell_type": "code", "execution_count": 8, "id": "step4-postprocess", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "\"Output" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# Compute fermion numbers\n", "N_mid_sim = get_number(raw_expvals_mid, num_lattice_point)\n", "N_sim = get_number(raw_expvals, num_lattice_point)\n", "N_diff_sim = calculate_difference(N_mid_sim, N_sim, num_lattice_point)\n", "\n", "# --- Reproduce Figure 3 style: Staggered Fermionic Occupation Number Dynamics ---\n", "fig, axes = plt.subplots(1, 3, figsize=(18, 5))\n", "\n", "# Convert to numpy arrays for plotting\n", "N_mid_arr = np.array(N_mid_sim)\n", "N_arr = np.array(N_sim)\n", "N_diff_arr = np.array(N_diff_sim)\n", "\n", "# Color scheme\n", "vmax = max(max(sublist) for sublist in N_arr)\n", "vmin = -vmax\n", "\n", "# Meson evolution\n", "norm1 = TwoSlopeNorm(vmin=vmin, vcenter=0, vmax=vmax)\n", "im0 = axes[0].imshow(\n", " N_mid_arr,\n", " aspect=\"auto\",\n", " origin=\"lower\",\n", " cmap=\"RdBu_r\",\n", " norm=norm1,\n", " extent=[0, num_lattice_point, 0.5, len(trotter_steps) + 0.5],\n", ")\n", "axes[0].set_xlabel(\"Lattice site $r$\", fontsize=12)\n", "axes[0].set_ylabel(\"Trotter step $t$\", fontsize=12)\n", "axes[0].set_title(\"$n_f(r,t)$ — Meson initial state\", fontsize=12)\n", "plt.colorbar(im0, ax=axes[0], label=\"$n_f(r,t)$\")\n", "\n", "# Vacuum (SCV) evolution\n", "im1 = axes[1].imshow(\n", " N_arr,\n", " aspect=\"auto\",\n", " origin=\"lower\",\n", " cmap=\"RdBu_r\",\n", " norm=norm1,\n", " extent=[0, num_lattice_point, 0.5, len(trotter_steps) + 0.5],\n", ")\n", "axes[1].set_xlabel(\"Lattice site $r$\", fontsize=12)\n", "axes[1].set_ylabel(\"Trotter step $t$\", fontsize=12)\n", "axes[1].set_title(\"$n_f(r,t)$ — Vacuum (SCV)\", fontsize=12)\n", "plt.colorbar(im1, ax=axes[1], label=\"$n_f(r,t)$\")\n", "\n", "# Differential: meson - vacuum\n", "norm2 = TwoSlopeNorm(vmin=vmin, vcenter=0, vmax=vmax)\n", "im2 = axes[2].imshow(\n", " N_diff_arr,\n", " aspect=\"auto\",\n", " origin=\"lower\",\n", " cmap=\"RdBu_r\",\n", " norm=norm2,\n", " extent=[0, num_lattice_point, 0.5, len(trotter_steps) + 0.5],\n", ")\n", "axes[2].set_xlabel(\"Lattice site $r$\", fontsize=12)\n", "axes[2].set_ylabel(\"Trotter step $t$\", fontsize=12)\n", "axes[2].set_title(\n", " \"Staggered Fermionic Occupation Number Dynamics\\n$|n_f^{\\\\mathrm{meson}} - n_f^{\\\\mathrm{vacuum}}|$\",\n", " fontsize=12,\n", ")\n", "plt.colorbar(im2, ax=axes[2], label=\"$n_f(r,t)$\")\n", "\n", "plt.suptitle(\n", " f\"Hadron propagation — {num_lattice_point}-site lattice (StatevectorEstimator)\",\n", " fontsize=14,\n", " y=1.02,\n", ")\n", "plt.tight_layout()\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "hardware-header", "metadata": {}, "source": [ "\n", "\n", "## 大規模なハードウェアの例\n", "\n", "ここでは、 IBM Quantum ハードウェア上で、30サイトの格子（60キュービット）へとスケールアップします。 このスケールでは、10トロッターステップの回路は、3400個以上の2量子ビットゲートと14,000個の単一量子ビットゲートで構成されています。\n", "\n" ] }, { "cell_type": "markdown", "id": "hardware-steps", "metadata": {}, "source": [ "\n", "\n", "### 手順 1～4（1つのコードブロックにまとめました）\n", "\n", "ハードウェア・ワークフローの主なポイント：\n", "\n", "* メソン回路および真空回路用の10ステップのトロッター法（ドリフトを最小限に抑えるためインターリーブ処理済み）\n", "* — による `optimization_level=1` トランスパイレーション — 回路レイアウトはすでにデバイスのトポロジー（線形チェーン）と同型であるため、配線によるSWAPは不要です。 このトランスパイラーは、ノイズの少ない物理量子ビットの連鎖を選択し、ゲートをネイティブゲートセットに分解することのみを目的として使用されます。\n", "* `EstimatorV2` TREXの読み出しエラーの低減とパウリ・トゥワーリングを用いて\n", "* `Batch` すべてのジョブを一括で送信するセッション\n", "\n" ] }, { "cell_type": "code", "execution_count": null, "id": "hardware-code", "metadata": {}, "outputs": [], "source": [ "# -------------------------Step 1: Define parameters & build circuits-------------------------\n", "\n", "from qiskit_ibm_runtime import QiskitRuntimeService\n", "from qiskit.transpiler.preset_passmanagers import generate_preset_pass_manager\n", "from qiskit_ibm_runtime import EstimatorV2, Batch\n", "from qiskit_ibm_runtime.options import (\n", " EstimatorOptions,\n", " ResilienceOptionsV2,\n", " TwirlingOptions,\n", " DynamicalDecouplingOptions,\n", ")\n", "\n", "service = QiskitRuntimeService()\n", "\n", "num_lattice_point_hw = 30\n", "num_qubits_hw = 2 * num_lattice_point_hw # 60 qubits\n", "c_hw = 0.15\n", "theta_hw = 0.01\n", "m_hw = 0.03\n", "trotter_steps_hw = range(1, 11) # 10 Trotter steps\n", "\n", "# Build meson and vacuum circuits\n", "circuits_mid_hw = [\n", " construct_circuit(\n", " num_lattice_point_hw,\n", " d,\n", " c_hw,\n", " theta_hw,\n", " m_hw,\n", " barriers=False,\n", " measurement=False,\n", " add_init_state=True,\n", " inverse_mid=True,\n", " )\n", " for d in trotter_steps_hw\n", "]\n", "\n", "circuits_hw = [\n", " construct_circuit(\n", " num_lattice_point_hw,\n", " d,\n", " c_hw,\n", " theta_hw,\n", " m_hw,\n", " barriers=False,\n", " measurement=False,\n", " add_init_state=True,\n", " inverse_mid=False,\n", " )\n", " for d in trotter_steps_hw\n", "]\n", "\n", "print(f\"Built {len(circuits_hw)} circuit pairs for {num_qubits_hw} qubits\")\n", "\n", "# -------------------------Step 2: Transpile for hardware-------------------------\n", "# The circuit topology is a linear chain, isomorphic to the device topology.\n", "# We use optimization_level=1 since no routing SWAPs are needed — the transpiler\n", "# only needs to select a low-noise qubit chain and decompose to native gates.\n", "\n", "backend = service.backend(\"ibm_boston\")\n", "\n", "layout = [\n", " 140,\n", " 141,\n", " 142,\n", " 143,\n", " 136,\n", " 123,\n", " 122,\n", " 121,\n", " 116,\n", " 101,\n", " 102,\n", " 103,\n", " 96,\n", " 83,\n", " 82,\n", " 81,\n", " 76,\n", " 61,\n", " 62,\n", " 63,\n", " 64,\n", " 65,\n", " 66,\n", " 67,\n", " 68,\n", " 69,\n", " 78,\n", " 89,\n", " 88,\n", " 87,\n", " 97,\n", " 107,\n", " 106,\n", " 105,\n", " 117,\n", " 125,\n", " 126,\n", " 127,\n", " 137,\n", " 147,\n", " 148,\n", " 149,\n", " 150,\n", " 151,\n", " 152,\n", " 153,\n", " 154,\n", " 155,\n", " 139,\n", " 135,\n", " 134,\n", " 133,\n", " 132,\n", " 131,\n", " 130,\n", " 129,\n", " 118,\n", " 109,\n", " 110,\n", " 111,\n", "]\n", "\n", "\n", "pm = generate_preset_pass_manager(\n", " optimization_level=1, backend=backend, initial_layout=layout\n", ")\n", "\n", "isa_circuits_mid = pm.run(circuits_mid_hw)\n", "isa_circuits = pm.run(circuits_hw)\n", "\n", "print(f\"Transpiled circuits. Example depth: {isa_circuits[0].depth()}\")\n", "\n", "# Define and layout-map observables\n", "observables_hw = [\n", " SparsePauliOp(\"I\" * i + \"Z\" + \"I\" * (num_qubits_hw - i - 1))\n", " for i in range(num_qubits_hw)\n", "]\n", "\n", "isa_observables_mid = [\n", " [obs.apply_layout(isa_circuits_mid[i].layout) for obs in observables_hw]\n", " for i in range(len(isa_circuits_mid))\n", "]\n", "isa_observables = [\n", " [obs.apply_layout(isa_circuits[i].layout) for obs in observables_hw]\n", " for i in range(len(isa_circuits))\n", "]\n", "\n", "# Build PUBs — interleave meson and vacuum for each Trotter step\n", "isa_pubs_mid = [\n", " (circ, obs) for circ, obs in zip(isa_circuits_mid, isa_observables_mid)\n", "]\n", "isa_pubs = [(circ, obs) for circ, obs in zip(isa_circuits, isa_observables)]\n", "\n", "pubs_to_execute = [\n", " [isa_pubs_mid[i], isa_pubs[i]] for i in range(len(isa_pubs))\n", "]\n", "\n", "# -------------------------Step 3: Execute on hardware-------------------------\n", "\n", "twirling_options = TwirlingOptions(\n", " enable_gates=True,\n", " enable_measure=True,\n", " shots_per_randomization=\"auto\",\n", " strategy=\"active-circuit\",\n", ")\n", "\n", "resilience_options = ResilienceOptionsV2(\n", " measure_mitigation=True, # TREX readout error mitigation\n", " zne_mitigation=False, # ZNE turned off\n", ")\n", "\n", "dd_options = DynamicalDecouplingOptions(\n", " enable=False # Circuit is sufficiently dense\n", ")\n", "\n", "options = EstimatorOptions(\n", " resilience=resilience_options,\n", " twirling=twirling_options,\n", " dynamical_decoupling=dd_options,\n", " default_shots=10_000,\n", ")\n", "\n", "ids = []\n", "with Batch(backend=backend) as batch:\n", " for idx, pub in enumerate(pubs_to_execute):\n", " print(f\"Submitting job for Trotter step {idx + 1}\")\n", " estimator = EstimatorV2(mode=batch, options=options)\n", " estimator.skip_transpilation = True\n", " job = estimator.run(pub)\n", " ids.append(job.job_id())\n", " batch_id = batch.session_id\n", "\n", "job_info = {\"ids\": ids, \"batch_id\": batch_id}\n", "print(f\"Submitted {len(ids)} jobs. Batch ID: {batch_id}\")" ] }, { "cell_type": "code", "execution_count": null, "id": "f03fb6e6-2ba7-49bb-b1f5-eb9b6eda993b", "metadata": {}, "outputs": [], "source": [ "print(ids)" ] }, { "cell_type": "code", "execution_count": 16, "id": "72d09009-e0d2-4bb0-9157-2e88a9d973ea", "metadata": {}, "outputs": [], "source": [ "# -------------------------Step 4: Post-process results-------------------------\n", "\n", "jobs = [service.job(job_id) for job_id in ids]\n", "results = [job.result() for job in jobs]\n", "\n", "# Extract expectation values (index 0 = meson, index 1 = vacuum)\n", "raw_expvals_mid_hw = [result[0].data.evs[::-1] for result in results]\n", "raw_expvals_hw = [result[1].data.evs[::-1] for result in results]\n", "\n", "# Compute fermion numbers and differential\n", "N_mid_hw = get_number(raw_expvals_mid_hw, num_lattice_point_hw)\n", "N_hw = get_number(raw_expvals_hw, num_lattice_point_hw)\n", "N_diff_hw = calculate_difference(N_mid_hw, N_hw, num_lattice_point_hw)" ] }, { "cell_type": "code", "execution_count": null, "id": "19ee420d", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "\"Output" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "N_diff_hw_arr = np.array(N_diff_hw)\n", "\n", "fig, ax = plt.subplots(figsize=(10, 6))\n", "vmax = np.abs(N_hw).max()\n", "vmin = -vmax\n", "norm = TwoSlopeNorm(vmin=vmin, vcenter=0, vmax=vmax)\n", "im = ax.imshow(\n", " N_diff_hw_arr,\n", " aspect=\"auto\",\n", " origin=\"lower\",\n", " cmap=\"RdBu_r\",\n", " norm=norm,\n", " extent=[0, 8, 0.5, len(trotter_steps_hw) + 0.5],\n", ")\n", "ax.set_xlabel(\"Lattice site $r$\", fontsize=13)\n", "ax.set_ylabel(\"Trotter step $t$\", fontsize=13)\n", "ax.set_title(\n", " \"Staggered Fermionic Occupation Number Dynamics\\nQuantum Simulation on IBM Hardware — 30-site lattice (60 qubits)\",\n", " fontsize=13,\n", ")\n", "cbar = plt.colorbar(im, ax=ax)\n", "cbar.set_label(\"$n_f(r,t)$\", fontsize=12)\n", "plt.tight_layout()\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "de8e2aa6", "metadata": {}, "source": [ "\n", "\n", "## パウリ伝播による古典的なベンチマーク\n", "\n", "パウリ伝播法（PPM）は、ハイゼンベルク図式において、測定された観測量を回路全体に逆伝播させることで、量子回路のノイズのない古典的シミュレーションを実現する。 クリフォード層（CNOT、H、S、Xゲート）の下では、パウリ演算子は、項の数を増やすことなく、他のパウリ演算子に写像される。 非クリフォード層（回路内の $R_z$ ゲート）は分岐を引き起こす可能性があり、最悪の場合、項の数が2倍になることもありますが、多くの分岐は係数が小さいため、切り捨てることができます。\n", "\n", "の [`pauli-prop`](https://github.com/Qiskit/pauli-prop) ワークフローは以下の通りです：\n", "\n", "1. `evolve_through_cliffords`を使用して、回路をクリフォード部分と非クリフォード部分**に分割する**。\n", "2. `atol``propagate_through_circuit`各観測量を、を用いて非クリフォード部分を通じて**伝播させ**、最大 `max_terms` 個のパウリ項まで保持し、係数が切り捨て閾値 未満の項は除外する。\n", "3. Qiskitに組み込まれているクリフォード演算のサポート機能を使用して、クリフォード演算を通じて結果を**導出します**。\n", "4. 対角のパウリ項（ $I$ および $Z$ のみを含む）の係数を合計することで、期待値を**求めます**。\n", "\n", "\n", "\n", "### 切り捨て閾値\n", "\n", "の `propagate_through_circuit` パラメータは `atol` 、小さなパウリ分岐がどの程度積極的に剪定されるかを制御します。 `1e-12`しきい値を非常に厳しく設定すると（例えば、）、ほぼすべての分岐が保持され、正確な結果が得られますが、回路の深さに応じてシミュレーション時間が急激に増加します。 [本論文](https://arxiv.org/abs/2602.18080)で取り上げられた120キュービットのシミュレーションは、デフォルト設定では約 8.5 時間かかりました。 `1e-3`閾値を引き上げると（例えば、 または `1e-6` まで）、係数がその値を下回る項が除外され、追跡対象となる項の数が劇的に減少し、計算速度が向上します。 その代償として、わずかで制御可能な近似誤差が生じますが、これは異なる閾値での結果を比較することで検証することができます。\n", "\n" ] }, { "cell_type": "code", "execution_count": 24, "id": "0ed2dd40", "metadata": {}, "outputs": [ { "name": "stdout", "output_type": "stream", "text": [ "PPM settings: atol=0.001, max_terms=66000\n", "Trotter step 1: 5.0 s\n", "Trotter step 2: 7.5 s\n", "Trotter step 3: 11.2 s\n", "Trotter step 4: 14.7 s\n", "Trotter step 5: 18.3 s\n", "Trotter step 6: 22.1 s\n", "Trotter step 7: 25.6 s\n", "Trotter step 8: 29.4 s\n", "Trotter step 9: 33.2 s\n", "Trotter step 10: 36.6 s\n", "\n", "Total PPM simulation time: 203.6 s\n", "Truncation threshold used: 0.001\n" ] } ], "source": [ "import time\n", "from pauli_prop import evolve_through_cliffords, propagate_through_circuit\n", "\n", "# ── PPM Configuration ──\n", "# Truncation threshold: controls the speed/accuracy trade-off.\n", "PPM_THRESHOLD = 1e-3\n", "\n", "# Maximum Pauli terms to track per observable (hard cap on memory/time)\n", "PPM_MAX_TERMS = 66_000\n", "\n", "print(f\"PPM settings: atol={PPM_THRESHOLD}, max_terms={PPM_MAX_TERMS}\")\n", "\n", "# We propagate each single-qubit Z observable through each circuit.\n", "# For PPM, we work with the un-transpiled circuits (ideal noiseless simulation).\n", "\n", "observables_pp = [\n", " SparsePauliOp(\"I\" * i + \"Z\" + \"I\" * (num_qubits_hw - i - 1))\n", " for i in range(num_qubits_hw)\n", "]\n", "\n", "\n", "def ppm_expectation_values(\n", " circuit, observables, max_terms=PPM_MAX_TERMS, atol=PPM_THRESHOLD\n", "):\n", " \"\"\"Compute expectation values of single-qubit Z observables\n", " via Pauli propagation.\n", "\n", " Args:\n", " circuit: The quantum circuit to simulate.\n", " observables: List of single-qubit Z observables.\n", " max_terms: Maximum number of Pauli terms to retain (hard cap).\n", " atol: Absolute tolerance — Pauli terms with coefficients below this\n", " value are discarded during propagation. Larger values give\n", " faster simulation at the cost of approximation accuracy.\n", " \"\"\"\n", " circuit = circuit.decompose([\"swap\"]) # decompose SWAPs into 3 CX gates\n", " cliff, non_cliff = evolve_through_cliffords(circuit)\n", "\n", " evs = []\n", " for obs in observables:\n", " evolved_obs = propagate_through_circuit(\n", " obs, non_cliff, max_terms=max_terms, atol=atol, frame=\"h\"\n", " )[0]\n", " evolved_obs.paulis = evolved_obs.paulis.evolve(cliff, frame=\"h\")\n", " diagonal_mask = ~evolved_obs.paulis.x.any(axis=1)\n", " ev = float(evolved_obs.coeffs[diagonal_mask].sum().real)\n", " evs.append(ev)\n", " return np.array(evs)\n", "\n", "\n", "# Run PPM for each Trotter step and record wall-clock time\n", "pp_expvals_mid = []\n", "pp_expvals = []\n", "pp_times = []\n", "\n", "for idx, d in enumerate(trotter_steps_hw):\n", " t_start = time.perf_counter()\n", "\n", " # Meson circuit\n", " evs_mid = ppm_expectation_values(circuits_mid_hw[idx], observables_pp)\n", "\n", " # Vacuum circuit\n", " evs_vac = ppm_expectation_values(circuits_hw[idx], observables_pp)\n", "\n", " elapsed = time.perf_counter() - t_start\n", " pp_times.append(elapsed)\n", "\n", " pp_expvals_mid.append(evs_mid[::-1])\n", " pp_expvals.append(evs_vac[::-1])\n", "\n", " print(f\"Trotter step {d:2d}: {elapsed:.1f} s\")\n", "\n", "print(f\"\\nTotal PPM simulation time: {sum(pp_times):.1f} s\")\n", "print(f\"Truncation threshold used: {PPM_THRESHOLD}\")" ] }, { "cell_type": "code", "execution_count": 25, "id": "pauli-prop-timing-plot", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "\"Output" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# --- PPM simulation time vs. Trotter steps ---\n", "fig, ax = plt.subplots(figsize=(8, 5))\n", "ax.plot(\n", " list(trotter_steps_hw),\n", " pp_times,\n", " \"o-\",\n", " color=\"tab:blue\",\n", " linewidth=2,\n", " markersize=6,\n", ")\n", "ax.set_xlabel(\"Trotter step\", fontsize=13)\n", "ax.set_ylabel(\"Wall-clock time (s)\", fontsize=13)\n", "ax.set_title(\n", " \"Pauli Propagation simulation time vs. Trotter steps\\n(30-site lattice, 60 qubits)\",\n", " fontsize=13,\n", ")\n", "ax.grid(True, alpha=0.3)\n", "plt.tight_layout()\n", "plt.show()" ] }, { "cell_type": "code", "execution_count": 37, "id": "pauli-prop-heatmap", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "\"Output" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "# --- PPM heatmap and comparison with hardware ---\n", "N_mid_pp = get_number(pp_expvals_mid, num_lattice_point_hw)\n", "N_pp = get_number(pp_expvals, num_lattice_point_hw)\n", "N_diff_pp = calculate_difference(N_mid_pp, N_pp, num_lattice_point_hw)\n", "\n", "N_diff_pp_arr = np.array(N_diff_pp)\n", "\n", "fig, axes = plt.subplots(1, 2, figsize=(18, 6))\n", "vmax = np.abs(N_hw).max()\n", "vmin = -vmax\n", "norm = TwoSlopeNorm(vmin=vmin, vcenter=0, vmax=vmax)\n", "\n", "# PPM result\n", "im0 = axes[0].imshow(\n", " N_diff_pp_arr,\n", " aspect=\"auto\",\n", " origin=\"lower\",\n", " cmap=\"RdBu_r\",\n", " norm=norm,\n", " extent=[0, num_lattice_point_hw, 0.5, len(trotter_steps_hw) + 0.5],\n", ")\n", "axes[0].set_xlabel(\"Lattice site $r$\", fontsize=12)\n", "axes[0].set_ylabel(\"Trotter step $t$\", fontsize=12)\n", "axes[0].set_title(\n", " \"Pauli Propagation\\n(classical noiseless simulation)\", fontsize=12\n", ")\n", "plt.colorbar(im0, ax=axes[0], label=\"$n_f(r,t)$\")\n", "\n", "# Hardware result\n", "im1 = axes[1].imshow(\n", " N_diff_hw_arr,\n", " aspect=\"auto\",\n", " origin=\"lower\",\n", " cmap=\"RdBu_r\",\n", " norm=norm,\n", " extent=[0, num_lattice_point_hw, 0.5, len(trotter_steps_hw) + 0.5],\n", ")\n", "axes[1].set_xlabel(\"Lattice site $r$\", fontsize=12)\n", "axes[1].set_ylabel(\"Trotter step $t$\", fontsize=12)\n", "axes[1].set_title(\n", " \"Quantum Simulation\\n(IBM Hardware, readout error mitigation only)\",\n", " fontsize=12,\n", ")\n", "plt.colorbar(im1, ax=axes[1], label=\"$n_f(r,t)$\")\n", "\n", "plt.suptitle(\n", " \"Staggered Fermionic Occupation Number Dynamics — 30-site lattice\",\n", " fontsize=14,\n", " y=1.02,\n", ")\n", "plt.tight_layout()\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "next-steps", "metadata": {}, "source": [ "\n", "\n", "## 次のステップ\n", "\n", "この作品に興味を持たれた方は、以下の資料もぜひご覧ください：\n", "\n", "\n", " * [Qiskit Estimator プリミティブのドキュメント](/docs/guides/get-started-with-estimator) — エラー軽減オプションの設定に関する詳細\n", " * [エラーの軽減および抑制手法](/docs/guides/error-mitigation-and-suppression-techniques) — TREX、ZNE、その他の軽減手法について学ぶ\n", " * [Qiskit Pauli Propagation (pauli-prop)](https://github.com/Qiskit/pauli-prop) — パウリ逆伝播を用いたRustによる高速化された古典シミュレーション\n", "\n", "\n", "\n", "\n", "## 参照\n", "\n", "\\[1] 原著論文：Ilčić, Majumdar, Mathew ほか「ノイズの多い量子プロセッサにおける堅牢かつコヒーレントな非アベル型ハドロンダイナミクスの観測」『 [arXiv:2602.18080](https://arxiv.org/abs/2602.18080) 』（2026年）\n", "\n" ] }, { "cell_type": "markdown", "metadata": {}, "id": "a1b8767d", "source": "© IBM Corp., 2017-2026" } ], "metadata": { "kernelspec": { "display_name": "Python 3", "language": "python", "name": "python3" }, "language_info": { "codemirror_mode": { "name": "ipython", "version": 3 }, "file_extension": ".py", "mimetype": "text/x-python", "name": "python", "nbconvert_exporter": "python", "pygments_lexer": "ipython3", "version": "3" }, "hours": 1, "qpuSeconds": 360 }, "nbformat": 4, "nbformat_minor": 5 }