こんにちは、さっそく組合せ最適化問題をやってみたいと思います。
今回はオンライン勉強会の題材として準備してみました。
参考にするのはだいぶ前に流行ったVWの論文です。これを改造してより正解率を高めます。
Traffic flow optimization using a quantum annealer
Florian Neukart, Gabriele Compostella, Christian Seidel, David von Dollen, Sheir Yarkoni, Bob Parney
https://arxiv.org/abs/1708.01625
英語版のチュートリアルを作りました。
https://blueqat.com/yuichiro_minato2/b2fd6995-b7a1-4f2f-926f-9ed2bfb4fa6f
https://github.com/Blueqat/blueqat-tutorials/blob/master/tutorial/316_trafficflow.ipynb
問題設定
シンプルな問題を作りました。交差点が９つ、道路が１２あります。0の地点から8の地点に向かいます。
車が２台あり、それぞれカーナビで経路が３つ提案されます。同じ通りをなるべく通らないように、
混雑を発生させないように二台のカーナビの経路を提案するような組合せを計算します。
交差点と通りに名前を付けました。通りにはsがついてます。
経路提案
経路の提案は比較的適当でいいですが、２台ともこうしました。
車１への経路提案
Q0 : s0,s3,s8,s11
Q1 : s2,s7,s10,s11
Q2 : s0,s1,s4,s9
車２への経路提案
Q3 : s0,s3,s8,s11
Q4 : s2,s7,s10,s11
Q5 : s0,s1,s4,s9
Q0からQ5は量子ビットが６つ割り振られます。経路が選ばれれば１、選ばれなければ０です。
最適化問題
今回は最適化問題を解きます。これは、量子コンピュータで、ある式の最小値を出してくれるQAOAという計算を使います。このQAOAを実行するには式が必要です。問題を式の形に落とすことで、自動的に量子ビットの値が0なのか1なのか計算してくれます。
混雑度の計算
まずは混雑度を計算します。計算は、上記の経路に登場する道路の登場回数を足し合わせて二乗します。これによりきれいに混雑度が提案できます。
ちょっと画像小さいですが、上記のようになります。計算してみてください。
制約
これだけではありませんもう一つ条件があります。実は上の計算をQAOAに入れると全部0が出てしまいます。なぜなら条件を付けてないので、すべて0の時に最小値0を取るからです。そこで、条件を付けます。それが、
車１は必ず一つ経路を選ぶ。
車２は必ず一つ経路を選ぶ。
というものです。これはこれまで習ってきた量子もつれとQAOAの解探索方法を使います。
mixer
答えの探し方はいくつかありますが、mixerを使います。今回は最初の状態を100100みたいなすでに一つ答えが選ばれている状態からスタートし、Q0からQ2の中で値を入れ替えながら計算、同様にQ3からQ5を入れ替えながら計算をします。
もともとQ0からQ2は100からスタートし、1が一つ0が二つの状態です。ここから例えばQ0とQ1を入れ替えても、010のように1と0の個数は変わりません。古典計算では2-optという手法になります。これを使います。
使い方はQAOAにmixerの設定で指定をします。デフォルト状態では単に量子ビット単体を0か1かにする操作を続けて探しているだけになります。今回はXYmixerを使って、入れ替えていきます。
入れ替え方は、Q0とQ1を入れ替える、Q1とQ2を入れ替える、Q2とQ0を入れ替えるようにまんべんなく入れ替えていきます。Q3からQ5も同様です。
QAOA
では、さっそくプログラムの形にします。
from blueqat import vqe, Circuit
from blueqat.pauli import X, Y
from blueqat.pauli import qubo_bit as q
hamiltonian = 2*q(0)*q(1)+2*q(0)*q(2)+8*q(0)*q(3)+2*q(0)*q(4)+2*q(0)*q(5)+4*q(0)+2*q(1)*q(3)+7*q(1)*q(4)+3*q(1)+2*q(2)*q(3)+8*q(2)*q(5)+4*q(2)+2*q(3)*q(4)+2*q(3)*q(5)+4*q(3)+3*q(4)+4*q(5)
step = 4
#mixer and init state
mixer = 0.5*X[0]*X[1] + 0.5*Y[0]*Y[1] + 0.5*X[1]*X[2] + 0.5*Y[1]*Y[2] + 0.5*X[2]*X[0] + 0.5*Y[2]*Y[0]
mixer += 0.5*X[3]*X[4] + 0.5*Y[3]*Y[4] + 0.5*X[4]*X[5] + 0.5*Y[4]*Y[5] + 0.5*X[5]*X[3] + 0.5*Y[5]*Y[3]
init = Circuit(6).x[0].x[3]
result = vqe.Vqe(vqe.QaoaAnsatz(hamiltonian, step, init, mixer)).run()
result.most_common(10)
最初の３行はツールの読み込みです。
次のhamiltonianはコスト関数を入れます。stepは精度です。高いほうが精度の高い計算が理論的にできます。
いよいよmixerです。XYmixerは、(XiXj + YiYj)/2をします。
ですので、上記の入れ替える操作をプログラムします。
mixerを使う前に量子ビットの値を100100にしておく操作が.x[0].x[3]の部分です。
これで準備が整ったので、計算します。
result = vqe.Vqe(vqe.QaoaAnsatz(hamiltonian, step, init, mixer)).run()
実行すると、
(((0, 1, 0, 0, 0, 1), 0.37291765538289845),
((0, 0, 1, 0, 1, 0), 0.3729176553828979),
((0, 0, 1, 1, 0, 0), 0.08145920498348241),
((1, 0, 0, 0, 0, 1), 0.0814592049834823),
((1, 0, 0, 1, 0, 0), 0.02375992120568354),
((0, 1, 0, 0, 1, 0), 0.022701464389553148),
((0, 0, 1, 0, 0, 1), 0.017478879274962108),
((0, 1, 0, 1, 0, 0), 0.013653007198510678),
((1, 0, 0, 0, 1, 0), 0.013653007198510647),
((0, 0, 0, 0, 1, 0), 1.1219347803421538e-31))
いかがでしょうか。010001と001010がでました。この組み合わせは一切車同士が同じ道を通らない混雑が最小の組合せとなりました。
実機
IonQやRigettiへの量子コンピュータへの投げ方は以下の通りです。途中全部の計算をするのはものすごいコストがかかりますので、最適化後の計算をお勧めします。以上です。
from bqcloud import load_api
from blueqat import Circuit
from bqcloud import Device
api = load_api()
Circuit, Device, Number of shots, Group name(Optional) for IonQ
task_ionq = api.execute(result.circuit, Device.IonQDevice, 10)
Circuit, Device, Number of shots, Group name(Optional) for Rigetti
task_rigetti = api.execute(result.circuit, Device.Aspen11, 10)