Quantum phase estimation is one of the most important subroutines in quantum computation. It serves as a central building block for many quantum algorithms and implements a measurement for essentially any Hermitian operator. Recall that a quantum computer initially only permits us to measure individual qubits. If we want to measure a more complex observable, such as the energy described by a Hamiltonian \(H\), we resort to quantum phase estimation. The routine prepares an eigenstate of the Hermitian operator in one register and stores the corresponding eigenvalue in a second register.

Quantum phase estimation is a discretization of von Neumann’s prescription to measure a Hermitian observable \(H = \sum_j E_j |\psi_j \rangle \langle \psi_j|\). The scheme that von Neumann envisioned is the following. We consider a quantum system that supports the observable \(H\), which we want to measure. We assume that we are only able to measure simpler observables, in our case single qubits, or as in the original setting the location of a single particle. It is therefore our goal to reduce the measurement of the complex observable \(H\) to a measurement of the simpler observable, e.g. the location. This simple observable is then referred to as the pointer. To map the complex observable on to the simpler one we’ll make use of a convenient observation from quantum mechanics. It is known that the momentum operator \(\hat{p}\) generate shifts for single particles.

That is, if we apply the unitary \(\exp(- i \hat{p} \lambda)\) to some wave packet \(\psi(x)\), then this wave packet will be shifted by \(\lambda\) in the positive direction.

The scheme now assumes that we can apply the unitary evolution
\(\exp(-i H \otimes \hat{p} t )\) to both the system and the
pointer register as illustrated in the following picture

This picture essentially describes von Neumann”s measurement
scheme. We now follow the steps and first adjoin an ancilla – the
pointer – which is a continuous quantum variable initialized in the
state \(|0\rangle\) (the origin), so that the system+pointer is
initialized in the state \(|\psi\rangle|0\rangle\), where
\(|\psi\rangle\) is the initial state of the system. Then we evolve
according to the new Hamiltonian \(K = H\otimes\hat{p}\) for a time
\(t\), so the evolution is given by

\(e^{-it H\otimes \hat{p}} = \sum_{j=1}^{2^N} |\psi_{j}\rangle\langle \psi_{j}|\otimes e^{-itE_j \hat{p}}.\)

We now observe the action of this measurement apparatus. Suppose that \(|\psi\rangle\) is an eigenstate \(|\psi_{j}\rangle\) of \(H\), we find that the system evolves to \(e^{-it H\otimes \hat{p}}|\psi_{j}\rangle|0\rangle = |\psi_{j}\rangle |x = tE_j\rangle.\) A measurement of the position of the pointer with sufficiently high accuracy will provide an approximation to \(E_j\).

To carry out the above operation efficiently on a quantum computer, we
discretize the pointer using \(r\) qubits, replacing the continuous
quantum variable with a \(2^r\)-dimensional space, where the
computational basis states \(|z\rangle\) of the pointer represent
the basis of momentum eigenstates of the original continuous quantum
variable. The label \(z\) is the binary representation of the integers
\(0\) through \(2^r-1\). In this representation, the discretization of
the momentum operator becomes

\(\hat{p} = \sum_{j=1}^r 2^{-j}
\frac{\mathbb{I}-\sigma^z_j}{2}.\)
xx

With this normalization \(\hat{p}|z\rangle =
\frac{z}{2^r}|z\rangle\). Now the discretized Hamiltonian \(K =
H\otimes \hat{p}\) is a sum of terms involving at most \(k+1\)
qubits, if \(H\) is a Hamiltonian involving terms of at most \(k\)
qubits. Thus we can simulate the dynamics of \(K\) using standard
methods. In terms of the momentum eigenbasis the initial (discretized)
state of the pointer is written \(| x=0\rangle =
\frac{1}{2^{r/2}}\sum_{z=0}^{2^r-1} |z\rangle\). This state can
be prepared efficiently on a quantum computer by first initializing
the qubits of the pointer in the state \(|0\rangle \cdots
|0\rangle\) and applying an (inverse) quantum Fourier
transform.
Since we have a very simple initial state, the Fourier transform can
be represented by a product of Hadamard matrices. The discretized
evolution of the system+pointer now can be written

\(e^{-it H\otimes \hat{p}}|\psi_{j}\rangle|x=0\rangle =
\frac{1}{2^{r/2}}\sum_{z=0}^{2^r-1} e^{-iE_j
zt/2^r}|\psi_{j}\rangle z\rangle.\)

Performing an inverse quantum Fourier transform on the pointer leaves
the system in the state \(|\psi_{j}\rangle\otimes|\phi\rangle\),
where

\(| \phi\rangle = \sum_{x=0}^{2^r-1} \left( \frac{1}{2^{r}}\sum_{z=0}^{2^r-1}e^{\frac{2\pi i}{2^r}\left(x-\frac{E_j t}{2\pi}\right)z} \right)|x\rangle,\)

which is strongly peaked near the location \(x = \lfloor
\frac{E_jt}{2\pi} \rfloor\). To ensure that there are no overflow
errors we need to choose \(t < \frac{2\pi}{\|H\|}\). (We assume
here, for simplicity, that \(H\geq 0\).) It is easy to see that
actually performing the simulation of \(K\) for \(t=1\) is a product of
\(r\) simulations of the evolution according to \(\frac{1}{2^{r}}
H\otimes \frac{\mathbb{I}-\sigma^z_k}{2}\) for \(1, 2, 2^2,
\ldots, 2^{r-1}\) units of time, respectively. This results in the
general circuit for quantum phase estimation:

In order to implement the full circuit on a quantum computer, we still need to decompose the controlled unitaries \(e^{-i H \frac{t}{2^k}}\) as well as the inverse quantum Fourier transform denoted by \(QFT^{-1}\) into our elementary gates.

The example below demonstrates quantum phase estimation for a toy
single-qubit Hamiltonian \(\sigma^x\) acting on qubit \(Q_2\). Qubit
\(Q_3\) serves as a pointer system. In this example the quantum
Fourier transform on the pointer system is equivalent to the Hadamard
gate \(H\) on \(Q_3\). The discretized evolution of the
system+pointer system is described by the CNOT gate. The final
measurement outcome on the pointer qubit \(Q_3\) is \(0\) or \(1\)
depending on whether \(Q_2\) is prepared in the \(+1\) or \(-1\)
eigenstate of \(\sigma^x\). In this example, qubit \(Q_2\) is
initialized in a state \(Z H|0\rangle\) which is \(-1\) eigenvector
\(\sigma^x\). Accordingly, the measurement outcome is \(1\).