CirqでDeustch-Jozsaのアルゴリズムを実装する

量子アルゴリズムの中で有名なアルゴリズムが3つあります。
1つ目は今回のDeustch-Jozsaのアルゴリズム
2つめはGroverの検索アルゴリズム。
3つめは皆さんご存じ素因数分解が高速で解けるというShorのアルゴリズム。
Deutch-Jozsaのアルゴリズムはこの中でも役に立たなさそうなやつNO.1なんですけど、どんな教科書にも必ず載っているのできっと大事なんでしょう！！これをPythonのパッケージのCirqってやつで実装していきます。CirqはWindowsだとコマンドプロンプトで

```
pip install cirq
```
やるとimportできるようになります。
#Deuch-Jozsaのアルゴリズムでできること

Deuch-Jozsaのアルゴリズムでは与えられた関数が「定数関数」もしくは「バランス関数」であるかを判定することができます。
ここで関数f(x)は
$$
入力：0\leq x\leq2^n-1\in N
$$
$$
出力：0 または 1
$$
のような関数です。
定数関数は任意のxに対して$f(x)=0$ となるような関数。
もしくは、任意のxに対して$f(x)=1$ となるような関数です。(xの範囲は入力の範囲）

バランス関数とは入力のxの丁度半分が$f(x)=0$となり、残りの半分が$f(x)=1$となるような関数です。つまりxをランダムに選ぶと$f(x)$は1/2の確率で0となり、1/2の確率で1となるような関数です。

ある関数が与えられ、それが「定数関数」か「バランス関数」かを当ててほしいと言われたら何回この関数を使えばよいでしょうか？古典コンピュータならx=0から1ずつ増やして関数を使っていけば最短で2回、最長で2^(n-1)+1回、関数を使わなければなりません。
n=4のとき
最短例:f(0)=1 f(1)=1 ->バランス関数
最長例:f(0)=0 f(1)=0,...,f(7)=0,f(8)=1 ->バランス関数

量子コンピュータを使うと関数を1回使うだけで与えられた関数が「定数関数」か「バランス関数」かを判定できるというのです。
#Deutch-Jozsaのアルゴリズム
ここから本題のアルゴリズムの解説に入ります。
まずはn個のqubitと1個の補助qubit（ancilla qubit)を用意します。(慣例的に初期量子ビットはすべて0になっています。

```math
|0^n\rangle |0\rangle　=|0\rangle...|0\rangle |1\rangle (0がnこ）
```
n+1ビット目（補助qubit)にXゲート(ビット反転するゲート）を作用させます
Xゲートは以下の働き

```math
 X|0\rangle=|1\rangle\quad X|1\rangle=|0\rangle
```
すると

```math
|0^n\rangle |1\rangle
```
次にn+1ビット全体にHadamardゲートを作用させます。
Hゲートは以下の働き

```math
H|0\rangle=\frac{1}{\sqrt{2}}(|0\rangle +|1\rangle)\\
H|1\rangle=\frac{1}{\sqrt{2}}(|0\rangle -|1\rangle)
``` 
すると

```math
\frac{1}{\sqrt{2^{n+1}}}\sum_{x_1=0}^{1}
\sum_{x_2=0}^{1}...\sum_{x_n=0}^{1}|x_1\rangle |x_2\rangle...|x_n\rangle (|0\rangle -|1\rangle)\\
=\frac{1}{\sqrt{2N}}\sum_{x=0}^{2^{n-1}}|x\rangle (|0\rangle -|1\rangle) \qquad ここでN=2^{n}とおいた
```
1行目から2行目への変形でビットの2進数表記をまとめて10進数表記にしていることに注意。
例

```math
|1\rangle |0\rangle |0\rangle |1\rangle =|9\rangle  
```

つぎに以下の性質を満たす $U_f$を考える

```math
U_f|x\rangle |y\rangle=|x\rangle |y\oplus f(x)\rangle\\
ここで\oplus は足してmod2をする記号
```
このような$U_f$ゲートを1回使うと関数$f(x)$を1回使っていることになる。問題なのはこのゲートをどうやって実装するかという点だが、後回しにする。この$U_f$ゲートをn+1ビット全体に作用させると

```math
\frac{1}{\sqrt{2N}}\sum_{x=0}^{2^{n-1}}|x\rangle (|0\oplus f(x)\rangle -|1\oplus f(x)\rangle) \\
=\frac{1}{\sqrt{2N}}\sum_{x=0}^{2^{n-1}}(-1)^{f(x)}|x\rangle (|0\rangle -|1\rangle) \qquad 
```
ここで補助キュビットは用済みなので、以下からn個のqubitのみで考えます

```math
\frac{1}{\sqrt{N}}\sum_{x=0}^{2^{n-1}}(-1)^{f(x)}|x\rangle
```
n個のqubitにそれぞれHadamardゲートを通すと

```math
\frac{1}{N}\sum_{x=0}^{2^{n-1}}\sum_{y=0}^{2^{n-1}}(-1)^{f(x)}(-1)^{\langle x,y\rangle}|y \rangle\\
ここでx,yを二進数でn桁で表し、i番目をx_i,y_iとして、\\
\langle x,y\rangle=\sum_{i=1}^{n}x_iy_i
```
(ここの変形がこのアルゴリズムの山場の1つです。以下の補題が理解できれば変形できます。頑張りましょう！！）
補題１

```math
任意のx(0\leq x\leq2^n-1\in N)に対して\\
\begin{eqnarray*}
H^{\otimes n}|x\rangle
&=& H|x_1\rangle H|x_2\rangle...H|x_n\rangle\\
&=& \frac{1}{\sqrt{2^n}}\sum_{y_1=0}^{1}(-1)^{x_1y_1}|y_1\rangle\sum_{y_2=0}^{1}(-1)^{x_2y_2}|y_2\rangle ...\sum_{y_n=0}^{1}(-1)^{x_ny_n}|y_n\rangle\\
&=&\frac{1}{\sqrt{2^n}}\sum_{y=0}^{2^n-1}(-1)^{(x_1y_1+x_2y_2+...+x_ny_n)}|y\rangle\\
&=&\frac{1}{\sqrt{2^n}}\sum_{y=0}^{2^n-1}(-1)^{\langle x,y\rangle}|y\rangle\\
\end{eqnarray*}
```

ここでnこのqubitを観測して$|0^n\rangle$ (つまり上の状態の$|y\rangle =|0\rangle$)が得られる確率を考えます。観測の確率はその振幅（係数）の二乗で与えられるので

```math
(|0^n\rangle が観測される確率)=\Bigg|\frac{1}{N}\sum_{x=0}^{2^{n-1}}(-1)^{f(x)}\Bigg|^2
```
ここで$f(x)$が定数関数だった場合、$\sum_{x=0}^{2^{n-1}}(-1)^{f(x)}=\pm2^n$ なので

```math
\begin{eqnarray*}
(|0^n\rangle が観測される確率)&=&\Bigg|\frac{1}{N}\sum_{x=0}^{2^{n-1}}(-1)^{f(x)}\Bigg|^2\\
&=&1
\end{eqnarray*}
```
つまり観測すると$|0^n\rangle $の状態のみで、それ以外が得られることはありません。

$f(x)$がバランス関数だった場合、$\sum_{x=0}^{2^{n-1}}(-1)^{f(x)}=0$なので

```math
\begin{eqnarray*}
(|0^n\rangle が観測される確率)&=&\Bigg|\frac{1}{N}\sum_{x=0}^{2^{n-1}}(-1)^{f(x)}\Bigg|^2\\
&=&0
\end{eqnarray*}
```
つまり観測すると$|0^n\rangle $の以外の状態が得られ、$|0^n\rangle $は絶対に観測されません。

まとめると、上記の操作を行い観測結果が

- $|0^n\rangle $　　なら$f(x)$は定数関数
- $|0^n\rangle $以外なら$f(x)$はバランス関数

このようにして、$f(x)$($U_f$ゲート)を一回使うだけで$f(x)$が定数関数かバランス関数かを判定することができます。
#アルゴリズムの概要
流れだけまとめていきます
(1) n+1ビット目にXゲートを作用する
(2) n+1ビット全体にHadamardゲートを作用させる
(3) n+1ビット全体に$U_f$ゲートを作用させる
(4) 最初のnビットにHadamardゲートを作用させる
(5) 最初のnビットを観測して
　　　$|0^n\rangle$　　なら定数関数
　　　$|0^n\rangle$以外ならバランス関数
　　と結論づける
#n=2で実装してみる
実装の肝となるのが手順３の$U_f$ゲートで、これをはじめから一般のnでやるのは厳しいのでまずはn=2で実装してみます。
手順1,2,4,5はそのままなので下に挙げます。

```python
import cirq
import numpy as np
import matplotlib.pyplot as plt

qubit = cirq.GridQubit(0,0)
n_qubits=2
circuit=cirq.Circuit()
qubits = [cirq.GridQubit(i,0) for i in range(n_qubits)]
ancilla = cirq.GridQubit(n_qubits,0)

def scheme1(circuit,qubits,ancilla):
    circuit.append([
    cirq.X(ancilla)
    ])
    return circuit

def scheme2(circuit,qubits,ancilla):
    circuit.append([
        cirq.H.on_each(qubits),
        cirq.H(ancilla)
    ])
    return circuit


def scheme4(circuit,qubits,ancilla):
    circuit.append([
        cirq.H.on_each(qubits),
    ])
    return circuit

def scheme5(circuit,qubits,ancilla):
    circuit.append([
        cirq.measure(*qubits,key="m")
    ])
    print("Circuit:")
    print(circuit)

    simulator=cirq.Simulator()
    result=simulator.run(circuit,repetitions=50000)
    print("Results:")
    print(result.histogram(key="m"))

    list = get_result_list(result)

    plt.figure(figsize=(4,3))
    plt.bar(list[:,0],list[:,1])
    return circuit

```
関数$f(x)$は$2^n$桁のビットで表し、左からi番目は$f(i-1)$を表します。
例
"1001"は$f(0)=1$,$f(1)=0$,$f(2)=0$,$f(3)=1$を表す。
以下を満たす$U_f$ゲートは

```math
U_f|x\rangle |y\rangle=|x\rangle |y\oplus f(x)\rangle\\
```
目標のビット以外の場所以外が変更されるのを防ぐため、エンタングルドしている$|x\rangle$を$f(x)=1$なら、Xゲートで$|1\rangle |1\rangle$に変えてからTOFFOLIゲートで$|y\rangle$をビット反転して、Xゲートで$|1\rangle |1\rangle$を$|x\rangle$に戻す、$f(x)=0$なら何もしないという風にしています。以下が手順3のコードです

```python
def seiri(i,n):
    X=str(bin(i)[2:])
    X=X.zfill(n)
    index=[]
    for i in range(n):
        if X[i]=="0":
            index.append(i)
    return index

def scheme3(circuit,qubits,ancilla):
    for i in range(len(bits)):
        
        p=seiri(i,n_qubits)    
        if  bits[i]=="1":
            for j in p:
                circuit.append([
                    cirq.X(qubits[j])
                ])
            circuit.append([
                cirq.TOFFOLI(qubits[0],qubits[1],ancilla)
            ])
            for j in p[::-1]: 
                circuit.append([
                    cirq.X(qubits[j])
                ])
    return circuit
```
実際にbits="1001"として実行してみると

```python
bits="1001"

circuit=scheme1(circuit,qubits,ancilla)
circuit=scheme2(circuit,qubits,ancilla)
circuit=scheme3(circuit,qubits,ancilla)
circuit=scheme4(circuit,qubits,ancilla)
circuit=scheme5(circuit,qubits,ancilla)
```
結果は
Circuit:
(0, 0): ───────H───X───────@───────X───X───@───X───H───M('m')───
                           │               │           │
(1, 0): ───────H───────X───@───X───────────@───────H───M────────
                           │               │
(2, 0): ───X───H───────────X───────────────X────────────────────
Results:
Counter({2: 50000})
![image.png](https://qiita-image-store.s3.amazonaws.com/0/319467/6d6632a6-ad46-f9e7-c9cc-5c056e15f929.png)

観測結果は$|0^2\rangle$ではないので$f(x)$はバランス関数であると推定でき、実際にそうなっています。
次にbits="1111"で試して結果は
Circuit:
(0, 0): ───────H───X───────@───────X───X───@───X───────@───────@───H───M('m')───
                           │               │           │       │       │
(1, 0): ───────H───────X───@───X───────────@───────X───@───X───@───H───M────────
                           │               │           │       │
(2, 0): ───X───H───────────X───────────────X───────────X───────X────────────────
Results:
Counter({0: 50000})
![image.png](https://qiita-image-store.s3.amazonaws.com/0/319467/45d9f89c-f9db-f90e-be22-1eff910c653f.png)
観測結果は$|0^2\rangle$のみなので$f(x)$は定数関数であると推定でき、実際にそうなっています。

ここで関数$f(x)$を定数関数でもバランス関数でもないものにしたらどうなるのでしょう？
bits="1011"とした結果は
Circuit:
(0, 0): ───────H───X───────@───────X───────@───────@───H───M('m')───
                           │               │       │       │
(1, 0): ───────H───────X───@───X───────X───@───X───@───H───M────────
                           │               │       │
(2, 0): ───X───H───────────X───────────────X───────X────────────────
Results:
Counter({2: 12651, 0: 12520, 3: 12511, 1: 12318})
![image.png](https://qiita-image-store.s3.amazonaws.com/0/319467/7560b32a-f8e8-550f-487e-174709fb41b5.png)
観測結果からは何も得られないことがわかります。（強いて言うならこの関数は定数関数でもバランス関数でもないことがわかります）

#終わりに
今回はn=2までだったが、$|x\rangle$を一時的に$|11\rangle$に変換してCNOTゲートやTOFFOLIゲートで目標を変換する技は簡単に一般化できそう。ただ$U_f$ゲートを使うのは1回だがその中でf(x)を何度も使ってるのはいいのだろうか...という気持ちです（$U_f$ゲートをうまく実装すればf(x)の参照が1回で済むようにできるのだろうか）

#参考文献

- 量子コンピュータ入門[第２版]　日本評論社
- 量子情報科学入門　共立出版