Quantum phase estimation

Kitaev Algorithm for Quantum Phase Estimation

In this post, we introduce the quantum phase estimation based on quantum Fourier transform. There is another way to implement quantum phase estimation via an interferometric scheme, usually called Kitaev algorithm.

As before, we assume that we can prepare an eigenstate $|\psi\rangle$ of a unitary operator $U$ ,

$U|\psi\rangle = e^{2\pi i \theta} |\psi\rangle ,$

where we take $\theta \in [-1/2,1/2)$ for convenience. Kitaev’s idea is to encode the phase information into the amplitudes of certain superposition states and extract $\theta$ using the interferometric scheme.

^{Fig 1: interferometric scheme for cosine part.}

The quantum circuit is shown in Figure 1.

Acting with a Hadamard gate on the control qubit gives

$\frac{1}{\sqrt{2}}(|0\rangle|\psi\rangle + |1\rangle|\psi\rangle).$

Applying the controlled- $U$ gate yields

$\frac{1}{\sqrt{2}}(|0\rangle|\psi\rangle + e^{2\pi i\theta}|1\rangle|\psi\rangle),$

so the phase is now encoded in the relative amplitude. A final Hadamard gate transforms the state to

$|\Phi_f^{c}\rangle = \frac{1}{\sqrt{2}}(|+\rangle|\psi\rangle + e^{2\pi i\theta}|-\rangle|\psi\rangle)$

$= \frac{1+e^{2\pi i \theta}}{2}|0\rangle|\psi\rangle + \frac{1-e^{2\pi i\theta}}{2}|1\rangle|\psi\rangle$

Thus the measurement probabilities are

$p(0) = \left|\frac{1+e^{2\pi i\theta}}{2}\right|^2 = \frac{1+\cos(2\pi\theta)}{2}$

$p(1) = \left|\frac{1-e^{2\pi i\theta}}{2}\right|^2= \frac{1-\cos(2\pi\theta)}{2}$

These measurements allow us to extract the cosine of the phase,

$\theta = \frac{1}{2\pi}\cos^{-1}(2p(0)-1)$

or equivalently,

$\theta = \frac{1}{2\pi}\cos^{-1}(1-2p(1))$

However, we know that $\cos(2\pi \theta)$ cannot distinguish between
$\theta$ and $-\theta$ .
It is therefore natural to incorporate the sine component $\sin(2\pi\theta)$ to resolve this ambiguity.
By first determining the sign of $\theta$ from the sine measurement, we can then output the final estimate as

$\theta = \pm \frac{1}{2\pi}\cos^{-1}(2p(0)-1),$

where the sign is chosen according to whether $\sin(2\pi\theta)$ is positive or negative.

^{Fig 2. interferometric scheme for sine part.}

To obtain the sine component, we slightly modify the circuit by inserting a phase gate $S$ on the control qubit prior to the controlled- $U$ operation, see Figure 2.

A parallel computation yields

$|\Phi_f^{s}\rangle = \frac{1+ie^{2\pi i\theta}}{2}|0\rangle|\psi\rangle+\frac{1-ie^{2\pi i\theta}}{2}|1\rangle|\psi\rangle.$

The measurement probabilities are now

$p(0) = \left|\frac{1+ie^{2\pi i\theta}}{2}\right|^2= \frac{1-\sin(2\pi\theta)}{2}, \quad p(1) = \left|\frac{1-ie^{2\pi i\theta}}{2}\right|^2= \frac{1+\sin(2\pi\theta)}{2}$

Therefore,

$\theta = \frac{1}{2\pi}\sin^{-1}(1-2p(0))$

or equivalently,

$\theta = \frac{1}{2\pi}\sin^{-1}(2p(1)-1)$

with $\sin^{-1}$ taking values in $[-\pi,\pi)$ .

Combining the estimates of $\cos(2\pi\theta)$ and $\sin(2\pi\theta)$ allows us to reconstruct $\theta$ uniquely in the range $[-1/2,1/2)$ .

Quantum Order-Finding Algorithm

The order-finding problem is a fundamental problem that plays a central role in several quantum algorithms; in particular, it underlies Shor’s celebrated factoring algorithm. Classically, order-finding is hard—no efficient classical algorithm is known—but quantum computers can solve it efficiently. Efficiently determining the order of an integer modulo $N$ allows integers to be factored in polynomial time, which would compromise RSA encryption, as we will discuss in the subsequent sections.

The order-finding problem is defined rigorously as follows:

Definition (Order-Finding Problem) Let $N \in \mathbf{Z}_+$ be a positive integer, and let $a \in \mathbf{Z}_N^{\times}$ be an integer coprime to $N$ (that is, the greatest common divisor $\gcd(a,N)$ of $a$ and $N$ equals one, $\gcd(a,N) = 1$ ). The notation $\mathbf{Z}_N^\times$ denotes the set of invertible elements in $\mathbf{Z}_N$ under multiplication—those integers that are coprime to $N$ . The order of $a$ modulo $N$ , denoted $\mathrm{ord}_N(a)$ or simply $r$ , is the smallest positive integer $r$ such that
$a^r \equiv 1 (\mathrm{mod}\, N)$
The order-finding problem is: given $a$ and $N$ , compute $r = \mathrm{ord}_N(a)$ .

At first glance, this definition may appear somewhat abstract.
Do not worry—shortly, we will unpack each of the mathematical ingredients and illustrate their meaning with concrete examples.

The quantum algorithm for order-finding exploits the periodic structure of the function

$f_{a,N}(x) \equiv a^x (\mathrm{mod}\, N) \ \ \ \ \ (1)$

which is periodic with period exactly $r$ . By preparing a quantum superposition over all possible exponents and applying the quantum phase estimation algorithm to the unitary operator

$U |y\rangle = |a y (\mathrm{mod}\, N)\rangle,$

we obtain an efficient quantum procedure that succeeds with high probability.

Preliminaries of Modular Arithmetic

Before presenting the complete quantum algorithm, let us first review a few essential results from number theory. This section may appear somewhat technical compared to others; if you prefer not to delve into the details, simply remember that for coprime integers $a$ and $N$ , there exists a positive order $r$ , and the function in Eq.(1) is periodic and turn to the next section.

An integer $N > 1$ is called a prime number if it has no nontrivial divisors; for example, $N = 2, 3, 5, 7, 11, \dots$ are prime numbers. Otherwise, $N$ is called a composite number. We use $a \mid b$ to denote that $b$ is divisible by $a$ , for instance $3 \mid 12$ , $2 \mid 6$ . Similarly, $a \nmid b$ means that $a$ does not divide $b$ , for example $3 \nmid 14$ .

Every integer $N > 1$ can be uniquely expressed as a product of prime powers:

$N = p_1^{e_1} \times p_2^{e_2} \times \cdots \times p_m^{e_m},$

where $p_1, \ldots, p_m$ are distinct prime numbers and $e_1, \ldots, e_m$ are positive integers. This result is known as the Fundamental Theorem of Arithmetic.

To obtain the prime factorization of a given number $N$ , one typically follows these steps:

Divide $N$ by $2$ repeatedly until the quotient is no longer divisible by $2$ .
Let $e_2$ be the number of times we divide by $2$ , and denote the resulting quotient by $N_2 = N / 2^{e_2}$ , for which $2 \nmid N_2$ .
For the quotient $N_2$ , proceed with the next prime, $3$ , and divide repeatedly until the quotient is no longer divisible by $3$ .
Let $e_3$ be the number of divisions, and denote the resulting quotient by $N_3 = N_2 / 3^{e_3}$ .
Continue in this manner with successive primes $5, 7, 11, \dots$ until the quotient becomes $1$ .

Exercise. Decompose $N = 84$ into its prime factors: $84 = 2^2 \times 3 \times 7$ .

The key operation we will use is modular arithmetic.
In everyday life we count linearly: 1, 2, 3, …, but in modular arithmetic we count on a circle.
Take a clock that shows only 12 hours. When the hand reaches 12, it jumps back to 0.
Mathematically, we say time is measured modulo 12.
This simple idea—wrapping around after a fixed modulus $N$ —is crucial in number theory and its applications in modern public-key cryptography.

Definition (Congruence) Fix an integer $N \ge 2$ , called the modulus.
Two integers $a$ and $b$ are congruent modulo $N$ if their difference is a multiple of $N$ :

$a \equiv b (\mathrm{mod}\, N) \quad \Longleftrightarrow \quad N \mid (a - b) \quad \Longleftrightarrow \quad \exists\, k \in \mathbf{Z} \text{ such that } a = b + kN.$

We read $a \equiv b (\mathrm{mod}\, N)$ as “ $a$ is congruent to $b$ modulo $N$ ”.

The congruence class of $a$ is usually represented by its smallest non-negative residue:

$a (\mathrm{mod}\, N) \in \mathbf{Z}_N := \{0, 1, 2, \dots, N-1\}$

It is also convenient to use $[a]$ to denote the congruence class whenever there is no ambiguity about the modulus.

Example. Let $N = 7$ . We have $\mathbf{Z}_7 :=\{0, 1, 2, \dots, 6\}$ . Then $1 \equiv 8 (\mathrm{mod}\, 7)$ and $2 \equiv 9 (\mathrm{mod}\, 7)$ . The congruence class $[0]$ consists of $\{0, 7, 14, 21, \dots\}$ , the congruence class $[1]$ consists of $\{1, 8, 15, 22, \dots\}$ , and $[2]$ consists of $\{2, 9, 16, 23, \dots\}$ , etc.

The order-finding problem considers the modular exponentiation for given integers $a$ and $N$ :

$a^x \equiv b (\mathrm{mod}\, N)$

The remainder $b$ can take values in $0, 1, \dots, N-1$ :

$a \equiv 0 (\mathrm{mod}\, N)$
$a \equiv 1 (\mathrm{mod}\, N)$
$\vdots$
$a \equiv N-1 (\mathrm{mod}\, N)$

However, not all values of $b = 0, 1, \dots, N-1$ are necessarily attainable.

For example, let $a = 2$ and $N = 6$ . We find that

$2^0 \equiv 1 (\mathrm{mod}\, 6)$
$2^1 \equiv 2 (\mathrm{mod}\, 6)$
$2^2 \equiv 4 (\mathrm{mod}\, 6)$
$2^3 \equiv 2 (\mathrm{mod}\, 6)$
$2^4 \equiv 4 (\mathrm{mod}\, 6)$
$2^5 \equiv 2 (\mathrm{mod}\, 6)$
$\vdots$

We observe a repeating pattern. In this case, the residues $0$ , $3$ , and $5$ never appear.

Note also that $a^0 \equiv 1 (\mathrm{mod}\, N)$ for any integers $a$ and $N$ . Hence, the order of $a$ modulo $N$ is defined as the smallest positive integer $r$ such that $a^r \equiv 1 (\mathrm{mod}\, N)$ . However, we must still ensure that such an order actually exists. In fact, if $\gcd(a, N) = d > 1$ , then no positive integer $r$ satisfies $a^r \equiv 1 (\mathrm{mod}\, N)$ . This can be seen in the previous example with $a = 2$ and $N = 6$ , where the sequence never returns to $1$ modulo $6$ .

Exercise. Let $N > 1$ and $\gcd(a, N) = d > 1$ . Then the only integer $r$ satisfying

$a^r \equiv 1 (\mathrm{mod}\, N)$

is $r = 0$ .

Proof. Suppose there exists an integer $r$ such that $a^r \equiv 1 (\mathrm{mod}\, N)$ .
Then $N$ divides $a^r - 1$ , so in particular

$d \mid (a^r - 1)$ .

Since $d \mid a$ , we also have $d \mid a^r$ , and therefore

$d \mid (a^r - 1) \implies d \mid 1$ ,

which contradicts $d > 1$ . Hence no such $r \neq 0$ exists.

For $r = 0$ , the convention $a^0 = 1$ gives

$a^0 = 1 \equiv 1 (\mathrm{mod}\, N)$

trivially. Thus $r = 0$ is the only solution ∎.

In Shor’s algorithm, we either assume $\gcd(a, N) = 1$ , or first compute the greatest common divisor classically, which already yields a nontrivial factor of $N$ if $d > 1$ .

Another useful concept we will use is Euler’s Totient function: $\varphi(N)$ counts the number of integers $k \in {1, 2, \dots, N-1}$ coprime to $N$ .
For example, $\varphi(1) = 1$ , $\varphi(2) = 1$ , $\varphi(3) = 2$ , $\varphi(4) = 2$ , etc.

For fixed $N$ , consider $\mathbf{Z}_N$ .
We call $a \in \mathbf{Z}_N$ invertible if there exists an integer $c$ such that $a c \equiv 1 (\mathrm{mod}\, N)$ .

Proposition (Bézout’s Identity). Let $a$ and $b$ be integers with greatest common divisor $d = \gcd(a,b)$ . Then there exist integers $x$ and $y$ such that

$a x + b y = d$ .

Moreover, the integers of the form $a z + b t$ are exactly the multiples of $d$ .

The set $\mathbf{Z}_N$ is not a group under multiplication modulo $N$ in general. We denote by $(\mathbf{Z}_N^{\times}, \times)$ the set of invertible elements, which indeed forms a group. By Bézout’s identity, $\mathbf{Z}_N^{\times}$ consists of all integers coprime to $N$ , and this group has order $\phi(N)$ , where $\phi$ denotes Euler’s totient function.

From a group-theoretic perspective, for any $a \in \mathbf{Z}_N^{\times}$ , there exists a positive integer $r$ such that

$a^r \equiv 1 (\mathrm{mod}\, N)$ .

Since $a$ must be invertible, it follows that $\gcd(a, N) = 1$ . In the multiplicative group $\mathbf{Z}_N^{\times}$ , the order of any element divides the order of the group, that is,

$r \mid |\mathbf{Z}_N^{\times}| = \phi(N)$ .

Theorem (Euler’s Theorem). If $\gcd(a, N) = 1$ , then

$a^{\phi(N)} \equiv 1 (\mathrm{mod}\, N)$

and hence the order of $a$ modulo $N$ divides $\phi(N)$ .

By Euler’s theorem, for integers $a$ and $N$ with $\gcd(a, N) = 1$ , there always exists a positive integer $r$ such that

$a^r \equiv 1 (\mathrm{mod}\, N).$

The order $r = \mathrm{ord}_N(a)$ always satisfies $1 \le r \le \phi(N)$ .

Since in the cyclic subgroup generated by $[a] \in \mathbf{Z}_N^{\times}$ we have

${[a^0], [a^1], [a^2], \dots, [a^{r-1}]}, \quad \text{with } [a^r] = [a^0] = [1],$

the sequence $a^n (\mathrm{mod}\, N)$ is periodic with period equal to the order $r$ :

$a^0 \equiv 1 (\mathrm{mod}\, N)$
$a^1 \equiv a (\mathrm{mod}\, N)$
$\dots$
$a^{r-1} \not\equiv 1 (\mathrm{mod}\, N)$
$a^r \equiv 1 (\mathrm{mod}\, N)$
$a^{r+1} \equiv a (\mathrm{mod}\, N)$
$a^{r+2} \equiv a^2 (\mathrm{mod}\, N)$
$a^{r+3} \equiv a^3 (\mathrm{mod}\, N)$

Thus, the sequence repeats with period $r$ .

Example. Set $a = 2$ and $N = 7$ . We have $f_{a,N}(x) = 2^x (\mathrm{mod}\, 7)$ as:

$2^0 \equiv 1 (\mathrm{mod}\, 7)$
$2^1 \equiv 2 (\mathrm{mod}\, 7)$
$2^2 \equiv 4 (\mathrm{mod}\, 7)$
$2^3 \equiv 8 \equiv 1 (\mathrm{mod}\, 7)$
$2^4 \equiv 16 \equiv 2 (\mathrm{mod}\, 7)$
$2^5 \equiv 32 \equiv 4 (\mathrm{mod}\, 7)$
$2^6 \equiv 64 \equiv 1 (\mathrm{mod}\, 7)$
$2^7 \equiv 128 \equiv 2 (\mathrm{mod}\, 7)$
$2^8 \equiv 256 \equiv 4 (\mathrm{mod}\, 7)$

We see that the order is $r = \mathrm{ord}_7(2) = 3$ .

Quantum Order-Finding Algorithm

The quantum order-finding algorithm can be divided into two main steps:

Encode the order into the phase of an eigenvalue of a suitable unitary operator, thereby reducing the order-finding problem to a quantum phase estimation problem.
Recover the order $r$ from the estimated phase $\varphi$ using the method of continued fractions.

Let us first examine how to encode the order into the phase of a unitary operator.

Fix integers $a$ and $N$ such that $\gcd(a, N) = 1$ .
Define the unitary operator $U_{a,N}$ acting on computational basis states $|y\rangle$ by

$U_{a,N} |y\rangle = | a y (\mathrm{mod}\, N) \rangle, \quad y = 0, 1, \dots, N-1$

We will see that the order can be mapped to the phase of an eigenvalue of $U_{a,N}$ .

Proposition. Let $N \ge 2$ and let $a$ be an integer coprime to $N$ .
The operator $U_{a,N}$ on $\mathbf{C}^N$ defined by

$U_{a,N} |y\rangle = | a y (\mathrm{mod}\, N) \rangle, \quad y = 0, 1, \dots, N-1$

is unitary.

Proof. Since $\gcd(a, N) = 1$ , multiplication by $a$ modulo $N$ defines a bijection on the set ${0, 1, \dots, N-1}$ .
Thus $U_{a,N}$ merely permutes the computational basis ${ |0\rangle, |1\rangle, \dots, |N-1\rangle }$ . Any permutation of an orthonormal basis extends to a unitary operator, so $U_{a,N}$ is unitary.

Alternatively, suppose $U_{a,N} |y\rangle = U_{a,N} |y'\rangle$ . Then

$a y \equiv a y' (\mathrm{mod}\, N) \implies a (y - y') \equiv 0 (\mathrm{mod}\, N) \implies N \mid a (y - y')$ .

Because $\gcd(a, N) = 1$ , we may cancel $a$ (more precisely, multiply both sides by the modular inverse of $a$ modulo $N$ ) to obtain

$N \mid (y - y')$ .

Since $0 \le y, y' < N$ , the only possibility is $y - y' = 0$ , so $y = y'$ .
Thus $U_{a,N}$ is injective on a basis of $\mathbf{C}^N$ .
Being an isometry that maps a basis to a basis, it is unitary ∎.

Proposition. Let $r$ be the order of $a$ modulo $N$ , i.e., the smallest positive integer such that

$a^r \equiv 1 (\mathrm{mod}\, N)$ .

For each $s = 0, 1, \dots, r-1$ , define the state

$|u_s\rangle := \frac{1}{\sqrt{r}} \sum_{k=0}^{r-1} e^{-2\pi i s k / r} \, | a^k (\mathrm{mod}\, N) \rangle.$

Then, for each $s = 0, 1, \dots, r-1$ , $|u_s\rangle$ is an eigenstate of $U_{a,N}$ with eigenvalue

$e^{2\pi i s / r}$ .

Proof. Apply $U_{a,N}$ to $|u_s\rangle$ :

$U_{a,N} |u_s\rangle = \frac{1}{\sqrt{r}} \sum_{k=0}^{r-1} e^{-2\pi i s k / r} \, U_{a,N} | a^k (\mathrm{mod}\, N) \rangle$

$\quad\; = \frac{1}{\sqrt{r}} \sum_{k=0}^{r-1} e^{-2\pi i s k / r} \, | a^{k+1} (\mathrm{mod}\, N) \rangle$

Reindex the sum by letting $k \mapsto k+1$ :

$\quad\; = \frac{1}{\sqrt{r}} \sum_{k=1}^{r} e^{-2\pi i s (k-1) / r} \, | a^k (\mathrm{mod}\, N) \rangle$

$\quad\; = \frac{1}{\sqrt{r}} \sum_{k=1}^{r} e^{-2\pi i s k / r} \, e^{2\pi i s / r} \, | a^k (\mathrm{mod}\, N) \rangle$

$\quad\; = e^{2\pi i s / r} \frac{1}{\sqrt{r}} \sum_{k=1}^{r} e^{-2\pi i s k / r} \, | a^k (\mathrm{mod}\, N) \rangle$

The sum now runs from $k = 1$ to $r$ , but the $k = r$ term is

$e^{-2\pi i s r / r} \, | a^r (\mathrm{mod}\, N) \rangle = e^{-2\pi i s} \, | 1 (\mathrm{mod}\, N) \rangle = | 1 (\mathrm{mod}\, N) \rangle,$

while the $k = 0$ term in the original definition of $|u_s\rangle$ is

$e^0 \, | a^0 (\mathrm{mod}\, N) \rangle = | 1 (\mathrm{mod}\, N) \rangle.$

Thus the missing $k = 0$ term is exactly replaced by the $k = r$ term, and we obtain

$U_{a,N} |u_s\rangle = e^{2\pi i s / r} \frac{1}{\sqrt{r}} \sum_{k=0}^{r-1} e^{-2\pi i s k / r} \, | a^k (\mathrm{mod}\, N) \rangle = e^{2\pi i s / r} |u_s\rangle.$

This holds for every $s = 0, 1, \dots, r-1$ ∎.

Now the procedure becomes “straightforward”. We apply quantum phase estimation (see this post) to the operators $U_a^{2^j}$ with the eigenstate $|u_s\rangle$ , which yields the phases

$\theta_s = \frac{s}{r}, \quad s = 0, 1, \dots, r-1$

The output of the phase estimation algorithm is thus

$0.\theta_0 \theta_1 \dots \theta_{n}$

represented as binary numbers. A classical algorithm called continued fractions can be applied to obtain the order $r$ .

A key application of the quantum Fourier transform (see this post) is the Quantum Phase Estimation algorithm. This algorithm allows one to estimate the phase (or eigenvalue) of a unitary operator and serves as a central component of many quantum algorithms, including Shor’s factoring algorithm and quantum simulation.

Consider a unitary matrix $U$ . Since $U^\dagger U = I$ , all eigenvalues of $U$ have the form $e^{2\pi i \theta}$ with $0 \leq \theta < 1$ . Expressing $\theta$ in binary gives

$\theta = 0.\theta_1 \theta_2 \theta_3 \cdots = \sum_{k=1}^\infty \theta_k \, 2^{-k},$

where each $\theta_k \in {0,1}$ . The task is to estimate $\theta$ to $n$ binary digits of precision.

More precisely, the quantum phase estimation algorithm aims to determine the phase $\theta$ in the eigenvalue equation

$U |\psi\rangle = e^{2\pi i \theta} |\psi\rangle$ ,

where $U$ is the given unitary operator, $|\psi\rangle$ is an eigenstate of $U$ that can be prepared as an input state, and $\theta$ is the phase to be determined. We further assume that controlled- $U^{2^l}$ operations can be implemented for arbitrary non-negative integer $l$ . The goal is to estimate $\theta$ to the desired number of digits with high probability.

Now suppose that we express $\theta$ in $n$ binary digits,

$\theta = 0.\theta_1\theta_2\cdots \theta_n$ .

Then, we have

$U^{2^l} \lvert \psi \rangle = e^{2 \pi i 0.\theta_1\theta_2\cdots\theta_n 2^l} \lvert \psi \rangle = e^{2 \pi i \theta_1\cdots\theta_l.\theta_{l+1}\cdots\theta_n} \lvert \psi \rangle$

where the integer part of $\theta_1\cdots \theta_l.\theta_{l+1}\cdots \theta_n$ can be omitted, as it contributes only a trivial phase factor of $1$ .

For the controlled unitary $C(U^{2^l})$ , its action on $|+\rangle|\psi\rangle$ (with $|+\rangle$ the control qubit) is:

$C(U^{2^l}) \frac{1}{\sqrt{2}} ( \lvert 0 \rangle + \lvert 1 \rangle ) \lvert \psi \rangle = \frac{1}{\sqrt{2}} \lvert 0 \rangle \otimes \lvert \psi \rangle + e^{2 \pi i 0.\theta_{l+1}\cdots \theta_n} \lvert 1 \rangle \otimes \lvert \psi \rangle$

$= \frac{1}{\sqrt{2}} (|0\rangle + e^{2\pi i 0.\theta_{l+1}\cdots \theta_n} |1\rangle) \otimes |\psi\rangle$

It is convenient to view this control qubit as representing the $(l+1)$-th digit of a binary number $k = k_1 k_2 \cdots k_n = k_1 \times 2^0 + k_2 \times 2^1 + \cdots + k_n \times 2^{n-1}$ .

With this notation, the above expression can be rewritten as

$C(U^{2^l}) \frac{1}{\sqrt{2}} (\lvert0\rangle + \lvert1\rangle)\lvert\psi\rangle = \frac{1}{\sqrt{2}} \sum_{k_{l+1}=0}^1 e^{2 \pi i \theta k_{l+1} 2^l} \lvert k_{l+1}\rangle \lvert\psi\rangle$

This expression for the first qubit is familiar from our earlier discussion of the quantum Fourier transform, and it naturally suggests that the quantum Fourier transform should be applied here. The eigenstate $\lvert \psi \rangle$ therefore remains unchanged and can be reused in subsequent steps.

Fig 1: quantum circuit for quantum phase estimation algorithm.

The quantum phase estimation algorithm uses two quantum registers as shown in Figure 1:

First register (phase record): Contains $n$ qubits, which will store the estimated phase.
Second register (eigenstate): Contains the eigenstate $\lvert \psi \rangle$ , which is an eigenstate of the unitary operator $U$ .

The quantum phase estimation algorithm consists of the following steps:

Prepare Initial State

The system starts with two registers:

The first register is initialized in the state $\lvert 0 \rangle^{\otimes n}$ .
The second register is initialized in the eigenstate $\lvert \psi \rangle$ of $U$ .

The initial state is: $\lvert 0 \rangle^{\otimes n} \otimes \lvert \psi \rangle$ .

Apply Hadamard Gates

Apply Hadamard gates to each qubit in the first register, thereby creating a uniform superposition over all possible computational basis states:

$\frac{1}{\sqrt{2^n}} \sum_{k=0}^{2^n-1} \lvert k \rangle \otimes \lvert \psi \rangle$

Here, $k = k_1 k_2 \cdots k_n$ denotes the binary representation of the integer $k$ , with the bits encoded from the bottom to the top qubit in the first register.

Apply Controlled- $U^{2^l}$ Operations

For each qubit $k_l$ in the first register, apply the controlled- $U^{2^l}$ operation.
This step entangles the first and second registers, yielding the state:

$\frac{1}{\sqrt{2^n}} \sum_{k=0}^{2^n-1} \lvert k \rangle \otimes U^k \lvert \psi \rangle$

Since $U \lvert \psi \rangle = e^{2 \pi i \theta} \lvert \psi \rangle$ , this becomes:

$\frac{1}{\sqrt{2^n}} \sum_{k=0}^{2^n-1} e^{2 \pi i k \theta} \lvert k \rangle \otimes \lvert \psi \rangle$

At this stage, the second register is no longer needed, as all phase information is now encoded in the first register.

Apply the Inverse Quantum Fourier Transform

Apply the inverse quantum Fourier transform ( $U_{\rm QFT}^\dagger$ ) on the first register to extract the phase information:

$\frac{1}{\sqrt{2^n}} \sum_{k=0}^{2^n-1} e^{2 \pi i k \theta} \lvert k \rangle$

Measure the First Register

Finally, measure the qubits in the first register. The measurement will yield the binary representation of the phase $\theta$ with an $n$ -bit approximation.
The output will be an $n$ -bit approximation of the phase $\theta$ , such that:

$\theta \approx 0.\theta_1 \theta_2 \dots \theta_n$

where $\theta_1, \theta_2, \dots, \theta_n$ are the binary digits of $\theta$ .