量子アニーリングのライブラリ PyQUBO を使ってみる #Python

約一年ぶりの投稿です。このブログの存在を忘れかけていました。
最近は機械学習やディープラーニングの延長、というわけでもないが量子アニーリングのモデルについてちょっと触ってみたので、せっかくなので備忘のためその導入部分についてつらつらと記述しておきます。書かないと絶対忘れますからね。そんなわけで今回の記事は基本的には自分のためのメモです。

量子コンピュータの 1 ビット

古典コンピュータのビット： 0 か 1 の状態
量子コンピュータのビット : 0 でも 1 でもありという状態。これに横の回転とも言える「位相」が加わる。

量子コンピュータの方式

演算方式	量子ゲート方式	量子アニーリング方式
用途	汎用	イジング問題に特化
量子ビット数	65	最近は 5,000 位まで増えた

量子を操る方法	超電導回路	イオン	半導体	光
量子の種類	電子	原子	電子	光子
エラー率	1% 以下	1% 以下	まだ高い	まだ高い
集積化	可能		可能
温度、設備	冷凍機	真空容器	冷凍機	室温で動作
量子ビット	不安定	安定

現在は超電導回路方式が主流。
冷凍機がめちゃめちゃでかい。

量子機械学習向けライブラリ

TensorFlow Quantum
https://www.tensorflow.org/quantum?hl=ja

古典コンピュータまたは Bristlecone プロセッサなどで動作。

PennyLane
https://github.com/PennyLaneAI/PennyLane

IBMQ, IonQ などの実機に接続して計算。

疑似量子コンピュータ

古典的なアニーリング法を GPU などで実装したもの。

量子アニーリングマシン D-Wave

2011 年頃いきなり登場。量子ゲートによる計算機ではなく量子アニーリング法による最適化計算に特化。
2018 年にアップグレードされた。これにより逆アニーリング（reverse annealing）と仮想グラフ（virtual graphs）が可能となった。

量子アニーリングの手順

1) 最小化したいハミルトニアンをイジングモデルの形で書く。

H_p = \sum_{i,j}J_{ij}\sigma_i^z\sigma_j^z + \sum_{i}h_i\sigma_i^z

2) 横磁場項を加えて時間とともに変換するハミルトニアンを構築。

IH = -\frac{A(s)}{2} ( \sum_i{{\hat\sigma_x^{(i)}}} )

FH = \frac{B(s)}{2} ( \sum_i{h_i \hat\sigma_z^{(i)}} + \sum_{i>j} J_{i,j} \hat\sigma_z^{(i)}\hat\sigma_z^{(j)} )

H_{ising} = IH + FH

3) 量子状態を 0-1 の重ね合わせ状態として用意し、ハミルトニアンを断熱的に変化させると、最終的に H_p を最小化する量子状態が観測される。

|\psi(t)> = T\int_0^tdt'e^{-iH(t')t'/h}|\psi(0)>

0 \leq t \leq_3 \tau

数式表現がいろいろおかしい気がするので突っ込みを募集します。

アニーリングマシンで組み合わせ最適化問題を解く手順

1) 整数計画問題として問題を定式化。
2) 制約式をペナルティ関数として目的関数に追加。
3) 整数をバイナリにエンコード。
4) ハミルトニアン数式の展開
5) 高次項の次数を削減。
6) 式を整理して係数から QUBO を作成。

すべての組み合わせ最適化問題は QUBO やイジングモデルとして定式化することが原理的には可能。

PyQUBO

PyQUBO
https://github.com/recruit-communications/pyqubo

C++ で実装された QUBO の Python 用ラッパーライブラリらしい。

$ wc -l site-packages/pyqubo/**/*.py
   14 site-packages/pyqubo/__init__.py
  800 site-packages/pyqubo/array.py
   19 site-packages/pyqubo/integer/__init__.py
   47 site-packages/pyqubo/integer/integer.py
   71 site-packages/pyqubo/integer/log_encoded_integer.py
   90 site-packages/pyqubo/integer/one_hot_enc_integer.py
  156 site-packages/pyqubo/integer/order_enc_integer.py
   68 site-packages/pyqubo/integer/unary_encoded_integer.py
  120 site-packages/pyqubo/logic.py
  145 site-packages/pyqubo/logical_constraint.py
   15 site-packages/pyqubo/package_info.py
    1 site-packages/pyqubo/utils/__init__.py
   42 site-packages/pyqubo/utils/asserts.py
  104 site-packages/pyqubo/utils/solver.py
 1692 合計

$ ls -la **/*cpp_pyqubo*
-rwxr-xr-x 1 ubuntu ubuntu 552816  8月 11 20:20 site-packages/cpp_pyqubo.cpython-37m-x86_64-linux-gnu.so

PyQUBO で QUBO を定式化し、コンパイルして辞書型に変換するまで

チュートリアルを見る。
https://openjij.github.io/OpenJijTutorial/build/html/ja/003-PyQUBO_2_OpenJij.html

# pyqubo.Array で変数を用意
x = Array.create(‘x’, shape=(N_VER,N_COLOR), vartype='BINARY')

# 第一項 (制約項)を定義。
H_A = Constraint(sum((1-sum(x[v,i] for i in range(1,N_COLOR)))**2 for v in range(N_VER)), label='HA')

# 第二項 (コスト、目的関数)を定義。
H_B = sum((-1+sum(x[v,i] for v in range (N_VER)))/2*(sum(x[v,i] for v in range (N_VER))) - sum(x[u,i]*x[v,i] for (u,v) in graph) for i in range (1,N_COLOR))

# ハミルトニアン全体を定義。
Q = H_A+H_B

# モデルをコンパイルしてソルバーに投げられるようにする。
model = Q.compile()
qubo, offset = model.to_qubo()

辞書型 QUBO を受け付けるソルバーに投げて最適化問題を解く

# D-Wave Ocean SDK の neal ライブラリを Simulated Annealing ソルバーとして利用
import neal
sampler = neal.SimulatedAnnealingSampler()
raw_solution = sampler.sample_qubo(qubo)

# 得られた結果を辞書型にデコードする。
decoded_sample = model.decode_sample(raw_solution.first.sample, vartype="BINARY")

# 解を表示するため.array(変数名, 要素番号)で要素の値を抽出する。
x_solution = {}
for i in range(N_VER):
    x_solution[i] = {}
    for j in range(1,N_COLOR):
        x_solution[i][j] = decoded_sample.array('x', (i, j))

# 制約が満たされているか確認。辞書が空なら OK 。
print(decoded_sample.constraints(only_broken=True))

今回はここまでです。