# Introduction

The quantum Fourier transform and quantum phase estimation provide the foundation for many quantum algorithms, including the quantum matrix inversion, which is extensively used in quantum machine learning. It is therefore worthwhile developing a good understanding of these building blocks before moving on to more complex algorithms.

We must emphasize that starting with this notebook, the algorithms presented are coherent quantum protocols. By that, we mean that the input and output of an algorithm is a quantum state that we do not have classical information about. The protocol itself might use measurements: in this sense, they are not fully coherent, since we gain some, but incomplete classical information about the quantum system. We might also perform post-selection, which means that a gate is controlled by the classical output of a measurement. In some cases, we entirely discard a calculation based on a measurement output.

Why does it matter that we begin and end up with quantum states? Can't we just use state preparation starting from classical data and then perform tomography on the final state? We could do that, but state preparation and tomography are resource-intensive, and they are likely to destroy any quantum advantage.

An additional problem is that the quantum Fourier transformation and other quantum algorithms similar in complexity require a very large number of gates on a large number of high-quality qubits. This is why the practical relevance of these algorithms is not immediate, but since they are the core of many quantum machine learning protocols, it is essential that we take a look at them.


```python
import numpy as np
from qiskit import QuantumCircuit, ClassicalRegister, QuantumRegister
from qiskit import execute
from qiskit import Aer
π = np.pi
```

# Quantum Fourier Transform

The quantum Fourier transform is a quantum algorithm for the discrete Fourier transform over the amplitudes of a wavefunction. The exposition here follows the introduction in [[1](#1)]. A similar approach can be found in the [Qiskit tutorials](https://github.com/Qiskit/qiskit-tutorials/blob/master/community/algorithms/shor_algorithm.ipynb).

The classical discrete Fourier transform acts on a vector $\vec{x}=\begin{bmatrix}x_0\\ \vdots\\ x_{N-1}\end{bmatrix}$ and maps it to the vector $\vec{y}=\begin{bmatrix}y_0\\ \vdots\\ y_{N-1}\end{bmatrix}$, where $y_k = \frac{1}{\sqrt{N}}\sum_{j=0}^{N-1}x_je^{ \boldsymbol{2\pi i} \frac{jk}{N}}$.

The quantum Fourier transform acts on an amplitude-encoded variant of this vector, the quantum state $|x\rangle=\sum_{i=0}^{N-1} x_i |i \rangle$ and maps it to the quantum state $|y\rangle=\sum_{k=0}^{N-1} y_k |k \rangle$, where
$y_k = \frac{1}{\sqrt{N}}\sum_{j=0}^{N-1}x_je^{\boldsymbol{2\pi i} \frac{jk}{N}}$. Since the transformed state is also in the superposition of the computational basis, in essence, only the amplitudes are transformed.

To derive a circuit for a power of two $N=2^n$, consider the transform acts on the state $| x \rangle = | x_1...x_n \rangle$ where $x_1$ is the most significant bit, unlike the usual convention followed in the course. We will rewrite a number $y$ in the fractional binary notation  as $j = 0.j_1...j_n = \sum_{k=1}^n j_k/2^k$. For example $0.8125_d = 0.1101_b = \sum_{k=1}^4 j_k/2^k = 1/2 + 1/4 + 0/8 + 1/16$, where the subscripts $d$ stands for decimal and $b$ for binary. The action of the unitary $U$ describing the transform can be expanded as
\begin{aligned}
U |x \rangle = U |x_1 x_2 \cdots x_n \rangle& = \frac{1}{\sqrt{N}} \sum_{k=0}^{2^n-1} e^{\boldsymbol{2\pi i} xk / 2^n} |k \rangle \\
& = \frac{1}{\sqrt{N}} \sum_{k_1=0}^{1}\ldots\sum_{k_n=0}^{1} e^{\boldsymbol{2\pi i} x\left(\sum_{l=1}^n k_l2^{-l}\right) } \vert k_1 ... k_n \rangle \\
& = \frac{1}{\sqrt{N}} \sum_{k_1=0}^{1}\ldots\sum_{k_n=0}^{1} \bigotimes_{l=1}^n e^{\boldsymbol{2\pi i} x k_l2^{-l}} | k_1 ... k_n \rangle \\
& = \frac{1}{\sqrt{N}} \bigotimes_{l=1}^n \sum_{k_l=0}^{1}  e^{\boldsymbol{2\pi i} x k_l2^{-l}} | k_1 ... k_n \rangle \\
& = \frac{1}{\sqrt{N}} \bigotimes_{l=1}^n  \left(|0\rangle + e^{\boldsymbol{2\pi i} x 2^{-l} } |1\rangle \right) \\
& = \frac{1}{\sqrt{N}} \left(|0\rangle + e^{\boldsymbol{2\pi i} 0.x_n} |1\rangle\right) \otimes...\otimes  \left(\vert0\rangle + e^{\boldsymbol{2\pi i} 0.x_1.x_2...x_{n-1}.x_n} |1\rangle\right) 
\end{aligned}

This form of the QFT is useful for deriving a circuit, since only the last qubit depends on the the
values of all the other input qubits. The remaining qubits depend less and less on the input qubits. The simple structure also allows to decompose the unitary as Hadamard gates and rotations. On three qubits, we can define the circuit as follows:


```python
q = QuantumRegister(3, 'q')
c = ClassicalRegister(1, 'c')
qft = QuantumCircuit(q, c)
qft.h(q[0])
qft.cu1(π/2, q[1], q[0])
qft.h(q[1])
qft.cu1(π/4, q[2], q[0])
qft.cu1(π/2, q[2], q[1])
qft.h(q[2])
```

We can plot it to make the qubit dependencies more apparent:


```python
from qiskit.tools.visualization import circuit_drawer
circuit_drawer(qft)
```

The conditional rotations dominate the complexity, which scales as $O(N^2)$.

# Quantum phase estimation

The goal of a quantum phase estimation algorithm is, given a unitary operator $U$ and an eigenvector $|\psi\rangle$ of $U$, to estimate $\theta$ in  $U|\psi \rangle =e^{2 i \pi \theta}|\psi \rangle$. Since $U$ is unitary, all of its eigenvalues have an absolute value of 1. By convention, $\theta$ is taken to be in $[0,1]$ and is called the *phase* of $U$ associated to $|\psi\rangle$.

The eigenvector $|\psi\rangle$ is encoded in one set of quantum registers. An additional set of $n$ qubits forms an ancilla register. At the end of the procedure, this ancilla register should contain an approximation of the [binary fraction](https://www.electronics-tutorials.ws/binary/binary-fractions.html) associated to $\theta$, with n-bits precision. A critical element is the ability to perform the controlled unitary $C-U^{2^k}$ -- it is usually assumed to be provided to the phase estimation protocol.

First, the uniform superposition is prepared in the ancilla register via the application of Hadamard gates $H$. These qubits will act as controls for the unitary operators at different time steps. Our goal is to create a superposition of $U$ as the unitary is applied for different durations. Since the eigenvalues are always situated on the complex unit circle, these differently evolved components in the superposition help reveal the eigenstructure. Given that the ancilla register we have a superposition of all possible time steps between $0$ and $2^{n-1}$, we will end up with a superposition of all possible evolutions to encode binary representations of the eigenvalues. At the end of this procedure, we have the state 
\begin{aligned}
& \frac {1}{2^{\frac {n}{2}}} (|0\rangle+{e^{2 i \pi \theta \cdot 2^{n-1}}}|1\rangle ) \otimes \cdots (|0\rangle+{e^{2 i \pi \theta \cdot 2^{1}}}|1\rangle ) \otimes (|0\rangle+{e^{2i \pi \theta \cdot 2^{0}}}|1\rangle ) = \\
&\frac {1}{2^{\frac {n}{2}}}\sum _{k=0}^{2^{n}-1}e^{2 i \pi \theta k}|k\rangle
\end{aligned}
in the ancilla. To write the ancilla in this form, we exploit that the controlled unitary operations when applied, introduce a global phase, and it is this global phase that we see in the ancilla. This phenomenon is also known as the phase kickback.

As a final step, we apply an inverse Fourier transform on the ancilla. Measuring out in the computational basis, we get the phase in the ancilla register:

\begin{align}
\frac {1}{2^{\frac {n}{2}}}\sum _{k=0}^{2^{n}-1}e^{2i \pi \theta k}|k\rangle \otimes | \psi \rangle \xrightarrow{\mathcal{QFT}_n^{-1}} | \tilde{\theta} \rangle \otimes | \psi \rangle
\end{align}

where $\tilde{\theta}$ is the n-bits approximation of the binary fraction representing $\theta$.

The circuit for phase estimation is the following:



As a toy example, let's take the $2\times 2$ unitary matrix $\begin{bmatrix}e^{0} & 0 \\0 & e^{i \pi}\end{bmatrix}=\begin{bmatrix}1 & 0 \\0 & -1\end{bmatrix}$, which has the eigenvectors $|0\rangle$ and $|1\rangle$, and phases $\theta_0=0$ and $\theta_1=\frac{1}{2}$. Therefore, the $C-U^{2^k}$ gate is a controlled-$Z$ gate for $k=0$ and the identity for $k\geq 1$.

Starting with $|\psi\rangle=|0\rangle$ in the main register, we prepare the superposition in the ancilla:


```python
q = QuantumRegister(3, 'q')
c = ClassicalRegister(2, 'c')

qpe = QuantumCircuit(q, c)
qpe.h(q[0])
qpe.h(q[1])
```

Next we perform the controlled unitary operations:


```python
# Controlled-U0
qpe.cz(q[1], q[2])
# Controlled-U1
# nothing: identity
```

We apply quantum inverse Fourier transformation to write the phase to the ancilla register:


```python
qpe.swap(q[0], q[1])
qpe.h(q[1])
qpe.cu1(-π / 2, q[0], q[1])
qpe.h(q[0])
qpe.swap(q[0], q[1])
```

We will get the result from the two first registers


```python
qpe.measure(q[0], c[0])
qpe.measure(q[1], c[1])
```

We can plot the circuit:


```python
circuit_drawer(qpe)
```

Let's now test our circuit:


```python
backend = Aer.get_backend('qasm_simulator')
job = execute(qpe, backend, shots=1000)
result = job.result()
result.get_counts(qpe)
```

As expected the result is $|2 \cdot \theta_0\rangle=|2\cdot 0\rangle=|00\rangle$

Let's now run the circuit for the eigenvector $|1\rangle$:


```python
qpe = QuantumCircuit(q, c)
qpe.h(q[0])
qpe.h(q[1])
qpe.x(q[2]) # create |1> in the main register

qpe.cz(q[1], q[2])

qpe.swap(q[0], q[1])
qpe.h(q[1])
qpe.cu1(-π / 2, q[0], q[1])
qpe.h(q[0])
qpe.swap(q[0], q[1])

qpe.measure(q[0], c[0])
qpe.measure(q[1], c[1])
```


```python
backend = Aer.get_backend('qasm_simulator')
job = execute(qpe, backend, shots=1000)
result = job.result()
result.get_counts(qpe)
```

The result should be $|10\rangle$. Indeed, $10 \rightarrow 1\cdot 2^{-1} + 0 \cdot 2^0=\frac{1}{2}=\theta_1$

# References
[1] M. Nielsen, I. Chuang. (2000). *Quantum Computation and Quantum Information*. Cambridge University Press. <a id='1'></a>

ceb28eac0d165fda097a705966928cabb887635f
ipynb
Jupyter Notebook
qiskit_version/13_Quantum_Phase_Estimation.ipynb
pedrorohde/qml-mooc
430a488b9fa16f89b36545bc61769c87997087b5
MIT
2020-03-29T15:57:03.000Z
2022-03-22T06:59:18.000Z
qiskit_version/13_Quantum_Phase_Estimation.ipynb
pedrorohde/qml-mooc
430a488b9fa16f89b36545bc61769c87997087b5
MIT
2020-07-13T19:31:15.000Z
2021-01-30T19:25:28.000Z
qiskit_version/13_Quantum_Phase_Estimation.ipynb
pedrorohde/qml-mooc
430a488b9fa16f89b36545bc61769c87997087b5
MIT
2020-02-20T09:53:42.000Z
2022-03-22T06:59:19.000Z
Qwen/Qwen-72B
1. YES
2. YES