量子フーリエ変換

はじめに

　先の回では「量子もつれ」という現象を紹介し、blueqatという名前のフレームワークを用いて量子もつれをシミュレートした。今回は、離散フーリエ変換を量子コンピュータ（量子回路）で行うアルゴリズム「量子フーリエ変換」を紹介し、先と同じようにblueqatでシミュレートする。

離散フーリエ変換

　離散値 $x_j\;(j=0,\cdots,N-1)$ の離散フーリエ変換は次式で定義される。

(1) $\begin{equation*} y_k=\frac{1}{\sqrt{N}}\sum_{j=0}^{N-1}x_j\exp{\left(i\frac{2\pi kj}{N}\right)} \end{equation*}$

ここで、 $k=0,\cdots,N-1$ であり、 $i$ は虚数単位である。いま

(2) $\begin{equation*} w_{kj}=\frac{1}{\sqrt{N}}\exp{\left(i\frac{2\pi kj}{N}\right)} \end{equation*}$

を定義すると、式(1)は

(3) $\begin{equation*} y_k=\sum_{j=0}^{N-1}w_{kj}x_j \end{equation*}$

と書くことができる。 $w_{kj}$ を成分に持つ行列 $W$ を導入すると、式(3)は

(4) $\begin{equation*} y=Wx \end{equation*}$

と書ける。ここで、 $y\in\mathbb{R}^N,x\in\mathbb{R}^N,W\in\mathbb{R}^{N\times N}$ である。行列 $W$ はユニタリー行列になっていることに注意する。すなわち

(5) $\begin{equation*} W^{\dagger}W=WW^{\dagger}=I \end{equation*}$

が成り立つ。ここで、 $W^{\dagger}$ は $W$ のエルミート共役（転置して複素共役をとる）を表す。 $I$ は単位行列である。ユニタリー性を満たさない演算は量子回路上に載せることができないので、離散フーリエ変換のユニタリー性（式(5)）は重要である。

振幅エンコーディング

　 $x$ や $y$ のようなベクトルをどのように量子コンピュータに載せるのかを考える。これまでのブログ（1,2）で1量子ビットと2量子ビットを扱った。任意の量子状態は、1量子ビットの場合は $|0\rangle,|1\rangle$ の2つの基底状態の重ね合わせで、2量子ビットの場合は $|00\rangle,|01\rangle,|10\rangle,|11\rangle$ の4つの基底状態の重ね合わせで表現されることを見た。ここで、0や1の並びは0以上の整数の2進数表記に相当し、2量子ビットの場合の00,01,10,11は、それぞれ、10進数表記で0,1,2,3に相当する。いま、4次元の単位ベクトル $x=(x_0,x_1,x_2,x_3)$ が与えられたとき、各成分を係数に持つ次の量子状態を考える。

(6) $\begin{eqnarray*} |x\rangle &\equiv& \sum_{j=0}^{3}x_j|j\rangle \\ &=&x_0|0\rangle+x_1|1\rangle+x_2|2\rangle+x_3|3\rangle \\ &=&x_0|00\rangle+x_1|01\rangle+x_2|10\rangle+x_3|11\rangle \end{eqnarray*}$

$x$ は単位ベクトルとしたから量子状態が満たすべき条件（規格化条件： $\sum_{j=0}^{3}|x_j|^2=1$ ）は満たされている。つまり、4次元ベクトルは2量子ビットを用いて上のように表すことができる。これを $n$ ビットに拡張する。

(7) $\begin{equation*} |x\rangle =\sum_{j=0}^{N-1}x_j|j\rangle \end{equation*}$

ただし、 $N=2^n$ である。式(7)は $N$ 次元ベクトルを表す量子状態である。 $n=3$ なら

(8) $\begin{eqnarray*} |x\rangle &=& \sum_{j=0}^{2^3-1}x_j|j\rangle \\ &=& x_0|0\rangle+x_1|1\rangle+\cdots+x_7|7\rangle \\ &=& x_0|000\rangle+x_0|001\rangle+\cdots+x_7|111\rangle \\ \end{eqnarray*}$

となり、8（ $=2^3$ ）次元ベクトルを表すことができる。このようにベクトルの成分を量子状態の係数（確率振幅）に埋め込むことを振幅エンコーディングと呼ぶ。どのようにしてこの状態を作り出すのかという話は少々複雑なので今回は言及しない。

量子フーリエ変換

　振幅エンコーディングにより、ベクトル $y$ を量子状態に変換する。

(9) $\begin{equation*} |y\rangle =\sum_{k=0}^{N-1}y_k|k\rangle \end{equation*}$

ここに式(3)を代入すると

(10) $\begin{eqnarray*} |y\rangle &=& \sum_{k=0}^{N-1}y_k|k\rangle \\ &=& \sum_{k=0}^{N-1}\left[\sum_{j=0}^{N-1}w_{kj}x_j\right]|k\rangle \\ &=& \sum_{j=0}^{N-1} x_j \left[\sum_{k=0}^{N-1}w_{kj}|k\rangle\right] \end{eqnarray*}$

を得る。ところで $x$ の量子状態は

(11) $\begin{equation*} |x\rangle =\sum_{j=0}^{N-1}x_j|j\rangle \end{equation*}$

であった。この2式を見比べると、 $|x\rangle$ から $|y\rangle$ に変換するには次の変換を行えばよいことが分かる。

(12) $\begin{equation*} |j\rangle \rightarrow \sum_{k=0}^{N-1}w_{kj}|k\rangle \end{equation*}$

すなわち、この変換を $F_n$ と書くことにすると

(13) $\begin{equation*} F_n|j\rangle =\sum_{k=0}^{N-1}w_{kj}|k\rangle \end{equation*}$

であるから

(14) $\begin{eqnarray*} F_n|x\rangle &=& F_n\sum_{j=0}^{N-1}x_j|j\rangle \\ &=& \sum_{j=0}^{N-1}x_j\left(F_n|j\rangle\right) \\ &=& \sum_{j=0}^{N-1}x_j\left[\sum_{k=0}^{N-1}w_{kj}|k\rangle\right] \\ &=& \sum_{k=0}^{N-1}\left[\sum_{j=0}^{N-1}w_{kj}x_j\right]|k\rangle \\ &=& \sum_{k=0}^{N-1}y_k|k\rangle \\ &=& |y\rangle \end{eqnarray*}$

となる。変換 $F_n$ を量子フーリエ変換（Quantum Fourier Transform：QFT）と呼ぶ。上式を見ると、 $F_n|x\rangle$ の結果からフーリエ係数 $y_k$ を求めるには $|k\rangle$ の係数を取り出す必要があることに注意しなければならない。これを行う手間についてはまた後で触れることにする。

QFTアルゴリズム（ステップ1）

　QFTを量子回路に載せる手順（QFTアルゴリズム）を説明する。まず最初に、QFTの式を量子回路に載せやすい形に変形する必要がある。ステップ1ではその変形の過程を示す。古典回路と同様に量子回路の場合も2進数の言葉で書き直す必要がある。

　式(12)の右辺に式(2)を代入し、 $N$ を $2^n$ に戻すと

(15) $\begin{equation*} \sum_{k=0}^{2^n-1}\frac{1}{\sqrt{2^n}}\exp{\left(i\frac{2\pi kj}{2^n}\right)} |k\rangle \end{equation*}$

となる。ところで、 $k$ の $n$ ビット2進数表記が $k_{0}k_{1}\cdots k_{n-1}$ となるとき（例えば、3の2進数表記は11なので $k_0=1,k_1=1$ ）、次式が成り立つ。

(16) $\begin{equation*} k=\sum_{l=0}^{n-1}2^{n-l-1} k_{l} \end{equation*}$

$k_l$ は0か1を取る変数である。状態 $|k\rangle$ を2進数で書き直すと

(17) $\begin{eqnarray*} |k\rangle &=& |k_{0}k_{1}\cdots k_{n-1}\rangle \\ &=& |k_{0}\rangle|k_{1}\rangle\cdots|k_{n-1}\rangle \end{eqnarray*}$

となる。表記を簡単にするため直積記号 $\otimes$ は省略した。さらに $k$ についての和は次のように書き直すことができる。

(18) $\begin{equation*} \sum_{k=0}^{2^n-1}\rightarrow\sum_{k_{0}=0}^{1}\;\sum_{k_{1}=0}^{1}\cdots\sum_{k_{n-1}=0}^{1} \end{equation*}$

以上のことを考慮すると式(15)は

(19) $\begin{eqnarray*} &&\sum_{k=0}^{2^n-1}\frac{1}{\sqrt{2^n}}\exp{\left(i\frac{2\pi kj}{2^n}\right)} |k\rangle \\ &=&\frac{1}{\sqrt{2^n}} \sum_{k_{0}=0}^{1}\;\sum_{k_{1}=0}^{1}\cdots\sum_{k_{n-1}=0}^{1}\exp{\left(i\frac{2\pi j\sum_{l=0}^{n-1}2^{n-l-1} k_l}{2^n}\right)} |k_{0}\rangle|k_{1}\rangle\cdots|k_{n-1}\rangle \\ &=&\frac{1}{\sqrt{2^n}} \sum_{k_{0}=0}^{1}\;\sum_{k_{1}=0}^{1}\cdots\sum_{k_{n-1}=0}^{1}\exp{\left(i2\pi j\sum_{l=0}^{n-1}2^{-l-1} k_l\right)} |k_{0}\rangle|k_{1}\rangle\cdots|k_{n-1}\rangle \\ &=&\frac{1}{\sqrt{2^n}} \prod_{l=0}^{n-1}\left[\sum_{k_{l}=0}^{1}\exp{\left(i2\pi j2^{-l-1} k_l\right)} |k_{l}\rangle\right] \\ &=&\frac{1}{\sqrt{2^n}} \prod_{l=0}^{n-1}\left[|0\rangle+\exp{\left(i2\pi j2^{-l-1}\right)} |1\rangle\right] \end{eqnarray*}$

となる。ここで、乗積を展開するときは $l=0$ の項が一番左側に来ると約束する。 $j$ についても2進数表記を行うと

(20) $\begin{equation*} j=\sum_{l=0}^{n-1}2^{n-l-1} j_l \end{equation*}$

従って

(21) $\begin{eqnarray*} j2^{-l-1} &=& 2^{-l-1}\left(2^{n-1}j_{0}+\cdots+2^{l+1}j_{n-l-2}+2^{l}j_{n-l-1}+\cdots+2^0j_{n-1}\right) \\ &=& 2^{n-l-2}j_{0}+\cdots+2^{0}j_{n-l-2}+2^{-1}j_{n-l-1}+\cdots+2^{-l-1}j_{n-1} \end{eqnarray*}$

となる。右辺の先頭の $n-l-1$ 個の和は正の整数である。これを整数 $M=2^{n-l-2}j_{0}+\cdots+2^{0}j_{n-l-2}$ に置き換え、式(21)を指数関数の肩に載せると

(22) $\begin{eqnarray*} \exp{\left(i2\pi j2^{-l-1}\right)} &=& \exp{\left(i2\pi\left[M+2^{-1}j_{n-l-1}+\cdots+2^{-l-1}j_{n-1}\right]\right)}　\\ &=& \exp{\left(i2\pi\left[2^{-1}j_{n-l-1}+\cdots+2^{-l-1}j_{n-1}\right]\right)} \end{eqnarray*}$

となる。ここで、 $\exp{\left(i2\pi M\right)}=1$ を用いた。浮動小数点の2進数表記を用いると

(23) $\begin{equation*} 2^{-1}j_{n-l-1}+\cdots+2^{-l-1}j_{n-1}=0.j_{n-l-1}j_{n-l}\cdots j_{n-1} \end{equation*}$

であるから、結局

(24) $\begin{equation*} \exp{\left(i2\pi j2^{-l-1}\right)}= \exp{\left(i2\pi\;0.j_{n-l-1}j_{n-l}\cdots j_{n-1}\right)} \end{equation*}$

を得る。以上の結果をまとめると、式(19)は最終的に次のように変形される。

(25) $\begin{equation*} &&\sum_{k=0}^{2^n-1}\frac{1}{\sqrt{2^n}}\exp{\left(i\frac{2\pi kj}{2^n}\right)} |k\rangle= \frac{1}{\sqrt{2^n}} \prod_{l=0}^{n-1}\left[|0\rangle+\exp{\left(i2\pi0.j_{n-l-1}j_{n-l}\cdots j_{n-1} \right)} |1\rangle\right] \end{equation*}$

改めて量子フーリエ変換を書くと

(26) $\begin{eqnarray*} F_n|j\rangle &=& \frac{1}{\sqrt{2^n}} \prod_{l=0}^{n-1}\left[|0\rangle+\exp{\left(i2\pi0.j_{n-l-1}j_{n-l}\cdots j_{n-1} \right)}|1\rangle\right] \\ &=& \frac{1}{\sqrt{2^n}}\left(|0\rangle+e^{i2\pi0.j_{n-1}}|1\rangle\right)\left(|0\rangle+e^{i2\pi0.j_{n-2}j_{n-1}}|1\rangle\right)\cdots\left(|0\rangle+e^{i2\pi 0.j_{0}j_{1}\cdots j_{n-1}}|1\rangle\right) \end{eqnarray*}$

となる。乗積を展開するときは $l=0$ の項が一番左側に来るように配置する。あとはこの式を量子回路で表現すれば良い。その手順を次に示す。

QFTアルゴリズム（ステップ2）

　いきなり $n$ 量子ビットの場合を考えると難しいので、最初に1量子ビット、2量子ビットの場合を考える。

QFTアルゴリズム：1量子ビットの場合

　 $n=1$ のときQFTは次式になる。

(27) $\begin{equation*} F_1|j\rangle &=& \frac{1}{\sqrt{2}}\left(|0\rangle+e^{i2\pi0.j_0}|1\rangle\right) \end{equation*}$

$j$ の2進数表記は $j_0$ なので

(28) $\begin{equation*} F_1|j_0\rangle = \frac{1}{\sqrt{2}}\left(|0\rangle+e^{i2\pi0.j_0}|1\rangle\right) \end{equation*}$

である。 $j_0=0,1$ を代入すると

(29) $\begin{eqnarray*} F_1|0\rangle &=& \frac{1}{\sqrt{2}}\left(|0\rangle+|1\rangle\right) \\ F_1|1\rangle &=& \frac{1}{\sqrt{2}}\left(|0\rangle-|1\rangle\right) \\ \end{eqnarray*}$

を得る。ここで、 $e^{i2\pi0.1}=e^{i2\pi/2}=-1$ を用いた。2進数 $0.1$ は10進数で表すと $1/2$ であることに注意する。式(29)は先のブログで紹介したアダマールゲート $H$ で実現できる変換である。すなわち、 $F_1=H$ である。 $n=1$ の場合のQFTを実現する量子回路は下図である。

アダマールゲート $H$ の行列表記は

(30) $\begin{equation*} H=\frac{1}{\sqrt{2}}\left( \begin{array}{cc} 1 & 1 \\ 1 & -1 \\ \end{array} \right) \end{equation*}$

である。 $H$ を $|0\rangle=(1,0)$ と $|1\rangle=(0,1)$ に作用させると式(29)が得られることを確認してほしい。

QFTアルゴリズム：2量子ビットの場合

　天下り的であるが先に2量子ビットのQFTを実現する量子回路を示す。

最初の状態は

(31) $\begin{equation*} |\psi_0\rangle\equiv|j_0j_1\rangle=|j_0\rangle|j_1\rangle \end{equation*}$

である。先頭の量子ビット $|j_0\rangle$ にアダマールゲート $H$ をかけると、 $n=1$ の場合のQFTの式を参考にして

(32) $\begin{eqnarray*} |\psi_0\rangle &\rightarrow& \left(H|j_0\rangle\right)|j_1\rangle \\ &=& \frac{1}{\sqrt{2}}\left(|0\rangle+e^{i2\pi0.j_0}|1\rangle\right)|j_1\rangle\equiv|\psi_1\rangle \end{eqnarray*}$

を得る。次に $j_1=1$ のときだけ先頭の状態 $|j_0\rangle$ に位相ゲート $R_1$ を作用させる（このとき量子ビット $|j_1\rangle$ は制御ビットと呼ばれる）。 $R_1$ は次式で定義される行列である。

(33) $\begin{equation*} R_1=\left( \begin{array}{cc} 1 & 0 \\ 0 & e^{i2\pi/4} \\ \end{array} \right) \end{equation*}$

この変換を丁寧に見ていく。 $j_1=0$ のときは何もしないので状態は $|\psi_1\rangle$ のままである。 $j_1=1$ のときは、量子ビット $|j_0\rangle$ に $R_1$ を作用させるので

(34) $\begin{eqnarray*} |\psi_1\rangle &\rightarrow& \left(R_1\frac{1} {\sqrt{2}}\left(|0\rangle+e^{i2\pi0.j_0}|1\rangle\right)\right)|j_1=1\rangle \\ &=&\frac{1} {\sqrt{2}}\left(|0\rangle+e^{i2\pi0.j_01}|1\rangle\right)|j_1=1\rangle \end{eqnarray*}$

となる。 $j_1=0,1$ をまとめて書くと

(35) $\begin{equation*} |\psi_1\rangle\rightarrow \frac{1} {\sqrt{2}}\left(|0\rangle+e^{i2\pi0.j_0j_1}|1\rangle\right)|j_1\rangle\equiv|\psi_2\rangle \end{equation*}$

と変化することになる。この状態の後ろの量子ビット $|j_1\rangle$ にアダマールゲート $H$ をかけると

(36) $\begin{eqnarray*} |\psi_2\rangle &\rightarrow& \frac{1}{\sqrt{2}}\left(|0\rangle+e^{i2\pi0.j_0j_1}|1\rangle\right) \left(H|j_1\rangle\right) \\ &\rightarrow& \frac{1}{\sqrt{2}}\left(|0\rangle+e^{i2\pi0.j_0j_1}|1\rangle\right) \frac{1}{\sqrt{2}}\left(|0\rangle+e^{i2\pi0.j_1}|1\rangle\right)\equiv|\psi_3\rangle \end{eqnarray*}$

最後に、2つの量子ビットの順番を交換するSWAPゲート $S$ を作用させる（今回はSWAPゲートの詳細は省略する）。

(37) $\begin{equation*} |\psi_3\rangle\rightarrow \frac{1}{\sqrt{2}}\left(|0\rangle+e^{i2\pi0.j_1}|1\rangle\right) \frac{1}{\sqrt{2}}\left(|0\rangle+e^{i2\pi0.j_0j_1}|1\rangle\right)\equiv|\psi_4\rangle \end{equation*}$

最終的に得られた状態 $|\psi_4\rangle$ は式(26)で $n=2$ としたものである。すなわち

(38) $\begin{equation*} F_2|\psi_0\rangle=|\psi_4\rangle \end{equation*}$

である。

QFTアルゴリズム： $n$ 量子ビットの場合

　 $n=2$ の量子回路を拡張して $n=3$ の量子回路を得る。

さらに、拡張すれば、一般の $n$ の場合の量子回路を得る。

ただし

(39) $\begin{equation*} R_k=\left( \begin{array}{cc} 1 & 0 \\ 0 & e^{i2\pi/2^{k+1}} \\ \end{array} \right) \end{equation*}$

である。 $n\geq3$ の場合、SWAPゲート $S$ によりビット列は逆順に並べ替えられる。以上で量子フーリエ変換を量子回路で表現できた。

古典アルゴリズム（FFT）との比較

　ここでは、QFTアルゴリズムと古典アルゴリズム（FFT）の計算量の見積もりを行う。まず最初にQFTアルゴリズムから考える。量子ゲートである $H$ と $R_k$ の個数の和は、回路内の横線上にある四角の数を数えればよいから

(40) $\begin{equation*} \sum_{k=1}^{n}k=\frac{1}{2}n(n+1) \end{equation*}$

である。図の $S$ で示したSWAPゲートについては今回はその詳細を説明しないが、 $3n/2$ 個のCNOTゲートから構成される（CNOTゲートは先のブログで紹介した）。従って、QFTアルゴリズムで使う総ゲート数は

(41) $\begin{equation*} \frac{1}{2}n(n+1)+\frac{3n}{2} \end{equation*}$

となる。これは $n$ が大きくなると $n^2$ のオーダーで大きくなる量 $\mathcal{O}\left(n^2\right)$ である。フーリエ係数の個数 $N=2^n$ に換算すると $\mathcal{O}\left((\log_2{N})^2\right)$ となる。これが、QFTアルゴリズムの計算量である。一方、FFTの計算量は、 $\mathcal{O}\left(N\log_2{N}\right)$ である。下図に見る通り、QFTの方が高速であることが分かる。

式(14)の後の文章で触れたとおり、QFTの場合、計算後に $N$ 個のフーリエ係数を取り出す必要がある。この作業は指数関数的な時間を必要とするため、単にフーリエ係数を知りたいだけならFFTの方が速いというオチが付く。

blueqatによる実装

　今回作成したソースコードはここにある。まず最初に2量子ビットの場合のQFTを考える。このときの量子回路は以下であった。

初期状態が $|j_0j_1\rangle=|00\rangle$ の場合のQFTを実現するコードは以下の通り（main_2.py）。

def qft_with_2qbit_00():
    # |j0j1>=|00>
    c = Circuit(2)

    # j0ビットにアダマールゲートをかける。
    c = c.h[0]

    # j1ビットを制御ビットとしてj0ビットにR1をかける。
    c = c.cphase(math.pi / 2)[1, 0]

    # j1ビットにアダマールゲートをかける。
    c = c.h[1]

    # 実行。
    r = c.run()

    print("> 2qbit00")
    for i, v in enumerate(r):
        b = "{:02b}".format(i)
        print("[{}]{}".format(b, v))

def qft_with_2qbit_00():

# |j0j1>=|00>

c = Circuit(2)

# j0ビットにアダマールゲートをかける。

c = c.h[0]

# j1ビットを制御ビットとしてj0ビットにR1をかける。

c = c.cphase(math.pi / 2)[1, 0]

# j1ビットにアダマールゲートをかける。

c = c.h[1]

# 実行。

r = c.run()

print("> 2qbit00")

for i, v in enumerate(r):

b = "{:02b}".format(i)

print("[{}]{}".format(b, v))

コード内のコメントで示した通り、量子ゲートをひとつひとつ再現しているだけである。ただし、SWAPゲートは省略した。出力は以下の通り。

> 2qbit00
[00](0.4999999999999999+0j)
[01](0.4999999999999999+0j)
[10](0.4999999999999999+0j)
[11](0.4999999999999999+0j)

> 2qbit00

[00](0.4999999999999999+0j)

[01](0.4999999999999999+0j)

[10](0.4999999999999999+0j)

[11](0.4999999999999999+0j)

各基底状態（ $|00\rangle,|01\rangle,|10\rangle,|11\rangle$ ）の係数（複素振幅）が出力される。これらは複素数である（Pythonでは虚数単位はjと表記される）。上の回路図に従って手計算を行うとすべての係数が0.5となるので（試してほしい）上の結果と一致している。

　次のコードは同じ2量子ビットのQFTであるが、初期状態が $|j_0j_1\rangle=|11\rangle$ の場合である。

def qft_with_2qbit_11():
    # |j0j1>=|00>
    c = Circuit(2)

    # 両ビットにXゲートをかける。Xゲートはビットを反転させる計算である。状態は|11>になる。
    c = c.x[0, 1]

    # j0ビットにアダマールゲートをかける。
    c = c.h[0]

    # j1ビットを制御ビットとしてj0ビットにR1をかける。
    c = c.cphase(math.pi / 2)[1, 0]

    # j1ビットにアダマールゲートをかける。
    c = c.h[1]

    # 実行。
    r = c.run()

    print("> 2qbit11")
    for i, v in enumerate(r):
        b = "{:02b}".format(i)
        print("[{}]{}".format(b, v))

def qft_with_2qbit_11():

# |j0j1>=|00>

c = Circuit(2)

# 両ビットにXゲートをかける。Xゲートはビットを反転させる計算である。状態は|11>になる。

c = c.x[0, 1]

# j0ビットにアダマールゲートをかける。

c = c.h[0]

# j1ビットを制御ビットとしてj0ビットにR1をかける。

c = c.cphase(math.pi / 2)[1, 0]

# j1ビットにアダマールゲートをかける。

c = c.h[1]

# 実行。

r = c.run()

print("> 2qbit11")

for i, v in enumerate(r):

b = "{:02b}".format(i)

print("[{}]{}".format(b, v))

6行目で $|j_0\rangle$ と $|j_1\rangle$ の両ビットに $X$ ゲートを作用させている。 $X$ ゲートはビットを反転させる量子ゲートである。

(42) $\begin{eqnarray*} X|0\rangle&=&|1\rangle \\ X|1\rangle&=&|0\rangle \end{eqnarray*}$

従って、6行目で状態は $|00\rangle$ から $|11\rangle$ に変換される。以降の処理は先のコードと同じである。出力は以下の通り。

> 2qbit11
[00](0.4999999999999999+0j)
[01](-3.0616169978683824e-17-0.4999999999999999j)
[10](-0.4999999999999999+0j)
[11](3.0616169978683824e-17+0.4999999999999999j)

> 2qbit11

[00](0.4999999999999999+0j)

[01](-3.0616169978683824e-17-0.4999999999999999j)

[10](-0.4999999999999999+0j)

[11](3.0616169978683824e-17+0.4999999999999999j)

ここまでのコードでは、 $F_2|j\rangle$ までの計算しかしていない。フーリエ係数まで求めるコードは以下の通りである。

def qft(x):
    ll = len(x)
    lb = format(ll - 1, 'b')
    b = len(lb)

    # 状態|0>に演算子Fを作用させる。
    F = Circuit()  # 初期状態は|0>
    for i in range(b):
        F.h[i]
        for j in range(1, b - i):
            F.cphase(math.pi / 2**(j))[i + j, i]

    w = []
    for j in range(ll):
        c = Circuit()
        b = format(j, 'b')
        bl = len(b)

        # Fが作用する状態|j>を作る。
        for a in range(bl):
            if b[bl - a - 1] == '1':
                c = c.x[len(lb) - 1 - a]

        # F|j>を計算する。
        # |j>を作る回路cの後ろにフーリエ変換する回路Fを結合する。
        d = c + F
        e = d.run()  # w_k,k=0,1,...
        w.append(e)

    y = []
    for k in range(ll):
        y_k = 0
        for j in range(ll):
            y_k += w[j][k] * x[j]
        y.append(y_k)
    return y

def qft(x):

ll = len(x)

lb = format(ll - 1, 'b')

b = len(lb)

# 状態|0>に演算子Fを作用させる。

F = Circuit() # 初期状態は|0>

for i in range(b):

F.h[i]

for j in range(1, b - i):

F.cphase(math.pi / 2**(j))[i + j, i]

w = []

for j in range(ll):

c = Circuit()

b = format(j, 'b')

bl = len(b)

# Fが作用する状態|j>を作る。

for a in range(bl):

if b[bl - a - 1] == '1':

c = c.x[len(lb) - 1 - a]

# F|j>を計算する。

# |j>を作る回路cの後ろにフーリエ変換する回路Fを結合する。

d = c + F

e = d.run() # w_k,k=0,1,...

w.append(e)

y = []

for k in range(ll):

y_k = 0

for j in range(ll):

y_k += w[j][k] * x[j]

y.append(y_k)

return y

このコードは一般の $n$ 量子ビットの場合のQFTである。入力がx（ $=(x_0,\cdots,x_{N-1})$ ）、出力がy（ $=(y_0,\cdots,y_{N-1})$ ）である。ここで、 $N=2^n$ である。7行目から11行目までのコードで量子回路Fを構築している。14行目から28行目までのコードで $|j\rangle=|j_0\cdots j_{n-1}\rangle$ の取り得るすべての状態cを再現し、その各々をFに接続している（26行目：d=c+F）。31行目から35行目までの計算でフーリエ係数y_kを計算している。 $n$ 量子ビットの場合の量子回路を再掲しておく。

上で実装した関数qftを用いて次式を離散フーリエ変換してみる。

(43) $\begin{equation*} y=5\sin{2x} \end{equation*}$

これを離散化するため、 $0$ から $2\pi$ までを $N=2^4=16$ 等分する。

(44) $\begin{equation*} x_j=5\sin{2\left(\frac{2\pi j}{16}\right)}, j=0,\cdots,15 \end{equation*}$

このとき、フーリエ係数は

(45) $\begin{equation*} y_k=\frac{1}{4}\sum_{j=0}^{N-1}x_j\exp{\left(i\frac{2\pi kj}{16}\right)} \end{equation*}$

となる。この結果をqftで用いて再現してみる。

def fun(j):
    return 5.0 * np.sin(2.0 * 2.0 * np.pi * j / N)


def normal_ft(k, xs):
    s = 0.0
    for j in range(N):
        s += xs[j] * np.exp(complex(0, 2.0 * np.pi * k * j / N))
    return s / np.sqrt(N)

def fun(j):

return 5.0 * np.sin(2.0 * 2.0 * np.pi * j / N)

def normal_ft(k, xs):

s = 0.0

for j in range(N):

s += xs[j] * np.exp(complex(0, 2.0 * np.pi * k * j / N))

return s / np.sqrt(N)

1行目の関数は式(44)を、5行目の関数は式(45)を表す。これらの式を用いて計算したフーリエ係数と、上で実装したqftによる結果とを比較する。

if __name__ == "__main__":
    xs = []
    for i in range(N):
        xs.append(fun(i))

    print("> Classical Fourier Transformation")
    for k in range(N):
        y = normal_ft(k, xs)
        v = np.abs(y)
        b = "{:04b}".format(k)
        print("[{}]{}".format(b, v))

    print()
    print("> Quantume Fourier Transformation")
    ys = qft(xs)
    for (i, y) in enumerate(ys):
        v = np.abs(y)
        b = "{:04b}".format(i)
        print("[{}]{}".format(b, v))

if __name__ == "__main__":

xs = []

for i in range(N):

xs.append(fun(i))

print("> Classical Fourier Transformation")

for k in range(N):

y = normal_ft(k, xs)

v = np.abs(y)

b = "{:04b}".format(k)

print("[{}]{}".format(b, v))

print()

print("> Quantume Fourier Transformation")

ys = qft(xs)

for (i, y) in enumerate(ys):

v = np.abs(y)

b = "{:04b}".format(i)

print("[{}]{}".format(b, v))

ここでは、フーリエ係数の大きさ（ $|y_k|,k=0,\cdots,15$ ）を出力している。結果は以下の通り。

> Classical Fourier Transformation
[0000]1.1102230246251565e-15
[0001]1.5895974606912448e-15
[0010]10.000000000000004
[0011]1.1430445635548515e-15
[0100]2.180278542193117e-15
[0101]2.835419185914153e-15
[0110]4.986125195590856e-15
[0111]7.795513010265399e-15
[1000]2.90236210530186e-15
[1001]8.799524318285616e-15
[1010]3.1401849173675502e-15
[1011]4.718447854656915e-15
[1100]8.253197966801011e-15
[1101]1.7997002884907162e-14
[1110]10.000000000000002
[1111]1.3011534395271244e-14

> Quantume Fourier Transformation
[0000]7.771561172376096e-16
[0001]1.7554167342883505e-15
[0010]10.0
[0011]1.4936523181711914e-15
[0100]5.992897688303891e-16
[0101]1.217454592664719e-15
[0110]1.8209011484034673e-15
[0111]1.710067842969272e-15
[1000]1.2212453270876722e-15
[1001]1.807312143953211e-15
[1010]2.144196810752504e-15
[1011]1.0990647210786425e-15
[1100]3.7300332558822236e-16
[1101]2.159233940399377e-15
[1110]9.999999999999998
[1111]1.738660300550457e-15

> Classical Fourier Transformation

[0000]1.1102230246251565e-15

[0001]1.5895974606912448e-15

[0010]10.000000000000004

[0011]1.1430445635548515e-15

[0100]2.180278542193117e-15

[0101]2.835419185914153e-15

[0110]4.986125195590856e-15

[0111]7.795513010265399e-15

[1000]2.90236210530186e-15

[1001]8.799524318285616e-15

[1010]3.1401849173675502e-15

[1011]4.718447854656915e-15

[1100]8.253197966801011e-15

[1101]1.7997002884907162e-14

[1110]10.000000000000002

[1111]1.3011534395271244e-14

> Quantume Fourier Transformation

[0000]7.771561172376096e-16

[0001]1.7554167342883505e-15

[0010]10.0

[0011]1.4936523181711914e-15

[0100]5.992897688303891e-16

[0101]1.217454592664719e-15

[0110]1.8209011484034673e-15

[0111]1.710067842969272e-15

[1000]1.2212453270876722e-15

[1001]1.807312143953211e-15

[1010]2.144196810752504e-15

[1011]1.0990647210786425e-15

[1100]3.7300332558822236e-16

[1101]2.159233940399377e-15

[1110]9.999999999999998

[1111]1.738660300550457e-15

どちらの結果も、0010(=2)番目と1110(=14)番目の項が10、その他の項は0となる。 $y_2$ と $y_{14}$ だけが値を持つのは元の式 $y=5\sin{2x}$ が周波数2の正弦波のためである（ $y_{14}$ は $y_{2}$ の複素共役である）。

まとめ

　今回は、量子コンピュータを用いた計算の具体例としてフーリエ変換を取り上げ、かなり詳しく解説した（かなり長い記事になってしまった）。量子コンピュータで計算するとは、アルゴリズムを量子回路で表現することである。次回は、回帰を量子コンピュータで計算する方法を説明する予定である。

参考文献

量子コンピューティング基本アルゴリズムから量子機械学習まで

量子コンピュータ入門[第2版]

量子フーリエ変換とは

blueqatで量子フーリエ変換

Kumada Seiya

仕事であろうとなかろうと勉強し続ける、その結果”中身”を知ったエンジニアになれる

初歩的な数学のアルゴリズムを考えてみよう

アジャイル導入について

量子フーリエ変換

はじめに

離散フーリエ変換

振幅エンコーディング

量子フーリエ変換

QFTアルゴリズム（ステップ1）

QFTアルゴリズム（ステップ2）

QFTアルゴリズム：1量子ビットの場合

QFTアルゴリズム：2量子ビットの場合

QFTアルゴリズム： $n$ 量子ビットの場合

古典アルゴリズム（FFT）との比較

blueqatによる実装

まとめ

参考文献

Kumada Seiya

最近の記事

顔認証（DeepFace vs GPT4o）

時系列解析再訪

πの新しい公式

中心極限定理

DBSCANを用いたクラスタリング

レイアウト解析

イジングモデルと巡回セールスマン問題

量子位相推定

量子もつれ

量子ビットとは何か？

量子テレポーテーション

量子振幅増幅

Google Vision APIでOCR

アーカイブ

カテゴリー