{ "cells": [ { "cell_type": "markdown", "id": "1379da5d-023f-42e7-8a3f-ea8c1e3a1d46", "metadata": {}, "source": [ "---\n", "title: \"サンプルベースのクリロフ量子対角化法\"\n", "description: \"SKQDでは、QKDの収束保証とSQDの速度を組み合わせるために、これまでの授業で扱った内容を統合する。\"\n", "---\n", "\n", "{/* cspell:ignore arxiv Neel probs eigenvals Hatano */}\n", "\n", "\n", "\n", "# サンプルベース・クリロフ量子対角化法（SKQD）\n", "\n" ] }, { "cell_type": "markdown", "id": "97cfd899-4196-4d24-bf16-60bdd3899600", "metadata": {}, "source": [ "サンプルベースのクリロフ量子対角化（SKQD）に関するこのレッスンでは、以前のメソッドで説明したメソッドを組み合わせています。 Qiskitパターンフレームワークを活用した一つの例で構成されている：\n", "\n", "* ステップ1：問題を量子回路と演算子にマップする\n", "* ステップ2：ターゲット・ハードウェアに最適化する\n", "* ステップ 3: `IBM Quantum` プリミティブを使用して実行する\n", "* ステップ4：後処理\n", "\n", "サンプルベースの量子対角化法における重要なステップは、部分空間の品質ベクトルを生成することである。 前のレッスンでは、LUCJアサッツを使って化学ハミルトニアンの部分空間ベクトルを生成した。 このレッスンでは、レッスン2で説明したように、量子クリロフ状態 [\\[1\\]](#references) を使います。 まず、時間発展操作を使って量子コンピューター上でクリロフ空間を作る方法を復習する。 それからサンプリングする。 システムのハミルトニアンをサンプリングされた部分空間に射影し、対角化して基底状態のエネルギーを推定する。 このアルゴリズムは、レッスン2で説明した仮定の下で、証明可能かつ効率的に基底状態に収束する。\n", "\n", "\n", "\n", "## 0. クリロフ空間\n", "\n", "次数のクリロフ空間 $\\mathcal{K}^r$ $r$ は、行列 $A$ の高次のべき乗 $r-1$ までを参照ベクトル $\\vert v \\rangle$ と掛け合わせることによって得られるベクトルによってスパンされる空間である。\n", "\n", "$$\n", "\\mathcal{K}^r = \\left\\{ \\vert v \\rangle, A \\vert v \\rangle, A^2 \\vert v \\rangle, ..., A^{r-1} \\vert v \\rangle \\right\\}\n", "$$\n", "\n", "行列 $A$ がハミルトニアン $H$ の場合、対応する空間は*べき乗クリロフ空間* $\\mathcal{K}_P$ と呼ばれる。 $A$ がハミルトニアン $U=e^{-iH(dt)}$ によって生成される時間発展演算子である場合、その空間はユニタリー・クリロフ空間 $\\mathcal{K}_U$ と呼ばれる。パワークリロフ部分空間は、 $H$ がユニタリー演算子ではないため、量子コンピュータ上で直接生成することはできません。 その代わりに、時間発展演算子 $U = e^{-iH(dt)}$、パワークリロフ空間と同様の[収束保証が](https://arxiv.org/abs/2110.07492)得られることを示すことができる。 $U$ のべき乗は、異なる時間ステップ $U^k = e^{-iH(k dt)}$ となり、 $k = 0, 1, 2, ..., (r-1)$.\n", "\n", "$$\n", "\\mathcal{K}_U^r = \\left\\{ \\vert \\psi \\rangle, U \\vert \\psi \\rangle, U^2 \\vert \\psi \\rangle, ..., U^{r-1} \\vert \\psi \\rangle \\right\\}\n", "$$\n", "\n" ] }, { "cell_type": "markdown", "id": "864a0ce7-d094-48ea-9cf0-2a7e38791815", "metadata": {}, "source": [ "\n", "\n", "## 1. 問題を量子回路と演算子に写像する\n", "\n" ] }, { "cell_type": "markdown", "id": "5888e3e8-51a1-426a-9023-c137aa9c6c83", "metadata": {}, "source": [ "このレッスンでは、 $L = 22$ サイトを持つ反強磁性XX-Z spin-1/2 鎖のハミルトニアンを周期境界条件付きで考える：\n", "\n", "$$\n", " H = \\sum_{i, j}^{N} J_{xy} (X_{i} X_{j} + Y_{i} Y_{j}) + Z_{i} Z_{j}\n", "$$\n", "\n" ] }, { "cell_type": "code", "execution_count": 1, "id": "f0d3e153-17da-4457-b1d2-2a35d31aea3a", "metadata": {}, "outputs": [], "source": [ "from qiskit.transpiler import CouplingMap\n", "from qiskit_addon_utils.problem_generators import generate_xyz_hamiltonian\n", "\n", "num_spins = 22\n", "coupling_map = CouplingMap.from_ring(num_spins)\n", "H_op = generate_xyz_hamiltonian(coupling_map, coupling_constants=(0.3, 0.3, 1.0))" ] }, { "cell_type": "markdown", "id": "b0660e0d-5521-4095-a810-949729421b57", "metadata": {}, "source": [ "クリロフ空間を構築するには、3つの主要な材料が必要だ：\n", "\n", "1. クリロフ次元( $r$ )と時間ステップ( $dt$ )の選択。\n", "2. ターゲット（地上）状態と多項式にオーバーラップする初期（参照）状態（ベクトル $\\vert v \\rangle$ 上記）。ターゲット状態はスパースである。 この多項式オーバーラップの要件は、量子位相推定アルゴリズムと同じである。\n", "3. 時間発展演算子 $U^{k}=e^{-iH(k * dt)}$ ( $k = 0, 1, 2, ..., r-1$ ).\n", "\n", "$r$ （および、 $dt$ ）の選択された値に対して、 $r$ 量子回路を作成し、それらからサンプルを採取する。 各量子回路は、参照状態の量子回路表現と、 $k$ 値の時間発展演算子を結合することで作成される。\n", "\n", "クリロフ次元が大きいほど、推定エネルギーの収束性が向上する。 このレッスンでは、収束の傾向を説明するために、ディメンションを $5$ に設定する。\n", "\n", "文献 [\\[2\\]](#references) は、KQDの十分に小さい時間ステップは $\\pi / \\vert \\vert H \\vert \\vert$、この値を過大評価するよりも過小評価する方が望ましいことを示した。 一方、 $dt$ を小さくしすぎると、クリロフ基底ベクトルのタイムステップごとの違いが小さくなるため、クリロフ部分空間のコンディショニングが悪くなる。 さらに、この $dt$ の選択はSKQDの収束にとって証明可能なほど適切であるが、このサンプリングに基づくコンテキストでは、 $dt$ の実際の最適な選択は現在進行中の研究テーマである。 このレッスンでは、 $dt = 0.15$ を設定する。\n", "\n", "クリロフ次元と時間ステップの他に、時間発展のためのトロッターステップ数を設定する必要がある。 ステップ数が少なすぎるとトロッター化の誤差が大きくなり、ステップ数が多すぎると回路が深くなる。 このレッスンでは、トロッターのステップ数を $6$ に設定する。\n", "\n" ] }, { "cell_type": "code", "execution_count": 2, "id": "52eb15e3-e06c-4820-bad1-80fb7c05feee", "metadata": {}, "outputs": [], "source": [ "# Set parameters for quantum Krylov algorithm\n", "krylov_dim = 5 # size of krylov subspace\n", "dt = 0.15\n", "num_trotter_steps = 6" ] }, { "cell_type": "markdown", "id": "a716ea3e-b3d6-42de-8ed5-331e6416e586", "metadata": {}, "source": [ "次に、基底状態と重なる参照状態（ $\\vert \\psi \\rangle$ ）を選ぶ必要がある。 このハミルトニアンでは、 1s と 0s $\\vert ...101...010...101 \\rangle$ を交互に繰り返す*ニール*状態を参照状態として使用する。\n", "\n" ] }, { "cell_type": "code", "execution_count": 3, "id": "53347df7-b43c-4859-a9ea-a3dea2a44184", "metadata": {}, "outputs": [], "source": [ "# Prep `Neel` state as the reference state for evolution\n", "from qiskit import QuantumCircuit\n", "\n", "qc_state_prep = QuantumCircuit(num_spins)\n", "for i in range(num_spins):\n", " if i % 2 == 0:\n", " qc_state_prep.x(i)" ] }, { "cell_type": "markdown", "id": "6952aaf2-c598-4ce9-8dc5-e63f0725de10", "metadata": {}, "source": [ "最後に、時間発展演算子を量子回路にマッピングする必要がある。 これはレッスン2で行ったが、ここではQiskitのメソッド、特に *synthesisという*メソッドを活用する。 量子ゲートを持つ量子回路に数学演算子を合成する方法はさまざまだ。 [Qiskitの合成モジュールには](/docs/api/qiskit/synthesis)、このような技法が数多く用意されています。 合成には [`LieTrotter`](/docs/api/qiskit/qiskit.synthesis.LieTrotter) 合成のためのアプローチ [\\[3\\]](#references) [\\[4\\]](#references) を使用する。\n", "\n" ] }, { "cell_type": "code", "execution_count": null, "id": "c82afe58-a90a-49a9-9dfe-c24b2b5d0a8c", "metadata": {}, "outputs": [], "source": [ "from qiskit.circuit import QuantumRegister\n", "from qiskit.circuit.library import PauliEvolutionGate\n", "from qiskit.synthesis import LieTrotter\n", "\n", "evol_gate = PauliEvolutionGate(\n", " H_op, time=(dt / num_trotter_steps), synthesis=LieTrotter(reps=num_trotter_steps)\n", ") # `U` operator\n", "\n", "qr = QuantumRegister(num_spins)\n", "qc_evol = QuantumCircuit(qr)\n", "qc_evol.append(evol_gate, qargs=qr)\n", "\n", "circuits = []\n", "for rep in range(krylov_dim):\n", " circ = qc_state_prep.copy()\n", "\n", " # Repeating the `U` operator to implement U^0, U^1, U^2, and so on, for power Krylov space\n", " for _ in range(rep):\n", " circ.compose(other=qc_evol, inplace=True)\n", "\n", " circ.measure_all()\n", " circuits.append(circ)" ] }, { "cell_type": "code", "execution_count": 8, "id": "7ffc3eaa-9563-4c89-9d7b-808a3a37b275", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "\"Output" ] }, "execution_count": 8, "metadata": {}, "output_type": "execute_result" } ], "source": [ "circuits[1].decompose().draw(\"mpl\", fold=-1)" ] }, { "cell_type": "code", "execution_count": 9, "id": "bdbc0881-06b0-49f2-8919-a376aa7f08e3", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "\"Output" ] }, "execution_count": 9, "metadata": {}, "output_type": "execute_result" } ], "source": [ "circuits[2].decompose().draw(\"mpl\", fold=-1)" ] }, { "cell_type": "markdown", "id": "9b43f693-b6d0-43cb-b982-3f8a3b6e185d", "metadata": {}, "source": [ "\n", "\n", "## 2. 対象ハードウェア向けに最適化\n", "\n" ] }, { "cell_type": "markdown", "id": "6ed5a22f-4bdc-4a8d-ad39-5307308d0cc4", "metadata": {}, "source": [ "回路ができたので、ターゲットとなるハードウェアに最適化することができる。 実用規模のQPUを選ぶ。\n", "\n" ] }, { "cell_type": "code", "execution_count": null, "id": "8944305a-0f29-4f97-a04f-43249cdc02ea", "metadata": {}, "outputs": [], "source": [ "import warnings\n", "\n", "from qiskit import generate_preset_pass_manager\n", "from qiskit_ibm_runtime import QiskitRuntimeService\n", "\n", "warnings.filterwarnings(\"ignore\")\n", "\n", "service = QiskitRuntimeService()\n", "# Use the least-busy backend or specify a quantum computer using the syntax commented out below.\n", "backend = service.least_busy(operational=True, simulator=False)\n", "# backend = service.backend(\"ibm_brisbane\")" ] }, { "cell_type": "markdown", "id": "3ed422ad-04f5-4f1f-a142-f08a0ac0a34b", "metadata": {}, "source": [ "次に、プリセット・パス・マネージャーを使って、回路をターゲット・バックエンドにトランスパイルする。\n", "\n" ] }, { "cell_type": "code", "execution_count": 11, "id": "8a72b655-ccee-4927-8bf2-b99d0b92723a", "metadata": {}, "outputs": [], "source": [ "pm = generate_preset_pass_manager(backend=backend, optimization_level=3)\n", "isa_circuits = pm.run(circuits=circuits)" ] }, { "cell_type": "markdown", "id": "6405ab01-5cec-4a02-8846-f9ba626d9cdc", "metadata": {}, "source": [ "\n", "\n", "## 3. 対象ハードウェア上で実行する\n", "\n" ] }, { "cell_type": "markdown", "id": "6aa67e8c-30db-4e11-9bf2-150c157ec49c", "metadata": {}, "source": [ "ハードウェア実行用に回路を最適化した後、ターゲット・ハードウェア上で回路を実行し、基底状態のエネルギー推定用のサンプルを収集する準備が整った。\n", "\n" ] }, { "cell_type": "code", "execution_count": 12, "id": "0613c95b-e696-4ddb-8f51-dc1714d21a0e", "metadata": {}, "outputs": [], "source": [ "from qiskit_ibm_runtime import SamplerV2 as Sampler\n", "\n", "sampler = Sampler(mode=backend)\n", "job = sampler.run(isa_circuits, shots=100_000) # Takes approximately 2m 58s of QPU time" ] }, { "cell_type": "code", "execution_count": 13, "id": "1c60e828-9359-46c6-84ef-88aecdb85b81", "metadata": {}, "outputs": [], "source": [ "counts_all = [job.result()[k].data.meas.get_counts() for k in range(krylov_dim)]" ] }, { "cell_type": "markdown", "id": "fa2fcfe0-66ab-4230-9262-9a46ecda034b", "metadata": {}, "source": [ "\n", "\n", "## 4. 結果の後処理\n", "\n" ] }, { "cell_type": "markdown", "id": "75321a43-0740-4d13-bdc4-350b44727b62", "metadata": {}, "source": [ "次に、クリロフ次元が増加した場合のカウントを累積的に集計する。 累積カウントを用いて、クリロフ次元が大きくなる部分空間にまたがり、収束挙動を解析する。\n", "\n" ] }, { "cell_type": "code", "execution_count": 14, "id": "9954125c-5ad8-40e8-affa-4a644354c51f", "metadata": {}, "outputs": [], "source": [ "from collections import Counter\n", "\n", "counts_cumulative = []\n", "for i in range(krylov_dim):\n", " counter = Counter()\n", " for d in counts_all[: i + 1]:\n", " counter.update(d)\n", "\n", " counts = dict(counter)\n", " counts_cumulative.append(counts)" ] }, { "cell_type": "markdown", "id": "14fcdfc2-b37e-4a9c-9794-7c1fa9a6ab1b", "metadata": {}, "source": [ "ハミルトニアンを投影し対角化するために、以下の機能を使用する。 [`qiskit-addon-sqd`](https://qiskit.github.io/qiskit-addon-sqd/). このアドオンには、パウリ文字列に基づくハミルトニアンを部分空間に投影し、 `SciPy` を使って固有値を解く機能があります。\n", "\n" ] }, { "cell_type": "code", "execution_count": 15, "id": "5ace841d-59c6-4ea1-9ceb-6f6436546093", "metadata": {}, "outputs": [], "source": [ "from qiskit_addon_sqd.counts import counts_to_arrays\n", "from qiskit_addon_sqd.qubit import solve_qubit" ] }, { "cell_type": "markdown", "id": "104cc0fd-4fa0-469b-bfb8-4089cea5145f", "metadata": {}, "source": [ "原理的には、部分空間をスパニングする前に、不正確なパターンを持つビット列をフィルタリングすることができる。 例えば、このレッスンで習う反強磁性ハミルトニアンの基底状態は、通常、「アップ」と「ダウン」のスピンの数が等しく、つまり、ビット列の「1」の数は、系内のビット（スピン）の総数のちょうど半分でなければならない。 以下の関数は、\"1 \"の数が正しくないビット列をカウントから除外する。\n", "\n" ] }, { "cell_type": "code", "execution_count": 37, "id": "81d728a1-6ab4-4664-a7ca-2017dae10e1a", "metadata": {}, "outputs": [], "source": [ "# Filters out bitstrings that do not have specified number (`num_ones`) of `1` bits.\n", "def postselect_counts(counts, num_ones):\n", " filtered_counts = {}\n", " for bitstring, freq in counts.items():\n", " if bitstring.count(\"1\") == num_ones:\n", " filtered_counts[bitstring] = freq\n", "\n", " return filtered_counts" ] }, { "cell_type": "markdown", "id": "d77772f5-fe1c-4a6c-abcb-c2c3c7d7293a", "metadata": {}, "source": [ "正しい数のアップ/ダウン電子を持つビット列を使い、部分空間をスパンし、クリロフ次元を増加させる固有値を計算する。 問題のサイズと利用可能な古典的リソースによっては、（ [SQDのレッスンに](/learning/courses/quantum-diagonalization-algorithms/sqd-implementation)似た）サブサンプリングを採用して、部分空間の次元を抑える必要があるかもしれない。 さらに、レッスン4と同様にコンフィギュレーション・リカバリーの概念を応用することもできる。 再構成された固有状態からサイトごとの電子占有率を計算し、その情報を使ってアップ/ダウン電子数が正しくないビット列を修正することができる。 これは、興味のある読者のための練習として残しておく。\n", "\n" ] }, { "cell_type": "code", "execution_count": 38, "id": "f0e0f56f-ec56-40b6-b8be-059be35d2ce5", "metadata": {}, "outputs": [], "source": [ "import numpy as np\n", "\n", "num_batches = 10\n", "rand_seed = 0\n", "scipy_kwargs = {\"k\": 2, \"which\": \"SA\"}\n", "\n", "ground_state_energies = []\n", "for idx, counts in enumerate(counts_cumulative):\n", " counts = postselect_counts(counts, num_ones=num_spins // 2)\n", " bitstring_matrix, probs = counts_to_arrays(counts=counts)\n", "\n", " eigenvals, eigenstates = solve_qubit(\n", " bitstring_matrix, H_op, verbose=False, **scipy_kwargs\n", " )\n", " gs_en = np.min(eigenvals)\n", " ground_state_energies.append(gs_en)" ] }, { "cell_type": "markdown", "id": "67b2a04b-4e70-45dc-bd21-2c4340c18377", "metadata": {}, "source": [ "次に、計算されたエネルギーをクリロフ次元の関数としてプロットし、厳密エネルギーと比較する。 正確なエネルギーは、古典的なブルートフォース法を用いて別途計算される。 推定される基底状態のエネルギーは、クリロフ空間の次元が大きくなるにつれて収束していくことがわかる。 $5$ のクリロフ次元は制限的であるが、それでも結果は印象的な収束を示しており、クリロフ次元が大きくなるにつれて改善することが期待される [\\[1\\]](#references)。\n", "\n" ] }, { "cell_type": "code", "execution_count": 39, "id": "debca96a-e0e3-4342-a755-8204dba50d74", "metadata": {}, "outputs": [ { "data": { "text/plain": [ "\"Output" ] }, "metadata": {}, "output_type": "display_data" } ], "source": [ "import matplotlib.pyplot as plt\n", "\n", "exact_gs_en = -23.934184\n", "plt.plot(\n", " range(1, krylov_dim + 1),\n", " ground_state_energies,\n", " color=\"blue\",\n", " linestyle=\"-.\",\n", " label=\"estimate\",\n", ")\n", "plt.plot(\n", " range(1, krylov_dim + 1),\n", " [exact_gs_en] * krylov_dim,\n", " color=\"red\",\n", " linestyle=\"-\",\n", " label=\"exact\",\n", ")\n", "plt.xticks(range(1, krylov_dim + 1), range(1, krylov_dim + 1))\n", "plt.legend()\n", "plt.xlabel(\"Krylov space dimension\")\n", "plt.ylabel(\"Energy\")\n", "plt.ylim([-24, -22.50])\n", "plt.title(\n", " \"Estimating Ground state energy with Sample-based Krylov Quantum Diagonalization\"\n", ")\n", "plt.show()" ] }, { "cell_type": "markdown", "id": "5865567b-cc32-48dd-b7f9-32ce64d7dd7c", "metadata": {}, "source": [ "\n", "\n", "#### 理解度チェック\n", "\n", "下の問題を読んで、答えを考え、三角形をクリックして解答を明らかにしよう。\n", "\n", "\n", " \n", " **回答:**\n", "\n", " クリロフ次元を上げる。 一般的には、シュート数を増やすこともできるが、上記の計算ではすでにかなり多くなっている。\n", " \n", "\n", " \n", " **回答:**\n", "\n", " 他にも有効な答えがあるかもしれないが、完全な答えには以下を含めること：\n", "\n", " (a) SKQDにはSQDにはない収束保証がある。 SQDでは、計算基底における基底状態のサポートと非常によく重なるようなアナザッツを推測するか、アナザッツのファミリーをサンプリングするために変分要素を計算に導入する必要があります。\n", "\n", " (b)SKQDは、ハダマルドテストによる行列要素の計算を回避できるため、QPUの所要時間が大幅に短縮される。\n", " \n", "\n", "\n" ] }, { "cell_type": "markdown", "id": "a4cf119d-65cf-456c-be0d-731bb7a6799a", "metadata": {}, "source": [ "\n", "\n", "## 5. 概要\n", "\n", "* クリロフ基底状態のサンプリングによる基底状態のエネルギー推定は、スピン系、凝縮系問題、格子ゲージ理論などの格子モデルに非常に適している。 VQEやヒューリスティック・アンサッツに基づくSQD（例えば、前のレッスンの化学問題）のように、変分アサッツで多くのパラメータを最適化する必要がないため、このアプローチはVQEよりもはるかに優れたスケールを持つ。\n", " * 回路の深さを低く抑えるには、フォールト・トレラント・ハードウェアに適した格子問題に取り組むのが賢明である。\n", "* SKQDでは、VQEのような量子測定の問題は発生しない。 推定すべきパウリ作用素のグループはない。\n", "* SKQDはノイズの多いサンプルに対しても頑健である。なぜなら、問題固有のポストセレクションルーチン（例えば、問題固有のパターンに合致しないビット列を除外する）を使用したり、従来の対角化のオーバーヘッド（すなわち、より大きな部分空間で対角化を行う）を許容したりすることで、ノイズの影響を効果的に除去できるからである。\n", "\n" ] }, { "cell_type": "markdown", "id": "c52bada7-1912-41f2-a770-3145a947de63", "metadata": {}, "source": [ "\n", "\n", "## 参照\n", "\n", "\\[1] ジェフェリー・ユーほか \"サンプルベースクリロフ対角化のための量子中心アルゴリズム\" (2025). [arxiv:quant-ph/2501.09702](https://arxiv.org/abs/2501.09702).\n", "\n", "\\[2] イーサン・エッペリー、リン・リン、中務裕司。 「量子部分空間対角化の理論」。 SIAM Journal on Matrix Analysis and Applications 43, 1263-1290 (2022).\n", "\n", "\\[2] N. 波多野とM. 鈴木, \"高次の指数積公式の発見\" (2005). [arXiv:math-ph/0506007](https://arxiv.org/abs/math-ph/0506007).\n", "\n", "\\[4] D. ベリー、G.アホカス、R.クリーブ、B. Sanders, \"Efficient quantum algorithms for simulating sparse Hamiltonians\" (2006). [arXiv:quant-ph/0508139](https://arxiv.org/abs/quant-ph/0508139).\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" } }, "nbformat": 4, "nbformat_minor": 4 }