Quibako Note

Questions unveiled invite brilliant and keen observation

和が条件を満たす2数の組をGroverのアルゴリズムで探す

量子コンピュータ

目的：Groverのアルゴリズムで遊ぶ。

Groverのアルゴリズムは、$ N $個の状態のうちから条件を満たすものを探索するアルゴリズム。 これを使って、和が条件「4以上の奇数」を満たす2数の組を検索する。 条件は、オラクルが簡単に作れそうなものを選んだ。

Groverのアルゴリズム

まず、Groverのアルゴリズムについてまとめる。

候補(検索対象)となる状態の集合を$ S $として、解(検索条件を満たすもの)となる状態の集合を$ \Omega $とする。 $ S $の要素の数を$ N $、$ \Omega $の要素の数を$ K $とする。 $ n $量子ビットの状態全体を候補とする場合は$ N=2 ^ n $ であるが、初期条件を絞りこめるなら、$ N<2 ^ n $でもOK。

次のように、候補となる状態すべてを同じ重み(同位相)で重ね合わせた状態を$ \ket{s} $、 解となる状態すべてを同じ重み(同位相)で重ね合わせた状態を$ \ket{\omega} $、 候補とはなったけど結果的に解ではなかった状態すべてを同じ重み(同位相)で重ね合わせた状態を$ \ket{\omega'} $とする。

$$ \begin{align*} \require{physics} \ket{s} &= \frac{1}{\sqrt{N}} \sum_{j\in S} \ket{j}, \\ \ket{\omega} &= \frac{1}{\sqrt{K}} \sum_{j\in \Omega} \ket{j}, \\ \ket{\omega'} &= \frac{1}{\sqrt{N-K}} \sum_{j\in S \backslash \Omega} \ket{j}. \end{align*} $$

このとき、$ \ket{s} $は、 $ \ket{\omega} $および$ \ket{\omega'} $の 2つの状態が張る(2次元)空間の状態になり、次のように書ける。

$$ \begin{align*} \ket{s} &= \sqrt{\frac{K}{N}} \frac{1}{\sqrt{K}} \sum_{j\in \Omega} \ket{j} + \sqrt{\frac{N-K}{N}} \frac{1}{\sqrt{N-K}} \sum_{j\in S\backslash \Omega} \ket{j} \\ &= \sin\theta \ket{\omega} + \cos\theta \ket{\omega'}, \\ \sin\theta &= \sqrt{\frac{K}{N}}, \\ \cos\theta &= \sqrt{\frac{N-K}{N}}. \end{align*} $$

状態$ \ket{s} $は、初期化のユニタリ演算子$ U _ s $を初期状態$ \ket{0 \cdots 0} $に作用させて構成する。

$$ \ket{s} = U_s \ket{0 \cdots 0} $$

状態$ \ket{s} $をスタートとして、 オラクル$ P _ \omega $で条件を満たす状態の符号を反転した後、 $ \ket{s} $で折り返す演算$ P _ s $を行うとGrover回転になる。

$ \ket{s} $で折り返す演算$ P _ s $(拡散演算子)は、 $ \ket{s} $を法線とする線で反転させ、全体に-1をかけても同様の結果になる。

$ \ket{s} $を法線とする線で反転させる演算を$ P _ {s _ \perp} $とし、 $ \ket{\omega'} $および$ \ket{\omega} $を基底とすると、

$$ \begin{align*} \ket{\omega'} &= \mqty[ 1 \\ 0] , & \ket{\omega} &= \mqty[ 0 \\ 1] , & \ket{s} &= \mqty[ \cos\theta \\ \sin\theta] \end{align*} $$
$$ \begin{align*} P_\omega &= \mqty[ 1 & 0 \\ 0 & -1 ] , \\ P_s &= 1 - 2 \ket{s}\bra{s} \\ &= \mqty[ 1 & 0 \\ 0 & 1 ] - 2 \mqty[ \cos\theta & \sin\theta] \mqty[ \cos\theta \\ \sin\theta] \\ &= \mqty[ 1-2\cos^2\theta & -2\sin\theta\cos\theta \\ -2\sin\theta\cos\theta & 1-2\sin^2\theta ] \end{align*} $$
$$ \begin{align*} P_{s_\perp} P_\omega &= \mqty[ 1-2\cos^2\theta & -2\sin\theta\cos\theta \\ -2\sin\theta\cos\theta & 1-2\sin^2\theta ] \mqty[ 1 & 0 \\ 0 & -1 ] \\ &= -\mqty[ 2\cos^2\theta-1 & -2\sin\theta\cos\theta \\ 2\sin\theta\cos\theta & 1-2\sin^2\theta ] \\ &= -\mqty[ \cos 2\theta & -\sin 2\theta \\ \sin 2\theta & \cos 2\theta ] \end{align*} $$

1回のGrover回転は$ 2\theta $回転なので、$ m $回繰り返すと次の状態が得られる。 $ P _ {s _ \perp} $を用いる場合に現れる位相因子$ (-1) ^ m $は、全体の状態の絶対位相に係る因子なので、ここでは無視してよい。

$$ \cos\qty((2m+1) \theta) \ket{\omega'} + \sin\qty((2m+1) \theta) \ket{\omega} $$

$ (2m+1) \theta = \pi/2 $、つまり次の条件を満たす$ m $の回数だけGrover回転を繰り返せば、条件を満たす状態$ \ket{\omega} $が得られる。

$$ \begin{align*} m & = \frac12 \qty(\frac{\pi}{2\theta} -1) \\ & \simeq \frac{\pi}{4} \sqrt{\frac{N}{K}} \end{align*} $$

ここで、$ P _ {s _ \perp} $は、初期化のユニタリ演算子$ U _ s $とマルチ制御$ Z $ゲート$ Z ^ {x _ 1 x _ 2 \cdots x _ {n-1}} _ n $とを使って、次のように構成できる。

$$ P_{s_\perp} = 1-2\ket{s}\bra{s} = U_s X^{\otimes n} Z^{x_1 x_2 \cdots x_{n-1}}_n X^{\otimes n} U_s^\dagger $$

実際、状態$ \ket{s} $に対して、

$$ \begin{align*} P_{s_\perp} \ket{s} &= U_s X^{\otimes n} Z^{x_1 x_2 \cdots x_{n-1}} _n X^{\otimes n} U_s^\dagger \ket{s} \\ &= U_s X^{\otimes n} Z^{x_1 x_2 \cdots x _{n-1}} _n X^{\otimes n} \ket{0\cdots 0} \\ &= U_s X^{\otimes n} Z^{x_1 x_2 \cdots x _{n-1}} _n \ket{1\cdots 1} \\ &= -U_s X^{\otimes n} \ket{1\cdots 1} \\ &= -U_s \ket{0\cdots 0} \\ &= - \ket{s} \end{align*} $$

また、状態$ \ket{s} $に直交する適当な状態$ \ket{s _ \perp} $について、 $ U _ s ^ \dagger\ket{s _ \perp} $は、$ \ket{0\cdots 0} $と次のように直交する。

$$ \begin{align*} 0 &= \braket{s|s_\perp} \\ &= \braket{0\cdots 0|U_s^\dagger|s_\perp} \end{align*} $$

このため、$ U _ s ^ \dagger\ket{s _ \perp} $は、 少なくとも1つのビットが1であるような状態、例えば$ \ket{0 \cdots 1 \cdots 0} $などの状態になる。

これを使って、

$$ \begin{align*} P_{s_\perp} \ket{s_\perp} &= U_s X^{\otimes n} \underbrace{ Z^{x_1 x_2 \cdots x_{n-1}}_n \underbrace{ X^{\otimes n} \underbrace{ U_s^\dagger \ket{s_\perp} }_{少なくとも1つのビットが1} }_{少なくとも1つのビットが0} }_{=X^{\otimes n} U_s^\dagger \ket{s_\perp}\ (ビット0があるので位相反転なし)} \\ &= \underbrace{U_s X^{\otimes n} X^{\otimes n} U_s^\dagger}_{=I} \ket{s_\perp} \\ &= \ket{s_\perp} \end{align*} $$

つまり、問題に応じた初期化のユニタリ演算子$ U _ s $とオラクル$ P _ \omega $があれば、 Groverのアルゴリズムの計算は作れる。

初期化演算子

検索対象とする2数の組は、2ビット整数として$ a={0,1,2,3} $および$ b={1,2,3} $とする。 ここから和$ a+b $が条件「4以上の奇数」を満たすような組を検索する。

$ a $については、2進数で書くと00,01,10,11で、2量子ビットの取りうる全状態に対応するので、 上位ビットおよび下位ビットに$ H $ゲートを描ければOK。

$ b $については、2進数で書くと01,10,11なので、 いったん上位ビットを0と1の重みが1:2(同位相)になるように重ね合わせ、 そこから上位ビットが1の状態10について下位ビットの異なる状態11, 10に同じ重み(同位相)で分ける。 これで00,11,10ができるので、下位ビットを反転すると01,10,11ができる。

10の振幅を11と10に分けるのは$ H $ゲート(制御$ H $ゲート)でできるので、 次のように1:2に分ける方を考える(振幅は $ : $)。

$$ \begin{align*} \sqrt{\frac{1}{3}} \ket{0} + \sqrt{\frac{2}{3}} \ket{1} \end{align*} $$

ブロッホ球を考えると、適当な大きさの比の重み(同位相)で重ね合わせるには、(0状態から)$ y $軸周りに回転させればよい。 回転角は、ほしい重みの比に応じて適当に決める。

$ y $軸周りの回転$ R _ y(\theta) $ゲートは、

$$ \begin{align*} R_y(\theta) \ket{0} &= \cos \frac{\theta}{2} \ket{0} + \sin \frac{\theta}{2} \ket{1} \end{align*} $$

なので、$ \cos \theta _ 3 /2= \sqrt{1/3} $となる$ \theta _ 3 $を使えば1:2に分けられることになる。

以上から、初期化は以下のようにゲートを作用させればよい。

オラクル

ここでは、補助ビットに和を計算したうえで、補助ビットの値に基づいて条件が満たされるかを確認する。

全加算器については、下記を参照。

そして、条件「4以上の奇数」は、和の3ビット整数に対しては、 最上位および最下位ビットが1であることを確認すればよい。

以上から、オラクルは以下のようにゲートを作用させればよい($ F $は全加算器)。

実装

Qiskit(v1.2.2)で実装(シミュレーション)する。

まず、初期化演算子と、拡散演算子を作る。

import numpy as np
from qiskit import QuantumCircuit, QuantumRegister

def initialize(to_gate=True):
  qra=QuantumRegister(2, 'a')
  qrb=QuantumRegister(2, 'b')
  qc=QuantumCircuit(qra,qrb)

  th3 = 2 * np.asin(np.sqrt(2./3))
  qc.h(qra)
  qc.ry(th3, qrb[1])
  qc.ch(qrb[1], qrb[0])
  qc.x(qrb[0])

  return qc.to_gate(label=r"$U_s$") if to_gate else qc

# draw circuit
initialize(to_gate=False).draw(output='mpl')

次に、拡散演算子$ P _ {s _ \perp} $は、

$$ P_{s_\perp} = U_s X^{\otimes n} Z^{x_1 x_2 \cdots x_{n-1}}_n X^{\otimes n} U_s^\dagger $$

なので、次のように定義。 qiskitにmulti-controlled-Zがない？ようなので、 multi-controlled-Xと$ Z=HXH $を使って作る。 (回路図では左から順に作用する(数式と逆)ので注意)

def diffuse(to_gate=True):
  qra=QuantumRegister(2, 'a')
  qrb=QuantumRegister(2, 'b')
  qc=QuantumCircuit(qra,qrb)

  qr = qra[:]+qrb[:]

  Us = initialize()
  Usdg = Us.inverse()
  Usdg.name = r"$U_s^\dagger$"

  qc.append(Usdg, qr)
  qc.x(qr)

  qc.h(qr[-1])
  qc.mcx(qr[:-1],qr[-1])
  qc.h(qr[-1])
    
  qc.x(qr)
  qc.append(Us, qr)

  return qc.to_gate(label=r"$P_{s_\perp}$") if to_gate else qc

# draw circuit
diffuse(to_gate=False).draw(output='mpl')

次に、オラクルを作る。 まずは、全加算器。

def full_adder(to_gate=True):
    
  qr=QuantumRegister(4)
  qc=QuantumCircuit(qr)

  qc.ccx(0,1,3)
  qc.cx(0,1)
  qc.ccx(1,2,3)
  qc.cx(1,2)
  qc.cx(0,1)

  return qc.to_gate(label=r"$F$") if to_gate else qc

# draw circuit
full_adder(to_gate=False).draw(output='mpl')

全加算器を使って、オラクルを作る。 このとき、答え(a=2,b=3とa=3,b=2)の情報を使っていない。

def oracle(to_gate=True):
  qra=QuantumRegister(2, 'a')
  qrb=QuantumRegister(2, 'b')
  qrs=QuantumRegister(3, 's')
  qrt=QuantumRegister(1, 'tgt')
  qc=QuantumCircuit(qra,qrb,qrs,qrt)

  fa = full_adder()
  fadg = fa.inverse()
  fadg.name = r"$F^\dagger$"

  qc.append(fa, [qra[0],qrb[0],qrs[0],qrs[1]])
  qc.append(fa, [qra[1],qrb[1],qrs[1],qrs[2]])
  qc.mcx([qrs[0], qrs[2]], qrt)
  # uncomputation
  qc.append(fadg, [qra[1],qrb[1],qrs[1],qrs[2]])
  qc.append(fadg, [qra[0],qrb[0],qrs[0],qrs[1]])

  return qc.to_gate(label=r"$P_\omega$") if to_gate else qc

# draw circuit
oracle(to_gate=False).draw(output='mpl')

これらから、Grover回転を作る。

def rotation(to_gate=True):
  qra=QuantumRegister(2, 'a')
  qrb=QuantumRegister(2, 'b')
  qrs=QuantumRegister(3, 's')
  qrt=QuantumRegister(1, 'tgt')
  qc=QuantumCircuit(qra,qrb,qrs,qrt)
   
  qc.append(oracle(), qra[:]+qrb[:]+qrs[:]+qrt[:])
  qc.append(diffuse(), qra[:]+qrb[:])
    
  return qc.to_gate(label=r"$G$") if to_gate else qc

# draw circuit
rotation(to_gate=False).draw(output='mpl')

ここで、$ a={0,1,2,3} $、$ b={1,2,3} $なので候補の数$ N $は12個、 条件「4以上の奇数」を満たすものの数$ K $はa=2,b=3とa=3,b=2の2つ。 このとき、$ \frac{\pi}{4} \sqrt{\frac{N}{K}} \simeq 1.9 $なので、 2回Grover回転させればOK。

複数回ゲートを作用させるのは、repeatメソッドでできる。

まとめると、以下のようになる。

from qiskit import ClassicalRegister

qra=QuantumRegister(2, 'a')
qrb=QuantumRegister(2, 'b')
qrs=QuantumRegister(3, 's')
qrt=QuantumRegister(1, 'tgt')
cra=ClassicalRegister(2, 'cl-a')
crb=ClassicalRegister(2, 'cl-b')

qc=QuantumCircuit(qra,qrb,qrs,qrt,cra,crb)

#initialization
qc.append(initialize(), qra[:]+qrb[:])

# make target qubit into |-> state 
qc.h(qrt)
qc.z(qrt)

# Grover rotation 
rot = rotation().repeat(2)
rot.name=r"$G^2$"
qc.append(rot, qra[:]+qrb[:]+qrs[:]+qrt[:])

# measurement
qc.measure(qra[:]+qrb[:], cra[:]+crb[:])

qc.draw(output='mpl')

実行してみる。

from qiskit_aer import AerSimulator
from qiskit import transpile
from qiskit.visualization import plot_histogram

backend = AerSimulator()
tqc = transpile(qc, backend)
job = backend.run(tqc, shots=1024)

rst = job.result()
counts = rst.get_counts()

# plot
plot_histogram(counts)

ということで、a=2,b=3 (11 10状態)とa=3,b=2 (10 11状態)を検索できた。 初期状態として作っていないb=0の状態は、ちゃんと現れていない(ノイズがないので)。