NumPyでゼロから作る量子シミュレータ｜状態ベクトルで量子を"自作"する

Qiskitのような量子ライブラリは便利ですが、「中で何が起きているか」はブラックボックスになりがちです。いちばんの理解の近道は、自分で作ってみること。この記事ではNumPyだけで、状態ベクトル型の量子シミュレータをゼロから組み立てます。100行に満たないコードで、重ね合わせ・ゲート・もつれ・測定がすべて動きます。

発想: 量子状態は「複素数のベクトル」

n個の量子ビットの状態は、2ⁿ個の複素数の並び（状態ベクトル） で完全に表せます。1量子ビットなら [α, β] の2個。αは「0が出る振幅」、βは「1が出る振幅」で、測定すると |α|²・|β|² の確率で 0/1 が出ます。

準備

pip install numpy

必要なのはこれだけ。量子ライブラリは使いません。

① 状態を用意する

nビットの初期状態 |00...0⟩ は、先頭だけ1・残り0のベクトルです。

import numpy as np

def zero_state(n):
    state = np.zeros(2**n, dtype=complex)
    state[0] = 1.0          # |00...0⟩
    return state

print(zero_state(1))   # [1.+0.j 0.+0.j] = |0⟩
print(zero_state(2))   # [1 0 0 0] = |00⟩

② ゲートは「行列」、適用は「掛け算」

量子ゲートはユニタリ行列です。1量子ビットの代表を定義します。

H = (1/np.sqrt(2)) * np.array([[1, 1],
                               [1, -1]], dtype=complex)   # アダマール
X = np.array([[0, 1], [1, 0]], dtype=complex)             # NOT
Z = np.array([[1, 0], [0, -1]], dtype=complex)
I = np.eye(2, dtype=complex)

問題は「nビットの状態」に「1ビットのゲート」をどう掛けるか。答えはテンソル積で、対象ビット以外を単位行列 I で埋めた大きな行列を作ることです。

def apply_1q(state, gate, target, n):
    # 各ビットにI、targetにだけgateを置いてテンソル積で全体行列を作る
    op = np.array([[1]], dtype=complex)
    for q in range(n):
        op = np.kron(op, gate if q == target else I)
    return op @ state

さっそく重ね合わせを作ってみます。

n = 1
s = zero_state(n)
s = apply_1q(s, H, target=0, n=n)
print(s)   # [0.707+0.j 0.707+0.j] = (|0⟩+|1⟩)/√2

アダマール1発で [0.707, 0.707]。0と1が半々の重ね合わせが、行列の掛け算だけで再現できました。

③ 2量子ビットゲート（CNOT）でもつれを作る

もつれには2ビットゲートが要ります。CNOTを直接 4×4 行列で定義します。

CNOT = np.array([[1,0,0,0],
                 [0,1,0,0],
                 [0,0,0,1],
                 [0,0,1,0]], dtype=complex)

def apply_2q_full(state, gate4):
    return gate4 @ state

n = 2
s = zero_state(n)
s = apply_1q(s, H, target=0, n=n)   # 0番目を重ね合わせ
s = apply_2q_full(s, CNOT)          # もつれさせる
print(np.round(s, 3))
# [0.707 0 0 0.707] = (|00⟩+|11⟩)/√2

|00⟩ と |11⟩ にだけ振幅が立ち、|01⟩・|10⟩ はゼロ。これが**ベル状態（もつれ）**です。片方が0なら必ずもう片方も0、というあの相関を、数字で確認できました。

④ 測定する

測定は「確率分布に従ってサンプリングし、状態を収縮させる」操作です。

def measure(state, shots=1000):
    probs = np.abs(state)**2                 # 各基底の確率
    n = int(np.log2(len(state)))
    outcomes = np.random.choice(len(state), size=shots, p=probs)
    counts = {}
    for o in outcomes:
        key = format(o, f"0{n}b")            # 例: 3 -> "11"
        counts[key] = counts.get(key, 0) + 1
    return counts

print(measure(s))   # {'00': 498, '11': 502} のように出る

{'00': 498, '11': 502} ——もつれた2ビットは「00」と「11」しか出ません。Qiskitで同じ実験をしたときと、まったく同じ結果です。あなたが今、量子コンピュータのシミュレータを自作したということです。

なぜ実機が必要なのか（自作の限界）

このシミュレータは美しく動きますが、状態ベクトルは量子ビット数nに対して 2ⁿ で爆発します。30量子ビットで約10億個の複素数、40量子ビットで家庭用PCのメモリを完全に超えます。「古典で真似できなくなる」その先に、本物の量子コンピュータの存在意義があります。

まとめ

• 量子状態は「2ⁿ個の複素数ベクトル」で完全に表せる
• ゲートは行列、適用は np.kron で全体化して掛け算するだけ
• 測定は |振幅|² を確率としたサンプリング
• 状態ベクトルは2ⁿで爆発する——だから実機に意味がある

もう少し詳しく（背景と理論）

状態ベクトル法は「量子状態を 2ⁿ 個の複素振幅で厳密に保持する」もっとも素直なシミュレーション手法です¹。メモリが指数的に増えるため、現実的に扱えるのはおおよそ 30〜40 量子ビット程度で²、この壁の先に「古典で真似できない」量子超越性の議論があります³。

なお、量子回路すべてが古典で難しいわけではありません。アダマール・位相・CNOT だけからなるクリフォード回路は、状態ベクトルを陽に持たなくても多項式時間で古典シミュレートできます（Gottesman–Knill 定理）⁴。「量子だから速い」のは、こうした古典的に効率シミュレート可能なクラスの外側を使うときに限られる、という点が重要です。

次の一歩 🌸

ライブラリ版で同じことをやるならQiskit入門、ゲートの中身をもっと知りたいなら量子ゲート入門へどうぞ。

Footnotes

1. 状態ベクトル法（Schrödinger 型シミュレーション）は n 量子ビットの状態を 2ⁿ 個の複素数で表す。倍精度複素数は1個16バイトなので、必要メモリは 16×2ⁿ バイト。 ↩

2. 具体的には30量子ビットで約16GB、40量子ビットで約16TB、45量子ビットで約0.5PBとなり、家庭用〜スパコン級で扱える上限がこの付近にある。テンソルネットワーク法など、回路の構造を利用してこれを緩和する手法も研究されている。 ↩

3. Arute, F. et al. (2019). "Quantum supremacy using a programmable superconducting processor." Nature, 574, 505–510. 53量子ビットで古典シミュレーションが現実的でない領域に到達したと主張した（量子超越性参照）。 ↩

4. Gottesman, D. (1998); Aaronson, S. & Gottesman, D. (2004). "Improved Simulation of Stabilizer Circuits." Phys. Rev. A, 70, 052328. クリフォードゲートのみの回路はスタビライザー形式で多項式時間シミュレートできる。普遍量子計算にはT ゲートなど非クリフォード要素が必要。 ↩