Order Finding — The Quantum Idea

You have seen quantum phase estimation in ch. 70-72. Feed it a unitary U and an eigenstate |u\rangle with U|u\rangle = e^{2\pi i \varphi}|u\rangle, and it spits out a binary approximation of \varphi. It is a universal machine for reading eigenvalue phases.

Order finding uses that machine on a very specific, very structured unitary: multiplication by a modulo N. The phase that comes out of phase estimation is s/r — a fraction whose denominator is exactly the number you are trying to find. One classical post-processing step (continued fractions) recovers r from s/r. That is the whole quantum idea.

This chapter builds the idea from the inside. You will see why modular multiplication is unitary, what its eigenstates actually look like, why the state |1\rangle is the magic target register that lets the algorithm work even though you do not know r (and therefore could not prepare any single eigenstate), how phase estimation samples a random s/r, and how continued fractions turns that sample into r. The next chapter (ch. 77) wires this subroutine into the full Shor pipeline — but first you need to understand what is happening inside the quantum black box.

What "order" means, and why you want it

Take N = 15 and a = 7. Compute the powers:

After four steps the sequence returns to 1. The order of 7 modulo 15 is r = 4. From here on the powers cycle: 7^5 = 7, 7^6 = 4, and so on.

The order r is the length of the orbit \{a^0, a^1, \ldots, a^{r-1}\} = \{1, a, a^2, \ldots, a^{r-1}\} of a under multiplication mod N. By Lagrange's theorem, r divides the size of (\mathbb{Z}/N)^\times, so r is bounded above, but finding it for an n-bit N with classical methods takes exponential time — no algorithm is known that does better than factoring N itself.

From the previous chapter (ch. 75), you know why this matters: if you can find the order r of a random a modulo N, you can factor N in polynomial time. The quantum subroutine for order finding is therefore the quantum subroutine for factoring.

The unitary — multiplication by a

Define an operator U_a on the Hilbert space spanned by the basis states |0\rangle, |1\rangle, \ldots, |N-1\rangle by

This is a linear operator (extend by linearity to superpositions). To use phase estimation you need to know two things: that U_a is unitary, and what its eigenvalues and eigenstates look like.

U_a is unitary. Consider its action on the subset \{1, 2, \ldots, N-1\} that is coprime to N. Multiplication by a is a bijection of this set onto itself — its inverse is multiplication by a^{-1} \bmod N, which exists because \gcd(a, N) = 1. A linear operator that maps a basis to a permutation of itself is unitary: it takes orthonormal basis vectors to orthonormal basis vectors, and it preserves inner products. So U_a is unitary on the coprime-to-N subspace.

Why only on the coprime subspace: if y shares a factor with N, the orbit y, ay, a^2 y, \ldots also shares that factor, so U_a permutes the coprime elements among themselves. The non-coprime elements form other orbits. For Shor's algorithm, you only ever need U_a to act on the coprime-to-N subspace; the standard construction of the quantum circuit for U_a leaves the non-coprime basis states alone.

The orbit of 1. Starting from |1\rangle, apply U_a repeatedly: |1\rangle \to |a\rangle \to |a^2\rangle \to \cdots \to |a^{r-1}\rangle \to |1\rangle. The orbit has exactly r elements — the order of a. These r basis states |1\rangle, |a\rangle, |a^2\rangle, \ldots, |a^{r-1}\rangle form an r-dimensional invariant subspace of U_a: they are closed under the action of U_a (and under U_a^{-1}, which is multiplication by a^{-1}).

On this r-dimensional subspace, U_a acts like a cyclic-shift operator: each basis vector shifts to the next, and the last wraps around to the first. That is the structure you are going to diagonalise.

The eigenstates — Fourier combinations over the orbit

A cyclic-shift operator on r states is diagonalised by the discrete Fourier basis on \mathbb{Z}/r. Explicitly: for each s \in \{0, 1, \ldots, r-1\}, define

These r states are eigenstates of U_a, and the eigenvalue of |u_s\rangle is e^{2\pi i\, s/r}.

Deriving the eigenvalue equation. Apply U_a to |u_s\rangle and see what comes out:

Why pull U_a inside the sum: U_a is linear, so it distributes over the superposition. The coefficients e^{-2\pi i s k/r} are just numbers, untouched by the operator.

Now U_a |a^k \bmod N\rangle = |a^{k+1} \bmod N\rangle — multiplying by a shifts the exponent. Substitute:

Re-index the sum. Let j = k + 1; as k runs from 0 to r-1, j runs from 1 to r. The term j = r has |a^r \bmod N\rangle = |1\rangle = |a^0\rangle (because r is the order), so it corresponds to j = 0 wrapped around. Rewriting:

Why re-indexing works cleanly: the orbit is a cycle of length r, so shifting the index by 1 and wrapping at r is the same as shifting the sum variable. The exponential at j = 0 is e^{-2\pi i s(-1)/r} = e^{2\pi i s/r}, exactly what j = r would have been had the index not wrapped.

Pull out the common factor e^{2\pi i s / r}:

That is the eigenvalue equation. |u_s\rangle is an eigenstate of U_a with eigenvalue e^{2\pi i\, s/r}. One derivation, r eigenstates — one for each s = 0, 1, \ldots, r-1.

The catch, and the trick — |1\rangle is a uniform superposition of eigenstates

To run phase estimation on U_a, you need an eigenstate |u_s\rangle loaded into the target register. But look at |u_s\rangle:

Preparing this requires knowing r and s. But r is the number you are trying to find! You cannot prepare any single eigenstate without already having the answer.

Here is where quantum mechanics saves you. Consider the computational basis state |1\rangle. Sum the eigenstates:

Why swap the sums: both sums are finite and over disjoint variables (s and k), so Fubini applies. Pulling the sum over s inside lets you collect the geometric series on s for each fixed k.

The inner sum over s is a geometric series. For k = 0, every term is e^{0} = 1, so the sum is r. For k \neq 0 (and k < r), the sum is \sum_{s=0}^{r-1} \omega^s where \omega = e^{-2\pi i k/r} is a non-trivial r-th root of unity; by the geometric-series identity, this equals (1 - \omega^r)/(1 - \omega) = 0. So the sum collapses:

Equivalently,

The computational-basis state |1\rangle is a uniform superposition of all r eigenstates of U_a. You do not need to know r to prepare it — you just load the number 1 into the register. The quantum magic is that when you run phase estimation on U_a with target |1\rangle, linearity carries every eigenstate component through the circuit in parallel, and the measurement samples one eigenvalue at random.

Why this is a breakthrough: the superposition trick completely evades the "I don't know r so I can't prepare an eigenstate" problem. The state |1\rangle is the simplest basis state to prepare, and it secretly contains all r eigenstates with equal amplitude. Phase estimation does not mind — it processes the superposition linearly and picks out one eigenvalue per run.

The algorithm — phase estimation plus continued fractions

Putting the pieces together, the full order-finding algorithm looks like this.

Input. A composite N and an a \in \{2, 3, \ldots, N-1\} with \gcd(a, N) = 1.

Output. The order r of a mod N.

Algorithm.

1. Choose ancilla width t. Use t = 2n + 1 ancilla qubits where n = \lceil \log_2 N \rceil. This ensures enough precision that continued fractions will recover r with high probability.
2. Prepare registers. t ancillas in |0\rangle^{\otimes t}, and an n-qubit target register in |1\rangle.
3. Run phase estimation on U_a. Apply Hadamards to the ancillas, then the ladder of controlled-U_a^{2^k} gates, then inverse QFT, then measure the ancillas. Because |1\rangle is a uniform sum of eigenstates |u_s\rangle, the measurement yields a t-bit integer y with

for some s \in \{0, 1, \ldots, r-1\} drawn uniformly at random. The approximation is accurate to within 1/2^{t+1} (by the QPE accuracy bound). 4. Apply continued fractions. Treat \widetilde\varphi = y/2^t as a rational number and compute its continued-fraction expansion. The convergents are successive rational approximations in lowest terms. One of them will be s/r with r < N; stop when you find a convergent whose denominator r' satisfies |y/2^t - s'/r'| < 1/(2 N^2). 5. Verify. Compute a^{r'} \bmod N. If it is 1, output r'; otherwise, return to step 2 and try again.

Each run produces a candidate r'. The probability that continued fractions succeeds on a single run is \Omega(1/\log\log r) \geq \Omega(1/\log\log N); repeat O(\log\log N) times to find r with probability arbitrarily close to 1.

Most of the total work happens inside phase estimation — specifically, inside the controlled-U_a^{2^k} gates. The classical post-processing (continued fractions) is essentially free: it runs in O(n^3) time on the measurement outcome.

Worked examples

Common confusions

• "|1\rangle is itself an eigenstate of U_a." No. Applying U_a to |1\rangle gives |a\rangle, a different computational-basis state. |1\rangle is not fixed by U_a — it is the first element of an orbit of length r. The eigenstates are the Fourier-over-the-orbit superpositions |u_s\rangle. The magic of |1\rangle is that it equals \tfrac{1}{\sqrt r}\sum_s |u_s\rangle, a uniform superposition of eigenstates.

• "Phase estimation returns r directly." No. Phase estimation returns \widetilde\varphi \approx s/r for a random s — a fraction, not an integer. Recovering r from s/r is the job of continued fractions, a classical post-processing step. If you measure y = 0 (which corresponds to s = 0), you learn nothing; try again.

• "The whole algorithm is quantum." No. The quantum part is phase estimation on U_a — which is itself dominated by the classical-looking subroutine of modular exponentiation implemented reversibly on a quantum register. Continued fractions is purely classical. The reduction from factoring to order finding is classical. The verification a^{r'} \bmod N \stackrel{?}{=} 1 is classical. The quantum computer is doing one specific thing: extracting s/r from the eigenvalue structure of U_a.

• "You need to know r to pick the right ancilla count t." You only need an upper bound, and r < N gives you one. Choose t = 2n + 1 where n = \lceil \log_2 N\rceil, and the precision is sufficient for any r up to N. This is why the ancilla register is "roughly twice the size of the data register" in every Shor-circuit diagram.

• "Continued fractions always works." It succeeds when \gcd(s, r) = 1, which happens for \varphi(r)/r fraction of s values, where \varphi is Euler's totient function. For typical r this is \Omega(1/\log\log r), and the expected number of repetitions to find a coprime s is logarithmic. If you get an s with \gcd(s, r) = d > 1, continued fractions returns r/d, a proper divisor of r — sometimes useful, often enough to compute the true r by taking the LCM of two or three runs.

Going deeper

Formal derivation of the |1\rangle eigenstate expansion

Starting from the eigenstate definitions, the sum-over-s identity is a two-line Fourier-over-\mathbb{Z}/r argument. Writing \omega = e^{-2\pi i/r} (a primitive r-th root of unity):

where the inner sum equals r when k \equiv 0 \pmod r and 0 otherwise — the orthogonality of characters of \mathbb{Z}/r. More generally, any computational-basis state |a^j\rangle on the orbit decomposes as

so running QPE on any orbit state returns s/r for random s — the choice of |1\rangle is just the easiest to prepare. (Proof: apply the inverse Fourier on the orbit to both sides of the |u_s\rangle definition.)

The continued-fractions algorithm, formally

Given a rational x \in (0, 1), the continued-fraction expansion is the unique sequence [a_0; a_1, a_2, \ldots] with a_0 \ge 0, a_i \ge 1 for i \ge 1, defined by x_0 = x, a_i = \lfloor x_i\rfloor, x_{i+1} = 1/(x_i - a_i) (stopping if x_i = a_i).

The i-th convergent p_i/q_i is the rational obtained by truncating the expansion at a_i. Convergents satisfy the Legendre theorem: if |x - p/q| < 1/(2q^2) for some fraction p/q in lowest terms, then p/q is a convergent of x.

Application to order finding. By the QPE accuracy bound, the measured \widetilde\varphi = y/2^t satisfies |\widetilde\varphi - s/r| \le 1/2^{t+1}. Choose t = 2n + 1; then 1/2^{t+1} \le 1/(2^{2n+2}) \le 1/(2 N^2) \le 1/(2 r^2) (since r \le N). By Legendre, s/r (in lowest terms) is a convergent of \widetilde\varphi. Compute all convergents — there are O(n) of them by a theorem of Khinchin — and find the one whose denominator r' satisfies a^{r'} \equiv 1 \pmod N. That is r (if \gcd(s, r) = 1) or a divisor of r (otherwise).

Failure modes and success probability

A single run fails to return r in two cases:

1. s = 0. The measured phase is 0/r = 0, carrying no information. Happens with probability 1/r.
2. \gcd(s, r) > 1. The continued-fractions convergent s/r reduces to s'/r' with r' < r a proper divisor. Happens with probability 1 - \varphi(r)/r, where \varphi is Euler's totient.

Combining, the probability of success per run is \varphi(r)/r, which for typical r is \Omega(1/\log\log r) (by Mertens' theorem; the worst case is when r is a primorial). Repeat O(\log\log N) times and take the LCM of the results; with probability > 1/2, the LCM is r. Total expected runs: O(\log\log N).

A careful analysis (Nielsen and Chuang §5.3.1) shows that in practice two independent runs suffice to pin down r with probability > 3/4: on average the \gcd-returned divisors of r span enough information that two LCMs hit r.

Modular exponentiation cost — the expensive part

The bottleneck of the quantum order-finding circuit is implementing the ladder of controlled-U_a^{2^k} gates for k = 0, 1, \ldots, t-1 = 2n. Each gate U_a^{2^k} is multiplication by a^{2^k} \bmod N, which is a specific permutation of the coprime-to-N residues, built from binary adders and multipliers.

A straightforward implementation of U_a^{2^k} uses O(n^2) Toffoli gates. With t = 2n + 1 ancillas, the total Toffoli count is O(n^3) — cubic in the bit-size of N. For 2048-bit RSA, that is \sim 2048^3 \approx 10^{10} Toffoli gates per run. Each Toffoli requires \sim 7 T gates in the standard decomposition, pushing the T-count to \sim 10^{11}. This is why Shor's algorithm is expensive even for "small" cryptographically relevant N — the constant factors hidden in O(n^3) are enormous.

Modern optimisations (Coppersmith, windowed arithmetic, Montgomery form, Gidney's recent work on modular multipliers) reduce the logical-Toffoli count by factors of 10 to 100, but the asymptotic remains O(n^3).

Kitaev's alternative — one-ancilla iterative order finding

Kitaev's 1995 construction (arXiv:quant-ph/9511026) replaces the t-ancilla QFT-based phase estimation with a single-ancilla iterative version. For order finding, this means extracting the bits of s/r one at a time, feeding each measurement outcome back as a classical control for the next round's controlled-U_a^{2^k}.

Advantages: one ancilla instead of 2n + 1, at the cost of O(n) circuit depth rounds with classical feedforward. Useful on NISQ devices where qubit count is scarce but coherence time per round is adequate.

Disadvantages: the classical post-processing to reconstruct s/r is more intricate, and round-to-round phase errors accumulate differently. For full-scale fault-tolerant Shor, the QFT-based version is standard.

Indian context — Aadhaar, UPI, and the post-quantum horizon

Order finding is the real quantum threat to RSA-based cryptography, which underpins much of India's digital infrastructure. Aadhaar authentication uses RSA signatures on biometric queries; UPI transactions sign payment messages with RSA and ECDSA; GSTIN filings, DigiLocker, eSign, and MCA21 all rely on RSA-based public-key infrastructure today. If a fault-tolerant quantum computer ran the order-finding subroutine on the modulus of any one of these signing keys, the cryptosystem would fall.

The timeline is years to decades away — current hardware is \sim 10^4 qubits short of breaking 2048-bit RSA, as shown by Gidney-Ekerå (2021). But the "harvest now, decrypt later" threat is immediate: an adversary can record encrypted traffic today and decrypt it in 15 years when the quantum computer exists. India's National Quantum Mission (2023, ₹6000 crore) has a working group on post-quantum cryptography migration precisely to get ahead of this; the CERT-In and MeitY advisories from 2024 onward recommend lattice-based signatures (Kyber, Dilithium) for new deployments. The next chapter (ch. 77) assembles order finding into the full Shor circuit and puts concrete numbers on the migration deadline.

Where this leads next

• The Shor circuit end to end — the full factoring pipeline: pick a, run order finding, recover factors, with complete resource estimates for 2048-bit RSA.
• Factoring reduces to order finding — the classical reduction that explains why order finding factors N.
• Quantum phase estimation — the subroutine that order finding invokes.
• Simon's algorithm — the (\mathbb{Z}/2)^n analogue: hidden period over the Boolean hypercube, the template Shor transposed to \mathbb{Z}/N.
• Continued fractions — the classical number-theory algorithm used for post-processing.

References

1. Peter Shor, Polynomial-time algorithms for prime factorization and discrete logarithms on a quantum computer (1994) — arXiv:quant-ph/9508027. The founding paper — order finding is §5.
2. Wikipedia, Shor's algorithm — includes a clean derivation of the eigenstates of modular multiplication.
3. Nielsen and Chuang, Quantum Computation and Quantum Information §5.3 — Cambridge University Press. Order finding with the continued-fractions analysis.
4. John Preskill, Lecture Notes on Quantum Computation Chapter 7 — theory.caltech.edu/~preskill/ph229. Order finding as phase estimation, with the |1\rangle = \sum_s |u_s\rangle/\sqrt r derivation.
5. Alexei Kitaev, Quantum measurements and the Abelian stabilizer problem (1995) — arXiv:quant-ph/9511026. The one-ancilla iterative formulation of phase estimation.
6. Artur Ekert and Richard Jozsa, Quantum computation and Shor's factoring algorithm (1996) — arXiv:quant-ph/9508027. A readable review of order finding and its place in the factoring algorithm.