単純な量子回路をゼロから書く

CHSH不等式の破れを確認するで、量子回路の基本的な書き方と実機での実行の仕方を学びました。そこで出てきたように、QCでは量子レジスタに特定の量子状態を実現することが計算を行うこととなります。そこで今回は、基本的なゲートを組み合わせて様々な量子状態を作る実習を行います。

\(\newcommand{\ket}[1]{|#1\rangle}\)

状態ベクトルシミュレータで量子状態を調べる

第一回の課題で利用したQCシミュレータqasm_simulatorは、実機と同様に量子回路の測定結果をヒストグラムとして返してくれるものでした。測定結果を返すということは、その出力からは量子状態の複素位相を含めた振幅は読み取れません。今回は作った量子状態をより詳細に調べるために、statevector_simulator(状態ベクトルシミュレータ)を利用します。状態ベクトルシミュレータは、回路の終状態におけるすべての計算基底の確率振幅、つまりその量子状態に関する最も完全な情報を返します。

シミュレータへのアクセスはqasm_simulator同様Aerからです。`

# まずは全てインポート
import numpy as np
import matplotlib.pyplot as plt
from qiskit import QuantumCircuit, Aer, execute
# utilsはこのワークブック独自のモジュール(インポートエラーが出る場合はPYTHONPATHを設定するか、sys.pathをいじってください)
from utils.show_state import show_state

print('notebook ready')
simulator = Aer.get_backend('statevector_simulator')
print(simulator.name())
statevector_simulator

シミュレートする回路を作ります。ここでは例として前半で登場した回路1を使います。測定以外のゲートをかけていきます。

circuit = QuantumCircuit(2)
circuit.h(0)
circuit.cx(0, 1)
circuit.ry(3. * np.pi / 4., 1)

circuit.draw('mpl')
_images/circuit_from_scratch_6_0.png

Caution

回路に測定操作(measure_all())を加えてしまうと、statevector_simulatorを利用して測定前の回路の量子状態を確認することができなくなります。

状態ベクトルシミュレーションを実行する際もexecute関数を使い、ジョブオブジェクトから結果を得ます。ただし、今回は計算結果からカウントではなく状態ベクトル(=量子振幅の配列)を得るのでget_countsではなくresult.data()['statevector']を参照します。

job = execute(circuit, simulator)
result = job.result()
statevector = result.data()['statevector']

print(type(statevector), statevector.dtype)
print(statevector)
<class 'numpy.ndarray'> complex128
[ 0.27059805+0.j -0.65328148+0.j  0.65328148+0.j  0.27059805+0.j]

状態ベクトルはnumpyndarray配列オブジェクトとして実装されています。データタイプは128ビットの複素数(実部と虚部それぞれ64ビット)です。配列のインデックスがそのまま二進数としてみなした計算基底の値に対応しています。

matplotlibのテキスト表示を利用して、得られた状態ベクトルをケット表示してみましょう。

state_string = r'|\psi\rangle ='
for idx in range(statevector.shape[0]):
    angle = np.arcsin(statevector[idx].real * np.sqrt(2.))
    norm_angle = np.round(angle / np.pi * 8.).astype(int)
    state_string += r' sin({} \pi / 8)|{}\rangle'.format('' if norm_angle == 1 else norm_angle, idx)
    if idx != statevector.shape[0] - 1:
        state_string += '+'

plt.figure(figsize=[10., 0.2])
plt.axis('off')
t = plt.text(0.5, 0.5, '${}$'.format(state_string), fontsize='xx-large')
_images/circuit_from_scratch_10_0.png

その他のゲート

CHSH不等式の破れを確認するで、よく使われる1量子ビットゲートと2量子ビット制御ゲートを紹介しましたが、他にも知っておくと便利なゲートがいくつかあります。

1、2量子ビット

ゲート名

対象量子ビット数

説明

Qiskitコード

\(\sqrt{X}\), SX

1

計算基底それぞれに対して、以下の変換をする。

\[\begin{split}\sqrt{X} \ket{0} = \frac{1}{2} \left[(1 + i) \ket{0} + (1 - i) \ket{1}\right] \\ \sqrt{X} \ket{1} = \frac{1}{2} \left[(1 - i) \ket{0} + (1 + i) \ket{1}\right]\end{split}\]

circuit.sx(i)

\(U_1\), \(P\)

1

パラメータ\(\phi\)を取り、\(\ket{1}\)に位相\(e^{i\phi}\)をかける。\(P(\phi)\)\(R_{z}(\phi)\)と等価なので、1量子ビットゲートとして利用する分にはどちらでも同じ結果が得られる。

circuit.p(phi, i)

\(U_3\)

1

パラメータ\(\theta, \phi, \lambda\)を取り、単一量子ビットの任意の状態を実現する。

\[\begin{split}U_3(\theta, \phi, \lambda)\ket{0} = \cos\left(\frac{\theta}{2}\right) \ket{0} + e^{i\phi}\sin\left(\frac{\theta}{2}\right) \ket{1} \\ U_3(\theta, \phi, \lambda)\ket{1} = -e^{i\lambda}\sin\left(\frac{\theta}{2}\right) \ket{0} + e^{i(\phi+\lambda)}\cos\left(\frac{\theta}{2}\right) \ket{1}\end{split}\]

circuit.u3(theta, phi, lambda, i)

\(C^i_j[U_1]\), \(C^i_j[P]\)

2

パラメータ\(\phi\)を取り、ビット\(i, j\)が1であるような計算基底の位相を\(\phi\)前進させる。

\[\begin{split}C^i_j[U_1(\phi)] \ket{00} = \ket{00} \\ C^i_j[U_1(\phi)] \ket{01} = \ket{01} \\ C^i_j[U_1(\phi)] \ket{10} = \ket{10} \\ C^i_j[U_1(\phi)] \ket{11} = e^{i\phi} \ket{11}\end{split}\]

circuit.cp(phi, i, j)

SWAP

2

二つの量子ビットの量子状態を交換する。CNOT三回で実装できる。

\[\begin{split}\begin{align} & C^0_1[X]C^1_0[X]C^0_1[X] (\alpha \ket{00} + \beta \ket{01} + \gamma \ket{10} + \delta \ket{11}) \\ = & C^0_1[X]C^1_0[X] (\alpha \ket{00} + \beta \ket{01} + \gamma \ket{11} + \delta \ket{10}) \\ = & C^0_1[X] (\alpha \ket{00} + \beta \ket{11} + \gamma \ket{01} + \delta \ket{10}) \\ = & \alpha \ket{00} + \beta \ket{10} + \gamma \ket{01} + \delta \ket{11} \end{align}\end{split}\]

circuit.swap(i, j)

3量子ビット以上

量子レジスタに対する全ての操作が1量子ビットゲートと2量子ビット制御ゲート(実際にはCNOTのみ)の組み合わせで表現できることはすでに触れましたが、実際に量子回路を書くときには3つ以上の量子ビットにかかるゲートを使用したほうがアルゴリズムをコンパクトに表現できます。そのようなゲートを含む回路を実機で実行するときには、実行直前に「トランスパイル」(transpilation)という過程を経て回路を変換します。3量子ビット以上のゲートは最終的には全て1量子ビットゲートとCNOTに自動的に分解されます。

3量子ビット以上にかかるゲートも基本は制御ゲートで、よく使われるのは制御ビットが複数あるような多重制御(multi-controlled)ゲート\(C^{ij\dots}_{k}[U]\)です。そのうち制御ビットが2つで標的ビット1つに\(X\)をかける\(C^{ij}_{k}[X]\)は特別多用され、Toffoliゲートとも呼ばれます。

多重制御ゲートでない3量子ビットゲートの例として、controlled SWAPというものも存在します。制御ビット1つで他の2つの量子ビットの間の量子状態交換を制御します。当然制御ビットを増やしたmulti-controlled SWAPというものも考えられます。

Qiskitで多重制御ゲートに対応するコードは一般にcircuit.mc?(parameters, [controls], target)という形を取ります。例えばToffoliゲートは

circuit.mcx([i, j], k)

で、\(C^{ij}_{k}[R_y]\)

circuit.mcry(theta, [i, j], k)

となります。

量子状態生成

それでは、これまでに登場したゲートを使って量子状態生成のエクササイズをしましょう。

回路の量子状態を表示するshow_stateという関数はこのワークブックのレポジトリからインポートしています。この関数は内部でstatevector_simulatorを利用して回路の状態ベクトルを読み出しています。

問題1: 1量子ビット、相対位相付き

問題

1量子ビットに対して状態

\[ \frac{1}{\sqrt{2}}\left(\ket{0} + i\ket{1}\right) \]

を作りなさい。

circuit = QuantumCircuit(1)
##################
### EDIT BELOW ###
##################
circuit.h(0)
# circuit.?
##################
### EDIT ABOVE ###
##################

show_state(circuit, amp_norm=(np.sqrt(0.5), r'\frac{1}{\sqrt{2}}'))

解答

正解は無限に存在しますが、上の回路ですでにアダマールゲートがかかっているので、最も手っ取り早いのは

circuit.h(0)
circuit.p(np.pi / 2., 0)

でしょう。自分で入力して確認してみてください。

また、上で触れたように、単一量子ビットの場合、\(U_3\)ゲートを使えば任意の状態を一つのゲート操作で生成できます。\(U_3\)を使った解答は

circuit.u3(np.pi / 2., np.pi / 2., 0., 0)

です。

問題2: ベル状態、相対位相付き

問題

2量子ビットに対して状態

\[ \frac{1}{\sqrt{2}}\left(\ket{0} + i\ket{3}\right) \]

を作りなさい。

circuit = QuantumCircuit(2)

##################
### EDIT BELOW ###
##################
# circuit.?
circuit.cx(0, 1)
##################
### EDIT ABOVE ###
##################

show_state(circuit, amp_norm=(np.sqrt(0.5), r'\frac{1}{\sqrt{2}}'))

解答

問題1の結果を利用すると簡単です。初期状態は両ビットとも\(\ket{0}\)なので、2ビットレジスタに対して問題1の操作をすると状態\(1/\sqrt{2}(\ket{0} + i\ket{1})\ket{0}\)ができます。ここでビット0から1にCNOTをかけると、\(\ket{1}\ket{0} \rightarrow \ket{1}\ket{1}\)となるので、望む状態が実現します。

circuit.h(0)
circuit.p(np.pi / 2., 0)
circuit.cx(0, 1)

問題3: GHZ状態

問題

3量子ビットに対して状態

\[ \frac{1}{\sqrt{2}} (\ket{0} + \ket{7}) \]

を作りなさい。

circuit = QuantumCircuit(3)

##################
### EDIT BELOW ###
##################
# circuit.?
circuit.mcx([0, 1], 2)
##################
### EDIT ABOVE ###
##################

show_state(circuit, amp_norm=(np.sqrt(0.5), r'\frac{1}{\sqrt{2}}'))

解答

Bell状態の生成と発想は同じですが、多重制御ゲート\(C^{01}_2[X]\)を利用します。

circuit.h(0)
circuit.cx(0, 1)
circuit.mcx([0, 1], 2)

問題4: Equal superposition

問題

一般の\(n\)量子ビットに対して状態

\[ \frac{1}{\sqrt{2^n}} \sum_{k=0}^{2^n-1} \ket{k} \]

を作る回路を考え、\(n=4\)のケースを実装しなさい。

num_qubits = 4

circuit = QuantumCircuit(num_qubits)

##################
### EDIT BELOW ###
##################
# circuit.?
##################
### EDIT ABOVE ###
##################

sqrt_2_to_n = np.sqrt(2 ** num_qubits)
show_state(circuit, amp_norm=(1. / sqrt_2_to_n, r'\frac{1}{' + str(int(sqrt_2_to_n)) + '}'))

解答

Equal superpositionと呼ばれるこの状態は、量子計算でおそらく一番多様される状態です。一見複雑に思えますが、\(n\)桁の二進数で\(0\)から\(2^n-1\)までの数を表現すると、全ての0/1の組み合わせが登場するということを考えると、

\[ \frac{1}{\sqrt{2^n}} \sum_{k=0}^{2^n-1} \ket{k} = \frac{1}{\sqrt{2}}(\ket{0} + \ket{1}) \otimes \frac{1}{\sqrt{2}}(\ket{0} + \ket{1}) \otimes \cdots \otimes \frac{1}{\sqrt{2}}(\ket{0} + \ket{1}) \]

と、\(n\)個の\(1/\sqrt{2}(\ket{0} + \ket{1})\)の積と等しいことがわかります。つまり各ビットに対してアダマールゲートをかけた状態にほかならず、

for i in range(4):
    circuit.h(i)

が正解です。

問題5: 特定の基底の符号を反転させる

問題

問題4の4ビットequal superposition状態において、基底\(\ket{5}\)の符号を反転させなさい。

num_qubits = 4

circuit = QuantumCircuit(num_qubits)

##################
### EDIT BELOW ###
##################
for i in range(num_qubits):
    circuit.h(i)
    
# circuit.?
##################
### EDIT ABOVE ###
##################

sqrt_2_to_n = np.sqrt(2 ** num_qubits)
show_state(circuit, amp_norm=(1. / sqrt_2_to_n, r'\frac{1}{' + str(int(sqrt_2_to_n)) + '}'))

解答

\(\ket{5}=\ket{0101}\)の符号を反転させるには、4つの量子ビットが0, 1, 0, 1という値を取っている状態のみに作用するゲートを使います。必然的に多重制御ゲートが用いられますが、残念ながら制御ゲートはビットの値が1であるときにしか働きません。そこで、まず\(\ket{0101}\)\(\ket{1111}\)に変換し、その上で\(C^{012}_{3}[P(\pi)]\)を使います1。符号を反転させたら、最初の操作の逆操作をして\(\ket{1111}\)\(\ket{0101}\)に変換します。

for i in range(4):
    circuit.h(i)

circuit.x(1)
circuit.x(3)
circuit.mcp(np.pi, [0, 1, 2], 3)
circuit.x(1) # Xの逆操作はX
circuit.x(3)

ビット1と3にかかっている\(X\)ゲートは当然他の基底にも影響するので、例えば\(\ket{0011}\)\(\ket{0110}\)となります。しかし、1と3に\(X\)をかけて\(\ket{1111}\)になるのは\(\ket{0101}\)だけなので、逆操作でちゃんと「後片付け」さえすれば、この一連の操作で\(\ket{5}\)だけをいじることができます。

問題6: Equal superpositionに位相を付ける

問題

一般の\(n\)量子ビットに対して状態

\[ \frac{1}{\sqrt{2^n}}\sum_{k=0}^{2^n-1} e^{2\pi i s k/2^n} \ket{k} \quad (s \in \mathbb{R}) \]

を作る回路を考え、\(n=6, s=2.5\)のケースを実装しなさい。

num_qubits = 6

circuit = QuantumCircuit(num_qubits)

s = 2.5

##################
### EDIT BELOW ###
##################
#for i in range(num_qubits):
#    circuit.h(i)
#
# circuit.?
##################
### EDIT ABOVE ###
##################

sqrt_2_to_n = np.sqrt(2 ** num_qubits)
show_state(circuit, amp_norm=(1. / sqrt_2_to_n, r'\frac{1}{' + str(int(sqrt_2_to_n)) + '}'), phase_norm=(2 * np.pi / (2 ** num_qubits), r'\frac{2 \pi i}{' + str(2 ** num_qubits) + '}'))

解答

この回路は次の問題の量子フーリエ変換と深く関係しますが、量子フーリエ変換では任意の状態\(\ket{s}\)(このとき\(s\)は整数)を上の状態に変換するのに対し、この問題では\(s\)は外から与えられているパラメータです。そのため回路もだいぶ簡単になり、実は制御ゲートなしで実装できます。

問題4と同じように、この状態も個々の量子ビットにおける\(\ket{0}\)\(\ket{1}\)の重ね合わせの積として表現できます。そのためには\(k\)の二進数表現\(k = \sum_{m=0}^{n-1} 2^m k_m \, (k_m=0,1)\)を用いて、\(\ket{k}\)の位相を

\[ \exp\left(2\pi i \frac{sk}{2^n}\right) = \exp \left(2\pi i \frac{s}{2^n} \sum_{m=0}^{n-1} 2^m k_m\right) = \prod_{m=0}^{n-1} \exp \left(2\pi i \frac{s}{2^{n-m}} k_m\right) \]

と分解します。一方

\[ \ket{k} = \ket{k_{n-1}} \otimes \ket{k_{n-2}} \otimes \cdots \ket{k_0} \]

なので、

\[ e^{2\pi i sk/2^n} \ket{k} = e^{2\pi i sk_{n-1}/2} \ket{k_{n-1}} \otimes \cdots \otimes e^{2\pi i sk_{1}/2^{n-1}} \ket{k_1} \otimes e^{2\pi i sk_{0}/2^n} \ket{k_{0}}. \]

最後に、\(k_m = 0\)ならば\(e^{2\pi i sk_{m}/2^n} = 1\)であることから、

\[ \frac{1}{\sqrt{2^n}}\sum_{k=0}^{2^n-1} e^{2\pi i sk/2^n} \ket{k} = \frac{1}{\sqrt{2}} \left(\ket{0} + e^{2\pi i s/2} \ket{1}\right) \otimes \cdots \otimes \frac{1}{\sqrt{2}} \left(\ket{0} + e^{2\pi i s/2^n} \ket{1}\right) \]

がわかります。

for m in range(num_qubits):
    circuit.h(m)
    circuit.p(2. * np.pi * s / (2 ** (num_qubits - m)), m)

が正解となります。

問題7: 量子フーリエ変換

問題

\(n\)量子ビットレジスタの状態\(\ket{j} \, (j \in \{0,1,\dots,2^n-1\})\)を以下のように変換する回路を考え、\(n=6, j=23\)のケースを実装しなさい。

\[ \ket{j} \rightarrow \frac{1}{\sqrt{2^n}}\sum_{k=0}^{2^n-1} e^{2\pi i jk/2^n} \ket{k} \]
num_qubits = 6

circuit = QuantumCircuit(num_qubits)

j = 23

# jの二進数表現で値が1になっているビットに対してXを作用させる -> 状態|j>を作る
for i in range(num_qubits):
    if ((j >> i) & 1) == 1:
        circuit.x(i)

##################
### EDIT BELOW ###
##################
# circuit.?
##################
### EDIT ABOVE ###
##################

sqrt_2_to_n = np.sqrt(2 ** num_qubits)
show_state(circuit, amp_norm=(1. / sqrt_2_to_n, r'\frac{1}{' + str(int(sqrt_2_to_n)) + '}'), phase_norm=(2 * np.pi / (2 ** num_qubits), r'\frac{2 \pi i}{' + str(2 ** num_qubits) + '}'))

この操作は量子フーリエ変換(Quantum Fourier transform, QFT)と呼ばれ、Shorの素因数分解を含め多くのアルゴリズムに応用されています[NC00]。現在知られている中で最も重要な量子サブルーチン(アルゴリズムの部品)と言っても過言ではないでしょう。

古典計算機で離散フーリエ変換\(j \rightarrow \{e^{2 \pi i j k/2^n}\}_k\)を計算するには\(\mathcal{O}(n2^n)\)回の演算が必要であることが知られています。一方、QFTは最も効率的な実装[HalesHallgren00]\(\mathcal{O}(n \log n)\)個のゲートしか利用しません(下の実装では\(\mathcal{O}(n^2)\))。つまり、QCは古典計算機に比べて指数関数的に早くフーリエ変換を実行できます。

解答

前の問題同様、基本は\(H\)ゲート\(P\)ゲートを使っていきますが、今回は\(\ket{j}\)という「入力」があるので、制御ゲートを利用する必要があります。

ビットn-1に着目しましょう。前回同様\(j_m \, (m=0,\dots,n-1, \, j_m=0,1)\)を使って\(j\)の二進数表現を\(j=\sum_{m=0}^{n-1} 2^m j_m\)としておきます。ビットn-1の初期状態は\(\ket{j_{n-1}}_{n-1}\)なので、アダマールゲートをかけると

\[ H\ket{j_{n-1}}_{n-1} = \frac{1}{\sqrt{2}} \left[\ket{0}_{n-1} + e^{2 \pi i \frac{j_{n-1}}{2}} \ket{1}_{n-1}\right] \]

です。さらに、\(j_{0}, \dots j_{n-2}\)の情報をこのビットに盛り込むためにビット0からn-2まででそれぞれ制御した\(C^m_{n-1}[P]\)を使います。ただし、かける位相は制御ビットごとに異なります。

\[\begin{split} \begin{align} & C^{n-2}_{n-1}\left[P\left(\frac{2^{n-2} \cdot 2 \pi}{2^n}\right)\right] \cdots C^{0}_{n-1}\left[P\left(\frac{2 \pi}{2^n}\right)\right] (H\ket{j_{n-1}}_{n-1}) \ket{j_{n-2}}_{n-2} \cdots \ket{j_0}_0 \\ = & \frac{1}{\sqrt{2}} \left[\ket{0}_{n-1} + \exp \left(2 \pi i \frac{\sum_{m=0}^{n-1} 2^{m} j_m}{2^n}\right) \ket{1}_{n-1}\right] \ket{j_{n-2}}_{n-2} \cdots \ket{j_0}_0 \end{align} \end{split}\]

次に、ビットn-2にも同様の操作をします。ただし、制御はビットn-3以下からのみです。

\[\begin{split} \begin{align} & C^{n-3}_{n-2}\left[P\left(\frac{2^{n-2} \cdot 2 \pi}{2^n}\right)\right] \cdots C^{0}_{n-2}\left[P\left(\frac{2 \cdot 2 \pi}{2^n}\right)\right] (H\ket{j_{n-2}}_{n-2}) \cdots \ket{j_0}_0 \\ = & \frac{1}{\sqrt{2}} \left[\ket{0}_{n-2} + \exp \left(2 \pi i \frac{2 \sum_{m=0}^{n-2} 2^{m} j_m}{2^n}\right) \ket{1}_{n-2}\right] \ket{j_{n-3}}_{n-3} \cdots \ket{j_0}_0 \end{align} \end{split}\]

これをビット0まで繰り返します。制御は常に位数の小さいビットからです。

\[\begin{split} \begin{align} & C^{n-4}_{n-3}\left[P\left(\frac{2^{n-2} \cdot 2 \pi}{2^n}\right)\right] \cdots C^{0}_{n-3}\left[P\left(\frac{2^2 \cdot 2 \pi}{2^n}\right)\right] (H\ket{j_{n-3}}_{n-3}) \cdots \ket{j_0}_0 \\ = & \frac{1}{\sqrt{2}} \left[\ket{0}_{n-3} + \exp \left(2 \pi i \frac{2^2 \sum_{m=0}^{n-3} 2^{m} j_m}{2^n}\right) \ket{1}_{n-3}\right] \ket{j_{n-4}}_{n-4} \cdots \ket{j_0}_0 \end{align} \end{split}\]
\[ \dots \]
\[\begin{split} \begin{align} & C^{0}_{1}\left[P\left(\frac{2^{n-2} \cdot 2 \pi}{2^n}\right)\right] (H\ket{j_{1}}_{1}) \ket{j_0}_0 \\ = & \frac{1}{\sqrt{2}} \left[\ket{0}_{1} + \exp \left(2 \pi i \frac{2^{n-2} \sum_{m=0}^{1} 2^{m} j_m}{2^n}\right) \ket{1}_{1}\right] \ket{j_0}_0 \end{align} \end{split}\]
\[\begin{split} \begin{align} & H\ket{j_0}_0 \\ = & \frac{1}{\sqrt{2}} \left[\ket{0}_0 + \exp \left(2 \pi i \frac{2^{n-1} \sum_{m=0}^{0} 2^{m} j_m}{2^n}\right) \ket{1}_{0}\right]. \end{align} \end{split}\]

ここまでの全ての操作を順に適用すると、

\[\begin{split} (H_0) (C^{0}_{1}[P]H_1) \cdots (C^{n-2}_{n-1}[P] \cdots C^{0}_{n-1}[P]H_{n-1}) \, \ket{j_{n-1}}_{n-1} \ket{j_{n-2}}_{n-2} \cdots \ket{j_{0}}_0 \\ = \frac{1}{\sqrt{2^n}} \left[\ket{0}_{n-1} + \exp \left(2 \pi i \frac{\sum_{m=0}^{n-1} 2^{m} j_m}{2^n}\right) \ket{1}_{n-1}\right] \cdots \left[\ket{0}_0 + \exp \left(2 \pi i \frac{2^{n-1} \sum_{m=0}^{0} 2^{m} j_m}{2^n}\right) \ket{1}_0\right]. \end{split}\]

\(\newcommand{tk}{\tilde{k}}\)

例によって全ての量子ビットに\(\ket{0}\)\(\ket{1}\)が現れるので、右辺は振幅\(\{c_k\}\)を用いて\(\sum_{k=0}^{2^n-1} c_k \ket{k}\)の形で表せます。後の利便性のために\(k\)の代わりに整数\(\tk\)とその二進数表現\(\tk_{l}\)を使って、

\[ c_{\tk} = \exp \left(\frac{2 \pi i}{2^n} \sum_{l=0}^{n-1}\sum_{m=0}^{l} 2^{n-1-l+m} j_m \tk_l \right). \]

二重和が厄介な形をしていますが、\(j_m, \tk_l \in {0, 1}\)なので、\(n-1-l+m >= n\)のとき和の中身が\(2^n\)の倍数となり、\(\exp (2 \pi i / 2^n \cdot a 2^n) = 1 \, \forall a \in \mathbb{N}\)であることを用いると、\(m\)についての和を\(n-1\)までに拡張できます。さらに\(k_l = \tk_{n-1 - l}\)と定義すると、

\[ c_{\tk} = \exp \left(\frac{2 \pi i}{2^n} \sum_{m=0}^{n-1} 2^m j_m \sum_{l=0}^{n-1} 2^l k_l \right). \]

つまり、ここまでの操作で得られる状態は

\[ \frac{1}{\sqrt{2^n}} \sum_{k_l=0,1} \exp \left(\frac{2 \pi i}{2^n} j \sum_{l=0}^{n-1} 2^l k_l \right) \ket{k_0}_{n-1} \cdots \ket{k_{n-1}}_0 \]

であり、求める状態に対してビット順序が逆転したものになっていることがわかります。したがって、最後にSWAPを使ってビット順序を逆転させれば量子フーリエ変換が完成します。

正解は

for itarg in range(num_qubits - 1, -1, -1):
    circuit.h(itarg)
    for ictrl in range(itarg - 1, -1, -1):
        power = ictrl - itarg - 1 + num_qubits
        circuit.cp((2 ** power) * 2. * np.pi / (2 ** num_qubits), ictrl, itarg)
        
for i in range(num_qubits // 2):
    circuit.swap(i, num_qubits - 1 - i)

です。

参考文献

NC00

Michael A Nielsen and Isaac L Chuang. Quantum Computation and Quantum Information, chapter 5, pages 216–247. Cambridge University Press, 2000.

HalesHallgren00

L. Hales and S. Hallgren. An improved quantum fourier transform algorithm and applications. In Proceedings 41st Annual Symposium on Foundations of Computer Science, volume, 515–525. 2000. doi:10.1109/SFCS.2000.892139.

1

(多重)制御\(P\)ゲートはどのビットを制御と呼んでも標的と呼んでも同じ作用をするので、\(C^{123}_0[P(\pi)]\)でもなんでも構いません。