めもめも

このブログに記載の内容は個人の見解であり、必ずしも所属組織の立場、戦略、意見を代表するものではありません。

量子フーリエ変換を「位相差に情報を埋め込む」という視点で理解してみる

enakai00.hatenablog.com

何の話かと言うと

上記の記事で説明した量子フーリエ変換は、本質的には、

　\displaystyle U( {\mid j_1\rangle}\otimes\cdots\otimes {\mid j_n\rangle} ) = \frac{1}{2^{n/2}} \left( {\mid 0\rangle} + e^{2\pi i[0.j_n]}{\mid 1\rangle}\right) \left({\mid 0\rangle}+e^{2\pi i[0.j_{n-1}j_n]}{\mid 1\rangle}\right) \cdots \left({\mid 0\rangle}+e^{2\pi i[0.j_1\cdots j_n]}{\mid 1\rangle}\right)　--- (1)

という線形変換（ユニタリ変換）を量子回路で実装しただけのことですが、「なぜ (1) の表式が量子回路と相性がよいのか？」という部分を踏み込んで考えてみます。

位相差情報の抽出機としての（逆）フーリエ変換

もともとの話の流れでは、数学的に知られていたフーリエ変換

　\displaystyle y_k := \frac{1}{\sqrt{N}}\sum_{j=0}^{N-1}x_j\exp\left(2\pi i\times\frac{jk}{N} \right)

を量子ビットで表現したら、なんとも不思議な (1) の表式が得られた・・・・という事になりますが、実は、(1) はアダマール変換 H の自然な拡張と見なすことができます。

どいういうことかと言うと、容易に分かるように、アダマール変換 H は次のように表すことができます。

　\displaystyle H : {\mid j\rangle} \longrightarrow \frac{1}{\sqrt{2}}\left({\mid 0\rangle} + (-1)^j {\mid 1\rangle}\right)

これは、1量子ビットの生の値を同じ量子ビットの位相差に変換する。つまり、位相差に情報を埋め込む演算とみなすことができます。

実は (1) で1量子ビットの場合を考えると、まさにこの変換に一致しています。つまり、(1) は「量子ビットの値を位相差に埋め込む」という操作を複数量子ビットに自然に拡張したものになっているのです。

そして、逆フーリエ変換はこの反対の操作、つまり、「位相差の情報を取り出す操作」と理解することができます。

こう考えると、位相推定のアルゴリズム（というかそれをやる下記の量子回路）は、かなりスッキリと理解することができます。

ユニタリ演算の固有値という「位相情報」をうまいこと、位相差の情報として埋め込んでおき、それを逆フーリエ変換で取り出しているだけにすぎません。

Phase Estimator の量子回路を解説してみる（その２）

前回 は量子フーリエ変換を実装した量子回路を紹介しました。今回は、これを用いて、ユニタリ演算の固有値（の位相）を推定する方法を説明します。

逆フーリエ変換

有名な事実ですが、前回説明したフーリエ変換

　\displaystyle y_k := \frac{1}{\sqrt{N}}\sum_{j=0}^{N-1}x_j\exp\left(2\pi i\times\frac{jk}{N} \right)　--- (1)

には、逆変換があります。次になります。

　\displaystyle x_j := \frac{1}{\sqrt{N}}\sum_{k=0}^{N-1}y_k\exp\left(-2\pi i\times\frac{jk}{N} \right)　--- (2)

(1) と (2) を見比べると指数関数の引数の符号が変わっているだけですが、これが逆変換になっていることは直接計算で確認できます。(2) の右辺に (1) を代入して計算してみます。

　\displaystyle \frac{1}{\sqrt{N}}\sum_{k=0}^{N-1}y_k\exp\left(-2\pi i\times\frac{jk}{N} \right) =\frac{1}{N}\sum_{k=0}^{N-1}\sum_{j'=0}^{N-1}x_{j'}\exp\left(2\pi i\times\frac{j'k}{N} \right)\exp\left(-2\pi i\times\frac{jk}{N} \right)

　　\displaystyle =\frac{1}{N}\sum_{j'=0}^{N-1}x_{j'}\sum_{k=0}^{N-1}\exp\left(2\pi i\times\frac{(j'-j)k}{N} \right)

ここで、j'=j の場合は、

　\displaystyle \sum_{k=0}^{N-1}\exp\left(2\pi i\times\frac{(j'-j)k}{N} \right) = N

となりますが、一方、j'\ne j の場合、この和は、複素平面上の単位円において中心角を等分した点の値を足し合わせたものになり、対称性から、

　\displaystyle \sum_{k=0}^{N-1}\exp\left(2\pi i\times\frac{(j'-j)k}{N} \right) = 0

となります。したがって、j' についての和は j'=j の部分だけが残って、(2) が成り立ちます。

したがって、前回のフーリエ変換の量子回路において、\displaystyle R_k = Z^{1/2^{k-1}} を \displaystyle R'_k = Z^{-1/2^{k-1}} に置き換えれば逆フーリエ変換

　\displaystyle U^{-1}\left(\frac{1}{2^{n/2}} \left( {\mid 0\rangle} + e^{2\pi i[0.j_n]}{\mid 1\rangle}\right) \left({\mid 0\rangle}+e^{2\pi i[0.j_{n-1}j_n]}{\mid 1\rangle}\right) \cdots \left({\mid 0\rangle}+e^{2\pi i[0.j_1\cdots j_n]}{\mid 1\rangle}\right)\right) ={\mid j_1\rangle}\otimes\cdots\otimes {\mid j_n\rangle} 　--- (3)

の量子回路が得られます。

位相推定

唐突ですが、(3) の左辺にある

　\displaystyle \frac{1}{2^{n/2}} \left( {\mid 0\rangle} + e^{2\pi i[0.j_n]}{\mid 1\rangle}\right) \left({\mid 0\rangle}+e^{2\pi i[0.j_{n-1}j_n]}{\mid 1\rangle}\right) \cdots \left({\mid 0\rangle}+e^{2\pi i[0.j_1\cdots j_n]}{\mid 1\rangle}\right)　--- (3)

という量子状態を作る方法を考えます。結論から言うと、次の量子回路で実現できます。

ここに、{\mid u\rangle} はユニタリ演算子 U の固有状態で固有値は e^{2\pi i\varphi}（\varphi = [0.j_1j_2\cdots j_n]）とします。つまり、

　\displaystyle U^{2^0}{\mid u\rangle} = e^{2\pi i [0.j_1j_2\cdots j_n]}{\mid u\rangle}

　\displaystyle U^{2^1}{\mid u\rangle} = e^{2\pi i [j_1.j_2\cdots j_n]}{\mid u\rangle} = e^{2\pi i [0.j_2\cdots j_n]}{\mid u\rangle}

　\displaystyle U^{2^2}{\mid u\rangle} = e^{2\pi i [j_1j_2.j_3\cdots j_n]}{\mid u\rangle} = e^{2\pi i [0.j_3\cdots j_n]}{\mid u\rangle}

などが成り立ちます。これより、上記の量子回路の出力が (3) に一致することがわかります。

さらにこの出力に逆フーリエ変換の量子回路を適用すると・・・最後の結果は、

　\displaystyle{\mid j_1\rangle}\otimes\cdots\otimes {\mid j_n\rangle}

になります。つまり、この状態を観測して、j_1,j_2,\cdots,j_n の値を読み取れば、固有値の位相 \varphi = [0.j_1j_2\cdots j_n] が分かることになります。

つまり、固有値が未知のユニタリ演算、および、固有状態のペアがあった時に、この未知の固有値を測定できることになります。これが、量子回路による位相推定アルゴリズムというわけです。

固有値が n 桁の二進数小数で厳密に表現できない場合は、当然ながらこの手法は近似値を与えることになります。

Cirq による実装例

github.com

Phase Estimator の量子回路を解説してみる（その１）

まずは、Phase Estimator の基礎となる、量子フーリエ変換の実装を解説します。

（離散）フーリエ変換とは？

形式的に言うと、N 成分の複素ベクトル \mathbf x = (x_0,\cdots,x_{N-1}) を N 成分の複素ベクトル \mathbf y = (y_0,\cdots,y_{N-1}) に変換する、次の計算式です。

　\displaystyle y_k := \frac{1}{\sqrt{N}}\sum_{j=0}^{N-1}x_j\exp\left(2\pi i\times\frac{jk}{N} \right)　--- (1)

このような変換を考えて何が嬉しいのかと言うと・・・、実はこの変換には、元の成分 (x_0,\cdots,x_{N-1}) を数列とみなした時に、この数列が持つ「周期性」を検出する機能があります。

論より証拠で具体例を見てみましょう。

まず、\mathbf x の成分がすべて同じ（言うなれば周期 0）の場合で、 x_j = 1\,(j=0,\cdots,15) とします。対応する y_k\,(k=0,\cdots,15) の方は、（複素数をグラフ表示するのは難しいので）実部を取り出してグラフに描きます。

これを見ると、y_k の方は、y_0 だけが値を持っています。つまり、y_0 は周期 0 （定数）を示すフラグなのです。

次に、\displaystyle x_j = \cos\left(\frac{2 \pi j}{16}\right) として同様のグラフを描きます。

この場合、x_j は j=0,\cdots,15 の範囲でちょうど1回振動する（この範囲を単位として）周期1の三角関数です。y_k の方は、y_1 と y_{15} だけが値を持っています。つまり、この2つが周期1のフラグと言えます。（連続的なフーリエ変換を知っている方は、y_{15} が値を持つのが気持ち悪いかもしれませんが、今は値を持つ部分が離散的なので、周期 1/2 と周期 1/30 の区別がつかないことが理由です。）

あとは大体予想がつきますが、\displaystyle x_j = \cos\left(\frac{4 \pi j}{16}\right) とすると次の結果になります。

x_j は j=0,\cdots,15 の範囲で2回振動する（この範囲を単位として）周期0.5の三角関数です。y_k の方は、y_2 と y_{14} だけが値を持っています。

この変換は、複数の周期を合成したものを検知することもできます。次の関数を考えてみましょう。

　\displaystyle x_j = \frac{1}{2}\left\{\cos\left(\frac{2 \pi j}{16}\right)+\cos\left(\frac{4 \pi j}{16}\right)\right\}

結果は次の通りです。

周期1と周期0.5に対応するフラグがどちらも反応していることがわかります。

このようにフーリエ変換を用いると、さまざまな関数の「周期性」を発見するのに役立てることができます。

・・・これ自体は量子回路と何の関係もありませんが、フーリエ変換を愚直に実行すると計算量は O(N^2) になります。つまり、N 項の級数を N 変数分計算する必要があります。より高速に計算する FFT （高速フーリエ変換）のアルゴリズムが知られていますが、この場合の計算量は O(N\log N) です。一方、量子回路を用いると、O\left((\log N)^2\right) の計算量でこれを実行することができます。つまり、フーリエ変換は、量子コンピューターを用いることで、古典コンピューターよりも高速に実行できるのです。

量子ビットによるフーリエ変換の表現

(1) のフーリエ変換はベクトルの成分を用いて表されていましたが、線形変換になっていることがわかります。したがって、この成分表示に対応する正規直交基底を {\mid 0\rangle},\cdots,{\mid N-1\rangle} として、それぞれの基底ベクトルの変換法則がわかれば十分です。具体的には、

　\displaystyle\mathbf x = \sum_{j=0}^{N-1}x_j{\mid j\rangle},\ \mathbf y = \sum_{k=0}^{N-1}y_k{\mid k\rangle}

として、フーリエ変換を \mathbf y = U(\mathbf x) と表すことにすると、次の関係が成り立ちます。

　\displaystyle U(\mathbf x) = \sum_{k=0}^{N-1}y_k{\mid k\rangle} = \sum_{k=0}^{N-1}\frac{1}{\sqrt{N}}\sum_{j=0}^{N-1}x_j\exp\left(2\pi i\times\frac{jk}{N} \right){\mid k\rangle} = \sum_{j=0}^{N-1}x_j\left\{\frac{1}{\sqrt{N}}\sum_{k=0}^{N-1}\exp\left(2\pi i\times\frac{jk}{N} \right){\mid k\rangle}\right\}

これより、基底ベクトルの変換法則は次になることがわかります。

　\displaystyle U({\mid j\rangle}) = \frac{1}{\sqrt{N}}\sum_{k=0}^{N-1}\exp\left(2\pi i\times\frac{jk}{N} \right){\mid k\rangle}　--- (2)

ここで、N 個ある基底ベクトルを n = \log N 個の量子ビットで表現します。これは単純に j=0,\cdots,N-1 を二進数表記に直して、各桁のビットを量子ビットに置き換えます。たとえば、N=4 として、次のような対応になります。

　 {\mid 0\rangle} \Leftrightarrow {\mid 00\rangle} = {\mid 0\rangle}\otimes {\mid 0\rangle}
　 {\mid 1\rangle} \Leftrightarrow {\mid 01\rangle} = {\mid 0\rangle}\otimes {\mid 1\rangle}
　 {\mid 2\rangle} \Leftrightarrow {\mid 10\rangle} = {\mid 1\rangle}\otimes {\mid 0\rangle}
　 {\mid 3\rangle} \Leftrightarrow {\mid 11\rangle} = {\mid 1\rangle}\otimes {\mid 1\rangle}

一般に、n 個の量子ビットによる基底ベクトルを {\mid k_1\rangle}\otimes\cdots\otimes {\mid k_n\rangle} と表すと、次の対応関係が成り立ちます。

　\displaystyle \sum_{k=0}^{N-1}f(k){\mid k\rangle} \Leftrightarrow \sum_{k_1=0}^1\cdots\sum_{k_n=0}^1f\left(\sum_{l=1}^nk_l2^{n-l}\right){\mid k_1\rangle}\otimes\cdots\otimes {\mid k_n\rangle}　--- (3)

ここに、 k = [k_1k_2\cdots k_n] は k の二進数表記（つまり \displaystyle k = \sum_{l=1}^n k_l 2^{n-l}）とします。

そして、この「量子ビット基底ベクトル」を用いると、(2) の変換は量子ビットごとの変換に置き換えることができます。\displaystyle j=[j_1j_2\cdots j_n]= \sum_{l=1}^n j_l 2^{n-l} とする時、(3) のルールを用いて (2) を書き直すと次のようになります。

　\displaystyle U( {\mid j_1\rangle}\otimes\cdots\otimes {\mid j_n\rangle}) = \frac{1}{2^{n/2}}\sum_{k_1=0}^1\cdots\sum_{k_n=0}^1\exp\left(2\pi i \times \frac{j k}{2^n}\right){\mid k_1\rangle}\otimes\cdots\otimes {\mid k_n\rangle}

　　\displaystyle = \frac{1}{2^{n/2}}\sum_{k_1=0}^1\cdots\sum_{k_n=0}^1\exp\left(2\pi i \times j \sum_{l=1}^n k_l 2^{-l}\right){\mid k_1\rangle}\otimes\cdots\otimes {\mid k_n\rangle}

　　\displaystyle = \frac{1}{2^{n/2}}\sum_{k_1=0}^1\cdots\sum_{k_n=0}^1\prod_{l=1}^n \exp\left(2\pi i \times j k_l 2^{-l}\right){\mid k_1\rangle}\otimes\cdots\otimes {\mid k_n\rangle}

　　\displaystyle = \frac{1}{2^{n/2}}\sum_{k_1=0}^1\cdots\sum_{k_n=0}^1\bigotimes_{l=1}^n \exp\left(2\pi i \times j k_l 2^{-l}\right){\mid k_l\rangle}

　　\displaystyle = \frac{1}{2^{n/2}}\bigotimes_{l=1}^n \left\{ {\mid 0\rangle} + \exp\left(2\pi i \times j 2^{-l}\right){\mid 1\rangle}\right\}　--- (4)

最後の \displaystyle \exp\left(2\pi i \times j 2^{-l}\right) の部分は、二進数表記 j = [j_1j_2\cdots j_n] を用いると、

　\displaystyle j2^{-l} = [j_1\cdots j_{n-l}]+[0.j_{n-l+1}\cdots j_n]

であることから、\displaystyle \exp\left(2\pi i\times [j_1\cdots j_{n-l}]\right) = 1 より、次のように簡単化できます。

　\displaystyle \exp\left(2\pi i \times j 2^{-l}\right) = \exp\left(2\pi i \times[0.j_{n-l+1}\cdots j_n]\right)

これを (4) に代入して、最終的に次の関係が得られます。

　\displaystyle U( {\mid j_1\rangle}\otimes\cdots\otimes {\mid j_n\rangle} ) = \frac{1}{2^{n/2}} \left( {\mid 0\rangle} + e^{2\pi i[0.j_n]}{\mid 1\rangle}\right) \left({\mid 0\rangle}+e^{2\pi i[0.j_{n-1}j_n]}{\mid 1\rangle}\right) \cdots \left({\mid 0\rangle}+e^{2\pi i[0.j_1\cdots j_n]}{\mid 1\rangle}\right)　--- (5)

ここまでは量子回路とは関係のない、ただのテンソル計算と見ることもできますが、実はこの最後の表式から、この変換はアダマールゲートH と位相シフトゲート

　\displaystyle R_k := \begin{pmatrix} 1 & 0 \\ 0 & e^{2\pi i / 2^k} \end{pmatrix}（k=2,\cdots,n)

を用いた量子回路で実装できることがわかります。上記の位相シフトは、\displaystyle R_k = Z^{1/2^{k-1}} として実装できます。（Z 演算を物理実装する制御パルスの照射時間を変化させます。）

また、二進数表記を用いると、位相シフトゲートは次のようにも表現できます。（これが量子回路実装のポイントになります。）

　 R_2 = \begin{pmatrix} 1 & 0 \\ 0 & e^{2\pi i[0.01]}\end{pmatrix},\ R_3 = \begin{pmatrix} 1 & 0 \\ 0 & e^{2\pi i[0.001]}\end{pmatrix},\cdots　--- (6)

量子回路によるフーリエ変換の実装

結論から言うと、次の量子回路で実装できます。

はじめに、{\mid j_1\rangle} に対する H \rightarrow CR_2 \rightarrow \cdots \rightarrow CR_n の流れを見ます。まず、

　\displaystyle j_1 = 0 \Rightarrow e^{ 2\pi i[0.j_1] } = 1
　\displaystyle j_1 = 1 \Rightarrow e^{ 2\pi i[0.j_1] } = -1

の関係に注意すると、

　\displaystyle H{\mid j_1\rangle} = \frac{1}{2^{1/2}}\left({\mid 0\rangle}+e^{2\pi i [0.j_1]}{\mid 1\rangle}\right)

が成り立ちます。次に、CR_2 を演算すると、(6) の関係に注意して次が成り立ちます。

　\displaystyle CR_2 H{\mid j_1\rangle}\otimes{\mid j_2\rangle} = \frac{1}{2^{1/2}} \left({\mid 0\rangle}+e^{2\pi i [0.j_1]}e^{2\pi i [0.0j_2]} {\mid 1\rangle}\right)\otimes{\mid j_2\rangle} =\frac{1}{2^{1/2}} \left({\mid 0\rangle}+e^{2\pi i [0.j_1j_2]}{\mid 1\rangle}\right)\otimes{\mid j_2\rangle}

以下同様にして、CR_3,\cdots,CR_n を演算すると、最初の量子ビット \displaystyle{\mid j_1\rangle} は次のように変化します。

　 \displaystyle{\mid j_1\rangle} \longrightarrow \displaystyle \frac{1}{2^{1/2}} \left({\mid 0\rangle}+e^{2\pi i [0.j_1j_2\cdots j_n]}{\mid 1\rangle}\right)

2つ目以降の量子ビットについても同様の処理が行われるので、最後の H まで演算したところで、全体の状態は次になります。

　\displaystyle \frac{1}{2^{n/2}} \left({\mid 0\rangle}+e^{2\pi i[0.j_1\cdots j_n]}{\mid 1\rangle}\right) \left({\mid 0\rangle}+e^{2\pi i[0.j_2\cdots j_n]}{\mid 1\rangle}\right) \cdots \left( {\mid 0\rangle} + e^{2\pi i[0.j_n]}{\mid 1\rangle}\right)

これは、(5) の量子ビットの順番を反転させたものになっているので、最後にスワップゲートを組み合わせて、量子ビットの順序を並び替えれば完成です。

Cirq による実装例

github.com


続きはこちら

enakai00.hatenablog.com