$\require{cancel} \newcommand{\Ket}[1]{\left|{#1}\right\rangle} \newcommand{\Bra}[1]{\left\langle{#1}\right|} \newcommand{\Braket}[1]{\left\langle{#1}\right\rangle} \newcommand{\Rsr}[1]{\frac{1}{\sqrt{#1}}} \newcommand{\RSR}[1]{1/\sqrt{#1}} \newcommand{\Verti}{\rvert} \newcommand{\HAT}[1]{\hat{\,#1~}} \DeclareMathOperator{\Tr}{Tr}$
Quantum Phase Estimation - How¶
First created in September 2018
Ref:
Quantum Phase Estimation but we are not following its illustration, as it is confusing.
What¶
Given a unitary $U$ and a state $\Ket\psi$, such that $U\Ket\psi=e^{2\pi i\theta}\Ket\psi$, the algorithm estimates the value of $\theta$. In other words, it is to find the eigenvalue of a unitary operator, given an eigenstate.
The Circuit¶
Here is the circuit from Wikipedia. The top line is the MSB.
Please note that as we are dealing with global phases, $U^{2^n}$ needs to be consistent among each other without losing relative phase.
For example, if $U=R_x(\pi/2),~~U^2=R_x(\pi)=\begin{bmatrix}0&-i\\-i&0\end{bmatrix}\ne X.$ So the usual way of taking $X$ as turning $\pi$ about $X$-axis, i.e. $R_x(\pi)$, cannot be applied here.
Also, during calculation, no "standardisation" (making the coefficient of $\Ket0$ real and non-negative) is to be applied. e.g. $-i\Ket0$ is different from $\Ket0$.
The $QFT_n^{-1}$ (Inverse Quantum Fourier Transform) on the right is as following.
When input $\Ket{x}=\Ket1^{\otimes n}$, $\displaystyle H\Ket1 =\Rsr2\left(\Ket0-\Ket1\right) =\Rsr2\left(\Ket0+e^{2\pi i/2}\Ket1\right) .~$ We can say that $H$ is like $R_1$. So $\displaystyle R_2\Ket1=\Rsr2\left(\Ket0+e^{2\pi i\cdot 2^{-2}}\Ket1\right) .$
Generally, $R_n$ gates in $QFT_n$ are phase shift by $2\pi\cdot 2^{-n}$. In $QFT_n^{-1}$, the phase shift is negative, i.e. $\displaystyle \exp\left(-2\pi i~[0.\overbrace{0..1}^n]\right).$ So $R_2$ is $S^\dagger$, $R_3$ is $T^\dagger$ and so on.
Implementation with $U=R_x(\pi)$¶
Detailed explanation can be found in Quantum Phase Estimation - Why.
$U=R_x(\pi).$
When $\Ket\psi=\Ket+,~~R_x(\pi)\Ket\psi=-i\Ket+,~~\theta_+=0.75=[0.11].$
When $\Ket\psi=\Ket-,~~R_x(\pi)\Ket\psi=i\Ket-,~~\theta_-=0.25=[0.01].$
# Initialisation
import sys
sys.path.append('../')
from qtol import *
Case $n=1,~~N=2^n=2.~~\Ket x$ initialised to ${1\over\sqrt2}\left(\Ket0+\Ket1\right),$
# Expecting half-and-half
qbNum = 2
q = QuantumRegister(qbNum)
c = ClassicalRegister(qbNum)
qc = QuantumCircuit(q, c)
# Preparation
# \Ket x into unbiased superposition
# \Ket\psi into q_1
qc.h(q)
# \Ket-: qc.z will measure 0.1
# \Ket+: qc.iden will measure 0.0
qc.z(q[1])
# QPE
qc.h(q[1])
qc.crz(np.pi, q[0], q[1])
qc.h(q[1])
# QFT^{-1}
qc.h(q[0])
show_me(qc, q, c, show_latex=True)
qc.measure(q[0], c[0])
job = execute(qc, backend = 'local_qasm_simulator', shots=shots)
result = job.result()
print(result.get_data())
data = result.get_counts(qc)
plot_histogram(data)
circuit_drawer(qc)
Case $n=2,~~N=2^n=4.~~\Ket x$ initialised to ${1\over\sqrt2}\left(\Ket0+\Ket1\right)\otimes{1\over\sqrt2}\left(\Ket0+\Ket1\right).$
# Expecting 11 for \Ket+ and 01 for \Ket-
qbNum = 3
q = QuantumRegister(qbNum)
c = ClassicalRegister(qbNum)
qc = QuantumCircuit(q, c)
# Preparation
# \Ket x into unbiased superposition
# \Ket\psi into q_2
qc.h(q)
# \Ket-: qc.z will measure 0.01
# \Ket+: qc.iden will measure 0.11
qc.z(q[2])
# QPE
qc.h(q[2])
qc.crz(np.pi, q[0], q[2])
##qc.h(q[2])
##qc.h(q[2])
qc.crz(np.pi*2.0, q[1], q[2])
qc.h(q[2])
#show_me(qc, q, c, show_latex=True)
# QFT^{-1}
qc.h(q[1])
qc.cu1(-np.pi*0.5, q[0], q[1])
qc.h(q[0])
show_me(qc, q, c, show_latex=True)
qc.measure(q[0], c[0])
qc.measure(q[1], c[1])
job = execute(qc, backend = 'local_qasm_simulator', shots=shots)
result = job.result()
print(result.get_data())
data = result.get_counts(qc)
plot_histogram(data)
circuit_drawer(qc)
IMPORTANT: By the virtue of QFT, the qubit order is reversed. So for a $10$ output, it means eigenvalue of $0.01$.
Case $n=3,~~N=2^n=8.~~\Ket x$ initialised to ${1\over\sqrt2}\left(\Ket0+\Ket1\right)\otimes {1\over\sqrt2}\left(\Ket0+\Ket1\right)\otimes {1\over\sqrt2}\left(\Ket0+\Ket1\right).$
# Expecting 11 for \Ket+ and 01 for \Ket-
qbNum = 4
q = QuantumRegister(qbNum)
c = ClassicalRegister(qbNum)
qc = QuantumCircuit(q, c)
# Preparation
# \Ket x into unbiased superposition
# \Ket\psi into q_2
qc.h(q)
# \Ket-: qc.z will measure 0.01
# \Ket+: qc.iden will measure 0.11
qc.z(q[3])
# QPE
qc.h(q[3])
qc.crz(np.pi, q[0], q[3])
##qc.h(q[3])
##qc.h(q[3])
qc.crz(np.pi*2.0, q[1], q[3])
##qc.h(q[3])
##qc.h(q[3])
qc.crz(np.pi*4.0, q[2], q[3])
qc.h(q[3])
#show_me(qc, q, c, show_latex=True)
# QFT^{-1}
qc.h(q[2])
qc.cu1(-np.pi*0.5, q[1], q[2])
qc.cu1(-np.pi*0.25, q[0], q[2])
qc.h(q[1])
qc.cu1(-np.pi*0.5, q[0], q[1])
qc.h(q[0])
show_me(qc, q, c, show_latex=True)
qc.measure(q[0], c[0])
qc.measure(q[1], c[1])
qc.measure(q[2], c[2])
job = execute(qc, backend = 'local_qasm_simulator', shots=shots)
result = job.result()
print(result.get_data())
data = result.get_counts(qc)
plot_histogram(data)
circuit_drawer(qc)
# Expecting 11 for \Ket+ and 01 for \Ket-
qbNum = 5
q = QuantumRegister(qbNum)
c = ClassicalRegister(qbNum)
qc = QuantumCircuit(q, c)
# Preparation
# \Ket x into unbiased superposition
# \Ket\psi into q_2
qc.h(q)
# \Ket-: qc.z will measure 0.01
# \Ket+: qc.iden will measure 0.11
qc.z(q[4])
# QPE
qc.h(q[4])
qc.crz(np.pi, q[0], q[4])
##qc.h(q[3])
##qc.h(q[3])
qc.crz(np.pi*2.0, q[1], q[4])
##qc.h(q[3])
##qc.h(q[3])
qc.crz(np.pi*4.0, q[2], q[4])
##qc.h(q[3])
##qc.h(q[3])
qc.crz(np.pi*8.0, q[3], q[4])
qc.h(q[3])
#show_me(qc, q, c, show_latex=True)
# QFT^{-1}
qc.h(q[3])
qc.cu1(-np.pi*0.5, q[2], q[3])
qc.cu1(-np.pi*0.25, q[1], q[3])
qc.cu1(-np.pi*0.125, q[0], q[3])
qc.h(q[2])
qc.cu1(-np.pi*0.5, q[1], q[2])
qc.cu1(-np.pi*0.25, q[0], q[2])
qc.h(q[1])
qc.cu1(-np.pi*0.5, q[0], q[1])
qc.h(q[0])
show_me(qc, q, c, show_latex=True)
qc.measure(q[0], c[0])
qc.measure(q[1], c[1])
qc.measure(q[2], c[2])
qc.measure(q[3], c[3])
job = execute(qc, backend = 'local_qasm_simulator', shots=shots)
result = job.result()
print(result.get_data())
data = result.get_counts(qc)
plot_histogram(data)
circuit_drawer(qc)
