Chapter 11

Fourier sampling, phase estimation, and periods

The powers of 2 modulo 21 run $1,2,4,8,16,11$ and then repeat, so their period is 6. A quantum circuit can sample information about such a period with a Fourier transform, but one run returns only a single integer. What does that integer reveal about the period, and what classical work is still left to do?

Before you start: Controlled gates and reversible evaluation of a function (Chapter 4); conditional states after measuring one register (Chapter 2); phases written by an oracle and read by a change of basis (Chapter 10). Finite geometric sums, and remainders modulo an integer. Notation guide.

What a Fourier measurement can reveal

Let $N=2^n$, and identify $n$-bit computational-basis strings with the integers $0,\ldots,N-1$. Our positive-sign Fourier transform is

$$F_N|x\rangle=\frac1{\sqrt N}\sum_{k=0}^{N-1}e^{2\pi ixk/N}|k\rangle.$$

$F_N$ is unitary, so it is a legitimate gate on $n$ qubits.

Why $F_N$ is unitary

To check unitarity, compare columns $x$ and $y$. Their overlap is $N^{-1}\sum_ke^{2\pi i(y-x)k/N}$. For $x=y$, all $N$ terms are one. Otherwise set $z=e^{2\pi i(y-x)/N}$: then $z\ne1$, $z^N=1$, and $\sum_{k=0}^{N-1}z^k=(1-z^N)/(1-z)=0$. The columns are orthonormal, so the matrix is unitary.

Its inverse $F_N^\dagger$ has the opposite sign in the exponent.

For an input $\sum_xc_x|x\rangle$, the output amplitude at $k$ is $N^{-1/2}\sum_xe^{2\pi ixk/N}c_x$. Measuring produces one integer drawn from the squared magnitudes of these amplitudes. To obtain a period from that integer, we must first prepare amplitudes whose Fourier probabilities depend on the period. The following two-register construction explains how.

From a function table to a periodic state

Use input register $A$ and output register $B$. For this first calculation, assume the period $r$ divides $N$ and $f(x)=f(x')$ exactly when $x,x'$ have the same remainder modulo $r$. These assumptions make the peaks exact; order-finding later will not require $r$ to divide $N$.

Prepare $A$ uniformly and $B$ as zero, then evaluate $f$ reversibly:

$$\frac1{\sqrt N}\sum_{x=0}^{N-1}|x\rangle_A|0\rangle_B \xrightarrow{U_f}\frac1{\sqrt N}\sum_{x=0}^{N-1}|x\rangle_A|f(x)\rangle_B.$$

Now measure $B$. If its outcome is $f(x_0)$, the projector onto that value removes every term with a different function value. Before renormalisation the remaining joint vector is

$$\frac1{\sqrt N}\sum_{j=0}^{L-1}|x_0+jr\rangle_A|f(x_0)\rangle_B,\qquad L=N/r.$$

It contains $L$ orthogonal terms of squared amplitude $1/N$, so this outcome has probability $L/N=1/r$. Dividing by its square root gives the conditional input state

$$|\psi_{x_0}\rangle_A=\frac1{\sqrt L}\sum_{j=0}^{L-1}|x_0+jr\rangle_A.$$

The set of integers with a fixed remainder $x_0$ is a coset of the subgroup of multiples of $r$ in the integers modulo $N$. “Coset state” here simply means this uniform superposition over that set. Measuring the output reveals the remainder class but leaves the individual inputs within that class indistinguishable, so their amplitudes remain coherent.

For $N=16,r=4$, outcome $f(1)$ leaves $(|1\rangle+|5\rangle+|9\rangle+|13\rangle)/2$. Its probability is $1/4$. The experiment does not select the particular input one; it selects all four inputs sharing the observed output.

Calculate the Fourier peaks

Apply $F_N$ to the conditional state. Its amplitude at $k$ is

$$a_k=\frac{e^{2\pi ix_0k/N}}{\sqrt{NL}}\sum_{j=0}^{L-1}e^{2\pi ijrk/N}.$$

Set $z=e^{2\pi irk/N}$. Since $L=N/r$, $z^L=e^{2\pi ik}=1$. If $z\ne1$, the geometric sum is zero. If $z=1$, every term is one and the sum is $L$. The latter condition requires $rk/N$ to be an integer, so the only possible outcomes are $k=mN/r$ for $m=0,\ldots,r-1$.

Each surviving amplitude has magnitude $L/\sqrt{NL}=1/\sqrt r$, giving probability $1/r$. The offset $x_0$ contributes a phase of magnitude one and cannot change those probabilities. Consequently the measured peak distribution is also unchanged if we discard $B$ without reading it: we average over different offsets with identical Fourier probabilities.

In the sixteen-dimensional example the outcomes are $0,4,8,12$, each with probability $1/4$. Dividing by $N$ gives $0,1/4,1/2,3/4$. The nonzero fractions have reduced denominators four, two, and four. Even with exact peaks, a single draw need not reveal the full period: zero supplies no denominator and $1/2$ reveals only a divisor.

Phase estimation when the phase is not exactly representable

Suppose a unitary has an eigenstate $|u\rangle$ with $U|u\rangle=e^{2\pi i\phi}|u\rangle$, where $0\le\phi<1$. A control register in $N^{-1/2}\sum_j|j\rangle$ can apply $U^j$ coherently. To see the circuit construction, write $j=\sum_\ell j_\ell2^\ell$: controlled-$U^{2^\ell}$ on each bit multiplies to $U^j$ on that branch. The resulting state is

$$\frac1{\sqrt N}\sum_{j=0}^{N-1}|j\rangle U^j|u\rangle =\frac1{\sqrt N}\sum_{j=0}^{N-1}e^{2\pi ij\phi}|j\rangle|u\rangle.$$

Applying $F_N^\dagger$ to the control gives amplitude $a_k=N^{-1}\sum_je^{2\pi ij(\phi-k/N)}$. If $\phi=m/N$, this sum is one for $k=m$ and zero otherwise, by the same cancellation used above. For example, $N=8$ and $\phi=3/8$ give binary outcome $011$ with certainty.

Otherwise, write $\delta=\phi-k/N$ and sum the series: $a_k=(1-e^{2\pi iN\delta})/[N(1-e^{2\pi i\delta})]$. The identity $|1-e^{it}|=2|\sin(t/2)|$ then yields

$$p(k)=\frac{\sin^2(\pi N\delta)}{N^2\sin^2(\pi\delta)}.$$

Values at removable zeroes are obtained by a limit. The probability concentrates near $k=N\phi$, but a measurement can land away from the nearest integer. We therefore obtain an estimate, not an exact general phase. Increasing $N$ improves the resolution; repeated attempts address unfavourable samples.

Which unitary contains the order?

Choose a known modulus $M\ge2$ and an integer $a$ coprime to it, meaning their only common positive divisor is one. The order $r$ is the smallest positive integer satisfying $a^r\equiv1\pmod M$. Define modular multiplication on basis labels $0\le y<M$ by

$$U_a|y\rangle=|ay\bmod M\rangle.$$

Here $ay\bmod M$ is the remainder when $ay$ is divided by $M$. Coprimality means $a$ has an inverse modulo $M$, so this map permutes the basis and is unitary. In a qubit register whose dimension exceeds $M$, define it as identity on unused basis labels $y\ge M$. Multiplication of valid labels stays within the valid subspace, so the extended map is still a permutation.

Starting at $|1\rangle$, repeated applications visit $|1\rangle,|a\bmod M\rangle,\ldots,|a^{r-1}\bmod M\rangle$ and then return. On this cycle, define

$$|u_s\rangle=\frac1{\sqrt r}\sum_{j=0}^{r-1}e^{-2\pi isj/r}|a^j\bmod M\rangle,\qquad s=0,\ldots,r-1.$$

Checking the eigenvalue, and why $|1\rangle$ is an equal superposition of eigenstates

Applying $U_a$ shifts the basis exponent $j$ to $j+1$. Relabel the new exponent as $\ell$; its coefficient becomes $e^{-2\pi is(\ell-1)/r}=e^{2\pi is/r}e^{-2\pi is\ell/r}$. Thus $U_a|u_s\rangle=e^{2\pi is/r}|u_s\rangle$.

For the second claim, sum the Fourier coefficients of $r^{-1/2}\sum_s|u_s\rangle$ over $s$: only $j=0$ survives, leaving $|1\rangle$.

These are eigenstates, $U_a|u_s\rangle=e^{2\pi is/r}|u_s\rangle$, so the eigenphases are the fractions $s/r$. We do not need to prepare an individual eigenstate: the easily prepared state $|1\rangle$ is an equal superposition of them, $|1\rangle=r^{-1/2}\sum_s|u_s\rangle$. Phase estimation therefore samples the phase distributions for the different $s$ values with equal weights. The orthogonal target eigenstates prevent cross terms between those distributions when we measure the control alone.

For $a=2,M=21$, the cycle is $1,2,4,8,16,11,1$, giving $r=6$. Phases such as $1/6$ and $1/3$ cannot be represented exactly with a power-of-two denominator. This is where rational approximation becomes necessary.

Recover a phase fraction using continued fractions

The measured rational number is $k/N$. Reducing it in lowest terms still gives a power-of-two denominator, so that operation alone cannot recover an odd factor of $r$. Instead, expand $k/N$ as a continued fraction, written $[c_0;c_1,c_2,\ldots]=c_0+1/(c_1+1/(c_2+\cdots))$. Truncating after successive entries gives rational approximations called convergents.

For example, $256=3\cdot85+1$ and $85=85\cdot1$, so $85/256=[0;3,85]$. Its convergents include $0$, $1/3$, and $85/256$. The odd denominator three appears in an intermediate approximation; it is not the reduced denominator of the measured fraction.

The relevant theorem is precise: if $p/q$ is in lowest terms and $|k/N-p/q|<1/(2q^2)$, then $p/q$ occurs among the convergents of the continued-fraction expansion of $k/N$. The proof is a number-theoretic result; see the linked continued-fraction reference. Here the desired $p/q$ is $s/r$ after reducing that unknown phase fraction. An error smaller than $1/(2r^2)$ is sufficient because its reduced denominator $q$ is at most $r$.

The approximation theorem and the register size are separate statements. On a nearest-integer phase-estimation outcome, $|k/N-s/r|\le1/(2N)$. Choosing $N\ge2r^2$ is a sufficient margin: that error is at most $1/(4r^2)$, strictly below the required bound. It is not a necessary condition of the continued-fraction theorem itself, nor does it guarantee every measured outcome is close enough. Since $r$ is unknown but $r<M$, choosing a power of two $N\ge2M^2$ provides this margin in advance.

Verify the denominator instead of trusting a sample

Continue with $a=2,M=21$ and choose $N=1024$, which exceeds $2M^2$. A possible near-$1/3$ sample is $k=341$. Division gives $341/1024=[0;3,341]$, so one candidate phase fraction is $1/3$. Its approximation error is

$$\left|\frac{341}{1024}-\frac13\right|=\frac1{3072}<\frac1{2\cdot6^2}.$$

The approximation is excellent, yet its denominator three is not the order: $2^3\bmod21=8\ne1$. This sample corresponds to $s=2,r=6$, whose fraction reduces to $1/3$. The information lost in that reduction cannot be restored from this sample alone.

A second possible sample near $1/6$ is $171/1024=[0;5,1,84,2]$. Which convergent to test is fixed by the register-size margin, provided the sample is accurate enough. If $N\ge2M^2$ and $k$ is the integer nearest to $Ns/r$ for some $s$, then the reduced form of $s/r$ is a convergent of $k/N$, and since $r<M$ it is the convergent with the largest denominator below $M$. The register size does not make every sample accurate: an outcome away from every peak, which occurs with small but nonzero probability, carries no such guarantee. Its candidate fails the power test, and the run is repeated. The convergents of $171/1024$ are $0,\ 1/5,\ 1/6,\ 85/509,\ 171/1024$; the last two have denominators above $M=21$, so the rule selects $[0;5,1]=1/6$ and candidate six, and $2^6\bmod21=1$. Checking its proper divisors one, two, and three excludes a smaller order in this example. In general a passing power test establishes a multiple of the order; further reduction or suitable candidate selection establishes the least one. Repetitions and classical verification are part of the algorithm.

For factoring, an even order can yield factors through greatest common divisors of $a^{r/2}\pm1$ with $M$. Here $\gcd(7,21)=7$ and $\gcd(9,21)=3$. Some choices yield odd orders or unhelpful square roots and must be retried. The reference treats those cases and the success analysis.

Account for the operations used

The exact QFT circuit uses $n$ Hadamards, $n(n-1)/2$ controlled phase rotations, and qubit-order reversal: $O(n^2)$ elementary operations, meaning an upper bound proportional to $n^2$ for large $n$. This count assumes those specified rotations are available; synthesis into a fixed alphabet adds approximation overhead.

Phase estimation also needs controlled powers. If $U$ is available only as a black box, constructing $U^{2^\ell}$ by repetition requires $2^\ell$ calls, for $N-1$ calls across the whole register. Modular multiplication has additional structure: $U_a^{2^\ell}$ is multiplication by the classically precomputed constant $a^{2^\ell}\bmod M$. Reversible arithmetic can implement each such multiplication efficiently in the number of bits of $M$. This arithmetic construction, rather than gate universality alone, supports the runtime advantage.

Where this leaves us

One Fourier sample reveals an integer $k$ with $k/N$ close to $s/r$ for some unknown $s$: a noisy fraction, not the period. When $r$ divides $N$ the peaks are exact; otherwise they spread, and a sample can miss. What remains is classical: find the candidate fraction by continued fractions, test the candidate order with one modular power, and repeat when a sample reduces to a divisor of $r$ or falls away from a peak. For factoring there is one more classical step and sometimes a retry. The quantum part runs efficiently because the Fourier transform takes $O(n^2)$ gates and efficient reversible modular arithmetic implements the controlled powers; the existence of a Fourier transform alone would not be enough.

All of this depends on structure: a periodic function whose period the Fourier transform can expose. The next chapter asks what is possible with no structure at all, only a routine that checks whether a candidate is the answer.

Problems

1. Fourier-transform $(|0\rangle+|4\rangle)/\sqrt2$ for $N=8$. Show the amplitudes and verify their probabilities sum to one.

Worked answer

Each of the two input terms contributes a Fourier amplitude $1/(\sqrt2\sqrt8)=1/4$. Keeping both contributions gives

$$a_k=\frac14\left(e^{2\pi i0k/8}+e^{2\pi i4k/8}\right)=\frac{1+(-1)^k}{4}.$$

For even $k$, the two contributions add to $1/2$; for odd $k$, they cancel. The probability list in integer order is $(1/4,0,1/4,0,1/4,0,1/4,0)$. There are four surviving outcomes, so the total is $4(1/4)=1$.

2. For modular multiplication by $2$ modulo $21$, phase estimation with $N=1024$ returns $k=341$. Compute the first nonzero convergent and test whether its denominator is the order. Explain what a sample $k=171$ adds.

Worked answer

First use integer division: $1024=3\cdot341+1$, then $341=341\cdot1$. Thus $341/1024=[0;3,341]$, whose first nonzero convergent is $1/3$. Reducing $341/1024$ directly would leave denominator $1024$, not three.

Test the candidate with $2^3=8\not\equiv1\pmod{21}$, so three fails. The actual phase can be $2/6=1/3$, which explains how a close phase estimate can lose a factor of the order.

For the other sample, the divisions begin $1024=5\cdot171+169$ and $171=1\cdot169+2$. The convergent $[0;5,1]=1/(5+1)=1/6$ gives candidate six. Now $2^6=64\equiv1\pmod{21}$, while powers at the proper divisors $1,2,3$ give $2,4,8$. Hence six is the least positive order. The second sample supplied a numerator coprime to that order.

The tempting wrong turn is to reduce $341/1024$ and read off a denominator. The fraction is already in lowest terms, with denominator $1024$; only the continued-fraction expansion exposes the nearby fraction with a small denominator.

References and further reading

Nielsen & Chuang, §5.1–5.4.1 (QFT, phase estimation, order-finding, factoring, and periodic functions). For the approximation theorem, see Ben Lynn’s notes on continued-fraction convergence. For the full algorithm and runtime argument, see Shor, Polynomial-Time Algorithms for Prime Factorization and Discrete Logarithms on a Quantum Computer.