量子コンピュータで組み合わせ最適化問題を解くアルゴリズム、QAOA(Quantum Approximate Optimization Algorithm)。
Qiskitなどのライブラリを使えば、数行のコードで「自動的に」回路を生成して解を求めることができます。しかし、真に量子アルゴリズムを理解し、独自の改良を加えるためには、ブラックボックスを脱ぎ捨てて「回路をフルスクラッチで自作する」工程(レベル3)の理解が不可欠です。
今回は、数式をいかにして物理的な「量子ゲート」へと翻訳していくのか、その核心を解説します。
1. ハミルトニアンの構築:問題のルールを定義する
まずは、グラフの「Max-Cut問題」を量子的な言葉であるハミルトニアン(H)に変換します。
以前は文字列(['I', 'I', 'Z', 'Z'])を使って定義しましたが、よりプログラミング的なアプローチとして、Numpy配列を用いた手法を紹介します。どちらの手法をとっても、最終的に生成される行列は同一です。
import numpy as np
from qiskit.quantum_info import Pauli, SparsePauliOp
n = 4 # 量子ビット数
edges = [(0, 1), (0, 2), (1, 2), (1, 3), (2, 3)] # グラフの接続
pauli_list = []
for i, j in edges:
z = np.zeros(n, dtype=bool)
x = np.zeros(n, dtype=bool)
# i番目とj番目にZゲートを適用するフラグを立てる
z[i] = True
z[j] = True
# Pauli((z, x)) で Zi Zj というパウリ演算子を生成
pauli_list.append((Pauli((z, x)).to_label(), 1.0))
H = SparsePauliOp.from_list(pauli_list)
print(f"構築されたハミルトニアン:\n{H}")
2. 初期状態の作成:全探索の準備
QAOAの第一歩は、すべての組み合わせを平等に探索できるよう、全量子ビットを「重ね合わせ状態」にすることです。
from qiskit import QuantumCircuit
qc = QuantumCircuit(n)
# 全ビットにアダマールゲート(H)を適用: |+> 状態の作成
for i in range(n):
qc.h(i)
3. コストレイヤーの物理実装
ここがQAOAの心臓部です。ハミルトニアンの各項(ZiZj)を、量子回路に教え込みます。
物理学的な時間発展の数式 e^(-iHt) を実現するために、「CNOTサンドイッチ」と呼ばれる手法を用います。
# gammas[0] は最適化パラメータ
qc.cx(0, 1) # 1. 二つのビットを絡ませる(エンタングル)
qc.rz(2 * gammas[0], 1) # 2. 位相を回転させスコアを計算
qc.cx(0, 1) # 3. 状態を元に戻す
なぜこの3行で数式が表せるのか?
数学的な定義を確認すると、QiskitのRz ゲートは以下の通りです。
$$RZ(\theta) = e^{-i \frac{\theta}{2} Z}$$ (lateXですみません。。。)
私たちが作りたいのは $e^{-i \gamma Z}$ です。そこで、コード内で角度を $2\gamma$ と指定することで、分母の2を打ち消しています。
RZ(2\gamma) = e^{-i \frac{2\gamma}{2} Z} = e^{-i \gamma Z}ということですね。
さらに、両端を Cxゲート で挟むことで、単一ビットの Z が、二つのビットの相互作用であるZiZj へと変換されるのです。
数式で言うと、CNOT_{i,j} \cdot RZ_j(2\gamma) \cdot CNOT_{i,j} = e^{-i \gamma Z_i Z_j}
4. ミキサーレイヤーによる解の探索
コストレイヤーで計算した「スコア」をもとに、別の可能性を探索(ミックス)するのがミキサーレイヤーです。量子アニーリングにおける「横磁場」の役割を果たします。
Python
# 全ての量子ビットを X軸回りに回転させる
for i in range(n):
qc.rx(2 * betas[0], i)
この $RX$ ゲートによる「揺さぶり」が、局所最適解(たまたま見つけた浅い谷)から抜け出し、真の最適解へと導くトンネル効果を生み出します。
5. 汎用エンジン(VQE)による最適化
最後に、自作した回路(Ansatz)を汎用最適化エンジンである VQE (Variational Quantum Eigensolver) に渡します。
Python
from qiskit_algorithms import VQE
from qiskit_algorithms.optimizers import COBYLA
from qiskit.primitives import Estimator
optimizer = COBYLA()
estimator = Estimator()
# QAOA専用クラスではなく、汎用のVQEに自作回路(qc)を投入
vqe = VQE(estimator=estimator, ansatz=qc, optimizer=optimizer)
result = vqe.compute_minimum_eigenvalue(H)
print(f"最適化されたパラメータ: {result.optimal_parameters}")
6. 全体のコード
import numpy as np
from qiskit import QuantumCircuit
from qiskit.circuit import ParameterVector
from qiskit.quantum_info import Pauli, SparsePauliOp, Statevector
from qiskit.primitives import StatevectorEstimator
from qiskit_algorithms import VQE
from qiskit_algorithms.optimizers import COBYLA
# 頂点は 0,1,2,3(=頂点1〜4)に対応
# 辺: (1,2),(1,3),(2,3),(2,4),(3,4)
# → インデックスでは: (0,1),(0,2),(1,2),(1,3),(2,3)
n = 4
edges = [(0, 1), (0, 2), (1, 2), (1, 3), (2, 3)]
pauli_list = []
coeffs = []
#Max-Cutをイジング化
for i, j in edges:
pauli_str_list = ["I"]*n
pauli_str_list[i] = "Z"
pauli_str_list[j] = "Z"
pauli_str_list.reverse()
pauli_str = "".join(pauli_str_list)
pauli_list.append((pauli_str, 1.0))
#ハミルトニアン定義
H = SparsePauliOp.from_list(pauli_list)
print("Hamiltonian H:")
print(H)
##自作QAOAの量子ゲート
p = 2
gammas = ParameterVector("gamma", p) # gamma[0], gamma[1]
betas = ParameterVector("beta", p) # beta[0], beta[1]
qc = QuantumCircuit(n)
# 0. 初期状態 |+>^⊗4
qc.h(0) # 頂点1
qc.h(1) # 頂点2
qc.h(2) # 頂点3
qc.h(3) # 頂点4
# ===== レイヤー1: Cost(γ0) + Mixer(β0) =====
# Cost layer (gamma[0])
qc.cx(0, 1)
qc.rz(2 * gammas[0], 1)
qc.cx(0, 1)
qc.cx(0, 2)
qc.rz(2 * gammas[0], 2)
qc.cx(0, 2)
qc.cx(1, 2)
qc.rz(2 * gammas[0], 2)
qc.cx(1, 2)
qc.cx(1, 3)
qc.rz(2 * gammas[0], 3)
qc.cx(1, 3)
qc.cx(2, 3)
qc.rz(2 * gammas[0], 3)
qc.cx(2, 3)
# Mixer layer (beta[0])
qc.rx(2 * betas[0], 0)
qc.rx(2 * betas[0], 1)
qc.rx(2 * betas[0], 2)
qc.rx(2 * betas[0], 3)
# ===== レイヤー2: Cost(γ1) + Mixer(β1) =====
# Cost layer (gamma[1])
qc.cx(0, 1)
qc.rz(2 * gammas[1], 1)
qc.cx(0, 1)
qc.cx(0, 2)
qc.rz(2 * gammas[1], 2)
qc.cx(0, 2)
qc.cx(1, 2)
qc.rz(2 * gammas[1], 2)
qc.cx(1, 2)
qc.cx(1, 3)
qc.rz(2 * gammas[1], 3)
qc.cx(1, 3)
qc.cx(2, 3)
qc.rz(2 * gammas[1], 3)
qc.cx(2, 3)
# Mixer layer (beta[1])
qc.rx(2 * betas[1], 0)
qc.rx(2 * betas[1], 1)
qc.rx(2 * betas[1], 2)
qc.rx(2 * betas[1], 3)
print("\nQAOA ansatz circuit (p=2):")
print(qc.draw("text"))
##擬似量子コンピュータで計算
estimator = StatevectorEstimator()
optimizer = COBYLA(maxiter=100)
vqe = VQE(
estimator=estimator,
ansatz=qc,
optimizer=optimizer,
)
result = vqe.compute_minimum_eigenvalue(H)
print("\n=== VQE (hand-made QAOA ansatz) result ===")
print("最小エネルギー:", result.eigenvalue)
print("最適パラメータ (gamma[0], gamma[1], beta[0], beta[1]):")
print(result.optimal_parameters)
# 1) 最適パラメータを回路に代入して状態ベクトルを作る
bound = qc.assign_parameters(result.optimal_parameters, inplace=False)
state = Statevector.from_instruction(bound)
# 2) 各計算基底状態の確率を計算
probs = state.probabilities() # 長さ 2^4 = 16 の配列
# 最も確率の高い基底状態のインデックス
best_index = int(np.argmax(probs))
# 3) インデックス → 各量子ビットのビット値(q0〜q3)に変換
# bits[0] が qubit0(=頂点1), bits[1] が qubit1(=頂点2) ... という対応
bits = [ (best_index >> i) & 1 for i in range(n) ] # i=0がq0
print("\n=== Classical solution from QAOA state ===")
print("最も確率の高い基底状態 index:", best_index)
print("各qubitのビット列 [q0,q1,q2,q3] :", bits)
# 4) このビット列を「頂点のグループ分け」に解釈
# ここでは 0 = グループA, 1 = グループB とする
group_A = [i+1 for i, b in enumerate(bits) if b == 0] # 頂点番号は1始まりにして出力
group_B = [i+1 for i, b in enumerate(bits) if b == 1]
print("グループA (bit=0): 頂点", group_A)
print("グループB (bit=1): 頂点", group_B)
# 5) この切り方で何本の辺がカットされているかも計算
def cut_size(bits, edges):
"""bits: 長さ4の 0/1 リスト, edges: (i,j) のリスト"""
cnt = 0
for i, j in edges:
if bits[i] != bits[j]:
cnt += 1
return cnt
cut = cut_size(bits, edges)
print("この切り方で切れている辺の本数 (cut size):", cut)
まとめ:レベル3をマスターする価値
今回の手法は、Qiskitの自動生成機能を使わず、一からゲートを並べました。
このメリットは、「アルゴリズムのカスタマイズ性」にあります。
- 特定の制約を守らせるための特殊なミキサー(XYミキサーなど)への変更
- ハードウェアの特性に合わせたゲート配置の最適化
これらは、ブラックボックスの中身を知っている「レベル3」の理解があって初めて可能になります。数式 $e^{-i \gamma Z_i Z_j}$ が、コード上の cx-rz-cx という物理操作に重なって見えたとき、あなたは量子ネイティブな開発者への一歩を踏み出したと言えるでしょう。