Shor算法是体现量子计算相较于经典计算优越性的一个典型算法。大数分解问题是目前应用最广泛的非对称加密算法RSA算法的核心，其保密性来源于经典计算机无法在有限的时间内将一个大数分解为两个质因数。经典计算机求解这个问题的复杂度为 exp(O(n^{1/3}(logn)^{2/3})) 。而量子计算求解这个问题的复杂度为 O(n^2logn (loglogn))

1 大数分解问题与周期寻找问题

问题描述（大数分解问题）：设N为一个大数，它唯二的质因数为p和q，即N=pq。现在我们的需求是，在知道N的情况下，求出p和q。

为了解决大数分解问题，我们需要先解决周期寻找问题，因为周期寻找问题是大数分解问题的关键步骤。我们先介绍周期寻找问题，随后说明它们的关系。

问题描述（周期寻找问题）：设a和N为正整数，a＜N，且a与N互素。则a模N的周期（或者称为 阶）是满足下面等式的最小正整数r： a^r = 1\pmod N 。（根据欧拉定理，这样的r一定存在。欧拉定理：设a和N为正整数，且互素，则有 a^{\phi(N)} = 1 \pmod N , \phi(N) < N ）

下面说明这两个问题的关系：

执行如下步骤：

1. （大数分解问题）随机选取小于N的正整数 m 要求 m 和 N 互素（这很显然，否则m就必须等于p或q，问题已经得到解决）。
2. 计算m模N的周期r，即 m^r =1 \pmod N （周期寻找问题）
3. 如果r为奇数，则重新选择m，计算新的r，直至得到的r为偶数。
4. 取 x = m^{r/2} \pmod N , 则 x^2 = 1 \pmod N ，既 x^2 -1 =0\pmod N
   既 （x-1)(x+1) = 0 \pmod N ，既 (x-1)(x+1) = kN=kpq ，k为正整数
   因为p，q为素数，所以有 （x-1)=ap，（x+1)=bq ，a,b为整数且，ab=k
5. 到这里我们看出来了，（x-1)和N具有共同因数p=gcd(x-1,N)，gcd表示最大公因数，（x+1)和N有共同因数q=gcd(x+1,N)，我们只需要使用辗转相除法就可以求出p和q了。
   可能存在平庸解：gcd(x-1,N)=1 ，gcd(x+1,N)=N，即此时 x+1是N的整数倍。此时重新选择m并重复上述操作。
   

可见，我们只需要解决了周期寻找问题，就能只需要几步就可以解决大数分解问题。于是我们现在的任务转换成求解周期寻找问题。

2 Shor算法

Shor的解决办法是对酉正变换 U|y\rangle =|ay\pmod N\rangle 使用量子相位估计。为此，我们得先了解下U的本征态度情况。为方便理解，举一个具体的例子，假设a=3，N=35，最开始时 |y\rangle=|1\rangle ,则有如下关系

U|1\rangle=|3\pmod {35} \rangle \\ U^2|1\rangle=|9\pmod {35} \rangle \\ U^3|1\rangle = |27 \pmod {35} \rangle\\ \dots \\ U^{r-1}|1\rangle = |12 \pmod {35} \rangle\\ U^r|1\rangle = |1\pmod {35} \rangle

所以我们可以从上述等式中构造出U的一个本征态（这里省略了记号：modN）：

|u_0\rangle = \frac{1}{12^{1/2}}(|1\rangle+|3\rangle+|9\rangle+\dots+|12\rangle)

因为U作用在其上：

U|u_0\rangle = \frac{1}{12^{1/2}}(U|1\rangle+U|3\rangle+U|9\rangle+\dots+U|12\rangle)\\=\frac{1}{12^{1/2}}(|3\rangle+|9\rangle+\dots+|12\rangle+|1\rangle) =|u_0\rangle

类似的对于一般都态，我们可以构建U的本征态：

|u_0\rangle = \frac{1}{r^{1/2}}\sum_{k=0}^{r-1}|a^k\pmod N\rangle

由于这个本征态的本征值为1，这样我们无法获得更多关于周期的信息。为此我们构建一个新的本征态，它在每个计算基态上增加了一个相位因子：

|u_1\rangle=\frac{1}{r^{1/2}}\sum_{k=0}^{r-1}e^{-2\pi ik/r}|a^k \pmod N\rangle

此时U作用在其上，结果为为：

U|u_1\rangle = e^{2\pi i /r}|u_1\rangle

此时我们发现周期r，出现在了相位因子中。

在我们的a=3，N=35的例子中计算过程如下：

类似的我们可以构建更多的本征态：

|u_s\rangle=\frac{1}{r^{1/2}}\sum_{k=0}^{r-1}e^{-2\pi isk/r}|a^k \pmod N\rangle

U|u_s\rangle = e^{2\pi is /r}|u_s\rangle

0\leq s \leq r-1

这样做的好处在于，我们将这些本征态相加之后，除了 |1\rangle 以外的计算基态全部抵消了，即：

\frac{1}{r^{1/2}}\sum^{r-1}_{s=0}|u_s\rangle = |1\rangle 。

这是因为 \sum^{r-1}_{s=0}e^{-2\pi i s/r} =0 ，这里r为偶数。

例子：取a=7，N=15，则r=4，按照上述定义的四个本征态相加，结果如下

由于 |1\rangle =\frac{1}{r^{1/2}}\sum^{r-1}_{s=0}|u_s\rangle 为U的本征态的叠加，且 U|u_s\rangle = e^{2\pi is /r}|u_s\rangle 我们直接对 |1\rangle 执行量子相位估计，我们将测量得到一个相位， \psi = s / r ，s是0到r-1中的随机一个整数。多次测量后就能得到r。

电路图如下：

参考资料：

Qiskit

再思可矣：Shor 算法-无门槛学习！