Editing Shor's algorithm (section)

=== Quantum order-finding subroutine ===
The goal of the quantum subroutine of Shor's algorithm is, given [[coprime integers]] <math> N </math> and <math> 1< a<N </math>, to find the [[Multiplicative order|order <math> r </math> of <math>a</math> modulo <math>N</math>]], which is the smallest positive integer such that <math>a^r \equiv 1 \pmod N</math>. To achieve this, Shor's algorithm uses a quantum circuit involving two registers. The second register uses <math> n </math> qubits, where <math> n </math> is the smallest integer such that <math> N\le 2^n </math>, i.e., <math> n = \left \lceil {\log_2N} \right \rceil </math>. The size of the first register determines how accurate of an approximation the circuit produces. It can be shown that using <math> 2n </math> qubits gives sufficient accuracy to find <math> r </math>. The exact quantum circuit depends on the parameters <math>a</math> and <math>N</math>, which define the problem. The following description of the algorithm uses [[bra–ket notation]] to denote quantum states, and <math>\otimes</math> to denote the [[tensor product]], rather than [[logical AND]].

The algorithm consists of two main steps:

# Use [[Quantum phase estimation algorithm|quantum phase estimation]] with unitary <math>U</math> representing the operation of multiplying by <math>a</math> (modulo <math>N</math>), and input state <math>|0\rangle^{\otimes 2 n}\otimes|1\rangle</math> (where the second register is <math>|1\rangle</math> made from <math>n</math> qubits). The eigenvalues of this <math>U</math> encode information about the period, and <math>|1\rangle</math> can be seen to be writable as a sum of its eigenvectors. Thanks to these properties, the quantum phase estimation stage gives as output a random integer of the form <math>\frac{j}{r} 2^{2n}</math> for random <math>j=0,1,...,r-1</math>.
# Use the [[simple continued fraction|continued fractions algorithm]] to extract the period <math>r</math> from the measurement outcomes obtained in the previous stage. This is a procedure to post-process (with a classical computer) the measurement data obtained from measuring the output quantum states, and retrieve the period.

The connection with quantum phase estimation was not discussed in the original formulation of Shor's algorithm,<ref name="siam" /> but was later proposed by Kitaev.<ref>{{cite arXiv |eprint=quant-ph/9511026 |last1=Kitaev |first1=A. Yu |date=1995 |title=Quantum measurements and the Abelian Stabilizer Problem }}</ref>

==== Quantum phase estimation ====
[[File:Shor's algorithm.svg|frame|Quantum subroutine in Shor's algorithm]]
In general the [[quantum phase estimation algorithm]], for any unitary <math>U</math> and eigenstate <math>|\psi\rangle</math> such that <math>U|\psi\rangle=e^{2\pi i\theta} |\psi\rangle</math>, sends input states <math>|0\rangle|\psi\rangle</math> to output states close to <math>|\phi\rangle|\psi\rangle</math>, where <math>\phi</math> is a superposition of integers close to <math>2^{2n} \theta</math>. In other words, it sends each eigenstate <math>|\psi_j\rangle</math> of <math>U</math> to a state containing information close to the associated eigenvalue. For the purposes of quantum order-finding, we employ this strategy using the unitary defined by the action
<math display="block">
 U|k\rangle = \begin{cases}
  |ak \pmod N\rangle & 0 \le k < N, \\
  |k\rangle & N \le k < 2^n.
\end{cases}</math>
The action of <math>U</math> on states <math>|k\rangle</math> with <math> N \leq k < 2^n </math> is not crucial to the functioning of the algorithm, but needs to be included to ensure that the overall transformation is a well-defined quantum gate. Implementing the circuit for quantum phase estimation with <math>U</math> requires being able to efficiently implement the gates <math> U^{2^j} </math>. This can be accomplished via [[modular exponentiation]], which is the slowest part of the algorithm.

The gate thus defined satisfies <math>U^r = I</math>, which immediately implies that its eigenvalues are the <math>r</math>-th [[Root of unity|roots of unity]] <math>\omega_r^k = e^{2\pi ik/r}</math>. Furthermore, each eigenvalue <math>\omega_r^j</math> has an eigenvector of the form <math display="inline">|\psi_j\rangle=r^{-1/2}\sum_{k=0}^{r-1}\omega_r^{-kj}|a^k\rangle </math>, and these eigenvectors are such that
<math display="block">\begin{align}
 \frac{1}{\sqrt{r}} \sum_{j = 0}^{r - 1} |\psi_j\rangle &= \frac{1}{r} \sum_{j = 0}^{r - 1} \sum_{k = 0}^{r - 1} \omega_r^{jk}|a^k\rangle \\
 &= |1\rangle + \frac{1}{r} \sum_{k = 1}^{r - 1} \left(\sum_{j = 0}^{r - 1} \omega_r^{jk} \right) |a^k\rangle =|1\rangle,
\end{align}</math>
where the last identity follows from the [[geometric series]] formula, which implies <math display="inline">\sum_{j = 0}^{r - 1} \omega_r^{jk} = 0</math>.

Using [[Quantum phase estimation algorithm|quantum phase estimation]] on an input state <math>|0\rangle^{\otimes 2 n}|\psi_j\rangle</math> would then return the integer <math>2^{2n} j/r</math> with high probability. More precisely, the quantum phase estimation circuit sends <math>|0\rangle^{\otimes 2 n}|\psi_j\rangle</math> to <math>|\phi_j\rangle|\psi_j\rangle</math> such that the resulting probability distribution <math>p_k \equiv|\langle k|\phi_j\rangle|^2</math> is peaked around <math>k=2^{2n} j/r</math>, with <math>p_{2^{2n}j/r} \ge 4/\pi^2 \approx 0.4053</math>. This probability can be made arbitrarily close to 1 using extra qubits.

Applying the above reasoning to the input <math>|0\rangle^{\otimes 2 n}|1\rangle</math>, quantum phase estimation thus results in the evolution
<math display="block">
 |0\rangle^{\otimes 2 n}|1\rangle =
  \frac{1}{\sqrt{r}} \sum_{j = 0}^{r - 1} |0\rangle^{\otimes 2 n} |\psi_j\rangle \to
  \frac{1}{\sqrt{r}} \sum_{j = 0}^{r - 1} |\phi_j\rangle|\psi_j\rangle.
</math>
Measuring the first register, we now have a balanced probability <math>1/r</math> to find each <math>|\phi_j\rangle</math>, each one giving an integer approximation to <math>2^{2 n} j/r</math>, which can be divided by <math>2^{2n}</math> to get a decimal approximation for <math>j/r</math>.

==== Continued-fraction algorithm to retrieve the period ====
Then, we apply the [[continued fraction|continued-fraction]] algorithm to find integers <math>b</math> and <math>c</math>, where <math>b/c</math> gives the best fraction approximation for the approximation measured from the circuit, for <math>b, c < N</math> and [[Coprime integers|coprime]] <math>b</math> and <math>c</math>. The number of qubits in the first register, <math>2n</math>, which determines the accuracy of the approximation, guarantees that
<math display="block">
 \frac{b}{c} = \frac{j}{r},
</math>
given the best approximation from the superposition of <math>|\phi_j\rangle</math> was measured<ref name="siam"/> (which can be made arbitrarily likely by using extra bits and truncating the output). However, while <math>b</math> and <math>c</math> are coprime, it may be the case that <math>j</math> and <math>r</math> are not coprime. Because of that, <math>b</math> and <math>c</math> may have lost some factors that were in <math>j</math> and <math>r</math>. This can be remedied by rerunning the quantum order-finding subroutine an arbitrary number of times, to produce a list of fraction approximations
<math display="block">
 \frac{b_1}{c_1}, \frac{b_2}{c_2}, \ldots, \frac{b_s}{c_s},
</math>
where <math>s</math> is the number of times the subroutine was run. Each <math>c_k</math> will have different factors taken out of it because the circuit will (likely) have measured multiple different possible values of <math>j</math>. To recover the actual <math>r</math> value, we can take the [[least common multiple]] of each <math>c_k</math>:
<math display="block">
 \operatorname{lcm}(c_1, c_2, \ldots, c_s).
</math>
The least common multiple will be the order <math>r</math> of the original integer <math>a</math> with high probability. In practice, a single run of the quantum order-finding subroutine is in general enough if more advanced post-processing is used.<ref name="Ekerå24">{{cite journal |last1=Ekerå |first1=Martin |title=On the Success Probability of Quantum Order Finding |journal=ACM Transactions on Quantum Computing |date=May 2024 |volume=5 |issue=2 |pages=1–40 |doi=10.1145/3655026 |doi-access=free |arxiv=2201.07791 }}</ref>

==== Choosing the size of the first register ====
Phase estimation requires choosing the size of the first register to determine the accuracy of the algorithm, and for the quantum subroutine of Shor's algorithm, <math>2n</math> qubits is sufficient to guarantee that the optimal bitstring measured from phase estimation (meaning the <math>|k\rangle</math> where <math display="inline">k / 2^{2n}</math> is the most accurate approximation of the phase from phase estimation) will allow the actual value of <math>r</math> to be recovered.

Each <math>|\phi_j\rangle</math> before measurement in Shor's algorithm represents a superposition of integers approximating <math>2^{2 n} j/r</math>. Let <math>|k\rangle</math> represent the most optimal integer in <math>|\phi_j\rangle</math>. The following theorem guarantees that the continued fractions algorithm will recover <math>j/r</math> from <math>k/2^{2 {n}}</math>:
{{Math theorem
| math_statement = If <math>j</math> and <math>r</math> are <math>n</math> bit integers, and
<math display="block">\left\vert \frac{j}{r} - \phi\right\vert \leq \frac{1}{2 r^2}</math>
then the continued fractions algorithm run on <math>\phi</math> will recover both <math display="inline">\frac{j}{\gcd(j,\; r)}</math> and <math display="inline">\frac{r}{\gcd(j,\; r)}</math>.
}}
<ref name=":0" /> As <math>k</math> is the optimal bitstring from phase estimation, <math>k/2^{2 {n}}</math> is accurate to <math>j/r</math> by <math>2n</math> bits. Thus,<math display="block">\left\vert\frac{j}{r} - \frac{k}{2^{2n}}\right\vert \leq \frac{1}{2^{2 {n} + 1}} \leq \frac{1}{2N^2} \leq \frac{1}{2r^2}</math>which implies that the continued fractions algorithm will recover <math>j</math> and <math>r</math> (or with their greatest common divisor taken out).