量子基底エンコーディングから振幅エンコーディングへ

# 基底エンコーディングから振幅エンコーディングへの変換

教科書：[量子コンピューティング　基本アルゴリズムから量子機械学習まで]
(https://www.amazon.co.jp/%E9%87%8F%E5%AD%90%E3%82%B3%E3%83%B3%E3%83%94%E3%83%A5%E3%83%BC%E3%83%86%E3%82%A3%E3%83%B3%E3%82%B0-%E5%9F%BA%E6%9C%AC%E3%82%A2%E3%83%AB%E3%82%B4%E3%83%AA%E3%82%BA%E3%83%A0%E3%81%8B%E3%82%89%E9%87%8F%E5%AD%90%E6%A9%9F%E6%A2%B0%E5%AD%A6%E7%BF%92%E3%81%BE%E3%81%A7-%E6%83%85%E5%A0%B1%E5%87%A6%E7%90%86%E5%AD%A6%E4%BC%9A%E5%87%BA%E7%89%88%E5%A7%94%E5%93%A1%E4%BC%9A/dp/4274226212)

の、pp.47-8　**2.7.2 QRAMを使った振幅エンコーディング**　がちょっと難しかったので、フォローアップしてみます。

# QRAMによる基底エンコーディング

まず基底エンコーディングをします。これはQRAMというもの（アルゴリズムと理解したほうが良いかもしれません）を使えば出来ます。
QRAMの入力状態は、

$ \frac{1}{\sqrt{N}} \sum_{k=1}^{N} |i \rangle _{a} \otimes | 000...0 \rangle \otimes | 000...0 \rangle $

です。$i$がデータのindexを意味しています。

QRAMは次のように作用し、index $i$から対応するデータ$v_{i}$を引っ張ってきて量子状態の基底へ埋め込みます。

$ \frac{1}{\sqrt{N}} \sum_{k=1}^{N} |i \rangle_{a} \otimes | \text{b} (v_{i}) \rangle \otimes | 000...0 \rangle $

ここで、データは2進数化されて$b(v_{i})$なるbit列として埋め込まれています。

# 振幅エンコーディング

上記のままではちょっと使いにくい場合もあるので、振幅エンコーディングと呼ばれる
次のような状態へ変形することを考えたいです。

$ \propto \frac{1}{\sqrt{N}} \sum_{k=1}^{N} v_{i}|i \rangle _{a} $

振幅エンコーディングでは、文字通り、データ$v_{i}$が量子状態の確率振幅に埋め込まれます。
基底エンコーディングはデジタル、振幅エンコーディングはアナログなので、
この操作をQuantum digital to analog convert (QDAC)とも言うみたいです。

まず、必要な　逆計算(uncompute)の概念を説明しておきます。

## 逆計算(uncompute)

量子ゲート計算は古典NANDゲートを作れるので、古典で可能なあらゆる論理演算、ひいては数値演算が可能です。
量子算術演算とも呼ばれます。

例えば
$ |x \rangle \otimes |0 \rangle \to |x \rangle \otimes |f(x) \rangle $
という操作がいつも可能です。$f$は任意の関数です。(ただし$x$も$f(x)$も2進数化されます）

結果的にはこうなるのですが、量子ゲートでは大抵の場合”途中式”を経由してここへ至ります。
計算の過程で、補助ビットを付け足す必要があったり、興味のない補助ビットが状態遷移を起こしたりします。
そのような　ゴミ状態(garbage)　が出てくることが通常です。

$ |x \rangle \otimes |0 \rangle \otimes |0 \rangle \to |x \rangle \otimes |f(x) \rangle \otimes | garbage \rangle $

しかしgarbageと$f(x)$がエンタングルしていない限りは、garbageに逆演算を施して最初の状態に巻き戻しても$f(x)$には影響がありません。

$ |x \rangle \otimes |f(x) \rangle \otimes | garbage \rangle \to |x \rangle \otimes |f(x) \rangle \otimes | 0 \rangle $

このようにすることで、ゴミ状態の発生は（理論的には）問題になりません。

## 逆演算を使ってarccos(v)を埋め込む。

基底エンコーディングの状態
$ \frac{1}{\sqrt{N}} \sum_{k=1}^{N} |i \rangle_{a} \otimes | \text{b} (v_{i}) \rangle \otimes | 000...0 \rangle $

に対して、量子算術演算＋逆演算　によって、

$ \frac{1}{\sqrt{N}} \sum_{k=1}^{N} |i \rangle_{a} \otimes | \text{b} (v_{i}) \rangle \otimes | \text{b}(\arccos(v_{i})) \rangle $

と必ず出来ます。$v_{i}$から、その$\arccos(v_{i})$を計算したわけです。
狙いとしては、基底エンコーディングから確率振幅へ変換する操作が$\cos()$という関数を経由するので、ここで逆関数を取っておいたのです。

$\arccos(v_{i})$ は2進数 $\text{b}(\arccos(v_{i}))$になっています。2進数はどのように定義するかというと、

1. arccosの値域は$0$から$\pi$までである。
2. $\pi$で規格化することで$0$から$1$までになる。
3. そのような数は $b_{0}/2+b_{1}/4+...+b_{m-1}/2^{m-1}$ という2進数展開ができる。

このような方法で2進数化されてます。
ちゃんと書くと、

$ 0 \leq \arccos(x) \leq \pi$

$ 0 \leq \frac{\arccos(x)}{\pi} \leq 1$

$ 0 \leq \frac{\arccos(x)}{\pi} \leq \frac{1}{2}+\frac{1}{4}+\frac{1}{8}+....$

$ \arccos(x)/\pi \fallingdotseq \sum_{k=0}^{m-1} 2^{-k-1}b_{k}  $

$ \text{b}(x) :=  b_{0}b_{1}...b_{m-1}$

となっています。
後のために両辺$\pi$かけると、

$ \arccos(x) \fallingdotseq \pi \sum_{k=0}^{m-1} 2^{-k-1}b_{k}  $

となります。

## 回転ゲート

振幅エンコーディングにするためには回転ゲート$RY$を使います。RYゲートは以下のようなゲートです。

$\text{Ry}(\theta) |0 \rangle = \cos(\theta/2)|0 \rangle + \sin(\theta/2)|1 \rangle$

いま、$\arccos(v_{i}$を代入してみると、

$\text{Ry}(2\arccos(v_{i})) |0 \rangle = v_{i}|0 \rangle + \sqrt{1-v_{i}^{2}}|1 \rangle$

となり、$RY$ゲートによって $\cos(\arccos())$が作られて振幅エンコーディングが出来そうな気配です。
しかし実際には$\arccos(v_{i})$は2進数として埋め込まれてしまっていますので、この操作は直接実現できません。

なんとかするために、回転ゲートの　回転　という性質に着目します。
$ \text{Ry}(\alpha+\beta) = \text{Ry}(\alpha)\text{Ry}(\beta)$

という積和分解が成り立ちます。

いま、$\arccos(v_{i})$を代入してみると、

$\text{Ry}(2\arccos(v_{i})) = \text{Ry}(\pi \sum_{k=0}^{m-1} 2^{-k}b_{k}) = \Pi_{k=0}^{m-1} \text{Ry}(\pi 2^{-k}b_{k}) $

となり、$\arccos(v_{i})$の2進数化$b_{0}b_{1}...b_{m-1}$　の各桁に応じた角度での回転操作に分解できます。[^1]
かつ、$b_{k}$が0であれば$0 [deg.]$回転させ、1であれば$\pi2^{-k} [deg.]$回転させればいいことになります。
このような bitの0,1に応じた操作(if文相当) は、回転ゲートを制御ゲート化すれば実現できるはずです。

## 後処理

制御回転ゲートをした後の状態は、

$ \frac{1}{\sqrt{N}} \sum_{k=1}^{N} |i \rangle_{a} \otimes | \text{b} (v_{i}) \rangle \otimes | \text{b}(\arccos(v_{i})) \rangle \otimes (v_{i}|0 \rangle + \sqrt{1-v_{i}^{2}}|1 \rangle)$

となります。補助ビットを測定し、$0$となることを祈ります。

$ \propto \frac{1}{\sqrt{N}} \sum_{k=1}^{N} |i \rangle_{a} \otimes | \text{b} (v_{i}) \rangle \otimes | \text{b}(\arccos(v_{i})) \rangle \otimes (v_{i}|0 \rangle)$

こうなります。最後に逆計算によって、要らないものを巻き戻して

$ \propto \frac{1}{\sqrt{N}} \sum_{k=1}^{N} |i \rangle_{a} \otimes | 000...0 \rangle \otimes | 000...0 \rangle \otimes (v_{i}|0\rangle)$

とします。一部の量子ビットに着目すれば、

$ \propto \frac{1}{\sqrt{N}} \sum_{k=1}^{N} v_{i}|i \rangle _{a} $

となっていますよね。これで振幅エンコーディングに変換することができました。[^2]

## かんたんなデモ

qiskitでシンプルな場合を見てみます。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/878664/b64a8011-99af-2e78-49c8-c09a5f696be2.png)

この回路では、$q_{0,1}$が$\text{b}(\arccos(v_{i}))=b_{0}b_{1}$に相当します。

データの例は$v_{0} = v_{1} = 1/\sqrt{2}$、即ち$\arccos(v_{i}) = \pi/4$ です。$b_{0}b_{1}=01$になります。

$q_{2}$が回転させられる量子ビットです。
$q_{3}$がデータのindexです。

本来はindex 0 のデータも1のデータも両方埋め込むことで完成しますが、今回はindex 0 の場合だけ埋め込んでいます。[^3]

$b_{0}b_{1}=01$　なので、$q_{0}$から伸びている制御Ryゲートは何もしません。
$q_{1}$から伸びている制御Ryゲートは、角度$\pi/4$の回転を起こします。（$Ry(2\theta)$が$\theta$の回転であることに注意)
これで振幅エンコーディングが出来ています。

**追記**
補助ビット$q_{2}$を測定して0が出るまで祈る必要があります。その後に逆計算をします。

最後に逆計算で要らないところを巻き戻します。
ここでは、$q_{1}$を巻き戻しています。

一応コードは以下の通り。

```python
from qiskit import QuantumCircuit
import math

circ = QuantumCircuit(4)

# prepare \sum |b(arccos(vi)>|0>|0>
# (|01>|0>|0>)
circ.x(1)
circ.barrier()

# controlled rotation (control:q_1,2 target:q_3)
pos_control = [0,1]
pos_target = 2
for k in pos_control:
  circ.cry(math.pi*(2**(-k)),k,pos_target)
circ.barrier()

# uncompite
circ.x(1)

```


## (追記) index 0 と 1 の2データを同時埋め込み

$v_{0} = 1/\sqrt{2}, v_{1} = -1/\sqrt{2}$
$\arccos(v_{0})/\pi = 01, \arccos(v_{1})/\pi = 11$

とします。以下の状態が得られればOKです。
$|b_{0}b_{1} \rangle |ancilla \rangle |index \rangle = \frac{1}{\sqrt{2}}(|01\rangle|0\rangle|0\rangle - |11\rangle|0\rangle|1\rangle)$

```python
circ = QuantumCircuit(4,1)
#v_0 = 1/sqrt(2), v_1 = -1/sqrt(2)
# prepare \sum |b0b1>|ancilla>|index>
# 1/sqrt(2)(|01>|0>|0> + |11>|0>|1>)
# superposition of index
circ.h(3)
circ.barrier()

# for index = |1>
circ.cx(3,0)
circ.cx(3,1)
# for index = |0>
circ.x(3)
circ.cx(3,1)
circ.x(3)
circ.barrier()

# controlled rotation (control:q_1,2 target:q_3)
pos_control = [0,1]
pos_target = 2
for k in pos_control:
  circ.cry(math.pi*(2**(-k)),k,pos_target)
circ.barrier()

result = execute(circ,backend=Aer.get_backend('statevector_simulator'),shots=1).result()
```

![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/878664/617542f0-cac4-34b3-6e85-dab856e3783b.png)


本当は補助ビットの射影測定も含めて実装したいのですが、qiskitだとなかなか難しいです。
射影測定前の状態ベクトルを見てみます。

```python
import math
def print_Zbasis_expression(statevector):
    nqubit = int(math.log2(np.size(statevector)))
    for i in range(2**nqubit):
        print('+({:.2f})*|{number:0{width}b}>'.format(statevector[i],number=i,width=nqubit),end='\n')

State_out = result.get_statevector()
print('format : p*|index>|ancilla>|b1b0>')
print_Zbasis_expression(State_out)
```

`format : p*|index>|ancilla>|b1b0>
+(0.00+0.00j)*|0000>
+(0.00+0.00j)*|0001>
+(0.50-0.00j)*|0010>
+(0.00+0.00j)*|0011>
+(0.00+0.00j)*|0100>
+(0.00+0.00j)*|0101>
+(0.50+0.00j)*|0110>
+(0.00+0.00j)*|0111>
+(0.00+0.00j)*|1000>
+(0.00+0.00j)*|1001>
+(0.00+0.00j)*|1010>
+(-0.50-0.00j)*|1011>
+(0.00+0.00j)*|1100>
+(0.00+0.00j)*|1101>
+(0.00+0.00j)*|1110>
+(0.50+0.00j)*|1111>`

qiskitではビット順序が左右反転してますが、左から2番目が補助ビットです。
補助ビットが0であるものを抜き出すと、

$ \propto \frac{1}{2}(|00\rangle \otimes |1 \rangle \otimes |0 \rangle) - \frac{1}{2}(|10\rangle \otimes |1 \rangle \otimes |1 \rangle)$

となっています。実際の射影測定ではこれが大きさ1となるよう勝手に規格化されるので、

$  \frac{1}{\sqrt{2}}(|00\rangle \otimes |1 \rangle \otimes |0 \rangle) - \frac{1}{\sqrt{2}}(|10\rangle \otimes |1 \rangle \otimes |1 \rangle)$

という量子状態になっているはずです。
これは、(qiskit実装によるビット順序が左右反転していることを除いて）用意したかった振幅エンコーディング状態になっています。

$|b_{0}b_{1} \rangle |ancilla \rangle |index \rangle = \frac{1}{\sqrt{2}}(|01\rangle|0\rangle|0\rangle - |11\rangle|0\rangle|1\rangle)$

なお射影測定後の逆計算は、（私の理解では）$q_{2}$のRyゲート以外の全てを巻き戻せば良い・・・のですが
$q_{3}$の$H$ゲートの手前でやめてください。$H$ゲートを巻き戻してしまうと、index 0 と　1　が混ざってしまいます。[^4]


#　結論

- Ryゲートと積和分解でデジタルアナログ変換(QDAC)するというアイデアは高い一般性があると思った。

- 式変形は大変だが、アイデアはそれほど難しいものではない。


[^1]: 要は、回転角度を$n\pi/2$,$m\pi/4$,...という「だんだん細かくなる」回転の合成で表現したということです。
[^2]: ちなみに線形方程式を解くHHLアルゴリズムでもこの方法が使われています。HHLでは、まずQPEで固有値を基底エンコーディングで書き出します。これを振幅エンコーディングに変換すると同時に逆数を取っておくことで、「固有値で割る」という操作を実現しています。固有値で割ることと逆行列を掛けることは本質的に同じなので、これで方程式が解けるというわけです。
[^3]: index 0のデータも1のデータも埋め込みたい場合は、indexと回転角度をエンタングルで埋め込む。またはindexの部分にHゲートを掛けて重ねわせにし、indexを制御ビットとする制御振幅エンコーディング操作をする。
[^4]: 巻き戻してってどうやるの？と思うのかもしれませんが、掛けてきたゲートを逆順に掛けていけば戻っていきます。