第五章 进阶主题

第五章 量子计算进阶专题 (Advanced Topics in Quantum Computing)

前四章已经构建了理解量子计算的完整基础:从线性代数与复数(第一章)、量子力学原理(第二章)、量子比特与量子门(第三章),到量子傅里叶变换的数学形式和Shor算法的核心思想(第四章)。然而,这些章节中仍有一些关键的技术细节和前沿专题尚未深入展开。本章将填补这些空白,从五个维度提升读者对量子计算的理解深度:我们将展示如何将抽象的量子傅里叶变换转化为具体的量子电路;我们将学习量子相位估计——读取量子系统本征值的指数级加速技术;我们将直面真实量子计算机最大的敌人——噪声与退相干;我们将探索NISQ(含噪中等规模量子)时代最重要的算法范式——变分量子算法;最后,我们还将介绍5.55.5哈密顿量模拟简介——量子计算在量子化学和材料科学中最有前景的应用之一。这五个专题分别对应量子算法的电路实现、核心应用、物理限制、近期实践和核心应用前景,共同构成从理论到实验的完整桥梁。


5.1 量子傅里叶变换的电路实现 (QFT Circuit Construction)

从DFT到QFT:回顾与动机

回顾1.8节中学习的离散傅里叶变换(DFT)。对于N=2nN = 2^n个复数x0,x1,,xN1x_0, x_1, \ldots, x_{N-1},DFT将它们变换为一组新的复数y0,y1,,yN1y_0, y_1, \ldots, y_{N-1}

yk=1Nj=0N1xjωNjk,ωN=e2πi/Ny_k = \frac{1}{\sqrt{N}}\sum_{j=0}^{N-1}x_j\omega_N^{jk}, \quad \omega_N = e^{2\pi i/N}

在量子计算中,**量子傅里叶变换(Quantum Fourier Transform, QFT)**是DFT的量子实现。它不是对经典数据数组做变换,而是对量子态的振幅做变换。给定一个nn量子比特的态:

ψ=j=0N1xjj\lvert\psi\rangle = \sum_{j=0}^{N-1}x_j\lvert j\rangle

QFT将其变换为:

UQFTψ=k=0N1ykkU_{\text{QFT}}\lvert\psi\rangle = \sum_{k=0}^{N-1}y_k\lvert k\rangle

其中yky_k正是上述DFT系数。用基矢表示,QFT的作用为:

UQFTj=1Nk=0N1ωNjkkU_{\text{QFT}}\lvert j\rangle = \frac{1}{\sqrt{N}}\sum_{k=0}^{N-1}\omega_N^{jk}\lvert k\rangle

这是量子计算中最重要的变换之一。从矩阵角度看,UQFTU_{\text{QFT}}是一个2n×2n2^n \times 2^n的酉矩阵,其(j,k)(j,k)元为ωNjk/N\omega_N^{jk}/\sqrt{N}——这正是1.8节中DFT矩阵的量子版本。但关键在于:经典计算机需要O(NlogN)=O(2nn)O(N\log N) = O(2^n n)次操作来计算DFT,而量子计算机只需要O(n2)O(n^2)个量子门就能实现QFT——这是指数级的加速。

二进制分数表示:QFT的因式分解

为了将QFT转化为电路,我们需要将上述求和式因式分解为单量子比特的张量积形式。这是QFT电路构造的核心技巧。

首先引入二进制分数表示:对于一个二进制小数0.j1j2jm0.j_1 j_2 \ldots j_m,其值为:

0.j1j2jm=j12+j24++jm2m0.j_1 j_2 \ldots j_m = \frac{j_1}{2} + \frac{j_2}{4} + \cdots + \frac{j_m}{2^m}

例如,0.1=1/20.1 = 1/20.01=1/40.01 = 1/40.11=3/40.11 = 3/4

现在,将整数jjnn位二进制表示为j=j1j2jnj = j_1 j_2 \ldots j_n,即j=j12n1+j22n2++jn20j = j_1 2^{n-1} + j_2 2^{n-2} + \cdots + j_n 2^0。QFT作用在基矢j=j1j2jn\lvert j\rangle = \lvert j_1 j_2 \ldots j_n\rangle上,经过一系列代数运算(展开指数、分离变量),可以得到惊人的因式分解形式:

UQFTj1j2jn=12n/2(0+e2πi0.jn1)(0+e2πi0.jn1jn1)(0+e2πi0.j1j2jn1)U_{\text{QFT}}\lvert j_1 j_2 \ldots j_n\rangle = \frac{1}{2^{n/2}}\left(\lvert0\rangle + e^{2\pi i \cdot 0.j_n}\lvert1\rangle\right) \otimes \left(\lvert0\rangle + e^{2\pi i \cdot 0.j_{n-1}j_n}\lvert1\rangle\right) \otimes \cdots \otimes \left(\lvert0\rangle + e^{2\pi i \cdot 0.j_1 j_2 \ldots j_n}\lvert1\rangle\right)

让我们验证这个公式的正确性。以n=1n=1为例:j=j1j = j_1N=2N = 2,则:

UQFTj1=12(0+e2πi0.j11)U_{\text{QFT}}\lvert j_1\rangle = \frac{1}{\sqrt{2}}\left(\lvert0\rangle + e^{2\pi i \cdot 0.j_1}\lvert1\rangle\right)

j1=0j_1 = 0时,e2πi0=1e^{2\pi i \cdot 0} = 1,得12(0+1)\frac{1}{\sqrt{2}}(\lvert0\rangle + \lvert1\rangle);当j1=1j_1 = 1时,e2πi1/2=eπi=1e^{2\pi i \cdot 1/2} = e^{\pi i} = -1,得12(01)\frac{1}{\sqrt{2}}(\lvert0\rangle - \lvert1\rangle)。这正是阿达马门HH的作用!所以单量子比特QFT就是HH门。

这个因式分解的重要性在于:输出态已经写成了nn个单量子比特态的张量积。这意味着我们可以用nn个并行的单比特操作来制备输出态的每个因子——这是构造高效量子电路的关键。

受控相位旋转门 CRkCR_k

观察因式分解中的第ll个因子:

0+e2πi0.jnl+1jn1\lvert0\rangle + e^{2\pi i \cdot 0.j_{n-l+1} \ldots j_n}\lvert1\rangle

相位e2πi0.jnl+1jne^{2\pi i \cdot 0.j_{n-l+1} \ldots j_n}可以展开为:

e2πi0.jnl+1e2πi0.0jnl+2e2πi0.00jne^{2\pi i \cdot 0.j_{n-l+1}} \cdot e^{2\pi i \cdot 0.0j_{n-l+2}} \cdot \cdots \cdot e^{2\pi i \cdot 0.0\ldots 0j_n}

每个指数项对应一个受控相位旋转。为此,定义受控相位旋转门(Controlled-Phase Rotation) CRkCR_k

CRk=(100001000010000e2πi/2k)CR_k = \begin{pmatrix}1&0&0&0\\0&1&0&0\\0&0&1&0\\0&0&0&e^{2\pi i/2^k}\end{pmatrix}

它作用于两个量子比特:当控制比特为1\lvert1\rangle时,给目标比特的1\lvert1\rangle分量添加相位e2πi/2ke^{2\pi i/2^k};当控制比特为0\lvert0\rangle时,目标比特不变。注意CR1CR_1就是受控ZZ门(CZ),因为e2πi/2=eπi=1e^{2\pi i/2} = e^{\pi i} = -1

CRkCR_k的矩阵在计算基{00,01,10,11}\{\lvert00\rangle, \lvert01\rangle, \lvert10\rangle, \lvert11\rangle\}下是对角的,这使其在实验上相对容易实现(许多物理平台天然支持对角的两比特门)。

QFT的级联电路构造

现在我们可以构造完整的QFT电路。QFT电路从第一个量子比特到第nn个量子比特依次操作:

第1个量子比特:施加HH门,然后依次施加CR2,CR3,,CRnCR_2, CR_3, \ldots, CR_n,控制比特分别为第2, 3, \ldots, nn个量子比特。

  • HH门产生0+e2πi0.j11\lvert0\rangle + e^{2\pi i \cdot 0.j_1}\lvert1\rangle
  • CR2CR_2(控制=比特2)添加相位e2πi0.0j2=e2πij2/4e^{2\pi i \cdot 0.0j_2} = e^{2\pi i j_2/4},将相位变为e2πi0.j1j2e^{2\pi i \cdot 0.j_1 j_2}
  • CR3CR_3(控制=比特3)添加相位e2πij3/8e^{2\pi i j_3/8},将相位变为e2πi0.j1j2j3e^{2\pi i \cdot 0.j_1 j_2 j_3}
  • 以此类推…

第2个量子比特:施加HH门,然后依次施加CR2,CR3,,CRn1CR_2, CR_3, \ldots, CR_{n-1},控制比特分别为第3, 4, \ldots, nn个量子比特。

nn个量子比特:仅施加HH门。

33量子比特QFT电路如下:

j₁: |j₁⟩ —H—•—————•———     → 输出比特 3
              |     |
j₂: |j₂⟩ ————R₂—H—•———     → 输出比特 2
                    |
j₃: |j₃⟩ ————R₃—R₂—H—     → 输出比特 1

其中RkR_k表示CRkCR_k门,控制比特在上方,目标比特在下方。图中R2R_2从比特2控制比特1,R3R_3从比特3控制比特1,依此类推。

更一般地,nn量子比特QFT电路需要:

  • nnHH
  • k=1n1k=n(n1)/2\sum_{k=1}^{n-1}k = n(n-1)/2个受控相位门

总门数为O(n2)O(n^2)。由于同一层的HH门可以并行执行(不同比特上的HH门互不干扰),而受控相位门涉及不同比特对,电路深度也是O(n)O(n)(如果所有两比特门可以并行执行),或在最保守估计下为O(n2)O(n^2)

反向排序问题与SWAP门修正

仔细观察因式分解的输出:

第1个因子第2个因子n个因子\text{第1个因子} \otimes \text{第2个因子} \otimes \cdots \otimes \text{第$n$个因子}

其中第ll个因子对应e2πi0.jnl+1jne^{2\pi i \cdot 0.j_{n-l+1}\ldots j_n},即它只依赖于输入比特jnl+1,,jnj_{n-l+1}, \ldots, j_n。这意味着:

  • 第1个输出因子(最高位)只依赖于输入的最低位jnj_n
  • nn个输出因子(最低位)依赖于所有输入比特j1,,jnj_1, \ldots, j_n

因此,QFT电路的输出比特顺序是反转的。如果我们希望输出按照k1k2knk_1 k_2 \ldots k_n的正常顺序排列,需要在电路末尾添加n/2\lfloor n/2 \rfloor个SWAP门来反转比特顺序。

SWAP门交换两个量子比特的状态,其矩阵为:

SWAP=(1000001001000001)\text{SWAP} = \begin{pmatrix}1&0&0&0\\0&0&1&0\\0&1&0&0\\0&0&0&1\end{pmatrix}

SWAP可以用3个CNOT门实现:CNOT(1→2)、CNOT(2→1)、CNOT(1→2)。因此完整的QFT电路(含比特顺序修正)需要O(n2)O(n^2)个基本门。

3量子比特QFT的完整数值例子

让我们以一个具体的3量子比特例子验证QFT电路的正确性。取输入态ψ=011\lvert\psi\rangle = \lvert011\rangle(即j=3j = 3)。

步骤1:初始态 ψ0=011213\lvert\psi_0\rangle = \lvert0\rangle_1 \otimes \lvert1\rangle_2 \otimes \lvert1\rangle_3

步骤2:对比特1施加HH ψ1=12(01+11)1213\lvert\psi_1\rangle = \frac{1}{\sqrt{2}}(\lvert0\rangle_1 + \lvert1\rangle_1) \otimes \lvert1\rangle_2 \otimes \lvert1\rangle_3

步骤3:CR2CR_2(控制=比特2,目标=比特1)

控制比特2为1\lvert1\rangle,所以目标比特1的1\lvert1\rangle分量获得相位e2πi/4=eπi/2=ie^{2\pi i/4} = e^{\pi i/2} = i

ψ2=12(01+i11)1213\lvert\psi_2\rangle = \frac{1}{\sqrt{2}}(\lvert0\rangle_1 + i\lvert1\rangle_1) \otimes \lvert1\rangle_2 \otimes \lvert1\rangle_3

步骤4:CR3CR_3(控制=比特3,目标=比特1)

控制比特3为1\lvert1\rangle,所以目标比特1的1\lvert1\rangle分量获得额外相位e2πi/8=eπi/4=1+i2e^{2\pi i/8} = e^{\pi i/4} = \frac{1+i}{\sqrt{2}}

ieπi/4=eπi/2eπi/4=e3πi/4=12+i2i \cdot e^{\pi i/4} = e^{\pi i/2} \cdot e^{\pi i/4} = e^{3\pi i/4} = -\frac{1}{\sqrt{2}} + \frac{i}{\sqrt{2}}

所以:

ψ3=12(01+e3πi/411)1213\lvert\psi_3\rangle = \frac{1}{\sqrt{2}}(\lvert0\rangle_1 + e^{3\pi i/4}\lvert1\rangle_1) \otimes \lvert1\rangle_2 \otimes \lvert1\rangle_3

这对应0+e2πi0.0111=0+e2πi3/81\lvert0\rangle + e^{2\pi i \cdot 0.011}\lvert1\rangle = \lvert0\rangle + e^{2\pi i \cdot 3/8}\lvert1\rangle

步骤5:对比特2施加HH

比特2目前为1\lvert1\rangleH1=12(01)H\lvert1\rangle = \frac{1}{\sqrt{2}}(\lvert0\rangle - \lvert1\rangle)

ψ4=12(01+e3πi/411)12(0212)13\lvert\psi_4\rangle = \frac{1}{\sqrt{2}}(\lvert0\rangle_1 + e^{3\pi i/4}\lvert1\rangle_1) \otimes \frac{1}{\sqrt{2}}(\lvert0\rangle_2 - \lvert1\rangle_2) \otimes \lvert1\rangle_3

步骤6:CR2CR_2(控制=比特3,目标=比特2)

控制比特3为1\lvert1\rangle,比特2的1\lvert1\rangle分量获得相位ii

ψ5=12(01+e3πi/411)12(02+i(1)12)13\lvert\psi_5\rangle = \frac{1}{\sqrt{2}}(\lvert0\rangle_1 + e^{3\pi i/4}\lvert1\rangle_1) \otimes \frac{1}{\sqrt{2}}(\lvert0\rangle_2 + i(-1)\lvert1\rangle_2) \otimes \lvert1\rangle_3

=12(01+e3πi/411)12(02+e3πi/212)13= \frac{1}{\sqrt{2}}(\lvert0\rangle_1 + e^{3\pi i/4}\lvert1\rangle_1) \otimes \frac{1}{\sqrt{2}}(\lvert0\rangle_2 + e^{3\pi i/2}\lvert1\rangle_2) \otimes \lvert1\rangle_3

这对应0+e2πi0.111\lvert0\rangle + e^{2\pi i \cdot 0.11}\lvert1\rangle

步骤7:对比特3施加HH

ψ6=12(01+e3πi/411)12(02+e3πi/212)12(0313)\lvert\psi_6\rangle = \frac{1}{\sqrt{2}}(\lvert0\rangle_1 + e^{3\pi i/4}\lvert1\rangle_1) \otimes \frac{1}{\sqrt{2}}(\lvert0\rangle_2 + e^{3\pi i/2}\lvert1\rangle_2) \otimes \frac{1}{\sqrt{2}}(\lvert0\rangle_3 - \lvert1\rangle_3)

这对应0+e2πi0.11\lvert0\rangle + e^{2\pi i \cdot 0.1}\lvert1\rangle

步骤8:SWAP修正比特顺序(比特1↔比特3)

最终输出(忽略归一化因子1/23/21/2^{3/2})为8个基矢的叠加,振幅由e2πi3k/8/8e^{2\pi i \cdot 3k/8}/\sqrt{8}给出,其中k=0,1,,7k = 0, 1, \ldots, 7。这正是DFT作用于基矢3\lvert3\rangle的结果。

近似QFT与深度优化

在实际量子计算机上,非常小的相位旋转(如CRnCR_nnn很大时,e2πi/2n1e^{2\pi i/2^n} \approx 1)难以精确实现。一个实用的优化是近似QFT(Approximate QFT):丢弃所有k>kmaxk > k_{\max}CRkCR_k门,其中kmax=O(logn)k_{\max} = O(\log n)。这种近似的误差可以被控制在任意小的范围内,但电路深度降低到O(nlogn)O(n \log n),在某些架构下甚至可以优化到O(n)O(n)

小结 (Summary): 量子傅里叶变换通过将输出态因式分解为单量子比特的张量积,实现了高效的电路构造。QFT电路由nnHH门和n(n1)/2n(n-1)/2个受控相位旋转门组成,总复杂度为O(n2)O(n^2)个门。二进制分数表示揭示了每个输出比特只依赖于输入比特的一个子集,这是级联电路结构的关键。输出比特顺序反转需要通过SWAP门修正。3量子比特的数值例子验证了电路的正确性。

与量子计算的连接 (Connection to Quantum Computing): QFT是Shor算法和量子相位估计的核心子程序。从抽象的DFT矩阵到具体的HH门+受控相位门的级联电路,是量子算法从数学描述到物理实现的关键一步。理解QFT的电路构造,是理解这些指数级加速算法如何在真实量子硬件上运行的基础。O(n2)O(n^2)的电路复杂度意味着QFT可以在中等规模的量子计算机上实现,使其成为近期量子计算中最实用的变换之一。


5.2 量子相位估计 (Quantum Phase Estimation)

问题陈述:读取量子系统的本征值

量子相位估计(Quantum Phase Estimation, QPE)是量子计算中最强大的子程序之一。它解决如下问题:

给定:一个酉算子UU和它的一个本征向量u\lvert u\rangle,满足Uu=e2πiθuU\lvert u\rangle = e^{2\pi i \theta}\lvert u\rangle目标:估计相位θ[0,1)\theta \in [0, 1)的数值。

这个问题看似抽象,但它包含了量子计算中许多核心算法的本质。例如,在Shor算法中,UU是实现模幂运算的酉算子,u\lvert u\rangle是对应的周期态,而θ\theta的精确值直接给出因数分解所需的关键信息。

经典计算机上,估计一个酉算子的本征相位没有直接的方法——我们通常需要计算UU的全部矩阵元并做对角化,复杂度为O(N3)O(N^3)N=2nN = 2^n)。量子相位估计则利用量子并行性和QFT,以仅O(n2)O(n^2)个门的复杂度直接”读取”相位——这是指数级加速。

算法电路结构:辅助寄存器 + 受控酉操作 + 逆QFT

QPE电路由三部分组成,使用两个量子寄存器:

  1. 第一寄存器(辅助/计数寄存器)tt个量子比特,初始化为0t\lvert0\rangle^{\otimes t}
  2. 第二寄存器(目标寄存器):存储本征向量u\lvert u\rangle(假设已制备好)

电路结构如下:

|0⟩ —H—•——————•——————•———H†—S†—...—     → 测量 → θ的二进制估计
|0⟩ —H—|  U¹   |  U²   |   —H†—...—
|0⟩ —H—|       |       |     ...
  ⋮    |       |       |
|0⟩ —H—•———————•———————•———————...—
       |       |       |
|u⟩ ——U¹——U²——U⁴——...——U^(2^{t-1})—     → |u⟩ (不变)

步骤1:创建叠加态

对所有tt个辅助比特施加HH门:

12tx=02t1xu\frac{1}{\sqrt{2^t}}\sum_{x=0}^{2^t-1}\lvert x\rangle \otimes \lvert u\rangle

步骤2:受控酉操作编码相位

依次施加受控U2kU^{2^k}操作,控制比特为辅助寄存器的第kk个比特,目标为u\lvert u\rangle。受控U2kU^{2^k}的作用是:当控制比特为1\lvert1\rangle时,对目标施加U2kU^{2^k};当控制比特为0\lvert0\rangle时,目标不变。

数学上,对于辅助比特处于x=xt1xt2x0\lvert x\rangle = \lvert x_{t-1} x_{t-2} \ldots x_0\rangle(其中xk{0,1}x_k \in \{0,1\}),受控操作将目标态变为:

Ux0U2x1U4x2U2t1xt1u=Uxu=e2πiθxuU^{x_0} \cdot U^{2x_1} \cdot U^{4x_2} \cdots U^{2^{t-1}x_{t-1}}\lvert u\rangle = U^x\lvert u\rangle = e^{2\pi i \theta x}\lvert u\rangle

因此,整个辅助-目标系统的态变为:

12tx=02t1e2πiθxxu\frac{1}{\sqrt{2^t}}\sum_{x=0}^{2^t-1}e^{2\pi i \theta x}\lvert x\rangle \otimes \lvert u\rangle

这是关键的一步:相位θ\theta被编码到了辅助寄存器的振幅中!注意这个态的形式——它与QFT的输出态非常相似。

步骤3:逆QFT读取相位

对上述态的辅助寄存器部分施加逆量子傅里叶变换(Inverse QFT, QFT^\dagger。QFT^\dagger是QFT的厄米共轭,它将频率域的表示转换回时域。

如果θ\theta恰好可以用tt位二进制精确表示,即θ=0.θ1θ2θt\theta = 0.\theta_1 \theta_2 \ldots \theta_t,则逆QFT将辅助寄存器精确地变为θ1θ2θt\lvert\theta_1 \theta_2 \ldots \theta_t\rangle。测量辅助寄存器,我们就直接读出了θ\theta的二进制表示!

更精确地说,如果θ=j/2t\theta = j/2^t(即tt位二进制小数),则:

QFT(12tx=02t1e2πijx/2tx)=j\text{QFT}^\dagger\left(\frac{1}{\sqrt{2^t}}\sum_{x=0}^{2^t-1}e^{2\pi i j x/2^t}\lvert x\rangle\right) = \lvert j\rangle

这正是QFT的定义性质(正交性)的直接推论。

精度分析与成功概率

θ\theta不是tt位二进制小数时(即不能被精确表示),QPE给出的是θ\theta的最佳tt位近似。设θ\theta的真实值与最近的tt位二进制小数b/2tb/2^t的距离为δ\deltaδ1/2t+1|\delta| \leq 1/2^{t+1})。

测量结果落在bb附近的概率满足:

p(测量结果b<2tm)112(2mt1)p(|\text{测量结果} - b| < 2^{t-m}) \geq 1 - \frac{1}{2(2^{m-t}-1)}

特别地,要获得nn位精度的估计(即误差小于1/2n1/2^n),需要t=n+O(1)t = n + O(1)个辅助比特。更精确地说,如果希望以至少1ϵ1-\epsilon的概率获得nn位精度,需要:

t=n+log2(2+12ϵ)t = n + \left\lceil\log_2\left(2 + \frac{1}{2\epsilon}\right)\right\rceil

例如,要以99%的成功率获得nn位精度,大约需要t=n+8t = n + 8个辅助比特。

与Shor算法的关系

Shor算法的核心是阶数发现(Order Finding):给定与NN互质的整数aa,找到最小的正整数rr使得ar1(modN)a^r \equiv 1 \pmod{N}。这个问题可以表述为相位估计:

定义酉算子UUUx=axmodNU\lvert x\rangle = \lvert ax \mod N\rangle(作用于log2N\lceil\log_2 N\rceil个量子比特)。这个UU的本征态为:

us=1rk=0r1e2πisk/rakmodN\lvert u_s\rangle = \frac{1}{\sqrt{r}}\sum_{k=0}^{r-1}e^{-2\pi i s k/r}\lvert a^k \mod N\rangle

对应的本征值为e2πis/re^{2\pi i s/r},其中s=0,1,,r1s = 0, 1, \ldots, r-1

通过QPE估计s/rs/r,然后用连分数展开s/rs/r的近似中提取rr。这就是Shor算法中QFT/QPE发挥核心作用的机制。QPE将抽象的数论问题转化为量子系统本征相位的物理测量,是数学与物理的深刻交汇。

实现挑战:受控 U2kU^{2^k} 操作

QPE的理论框架优美而简洁,但其实验实现面临重大挑战:

  1. 制备本征态u\lvert u\rangle:通常我们不知道u\lvert u\rangle的具体形式。在实际中,我们制备1\lvert1\rangle(或其他简单态),它可以表示为所有本征态的叠加:1=sus/r\lvert1\rangle = \sum_s \lvert u_s\rangle/\sqrt{r}。QPE同时对多个本征相位进行估计,测量结果以一定概率给出某个s/rs/r的值——这足以完成阶数发现。

  2. 实现受控U2kU^{2^k}:这是最大的实验挑战。对于Shor算法中的模幂运算U2kx=a2kxmodNU^{2^k}\lvert x\rangle = \lvert a^{2^k}x \mod N\rangle,需要实现O((logN)3)O((\log N)^3)个基本门。虽然理论上可行,但当前量子硬件的门保真度和相干时间限制了可执行的U2kU^{2^k}规模。

  3. 量子比特数量:Shor算法分解一个LL位的整数,需要约2L2L个量子比特用于辅助寄存器,加上O(L)O(L)个量子比特用于目标寄存器。对于2048位RSA(L=2048L = 2048),需要约4000-5000个量子比特——这远超当前最先进的量子计算机(~1000量子比特)。

小结 (Summary): 量子相位估计通过创建辅助寄存器的叠加态、施加受控U2kU^{2^k}操作将相位编码到振幅中,再用逆QFT将相位信息转换为可测量的二进制表示,实现了对酉算子本征相位的指数级加速估计。tt个辅助比特提供约tO(1)t - O(1)位的精度估计。Shor算法中的阶数发现是QPE的直接应用。受控U2kU^{2^k}的实现是当前实验的主要瓶颈。

与量子计算的连接 (Connection to Quantum Computing): QPE是量子计算”读取”量子系统内禀信息的通用工具。它不仅用于Shor算法,还广泛应用于量子化学(计算分子能级)、量子模拟(读取系统演化频率)和量子机器学习等领域。理解QPE的电路结构和精度分析,是理解量子计算机如何从抽象的量子态中提取有用经典信息的关键。QPE展示了量子算法设计的典型范式:利用叠加态编码问题、用量子干涉提取答案。


5.3 噪声模型与退相干 (Noise Models & Decoherence)

开放量子系统:为什么噪声不可避免

在前面的章节中,我们讨论的都是封闭量子系统(Closed Quantum System)——一个完全孤立于外界环境、仅受精确控制的酉演化支配的理想系统。但在现实中,没有任何量子系统是完全孤立的。**开放量子系统(Open Quantum System)**与环境持续耦合,导致量子信息的退化和丢失。

这种耦合为什么不可避免?因为量子比特是真实的物理系统(电子自旋、超导电路、离子能级等),它们必然与周围的电磁场、晶格振动(声子)、热辐射等环境自由度相互作用。这些环境自由度构成了一个巨大的”热库”(Heat Bath),其希尔伯特空间的维度通常是量子比特的指数倍。量子比特与环境之间的纠缠,使得量子比特的纯态演化为我们将在2.6节中学习的混态——这正是退相干(Decoherence)的物理起源。

理解噪声不是可选的附加知识,而是评估任何量子计算方案实际可行性的基础。当前所有量子计算机都在”NISQ”(含噪中等规模量子)时代运行,噪声是限制计算规模和深度的首要因素。

算子和表示(Kraus表示)

对于开放量子系统,演化不再是纯态到纯态的酉映射,而是密度矩阵到密度矩阵的完全正定保迹映射(CPTP Map)。这种映射可以用算子和表示(Operator-Sum Representation),也称Kraus表示来刻画:

ρE(ρ)=kEkρEk\rho \to \mathcal{E}(\rho) = \sum_k E_k \rho E_k^\dagger

其中{Ek}\{E_k\}称为Kraus算子(Kraus Operators),满足完备性关系

kEkEk=I\sum_k E_k^\dagger E_k = I

这个条件保证了映射保迹:Tr(E(ρ))=Tr(ρ)=1\text{Tr}(\mathcal{E}(\rho)) = \text{Tr}(\rho) = 1

Kraus表示的物理意义是:环境对量子系统的作用可以看作一系列”量子操作”的统计混合。每个EkE_k对应环境的一种可能”响应”,而整个演化是这些响应的概率加权平均。这与经典噪声的信道模型完全类似——只是这里的”信道”作用在密度矩阵上,而非经典概率分布上。

如果只有一个Kraus算子E0=UE_0 = U(酉矩阵),则退化为封闭系统的酉演化:E(ρ)=UρU\mathcal{E}(\rho) = U\rho U^\dagger。因此,Kraus表示是酉演化的自然推广,统一描述了从纯酉演化到完全退极化的全部可能情况。

退极化信道(Depolarizing Channel)

**退极化信道(Depolarizing Channel)**是最简单的噪声模型之一,描述以概率pp将量子态完全随机化的过程:

E(ρ)=(1p)ρ+p3(XρX+YρY+ZρZ)\mathcal{E}(\rho) = (1-p)\rho + \frac{p}{3}(X\rho X + Y\rho Y + Z\rho Z)

对应的Kraus算子为:

E0=1pI,E1=p3X,E2=p3Y,E3=p3ZE_0 = \sqrt{1-p}\, I, \quad E_1 = \sqrt{\frac{p}{3}}\, X, \quad E_2 = \sqrt{\frac{p}{3}}\, Y, \quad E_3 = \sqrt{\frac{p}{3}}\, Z

验证完备性:E0E0+E1E1+E2E2+E3E3=(1p)I+p3(I+I+I)=IE_0^\dagger E_0 + E_1^\dagger E_1 + E_2^\dagger E_2 + E_3^\dagger E_3 = (1-p)I + \frac{p}{3}(I + I + I) = I

物理意义:以概率1p1-p,量子态不受影响;以概率p/3p/3分别施加XXYYZZ错误(比特翻转、相位翻转或两者兼有)。当p=1p = 1时,任意输入态都变为完全混态I/2I/2——信息完全丢失。

在布洛赫球上,退极化信道将所有点均匀地向球心收缩:

r(1p)r\vec{r} \to (1-p)\vec{r}

即布洛赫向量的长度缩小为原来的1p1-p倍,但方向不变。当p=1p = 1时,整个球面坍缩到球心r=0\vec{r} = \vec{0}

比特翻转信道(Bit-Flip Channel)

**比特翻转信道(Bit-Flip Channel)**以概率pp翻转量子比特(01\lvert0\rangle \leftrightarrow \lvert1\rangle),对应经典的”位错误”:

E(ρ)=(1p)ρ+pXρX\mathcal{E}(\rho) = (1-p)\rho + p X\rho X

Kraus算子:E0=1pIE_0 = \sqrt{1-p}\, IE1=pXE_1 = \sqrt{p}\, X

在布洛赫球上,比特翻转信道的效果是:

(x,y,z)(x,(12p)y,(12p)z)(x, y, z) \to (x, (1-2p)y, (1-2p)z)

不对,让我重新计算。设ρ=12(I+xX+yY+zZ)\rho = \frac{1}{2}(I + xX + yY + zZ),则:

E(ρ)=(1p)ρ+pXρX=12(I+xX+(12p)yY+(12p)zZ)\mathcal{E}(\rho) = (1-p)\rho + p X\rho X = \frac{1}{2}\left(I + xX + (1-2p)yY + (1-2p)zZ\right)

所以布洛赫坐标变换为:

xx,y(12p)y,z(12p)zx \to x, \quad y \to (1-2p)y, \quad z \to (1-2p)z

xx坐标保持不变,yyzz坐标按因子12p1-2p缩放。当p=1/2p = 1/2时,yyzz分量完全消失,态变为ρ=12(I+xX)\rho = \frac{1}{2}(I + xX)——如果初始x=0x = 0,则变为完全混态。比特翻转信道在xx轴上保持信息,但在yy-zz平面上压缩态。

相位翻转信道(Phase-Flip / Dephasing Channel)

相位翻转信道(Phase-Flip Channel),也称退相位信道(Dephasing Channel),以概率pp施加ZZ门(翻转1\lvert1\rangle的相位):

E(ρ)=(1p)ρ+pZρZ\mathcal{E}(\rho) = (1-p)\rho + p Z\rho Z

Kraus算子:E0=1pIE_0 = \sqrt{1-p}\, IE1=pZE_1 = \sqrt{p}\, Z

布洛赫球上的效果:

E(ρ)=12(I+(12p)xX+(12p)yY+zZ)\mathcal{E}(\rho) = \frac{1}{2}\left(I + (1-2p)xX + (1-2p)yY + zZ\right)

即:

x(12p)x,y(12p)y,zzx \to (1-2p)x, \quad y \to (1-2p)y, \quad z \to z

**这是退相干的核心机制!**相位翻转信道在xx-yy平面上压缩态(衰减相干项),但保持zz坐标(能量本征态的占据概率)不变。当p=1/2p = 1/2时,xxyy分量消失,态变为经典的概率混合ρ=1+z200+1z211\rho = \frac{1+z}{2}\lvert0\rangle\langle0\rvert + \frac{1-z}{2}\lvert1\rangle\langle1\rvert——量子相干性完全丢失,只剩下经典概率信息。

退相位信道的重要性在于:它是大多数物理系统中主导的噪声机制。由于能量守恒的约束,系统与环境交换能量的过程(导致zz变化)通常比纯相位信息丢失(xxyy衰减)慢得多。因此,T2T_2(相位相干时间)通常比T1T_1(能量弛豫时间)短——这是实验量子计算中最常见的限制因素。

幅值阻尼信道(Amplitude Damping)

**幅值阻尼信道(Amplitude Damping Channel)**描述能量从量子系统耗散到环境中的过程,例如激发态自发辐射衰变到基态:

E(ρ)=E0ρE0+E1ρE1\mathcal{E}(\rho) = E_0 \rho E_0^\dagger + E_1 \rho E_1^\dagger

其中:

E0=(1001γ),E1=(0γ00)E_0 = \begin{pmatrix}1&0\\0&\sqrt{1-\gamma}\end{pmatrix}, \quad E_1 = \begin{pmatrix}0&\sqrt{\gamma}\\0&0\end{pmatrix}

参数γ[0,1]\gamma \in [0, 1]是能量耗散的概率。E1E_11\lvert1\rangle(激发态)变为0\lvert0\rangle(基态),对应能量ω\hbar\omega释放到环境中。E0E_0缩减1\lvert1\rangle的振幅,保持0\lvert0\rangle不变。

验证完备性:

E0E0+E1E1=(1001γ)+(000γ)=IE_0^\dagger E_0 + E_1^\dagger E_1 = \begin{pmatrix}1&0\\0&1-\gamma\end{pmatrix} + \begin{pmatrix}0&0\\0&\gamma\end{pmatrix} = I

在布洛赫球上,幅值阻尼将态向北极(0\lvert0\rangle,基态)吸引:

ρ=12(1+zxiyx+iy1z)12(1+z+γ(1z)1γ(xiy)1γ(x+iy)(1γ)(1z))\rho = \frac{1}{2}\begin{pmatrix}1+z&x-iy\\x+iy&1-z\end{pmatrix} \to \frac{1}{2}\begin{pmatrix}1+z+\gamma(1-z)&\sqrt{1-\gamma}(x-iy)\\\sqrt{1-\gamma}(x+iy)&(1-\gamma)(1-z)\end{pmatrix}

坐标变换为:

x1γx,y1γy,zz+γ(1z)=1(1γ)(1z)x \to \sqrt{1-\gamma}\, x, \quad y \to \sqrt{1-\gamma}\, y, \quad z \to z + \gamma(1-z) = 1 - (1-\gamma)(1-z)

γ1\gamma \to 1时,任意态都趋近于0\lvert0\rangle。当γ1\gamma \ll 1时,系统经历指数衰减et/T1e^{-t/T_1},其中T1T_1能量弛豫时间(Energy Relaxation Time)

5.3.4 量子误差缓解 (Quantum Error Mitigation)

量子纠错(4.4节)通过冗余编码主动纠正错误,但需要大量物理量子比特和远低于阈值的物理错误率。在NISQ时代(含噪中等规模量子),硬件尚未达到纠错阈值,于是出现了另一条技术路线——量子误差缓解 (Quantum Error Mitigation, QEM)。与QEC不同,QEM不纠正错误,而是在经典后处理中移除噪声的影响,以无偏方式估计无噪声的期望值。QEM的核心观察是:我们虽然无法在量子硬件上消除噪声,但可以通过巧妙地改变噪声并测量其影响,在经典计算机上”减去”噪声贡献。

零噪声外推 (Zero-Noise Extrapolation, ZNE)

  • 核心思想:在不同噪声水平下测量同一算符的期望值,然后外推到零噪声极限
  • 如何控制噪声水平:通过脉冲拉伸(pulse stretching,在超导平台上延长门操作时间)或门折叠(gate folding,插入恒等门的分解来增加电路深度)来系统性地放大噪声
  • 拟合模型:O(λ)=O0+a1λ+a2λ2+\langle O \rangle(\lambda) = \langle O \rangle_0 + a_1\lambda + a_2\lambda^2 + \cdots,其中 λ\lambda 是噪声放大因子
  • Richardson外推或指数拟合
  • 实验验证:ZNE技术在VQE中可将能量估计误差降低5-10倍(IBM、Rigetti 2020-2022)

概率误差消除 (Probabilistic Error Cancellation, PEC)

  • 核心思想:将噪声量子门表示为理想门与噪声通道的组合,通过蒙特卡洛采样反转噪声
  • 实现方式:将每个噪声门分解为理想门的线性组合,以”准概率”方式采样,乘以符号因子 γ\gamma 来修正期望值
  • 代价:采样方差随 γ2\gamma^2 增长,γ\gamma 随电路深度指数级增长——因此PEC仅适用于浅层电路
  • 适用范围:与ZNE相比,PEC在低深度下提供更精确的修正但成本更高

虚拟态蒸馏 (Virtual Distillation)

  • 较新的方法(2021-2022),通过制备 MM 份副本并测量交换算符的期望值来”蒸馏”无噪声态
  • 不需要额外的辅助比特,但需要 SWAP\text{SWAP} 门网络
  • 适用于中等深度电路,是一种介于QEC和QEM之间的方法

测量误差修正

  • 最简单的QEM形式:通过制备已知基态(如 0n\lvert0\rangle^{\otimes n}1n\lvert1\rangle^{\otimes n})并测量,构建测量误差响应矩阵 AA,然后通过 Pcorrected=A1PrawP_{\text{corrected}} = A^{-1} P_{\text{raw}} 修正所有后续测量
  • 复杂度:O(2n)O(2^n) 的完整矩阵修正,或使用张量积分解近似

QEC vs QEM 对比

维度量子纠错 (QEC)量子误差缓解 (QEM)
目标主动纠正错误后处理移除噪声偏差
量子开销大量辅助比特(d2:1d^2:1额外电路执行(2-10x深度)
经典开销实时解码器蒙特卡洛/外推拟合
误差缩放指数级(低于阈值时)多项式级(外推次数)
适用期容错时代(2030+)NISQ时代(当前)
已展示低于阈值(Willow 10510^{-5}2-10x 精度提升(IBM/Google)
T1T_1T2T_2 时间

实验量子计算中,噪声用两个特征时间刻画:

T1T_1(能量弛豫时间,Longitudinal Relaxation Time):描述能量从激发态泄漏到基态的速率。如果t=0t = 0时系统处于1\lvert1\rangle,则zz坐标的演化满足:

z(t)=2et/T11z(t) = 2e^{-t/T_1} - 1

即激发态占据概率按et/T1e^{-t/T_1}衰减。T1T_1越大,量子比特保持能量的能力越强。

T2T_2(退相干时间,Transverse Relaxation Time):描述量子叠加态(xx-yy平面上的态)失去相干性的速率。xxyy坐标的衰减满足:

x(t)=x(0)et/T2,y(t)=y(0)et/T2x(t) = x(0)\, e^{-t/T_2}, \quad y(t) = y(0)\, e^{-t/T_2}

T2T_2T1T_1的关系为:

1T2=12T1+1Tϕ\frac{1}{T_2} = \frac{1}{2T_1} + \frac{1}{T_\phi}

其中TϕT_\phi纯退相位时间,描述纯相位噪声(能量守恒的相位随机化)的速率。通常TϕT1T_\phi \ll T_1,因此T2T1T_2 \ll T_1——相位相干比能量弛豫快得多地丢失。

当前最先进的超导量子比特:T1100500μsT_1 \sim 100\text{--}500\,\mu\text{s}T250200μsT_2 \sim 50\text{--}200\,\mu\text{s}。这意味着量子计算必须在远小于100μs\sim 100\,\mu\text{s}的时间内完成,否则噪声将淹没量子信号。由于典型单量子比特门时间约1050ns\sim 10\text{--}50\,\text{ns},当前系统可以执行约10310410^3\text{--}10^4个门操作——这定义了NISQ时代的计算边界。

噪声对量子纠错的启示

上述噪声模型的分析揭示了一个关键事实:单量子比特的错误可以看作是IIXXYYZZ的统计组合。这暗示了一种对抗噪声的策略:如果我们能够检测并纠正XXYYZZ错误,就能保护量子信息。

然而,量子纠错的挑战比经典纠错大得多,原因有三:

  1. 连续错误:量子错误是连续的(布洛赫球上的任意旋转),而经典错误是离散的(比特翻转)。虽然Kraus表示将连续演化分解为离散操作的混合,但纠错电路本身也受噪声影响。

  2. 不可克隆定理:我们不能复制量子信息来保护它,因此必须采用更精巧的编码策略(如将1个逻辑量子比特编码到多个物理量子比特中)。

  3. 测量破坏相干性:为了检测错误,我们需要测量某些可观测量,但测量会导致坍缩。量子纠错码(如Steane码、Shor码、表面码)通过测量**稳定子(Stabilizer)**而不是量子比特本身来绕过这一限制——稳定子测量揭示错误信息而不破坏编码的量子态。

**表面码(Surface Code)**是目前最有前景的量子纠错方案。它使用二维晶格上的物理量子比特编码一个逻辑量子比特,通过局域的稳定子测量检测XXZZ错误。表面码的阈值定理表明:如果单量子比特的错误率低于约1%1\%,逻辑错误率可以任意降低——这是容错量子计算的理论基础。

小结 (Summary): 开放量子系统与环境耦合导致信息丢失,用Kraus算子和表示统一描述。退极化信道以概率pp完全随机化量子态;比特翻转信道保留xx坐标但压缩yyzz;相位翻转(退相位)信道是退相干的核心机制,衰减xxyy但保留zz;幅值阻尼描述能量耗散,将态向基态吸引。T1T_1T2T_2时间刻画噪声强度,当前硬件限制NISQ时代电路深度约10310410^3\text{--}10^4门。量子纠错通过多比特编码和稳定子测量对抗噪声,表面码是最有前景的容错方案;在NISQ时代,量子误差缓解(ZNE、PEC、虚拟态蒸馏)通过经典后处理以无偏方式估计无噪声期望值,为当前硬件提供了另一种噪声应对策略,与面向容错时代的QEC形成互补。

与量子计算的连接 (Connection to Quantum Computing): 噪声是量子计算从理论走向实践的最大障碍。理解噪声模型不仅帮助我们评估量子算法的实际可行性,还指导量子纠错码的设计。NISQ时代的变分量子算法(5.4节)之所以被寄予厚望,正是因为它们对噪声的鲁棒性优于需要深度电路的算法(如Shor算法)。从退相位信道到T2T_2时间,从Kraus表示到表面码,噪声理论架起了抽象量子算法与真实物理硬件之间的桥梁。


5.4 变分量子算法 (Variational Quantum Algorithms)

量子-经典混合架构

我们已学习的量子算法(Shor、Grover、QPE)都假设理想的容错量子计算机——数百万量子比特、极低的错误率、足够长的相干时间。但当前和近期可预见的量子硬件远未达到这一标准。**NISQ(Noisy Intermediate-Scale Quantum)**时代的设备有501000\sim 50\text{--}1000个量子比特、有限相干时间和显著的门错误率。在这种限制下,需要新的算法范式。

**变分量子算法(Variational Quantum Algorithms, VQA)**应运而生。它的核心思想是:将量子计算机用作”参数化态制备器”,将困难的优化问题交给经典计算机处理。

混合架构的工作流程为:

  1. 经典优化器选择一组参数θ\vec{\theta}
  2. 量子计算机用参数化电路U(θ)U(\vec{\theta})制备试探态ψ(θ)\lvert\psi(\vec{\theta})\rangle
  3. 量子计算机测量该态的某些可观测量,得到能量(或其他代价函数)的估计值E(θ)E(\vec{\theta})
  4. 经典计算机根据E(θ)E(\vec{\theta})更新参数θθ\vec{\theta} \to \vec{\theta}'
  5. 重复步骤2—4,直到收敛

这个循环利用了量子计算机在态制备和测量方面的优势(指数级大希尔伯特空间中的操作),同时避开了深度量子电路的噪声积累——参数化电路通常很浅(O(1)O(1)O(logn)O(\log n)深度)。经典优化器处理参数更新,利用成熟的经典优化算法(梯度下降、Adam、L-BFGS等)。

变分原理:能量上界

变分量子算法的理论基础是变分原理(Variational Principle):对于任意量子态ψ\lvert\psi\rangle,哈密顿量HH的期望值满足:

ψHψE0\langle\psi|H|\psi\rangle \geq E_0

其中E0E_0HH的基态能量(最低本征值)。等号当且仅当ψ\lvert\psi\rangle是基态时成立。

这意味着:如果我们能够制备一个接近基态的试探态ψ(θ)\lvert\psi(\vec{\theta})\rangle,其能量期望值就给出了基态能量的上界。通过最小化Hθ\langle H \rangle_{\vec{\theta}},我们同时获得了基态能量的估计和近似基态本身。

这个原理在量子化学中尤为重要:分子的电子结构由薛定谔方程Hψ=EψH\lvert\psi\rangle = E\lvert\psi\rangle描述,但精确求解对于多电子系统(超过约20个电子)在经典计算机上是不可能的。变分量子算法提供了用NISQ设备计算分子基态能量的可能路径。

VQE:变分量子本征求解器

**变分量子本征求解器(Variational Quantum Eigensolver, VQE)**是VQA家族中最著名的算法,专门用于寻找量子系统的基态能量和基态波函数。

算法流程

1. 参数化量子电路(Ansatz)

选择一个参数化电路U(θ)U(\vec{\theta}),将初始态0n\lvert0\rangle^{\otimes n}变换为试探态:

ψ(θ)=U(θ)0n\lvert\psi(\vec{\theta})\rangle = U(\vec{\theta})\lvert0\rangle^{\otimes n}

常见的Ansatz包括:

  • 硬件高效Ansatz(Hardware-Efficient Ansatz, HEA):使用与目标量子硬件原生门集匹配的简单门层(如单比特旋转+纠缠门)。深度浅,但可能难以表示复杂的基态。
  • UCCSD(Unitary Coupled Cluster with Single and Double excitations):基于量子化学的激发算子构建,物理动机强,但电路深度较大。

2. 能量测量

将哈密顿量HH分解为泡利算子的线性组合(通过Jordan-Wigner变换Bravyi-Kitaev变换):

H=kckPk,Pk{I,X,Y,Z}nH = \sum_k c_k P_k, \quad P_k \in \{I, X, Y, Z\}^{\otimes n}

能量的期望值为:

Hθ=kckψ(θ)Pkψ(θ)\langle H \rangle_{\vec{\theta}} = \sum_k c_k \langle\psi(\vec{\theta})|P_k|\psi(\vec{\theta})\rangle

每个Pk\langle P_k \rangle可以通过适当的基变换后在计算基下测量得到。例如,测量X\langle X \rangle时,先对所有相关比特施加HH门,再测量ZZ

3. 经典优化

使用经典优化器最小化Hθ\langle H \rangle_{\vec{\theta}}。常用的优化器包括:

  • 梯度下降:需要估计梯度H/θi\partial \langle H \rangle / \partial \theta_i
  • SPSA(Simultaneous Perturbation Stochastic Approximation):只用2次函数评估估计所有梯度分量,适合噪声环境
  • L-BFGS:拟牛顿法,收敛快但对噪声敏感

梯度估计的量子方法

参数偏移规则(Parameter Shift Rule)允许精确计算梯度:

Hθi=12(Hθi+π/2Hθiπ/2)\frac{\partial \langle H \rangle}{\partial \theta_i} = \frac{1}{2}\left(\langle H \rangle_{\theta_i + \pi/2} - \langle H \rangle_{\theta_i - \pi/2}\right)

这只需要在θi±π/2\theta_i \pm \pi/2处分别评估能量,无需数值微分。

应用:量子化学

VQE最重要的应用是计算分子的基态能量。例如,计算氢分子H2H_2的基态能量需要4个量子比特(4个自旋轨道),VQE可以在当前的NISQ设备上实现。对于更大的分子(如N2N_2H2OH_2O),需要更多量子比特和更深的电路,但仍是近期最有希望的量子优势应用场景之一。

QAOA:量子近似优化算法

量子近似优化算法(Quantum Approximate Optimization Algorithm, QAOA)由Farhi、Goldstone和Gutmann于2014年提出,专门用于解决组合优化问题

问题设定

组合优化问题可以表述为:最小化一个定义在nn位二进制串上的代价函数C(z)C(z),其中z{0,1}nz \in \{0,1\}^n。例如,MaxCut问题:给定一个图G=(V,E)G = (V, E),将顶点分成两组,使得两组之间的边数最大。

将代价函数映射为问题哈密顿量(Problem Hamiltonian)CC

C=zC(z)zzC = \sum_{z} C(z) \lvert z\rangle\langle z\rvert

基态zopt\lvert z_{\text{opt}}\rangle对应最优解。但直接寻找基态是困难的(NP-hard)。

QAOA电路结构

QAOA使用pp层交替的酉演化:

ψ(γ,β)=eiβpBeiγpCeiβ1Beiγ1C+n\lvert\psi(\vec{\gamma}, \vec{\beta})\rangle = e^{-i\beta_p B}e^{-i\gamma_p C} \cdots e^{-i\beta_1 B}e^{-i\gamma_1 C} \lvert+\rangle^{\otimes n}

其中:

  • 混合哈密顿量(Mixer Hamiltonian)B=j=1nXjB = \sum_{j=1}^{n} X_j,即所有比特上的XX门之和。它的作用是驱动系统在解空间中进行”随机游走”。
  • 问题哈密顿量C=j,kE12(IZjZk)C = \sum_{\langle j,k\rangle \in E} \frac{1}{2}(I - Z_j Z_k)(以MaxCut为例)
  • 初始态+n\lvert+\rangle^{\otimes n},所有比特的等幅叠加

参数γ=(γ1,,γp)\vec{\gamma} = (\gamma_1, \ldots, \gamma_p)β=(β1,,βp)\vec{\beta} = (\beta_1, \ldots, \beta_p)通过经典优化器调整,以最小化期望值C\langle C \rangle

p=1p=1的QAOA电路示例(3比特MaxCut)

|0⟩ —H—Rz(γ₁)—•——————Rx(β₁)—
               | (ZZ)
|0⟩ —H—Rz(γ₁)—⊕—•—————Rx(β₁)—
                  | (ZZ)
|0⟩ —H—Rz(γ₁)———⊕——Rx(β₁)—

其中Rz(γ)Rz(\gamma)来自eiγCe^{-i\gamma C}的对角演化,ZZZZ门(受控ZZ门的变形)实现eiγZjZk/2e^{-i\gamma Z_j Z_k/2}

应用

QAOA已被应用于MaxCut、图着色、旅行商问题(TSP)、投资组合优化等问题。虽然对于一般问题,QAOA的近似比是否能超越经典算法仍是开放问题,但对于某些特定问题实例,QAOA已展示出优于经典贪心算法的性能。

变分算法的挑战

尽管变分量子算法为NISQ时代提供了实用路径,但它们面临几个根本性挑战:

1. 贫瘠高原(Barren Plateaus)

这是变分量子算法最严峻的理论挑战。研究表明,对于深度为O(poly(n))O(\text{poly}(n))的随机参数化电路,代价函数的梯度H/θi\partial \langle H \rangle / \partial \theta_i的方差随量子比特数nn指数级减小:

Var(Hθi)12n\text{Var}\left(\frac{\partial \langle H \rangle}{\partial \theta_i}\right) \sim \frac{1}{2^n}

这意味着当电路规模增大时,梯度变得极其微小且被噪声淹没——优化器无法找到下降方向。这种现象被称为贫瘠高原(Barren Plateaus),因为代价函数 landscape 像贫瘠的高原一样平坦,没有明显的梯度信息。

缓解策略包括:

  • 使用局部代价函数(而非全局代价函数)
  • 设计问题特定的浅层Ansatz
  • 采用分层优化策略

2. 测量采样开销

估计能量期望值H=kckPk\langle H \rangle = \sum_k c_k \langle P_k \rangle需要测量每一项Pk\langle P_k \rangle。每个Pk\langle P_k \rangle的估计精度为ϵ\epsilon需要O(1/ϵ2)O(1/\epsilon^2)次测量(由中心极限定理)。如果哈密顿量有MM项,则总采样复杂度为O(M/ϵ2)O(M/\epsilon^2)。对于量子化学问题,MM可以非常大(O(n4)O(n^4)),导致巨大的测量开销。

缓解策略包括:

  • 分组对易(Grouping Commuting Pauli Strings):同时测量对易的泡利算子串
  • 经典阴影(Classical Shadows):用随机测量一次性获取多个可观测量的信息

3. 电路深度与表达能力的权衡

更深的电路可以表示更复杂的量子态(表达能力更强),但也更容易受到噪声影响(保真度更低)。Ansatz的设计需要在表达能力和噪声鲁棒性之间找到平衡。硬件高效Ansatz深度浅、噪声鲁棒,但可能无法近似目标基态;UCCSD等化学启发Ansatz表达能力强,但电路深度可能超过NISQ设备的限制。

变分算法与容错量子计算的关系

变分量子算法不是容错量子计算的替代品,而是过渡期(NISQ时代)的核心算法范式。它们的关系可以概括如下:

特性NISQ变分算法容错量子算法
电路深度O(1)O(1)O(logn)O(\log n)O(poly(n))O(\text{poly}(n))
量子比特数10210310^2\text{--}10^3106+10^6+
错误率容忍10210^{-2}10410^{-4}或更低
经典优化必需不需要
应用场景量子化学、优化因子分解、大分子模拟

VQE在量子化学中的应用被认为是最有希望在NISQ设备上实现量子实用优势的领域。虽然经典方法(如密度矩阵重整化群DMRG、耦合簇CCSD(T))对于小分子非常精确,但对于强关联系统(如过渡金属催化剂、高温超导体的Hubbard模型),经典方法失效,而VQE可能提供新的计算路径。

从更广阔的视角看,变分量子算法代表了量子计算的一种”实用主义”转向:与其等待完美的容错量子计算机,不如利用现有(不完美)的量子硬件,结合经典计算,解决有实际价值的科学和工程问题。这种混合范式可能是量子计算在未来十年内产生实际影响的主要途径。

小结 (Summary): 变分量子算法采用量子-经典混合架构,将参数化量子电路用于态制备,经典优化器用于参数更新。VQE利用变分原理ψHψE0\langle\psi|H|\psi\rangle \geq E_0估计基态能量,在量子化学中有重要应用。QAOA通过交替应用问题哈密顿量和混合哈密顿量解决组合优化问题。变分算法面临三大挑战:贫瘠高原(梯度指数级消失)、测量采样开销O(1/ϵ2)O(1/\epsilon^2)、电路深度与表达能力的权衡。它们是NISQ时代最有希望的实用量子计算范式。

与量子计算的连接 (Connection to Quantum Computing): 变分量子算法架起了当前NISQ硬件与未来容错量子计算机之间的桥梁。它们展示了量子计算如何在噪声和规模限制下产生实际价值。从VQE的量子化学应用到QAOA的组合优化,变分算法将抽象的量子力学原理转化为解决实际问题的计算工具。理解变分算法的原理、优势和局限,是评估量子计算在近期内能否产生”量子优势”的关键。无论未来量子硬件如何发展,变分优化的思想——用量子态制备和经典优化的协同来解决复杂问题——都将继续在量子计算中扮演重要角色。


5.5 哈密顿量模拟简介 (Hamiltonian Simulation)

问题陈述:给定一个哈密顿量 HH(描述量子系统的总能量),计算时间演化算子 U(t)=eiHtU(t) = e^{-iHt} 并模拟其作用于初始态的效果。这是量子计算最有前景的应用之一——用于量子化学(分子基态计算)、材料科学(电子结构)、高能物理(晶格规范理论)。

为什么经典模拟困难HH 作用于 nn 个量子比特的系统,其矩阵大小为 2n×2n2^n \times 2^n——指数级。但 HH 通常是稀疏的(如量子化学中的 H=tpqapaq+vpqrsapaqarasH = \sum t_{pq} a_p^\dagger a_q + \sum v_{pqrs} a_p^\dagger a_q^\dagger a_r a_s),只涉及少数几项。

Trotter分解(Trotter-Suzuki分解——最基本的方法):

  • H=j=1mHjH = \sum_{j=1}^m H_j,且每个 HjH_j 容易模拟(如仅作用于少数比特的Pauli串),则 eiHt(j=1meiHjt/N)Ne^{-iHt} \approx \left(\prod_{j=1}^m e^{-iH_j t/N}\right)^N 其中 NN 是 Trotter 步数,误差 O(t2/N)O(t^2/N)
  • 物理直觉:将时间 tt 切分为 NN 个小段 t/Nt/N,在每个小段内近似认为各 HjH_j 对易
  • 复杂度:O(m2t2/ϵ)O(m^2 t^2/\epsilon) 达到精度 ϵ\epsilon

高阶 Trotter 方法

  • 二阶分解:ei(A+B)ΔteiAΔt/2eiBΔteiAΔt/2e^{-i(A+B)\Delta t} \approx e^{-iA\Delta t/2} e^{-iB\Delta t} e^{-iA\Delta t/2}
  • 四阶分解:通过嵌套二阶分解,误差从 O(Δt2)O(\Delta t^2) 降至 O(Δt4)O(\Delta t^4)

后 Trotter 方法(近年突破):

  • 量子信号处理 (Quantum Signal Processing):通过量子奇异值变换(QSVT)实现最优的 O(t+log(1/ϵ))O(t + \log(1/\epsilon)) 复杂度
  • 泰勒级数法:将 eiHte^{-iHt} 展开为泰勒级数,通过线性组合酉算子(LCU)实现
  • 量子行走法:将哈密顿量模拟转化为量子行走问题

应用

  • 量子化学:模拟分子 H2,LiH,Fe2S2H_2, LiH, Fe_2S_2 的电子结构(已经在小规模上实验验证)
  • 凝聚态物理:Hubbard模型、t-J模型的基态和动力学
  • 高能物理:Schwinger模型、晶格QED的模拟

当前局限

  • 实用有用的量子化学模拟需要 106-10810^6\text{-}10^8 个量子门,远超当前硬件能力
  • 纠错后的资源开销仍然巨大(可能需要数千逻辑量子比特运行数千小时)
  • Trotter误差与门保真度的耦合尚需深入研究

小结 (Summary): 哈密顿量模拟旨在计算时间演化算子 eiHte^{-iHt},是量子计算在量子化学、材料科学和高能物理中最有前景的应用。经典模拟因指数级矩阵规模而困难,但哈密顿量通常是稀疏的。Trotter分解将总哈密顿量拆分为易模拟项的乘积,是最基本的方法;高阶Trotter方法通过对称化分解降低误差。近年突破包括量子信号处理(QSVT)、泰勒级数法和量子行走法,实现了更优的复杂度。当前局限在于实用规模的模拟需要远超现有硬件能力的门数,纠错后的资源开销仍然巨大。

与量子计算的连接 (Connection to Quantum Computing): 哈密顿量模拟是量子计算”杀手级应用”的最有力候选之一。它直接利用量子系统的自然演化来模拟其他量子系统,避免了经典计算机指数级存储和计算的瓶颈。从Trotter分解到量子信号处理,哈密顿量模拟的发展展示了量子算法如何从基本物理直觉出发,通过数学深化实现复杂度突破。它与变分量子算法(如VQE)紧密相连:VQE用于寻找基态能量,而哈密顿量模拟用于研究系统动力学。理解哈密顿量模拟的原理和局限,是评估量子计算能否在量子化学和材料科学中产生实际优势的关键。


第五章总结

本章从四个关键维度补全了量子计算进阶专题的知识拼图:

  • 5.3.4 量子误差缓解(ZNE、PEC、虚拟态蒸馏)填补了NISQ时代”如何在纠错不可用的情况下运行算法”的核心缺口
  • 5.1-5.2 QFT电路实现与QPE构成了Shor算法和量子模拟的算法引擎,是量子指数级加速的技术核心
  • 5.3 噪声模型(Kraus算子、退极化/比特翻转/相位翻转/幅值阻尼信道)为理解真实量子硬件的物理限制提供了语言
  • 5.4 变分量子算法(VQE、QAOA)展示了量子-经典混合范式在NISQ时代的实用价值
  • 5.5 哈密顿量模拟作为量子计算的”杀手级应用”候选,连接了算法理论与实际应用场景

这五个专题的共同主题是:量子计算从理论走向实践所面临的核心挑战——噪声、退相干、纠错开销和算法设计约束。掌握这些专题,读者现在不仅能理解”量子计算如何工作”,更能评估”量子计算何时、以何种方式产生实际价值”。


附录

量子算法详解:从原理到电路实现

补充材料 — 量子计算前置教程(第三部分 3.6 节的深化与扩展)

本文档为教程第三部分 3.6 节(量子算法导引)中概述的五种核心量子算法提供完整的 数学推导、电路构造、复杂度证明和工作示例。读者应已完成教程 Part 1–3 的学习, 熟悉复数、线性代数、量子比特、量子门、量子电路以及量子傅里叶变换(QFT)的基础概念。

Part 5.1(QFT 电路实现)和 Part 5.2(量子相位估计 QPE)的内容将被 Shor 算法直接引用。


目录

  1. Deutsch-Jozsa 算法
  2. Bernstein-Vazirani 算法(附加)
  3. Simon 算法(附加)
  4. Grover 搜索算法
  5. Shor 因子分解算法
  6. 算法复杂度对比总表
  7. 参考文献与延伸阅读

1. Deutsch-Jozsa 算法

1.1 问题定义与经典复杂度

问题:给定一个布尔函数 f:{0,1}n{0,1}f: \{0,1\}^n \to \{0,1\},承诺 ff 要么是常数函数(对所有 x{0,1}nx \in \{0,1\}^nf(x)f(x) 取值相同),要么是平衡函数(恰好对一半输入输出 00,一半输出 11)。判定 ff 是常数函数还是平衡函数。

经典复杂度分析(最坏情况):

  • 确定性算法:在最坏情况下需要 2n1+12^{n-1}+1 次查询。因为即使前 2n12^{n-1} 次查询都得到了相同的结果,仍然不能判定 ff 是常数函数——可能剩下的 2n12^{n-1} 次全部是相反结果,使得 ff 恰好平衡。只有查询到第 2n1+12^{n-1}+1 次并得到相同结果时,才能确定 ff 是常数函数。因此确定性查询复杂度为 Θ(2n)\Theta(2^n)
  • 随机化算法:如果接受概率错误,可以更高效。但最坏情况下仍然需要指数级查询。

量子复杂度:仅需 1 次 查询。这是量子计算最早展示指数级加速的算法(尽管针对的是一个人为构造的问题)。


1.2 Oracle(预言机)的构造

量子算法中,函数 ff 不是作为”黑箱”被动查询的,而是通过一个量子 oracle(量子预言机) 实现的酉算子。oracle 是量子电路的基本组件:它将函数 ff 编码为可逆的酉变换。

标准的 Deutsch-Jozsa oracle 实现采用相位 oracle(Phase Oracle) 形式:

Uf:xy    xyf(x)U_f: \lvert x\rangle \lvert y\rangle \;\longmapsto\; \lvert x\rangle \lvert y \oplus f(x)\rangle

其中 x\lvert x\ranglenn 比特输入寄存器,y\lvert y\rangle 是单比特输出寄存器,\oplus 表示模 2 加法(XOR)。

矩阵表示:UfU_f 在计算基 {xy}\{\lvert x\rangle\lvert y\rangle\} 下是一个对角矩阵加上交换操作。具体地,对每个 xx

  • f(x)=0f(x)=0Ufxy=xyU_f\lvert x\rangle\lvert y\rangle = \lvert x\rangle\lvert y\rangle(不变)
  • f(x)=1f(x)=1Ufx0=x1U_f\lvert x\rangle\lvert 0\rangle = \lvert x\rangle\lvert 1\rangleUfx1=x0U_f\lvert x\rangle\lvert 1\rangle = \lvert x\rangle\lvert 0\rangle(在输出比特上施加 XX 门)

UfU_f 的酉性验证:Uf2=IU_f^2 = I,因为两次 XOR 操作恢复原态。UfU_f 也是一个置换矩阵(每行每列恰好一个 1),因此显然是酉矩阵。


1.3 相位反冲(Phase Kickback)机制

Phase Kickback 是 Deutsch-Jozsa 算法乃至许多量子算法的核心技巧。它的关键洞察是:将 oracle 的目标比特制备为 \lvert-\rangle 态,可以将函数值 f(x)f(x) “反冲”到输入寄存器的相位中

推导

将输出寄存器初始化为 1\lvert1\rangle 并施加 HH 门,得到:

y=H1=12(01)=\lvert y\rangle = H\lvert1\rangle = \frac{1}{\sqrt{2}}(\lvert0\rangle - \lvert1\rangle) = \lvert-\rangle

现在考察 UfU_f 作用在 x\lvert x\rangle\lvert-\rangle 上的效果:

Ufx=x12(0f(x)1f(x))U_f\lvert x\rangle\lvert-\rangle = \lvert x\rangle \cdot \frac{1}{\sqrt{2}}\big(\lvert0 \oplus f(x)\rangle - \lvert1 \oplus f(x)\rangle\big)

分两种情况讨论:

情况 Af(x)=0f(x) = 012(0010)=12(01)=\frac{1}{\sqrt{2}}(\lvert0 \oplus 0\rangle - \lvert1 \oplus 0\rangle) = \frac{1}{\sqrt{2}}(\lvert0\rangle - \lvert1\rangle) = \lvert-\rangle

情况 Bf(x)=1f(x) = 112(0111)=12(10)=\frac{1}{\sqrt{2}}(\lvert0 \oplus 1\rangle - \lvert1 \oplus 1\rangle) = \frac{1}{\sqrt{2}}(\lvert1\rangle - \lvert0\rangle) = -\lvert-\rangle

综合两种情况:

Ufx=(1)f(x)xU_f\lvert x\rangle\lvert-\rangle = (-1)^{f(x)}\lvert x\rangle\lvert-\rangle

这就是相位反冲:函数值 f(x)f(x)(1)f(x)(-1)^{f(x)} 的形式出现在输入寄存器 x\lvert x\rangle 的全局相位上,而输出寄存器 \lvert-\rangle 完全不变(可以丢弃)。换句话说,oracle 的作用等价于:

Ufphase:x    (1)f(x)xU_f^{\text{phase}}: \lvert x\rangle \;\longmapsto\; (-1)^{f(x)}\lvert x\rangle

此时 oracle 已经退化为一个对角酉矩阵 Ufphase=diag((1)f(00),(1)f(01),,(1)f(11))U_f^{\text{phase}} = \text{diag}((-1)^{f(0\cdots0)}, (-1)^{f(0\cdots1)}, \ldots, (-1)^{f(1\cdots1)})

物理直觉:输出寄存器 \lvert-\rangle 扮演了”相位参考”的角色。当 oracle 试图翻转 \lvert-\rangle 时(f(x)=1f(x)=1 的情况),由于 \lvert-\rangle 的两个分量 0\lvert0\rangle1\lvert1\rangle 被反相对称地翻转,整体获得 1-1 全局相位。这种相位因门的可逆性”反冲”到输入寄存器上。


1.4 Deutsch-Jozsa 算法完整描述

电路图(n 比特,标准记法)
|0⟩^⊗n —H^⊗n—•—H^⊗n—[M]
               |
|1⟩ —————H—— ⊕ ——————

其中:

  • 上层 nn 条线是输入寄存器,初始化为 0n\lvert0\rangle^{\otimes n}
  • 下层是辅助比特(输出寄存器),初始化为 1\lvert1\rangle
  • \bullet\oplus 之间的连线表示 UfU_f oracle(受控操作集合)
  • [M][M] 表示计算基测量
  • HnH^{\otimes n} 表示 nn 个并行 HH
逐步推导

步骤 0 — 初始化: ψ0=0n1\lvert\psi_0\rangle = \lvert0\rangle^{\otimes n} \otimes \lvert1\rangle

步骤 1 — 对辅助比特施加 HH 门: ψ1=0n\lvert\psi_1\rangle = \lvert0\rangle^{\otimes n} \otimes \lvert-\rangle

步骤 2 — 对所有 n+1n+1 个量子比特施加 HnIH^{\otimes n} \otimes I(即对输入寄存器做 nn 个并行 HH 门):

回顾 Hn0nH^{\otimes n}\lvert0\rangle^{\otimes n} 的作用: Hn0n=12nx=02n1xH^{\otimes n}\lvert0\rangle^{\otimes n} = \frac{1}{\sqrt{2^n}}\sum_{x=0}^{2^n-1}\lvert x\rangle

这是因为 H0=12(0+1)H\lvert0\rangle = \frac{1}{\sqrt{2}}(\lvert0\rangle + \lvert1\rangle),张量积展开后得到所有 2n2^n 个计算基矢的等幅叠加。因此:

ψ2=12nx=02n1x\lvert\psi_2\rangle = \frac{1}{\sqrt{2^n}}\sum_{x=0}^{2^n-1}\lvert x\rangle \otimes \lvert-\rangle

步骤 3 — 施加 oracle UfU_f

利用相位反冲:

ψ3=12nx=02n1(1)f(x)x\lvert\psi_3\rangle = \frac{1}{\sqrt{2^n}}\sum_{x=0}^{2^n-1}(-1)^{f(x)}\lvert x\rangle \otimes \lvert-\rangle

这是关键步骤:一次 oracle 调用同时标记了所有 2n2^n 个输入的函数值

步骤 4 — 对输入寄存器施加第二次 HnH^{\otimes n}

我们需要计算 HnxH^{\otimes n}\lvert x\rangle 的显式形式。回忆 H0=+H\lvert 0\rangle = \lvert+\rangleH1=H\lvert 1\rangle = \lvert-\rangle。对于 nn 比特,对单个基矢 x=x1x2xn\lvert x\rangle = \lvert x_1x_2\ldots x_n\rangle

Hnx=12nz=02n1(1)xzzH^{\otimes n}\lvert x\rangle = \frac{1}{\sqrt{2^n}}\sum_{z=0}^{2^n-1}(-1)^{x \cdot z}\lvert z\rangle

其中 xz=i=1nxizi(mod2)x \cdot z = \sum_{i=1}^n x_i z_i \pmod{2} 是比特内积(模 2 加法)。

这个公式的证明:Hxi=12(0+(1)xi1)H\lvert x_i\rangle = \frac{1}{\sqrt{2}}(\lvert0\rangle + (-1)^{x_i}\lvert1\rangle),因此:

Hnx1xn=i=1n12(0+(1)xi1)=12nz1,,zn{0,1}(1)ixiziz1zn=12nz=02n1(1)xzzH^{\otimes n}\lvert x_1\ldots x_n\rangle = \bigotimes_{i=1}^n\frac{1}{\sqrt{2}}(\lvert0\rangle + (-1)^{x_i}\lvert1\rangle) = \frac{1}{\sqrt{2^n}}\sum_{z_1,\ldots,z_n\in\{0,1\}}(-1)^{\sum_i x_i z_i}\lvert z_1\ldots z_n\rangle = \frac{1}{\sqrt{2^n}}\sum_{z=0}^{2^n-1}(-1)^{x \cdot z}\lvert z\rangle

现在将第二次 HnH^{\otimes n} 作用在 ψ3\lvert\psi_3\rangle 的输入寄存器上:

ψ4=12nx=02n1(1)f(x)(Hnx)\lvert\psi_4\rangle = \frac{1}{\sqrt{2^n}}\sum_{x=0}^{2^n-1}(-1)^{f(x)}\left(H^{\otimes n}\lvert x\rangle\right) \otimes \lvert-\rangle

=12nx=02n1(1)f(x)z=02n1(1)xzz= \frac{1}{2^n}\sum_{x=0}^{2^n-1}(-1)^{f(x)}\sum_{z=0}^{2^n-1}(-1)^{x\cdot z}\lvert z\rangle \otimes \lvert-\rangle

=z=02n1(12nx=02n1(1)f(x)+xz)振幅 αzz= \sum_{z=0}^{2^n-1}\underbrace{\left(\frac{1}{2^n}\sum_{x=0}^{2^n-1}(-1)^{f(x)+x\cdot z}\right)}_{\text{振幅 } \alpha_z}\lvert z\rangle \otimes \lvert-\rangle

步骤 5 — 测量输入寄存器:

测量 z\lvert z\rangle 得到结果 z{0,1}nz \in \{0,1\}^n。关键问题:我们能否从测量结果判断 ff 是常数还是平衡?


1.5 正确性证明

定理:在上述算法中,如果 ff 是常数函数,则以概率 11 测得 z=000z = 00\ldots 0;如果 ff 是平衡函数,则测得 z=000z = 00\ldots 0 的概率为 00

证明

考虑 000\lvert 00\ldots 0\rangle 的振幅 α0\alpha_0(对应 z=0z = 0):

α0=12nx=02n1(1)f(x)(1)x0=12nx=02n1(1)f(x)\alpha_0 = \frac{1}{2^n}\sum_{x=0}^{2^n-1}(-1)^{f(x)}(-1)^{x\cdot 0} = \frac{1}{2^n}\sum_{x=0}^{2^n-1}(-1)^{f(x)}

因为 x0=0x \cdot 0 = 0 对所有 xx 成立。

情况 1ff 是常数函数:

  • 若对所有 xxf(x)=0f(x) = 0,则 (1)f(x)=1(-1)^{f(x)} = 1α0=12n2n=1\alpha_0 = \frac{1}{2^n} \cdot 2^n = 1
  • 若对所有 xxf(x)=1f(x) = 1,则 (1)f(x)=1(-1)^{f(x)} = -1α0=12n(2n)=1\alpha_0 = \frac{1}{2^n} \cdot (-2^n) = -1

两种情况下 α02=1|\alpha_0|^2 = 1,因此测量以概率 1 得到 00000\ldots 0

情况 2ff 是平衡函数: 恰好一半输入的 f(x)=0f(x)=0,一半输入的 f(x)=1f(x)=1。因此:

x=02n1(1)f(x)=2n2(+1)+2n2(1)=0\sum_{x=0}^{2^n-1}(-1)^{f(x)} = \frac{2^n}{2} \cdot (+1) + \frac{2^n}{2} \cdot (-1) = 0

所以 α0=0\alpha_0 = 0,测量得到 00000\ldots 0 的概率为零。

推论:若测得 z0z \neq 0,则 ff 必为平衡函数;若测得 z=0z = 0,则 ff 必为常数函数。单次查询、确定性的判定!\blacksquare


1.6 工作示例:n=2n=2(Deutsch 算法)

n=2n=2 时,输入为 x1x2\lvert x_1x_2\ranglex{0,1,2,3}x \in \{0,1,2,3\}

示例 A:常数函数 f(x)=0f(x) = 0
xxf(x)f(x)(1)f(x)(-1)^{f(x)}
000+1
010+1
100+1
110+1

步骤 3 后的态: ψ3=12(+00+01+10+11)\lvert\psi_3\rangle = \frac{1}{2}\big(+\lvert00\rangle + \lvert01\rangle + \lvert10\rangle + \lvert11\rangle\big)\lvert-\rangle

步骤 4(第二次 H2H^{\otimes 2}):

H200=12(00+01+10+11)H^{\otimes 2}\lvert00\rangle = \frac{1}{2}(\lvert00\rangle + \lvert01\rangle + \lvert10\rangle + \lvert11\rangle) H201=12(0001+1011)H^{\otimes 2}\lvert01\rangle = \frac{1}{2}(\lvert00\rangle - \lvert01\rangle + \lvert10\rangle - \lvert11\rangle) H210=12(00+011011)H^{\otimes 2}\lvert10\rangle = \frac{1}{2}(\lvert00\rangle + \lvert01\rangle - \lvert10\rangle - \lvert11\rangle) H211=12(000110+11)H^{\otimes 2}\lvert11\rangle = \frac{1}{2}(\lvert00\rangle - \lvert01\rangle - \lvert10\rangle + \lvert11\rangle)

ψ4=12(H200+H201+H210+H211)\lvert\psi_4\rangle = \frac{1}{2}\big(H^{\otimes 2}\lvert00\rangle + H^{\otimes 2}\lvert01\rangle + H^{\otimes 2}\lvert10\rangle + H^{\otimes 2}\lvert11\rangle\big)\lvert-\rangle

合并同类项,00\lvert00\rangle 的系数为: 12(12+12+12+12)=1\frac{1}{2}\left(\frac{1}{2}+\frac{1}{2}+\frac{1}{2}+\frac{1}{2}\right) = 1

其他 z0\lvert z \neq 0\rangle 的系数为 0。因此 ψ4=00\lvert\psi_4\rangle = \lvert00\rangle\lvert-\rangle,测量必然得到 0000

示例 B:平衡函数 f(x)=x1f(x) = x_1(第一个比特的值)
xxf(x)f(x)(1)f(x)(-1)^{f(x)}
000+1
010+1
101-1
111-1

步骤 3 后的态: ψ3=12(+00+011011)\lvert\psi_3\rangle = \frac{1}{2}\big(+\lvert00\rangle + \lvert01\rangle - \lvert10\rangle - \lvert11\rangle\big)\lvert-\rangle

步骤 4: ψ4=12(H200+H201H210H211)\lvert\psi_4\rangle = \frac{1}{2}\big(H^{\otimes 2}\lvert00\rangle + H^{\otimes 2}\lvert01\rangle - H^{\otimes 2}\lvert10\rangle - H^{\otimes 2}\lvert11\rangle\big)\lvert-\rangle

00\lvert00\rangle 的系数: 12(12+121212)=0\frac{1}{2}\left(\frac{1}{2}+\frac{1}{2}-\frac{1}{2}-\frac{1}{2}\right) = 0

10\lvert10\rangle 的系数: 12(12+12+12+12)=1\frac{1}{2}\left(\frac{1}{2}+\frac{1}{2}+\frac{1}{2}+\frac{1}{2}\right) = 1

因此 ψ4=10\lvert\psi_4\rangle = \lvert10\rangle\lvert-\rangle,测量得到 1010(非零),判定 ff 为平衡函数。正确!

示例 C:平衡函数 f(x)=x1x2f(x) = x_1 \oplus x_2(XOR)
xxf(x)f(x)(1)f(x)(-1)^{f(x)}
000+1
011-1
101-1
110+1

步骤 3: ψ3=12(+000110+11)\lvert\psi_3\rangle = \frac{1}{2}\big(+\lvert00\rangle - \lvert01\rangle - \lvert10\rangle + \lvert11\rangle\big)\lvert-\rangle

步骤 4,00\lvert00\rangle 的系数: 12(121212+12)=0\frac{1}{2}\left(\frac{1}{2}-\frac{1}{2}-\frac{1}{2}+\frac{1}{2}\right) = 0

测量必然得到非零结果。注意不同平衡函数会产生不同的 zz 模式,但算法只需要检查是否全零即可。


1.7 工作示例:n=3n=3 电路

|0⟩ —H—•—H—[M]—
|0⟩ —H—•—H—[M]—
|0⟩ —H—•—H—[M]—
        |
|1⟩ —H—⊕——————

输入空间 {0,1}3\{0,1\}^3 共 8 个元素。

常数函数示例 f(x)=1f(x)=1ψ3=18x=07x\lvert\psi_3\rangle = \frac{-1}{\sqrt{8}}\sum_{x=0}^{7}\lvert x\rangle\lvert-\rangle

ψ4=000\lvert\psi_4\rangle = -\lvert000\rangle\lvert-\rangle

测得的概率分布:P(000)=1P(000) = 1P(其他)=0P(\text{其他}) = 0

平衡函数示例 f(x)=x1x2x3f(x) = x_1 \land x_2 \land x_3(三个比特的 AND,仅当 x=111x=111f=1f=1,其他情况 f=0f=0):

这个函数是平衡的吗?不,AND 函数对 8 个输入中仅 1 个输出 1,对 7 个输出 0——它既不常数也不平衡,因此不在 Deutsch-Jozsa 问题的承诺范围内。Deutsch-Jozsa 算法只对满足”常数或平衡”承诺的函数有效。

合法的平衡函数示例:f(x)=x1f(x) = x_1(仅取决于第一个比特)。这时 8 个输入中一半(x1=0x_1=0 的 4 个)输出 0,一半(x1=1x_1=1 的 4 个)输出 1。


1.8 复杂度分析

量子复杂度

  • 量子门数:O(n)O(n) 个单比特门(2n+12n+1HH 门)+ 1 次 oracle 调用
  • 电路深度:O(n)O(n)HH 门可以并行执行)
  • 总复杂度:O(n)O(n)

经典复杂度(确定性):

  • 最坏情况:2n1+12^{n-1}+1 次查询
  • 指数级差距!对 n=100n=100,经典需要  2996×1029~2^{99} \approx 6 \times 10^{29} 次查询(不可行),量子只需约 200 个门。

重要说明:Deutsch-Jozsa 算法解决的问题具有”常数或平衡”的承诺结构,且函数是 promise problem(承诺问题)而非一般的判定问题。它不是一个”实用”的算法,但它是量子算法设计思想的完美教学案例:叠加 → 并行评估 → 干涉 → 提取全局信息。


2. Bernstein-Vazirani 算法

2.1 问题定义

问题:给定一个函数 fs:{0,1}n{0,1}f_s: \{0,1\}^n \to \{0,1\},其形式为 fs(x)=sx(mod2)f_s(x) = s \cdot x \pmod{2},其中 s{0,1}ns \in \{0,1\}^n 是隐藏的比特串,sx=i=1nsixis \cdot x = \sum_{i=1}^n s_i x_i 是比特内积。求 ss

经典复杂度:逐个查询每个比特。设置 x=1000x = 100\ldots 0 得到 fs(x)=s1f_s(x) = s_1;设置 x=0100x = 010\ldots 0 得到 s2s_2;以此类推。总共需要 nn 次查询。

量子复杂度:仅需 1 次 查询。


2.2 算法电路

Bernstein-Vazirani 算法与 Deutsch-Jozsa 算法几乎完全相同,唯一区别在于测量结果的解释不同。

|0⟩^⊗n —H^⊗n—•—H^⊗n—[M] → s (直接读出!)
               |
|1⟩ —————H—— ⊕ ——————

其中 oracle 实现 Ufs:xyxy(sx)U_{f_s}: \lvert x\rangle\lvert y\rangle \to \lvert x\rangle\lvert y \oplus (s \cdot x)\rangle

2.3 完整的数学推导

步骤 0–3 与 Deutsch-Jozsa 完全相同。oracle 调用后的态为:

ψ3=12nx=02n1(1)sxx\lvert\psi_3\rangle = \frac{1}{\sqrt{2^n}}\sum_{x=0}^{2^n-1}(-1)^{s \cdot x}\lvert x\rangle\lvert-\rangle

步骤 4 — 施加第二次 HnH^{\otimes n}

回忆 Hnx=12nz=02n1(1)xzzH^{\otimes n}\lvert x\rangle = \frac{1}{\sqrt{2^n}}\sum_{z=0}^{2^n-1}(-1)^{x \cdot z}\lvert z\rangle。因此:

ψ4=12nx=02n1(1)sx12nz=02n1(1)xzz\lvert\psi_4\rangle = \frac{1}{\sqrt{2^n}}\sum_{x=0}^{2^n-1}(-1)^{s\cdot x}\frac{1}{\sqrt{2^n}}\sum_{z=0}^{2^n-1}(-1)^{x \cdot z}\lvert z\rangle\lvert-\rangle

=12nz=02n1(x=02n1(1)(sz)x)z= \frac{1}{2^n}\sum_{z=0}^{2^n-1}\left(\sum_{x=0}^{2^n-1}(-1)^{(s\oplus z)\cdot x}\right)\lvert z\rangle\lvert-\rangle

其中 szs \oplus z 是逐比特 XOR。注意求和 x(1)tx\sum_{x}(-1)^{t\cdot x} 当且仅当 t=0t = 0 时值为 2n2^n(因为 (1)0+(1)0+=2n(-1)^0 + (-1)^0 + \cdots = 2^n),否则为 00

因此:

ψ4=s\lvert\psi_4\rangle = \lvert s\rangle\lvert-\rangle

步骤 5 — 测量输入寄存器,直接得到 ss 的每一位!

直观解释:量子并行性一次评估了所有输入 xx,通过干涉构造了隐藏比特串 ss 的傅里叶变换形式,第二次 HnH^{\otimes n} 执行了逆傅里叶变换,将 ss 的信息聚焦到单一量子态上。

2.4 工作示例

n=3n=3s=101s = 101(即 f(x)=x1+x3(mod2)f(x) = x_1 + x_3 \pmod{2})。

oracle 调用后的态: ψ3=18x{0,1}3(1)x1+x3x\lvert\psi_3\rangle = \frac{1}{\sqrt{8}}\sum_{x\in\{0,1\}^3}(-1)^{x_1 + x_3}\lvert x\rangle\lvert-\rangle

具体展开(仅列出 (1)f(x)=1(-1)^{f(x)} = -1 的项,即 x1+x3x_1+x_3 为奇数的输入):

xx(1)f(x)(-1)^{f(x)}
001-1
010+1
011-1
100-1
101+1
110-1
111+1

第二次 H3H^{\otimes 3} 后,所有振幅除了 101\lvert101\rangle 外都干涉相消,测量得到 101101


2.5 Bernstein-Vazirani 与 Deutsch-Jozsa 的关系

  • BV 是 DJ 的”参数化”版本:DJ 判断 ff 是常数还是平衡,BV 找出 ff 的隐藏参数 ss
  • BV 的 oracle 结构更具体(线性函数而不是任意常数/平衡函数)
  • BV 展示了量子傅里叶采样的思想:通过傅里叶变换将隐藏结构转化为可测量的尖峰
  • 两种算法共享同一个电路,但 BV 提供了更强的结果(找出 ss,而不仅仅是分类)

复杂度对比

算法经典查询量子查询加速比
Deutsch-Jozsa2n1+12^{n-1}+11指数级
Bernstein-Vaziraninn1nn

BV 的加速是”线性”的(nn 倍),而不是指数的。但作为量子算法设计模板,它与 DJ 一样具有重要的教学意义。


3. Simon 算法

3.1 问题定义

问题:给定一个函数 f:{0,1}n{0,1}nf: \{0,1\}^n \to \{0,1\}^{n}(输出也是 nn 比特),承诺存在一个隐藏的非零比特串 s{0,1}ns \in \{0,1\}^n,使得对于所有 x,y{0,1}nx, y \in \{0,1\}^n

f(x)=f(y)    y=xsf(x) = f(y) \iff y = x \oplus s

ff二对一的,且冲突对恰好相差 ss。求 ss

几何理解:输入空间 {0,1}n\{0,1\}^n 被划分为 2n12^{n-1} 个”配对”{x,xs}\{x, x\oplus s\},每个配对的函数值相同。寻找决定这种配对结构的隐藏周期 ss

经典复杂度:最坏情况下需要 Θ(2n/2)\Theta(2^{n/2}) 次查询(由生日悖论给出)。

量子复杂度O(n)O(n) 次查询,指数级加速!


3.2 Simon 算法的电路

Simon 算法的核心思想与 DJ/BV 有深刻的相似性:叠加 → oracle → Hadamard → 测量。但 Simon 需要多次运行来收集线性方程。

单次运行电路

|0⟩^⊗n —H^⊗n—•—H^⊗n—[M] → 随机 z
               |
|0⟩^⊗n ——————⊕—————      → 丢弃(测量用于验证)

其中 oracle Uf:xyxyf(x)U_f: \lvert x\rangle\lvert y\rangle \to \lvert x\rangle\lvert y \oplus f(x)\rangle

3.3 逐步推导

步骤 0ψ0=0n0n\lvert\psi_0\rangle = \lvert0\rangle^{\otimes n} \otimes \lvert0\rangle^{\otimes n}

步骤 1 — 对第一寄存器施加 HnH^{\otimes n}ψ1=12nx=02n1x0n\lvert\psi_1\rangle = \frac{1}{\sqrt{2^n}}\sum_{x=0}^{2^n-1}\lvert x\rangle \otimes \lvert0\rangle^{\otimes n}

步骤 2 — 施加 oracle: ψ2=12nx=02n1xf(x)\lvert\psi_2\rangle = \frac{1}{\sqrt{2^n}}\sum_{x=0}^{2^n-1}\lvert x\rangle\lvert f(x)\rangle

步骤 3 — 对第一寄存器施加第二次 HnH^{\otimes n}

ψ3=12nx=02n1(12nz=02n1(1)xzz)f(x)\lvert\psi_3\rangle = \frac{1}{\sqrt{2^n}}\sum_{x=0}^{2^n-1}\left(\frac{1}{\sqrt{2^n}}\sum_{z=0}^{2^n-1}(-1)^{x\cdot z}\lvert z\rangle\right)\lvert f(x)\rangle

=12nz=02n1x=02n1(1)xzzf(x)= \frac{1}{2^n}\sum_{z=0}^{2^n-1}\sum_{x=0}^{2^n-1}(-1)^{x\cdot z}\lvert z\rangle\lvert f(x)\rangle

步骤 4 — 测量第一寄存器。测得特定 zz 的概率为:

P(z)=12nx=02n1(1)xzf(x)2P(z) = \left\|\frac{1}{2^n}\sum_{x=0}^{2^n-1}(-1)^{x\cdot z}\lvert f(x)\rangle\right\|^2

由于 ff 是二对一的,每个 f(x)f(x) 对应两个 xx 值:x0x_0x0sx_0 \oplus s。因此:

P(z)=12n每对{x,xs}(1)xzf(x)2P(z) = \left\|\frac{1}{2^n}\sum_{\text{每对}\{x,x\oplus s\}}(-1)^{x\cdot z}\lvert f(x)\rangle\right\|^2

对于任意”配对”{x,xs}\{x, x\oplus s\},其贡献为:

12n[(1)xz+(1)(xs)z]f(x)=12n(1)xz[1+(1)sz]f(x)\frac{1}{2^n}\left[(-1)^{x\cdot z} + (-1)^{(x\oplus s)\cdot z}\right]\lvert f(x)\rangle = \frac{1}{2^n}(-1)^{x\cdot z}\left[1 + (-1)^{s\cdot z}\right]\lvert f(x)\rangle

关键

  • sz=1s\cdot z = 1,则 1+(1)=01 + (-1) = 0,该配对的贡献为零
  • sz=0s\cdot z = 0,则 1+1=21 + 1 = 2,该配对的贡献非零

因此,只有当 sz=0(mod2)s \cdot z = 0 \pmod{2} 时,P(z)>0P(z) > 0;否则 P(z)=0P(z) = 0

结论:每次运行 Simon 算法,测得一个随机的 zz 满足 sz=0s \cdot z = 0。特别地,P(z=0)=2n+1P(z=0) = 2^{-n+1}P(z=s)=2n+1P(z=s) = 2^{-n+1}

3.4 从测量结果恢复 ss

单次运行给出一个 zz 满足 sz=0s \cdot z = 0。这是关于 ss 的一个线性方程。运行 O(n)O(n) 次,收集 n1n-1线性独立的 zz 值:

z1s=0z_1 \cdot s = 0 z2s=0z_2 \cdot s = 0 \vdots zn1s=0z_{n-1} \cdot s = 0

这构成一个 GF(2)GF(2) 上的线性方程组。解出非零解 ss(以及平凡解 s=0s=0)。由于 s0s \neq 0 是承诺条件,我们取非零解。

所需运行次数:在 GF(2)nGF(2)^n 中,随机均匀选取 mmnn 维向量,它们张成 (n1)(n-1) 维子空间(排除 ss 的正交补)的概率为:

P(成功)=k=1n1(12k1m)P(\text{成功}) = \prod_{k=1}^{n-1}(1-2^{k-1-m})

m=n+O(1)m = n + O(1) 时成功概率趋近于 1。经典后处理在 O(n3)O(n^3) 内完成(高斯消元)。

3.5 工作示例:n=2n=2s=11s = 11

函数定义(承诺示例):

xxf(x)f(x)
0001
0110
1010
1101

验证:f(00)=f(11)f(00)=f(11)f(01)=f(10)f(01)=f(10),且 0011=1100\oplus 11 = 110110=1101\oplus 10 = 11。所以 s=11s=11

运行 1:假设测得 z=11z = 11。方程:11s=s1+s2=0(mod2)11 \cdot s = s_1 + s_2 = 0 \pmod{2}

运行 2:假设测得 z=10z = 10。方程:10s=s1=0(mod2)10 \cdot s = s_1 = 0 \pmod{2}

s1=0s_1 = 0s1+s2=0s_1 + s_2 = 0s2=0s_2 = 0。但这给出 s=00s=00,与承诺矛盾(s0s\neq 0)。怎么回事?我们收集的方程不独立!实际上第二个方程 s1=0s_1=0 已经推出 s1=0s_1=0,再由第一个方程得 s2=0s_2=0,但我们知道 s=11s=11

这揭示了一个关键点:Simon 算法的测量结果 zz 的分布是均匀的(在满足 sz=0s\cdot z=0zz 中)。可能恰好 zz 都属于 {00,11}\{00, 11\} 这个 1 维子空间——这种情况我们只能得到 s1+s2=0s_1 + s_2 = 0 这一个方程,需要更多运行。

运行 3:测得 z=01z = 01。方程:01s=s2=001 \cdot s = s_2 = 0

现在 s1=0,s2=0s_1=0, s_2=0 还是 0000。仍然矛盾!这说明测量的 zz 全部来自 ss 的正交补空间 {00}\{00\}——如果每次测量都得到 0000(概率虽小但非零),我们需要更多运行。

实际上,我们需要至少 n1=1n-1=1 个独立的非零方程来约束 ssss 的正交补空间是 n1n-1 维的,测量均匀地从中采样。找到 n1n-1 个线性独立向量的期望运行次数为 O(n)O(n)。一旦我们找到非零 ss,就完成了。

3.6 复杂度分析

量子复杂度

  • 每次运行:O(n)O(n) 个门 + 1 次 oracle 调用
  • 运行次数:O(n)O(n) 次(期望)
  • 总 oracle 调用:O(n)O(n)

经典复杂度

  • 确定性:Θ(2n/2)\Theta(2^{n/2}) 次函数求值(通过生日悖论找到冲突对)
  • 随机化:同样需要 Ω(2n/2)\Omega(2^{n/2})

指数加速: Simon 算法是第一个展示指数级量子加速的算法(比 Deutsch-Jozsa 的”承诺指数”更令人信服,因为问题更自然)。它直接启发了 Shor 算法——Simon 的”隐藏子群”框架直接推广到有限阿贝尔群,Shor 的周期发现就是 Simon 在群 Z\mathbb{Z} 上的推广。

3.7 Simon 问题与隐藏子群问题

Simon 问题是隐藏子群问题 (Hidden Subgroup Problem, HSP) 在群 (Z2)n(\mathbb{Z}_2)^n 上的一个实例:

  • G=(Z2)nG = (\mathbb{Z}_2)^n
  • 隐藏子群 H={0,s}H = \{0, s\}(二阶子群)
  • 函数 ffHH 的陪集上为常数
  • 寻找 HH 的生成元 ss

Shor 因子分解算法中的周期发现问题对应 HSP 在群 Z\mathbb{Z} 上的实例。这一统一框架揭示了 Simon 与 Shor 之间的深刻联系。


4. Grover 搜索算法

4.1 问题定义

问题:在 N=2nN = 2^n 个无结构数据项中,找到满足某个条件的”标记”项。假设存在唯一的标记项 xx^*,可以通过 oracle 查询来检验。

形式化:存在 oracle 函数 f:{0,1}n{0,1}f: \{0,1\}^n \to \{0,1\},其中 f(x)=1f(x) = 1 当且仅当 x=xx = x^*(标记项),否则 f(x)=0f(x) = 0。通过查询 ff 找到 xx^*

经典复杂度:顺序搜索需要平均 N/2N/2 次查询,最坏 NN 次。

量子复杂度O(N)O(\sqrt{N}) 次查询(约 π4N\frac{\pi}{4}\sqrt{N} 次迭代),二次加速。


4.2 Oracle 构造

与 Deutsch-Jozsa 类似,Grover oracle 通过相位反冲实现:

Ufx=(1)f(x)xU_f\lvert x\rangle\lvert-\rangle = (-1)^{f(x)}\lvert x\rangle\lvert-\rangle

对于标记项 xx^*f(x)=1f(x^*) = 1,因此施加 1-1 相位;对其他 xx,相位保持不变。

在电路层面,可以将 Grover oracle 写作:

Uoracle=I2xxU_{\text{oracle}} = I - 2\lvert x^*\rangle\langle x^*\rvert

即对标记态 x\lvert x^*\rangle 施加相位翻转 (1)(-1),而保持所有其他态不变。

验证:(I2xx)x=x2x=x(I - 2\lvert x^*\rangle\langle x^*\rvert)\lvert x^*\rangle = \lvert x^*\rangle - 2\lvert x^*\rangle = -\lvert x^*\rangle;而对 xx\lvert x \neq x^*\rangle(I2xx)x=x(I - 2\lvert x^*\rangle\langle x^*\rvert)\lvert x\rangle = \lvert x\rangle


4.3 扩散算子(Diffusion Operator)的数学推导

Grover 算法的第二个关键组件是扩散算子(也称”反演关于均值的操作”,inversion about the mean):

Udiff=2ψψIU_{\text{diff}} = 2\lvert \psi\rangle\langle\psi\rvert - I

其中 ψ=Hn0n=1Nx=0N1x\lvert\psi\rangle = H^{\otimes n}\lvert0\rangle^{\otimes n} = \frac{1}{\sqrt{N}}\sum_{x=0}^{N-1}\lvert x\rangle 是所有基矢的等幅叠加。

为什么叫”反演关于均值”?

axa_xx\lvert x\rangle 的振幅,UdiffU_{\text{diff}} 作用于态 xaxx\sum_x a_x\lvert x\rangle 后:

Udiffxaxx=(2ψψI)xaxxU_{\text{diff}}\sum_x a_x\lvert x\rangle = (2\lvert\psi\rangle\langle\psi\rvert - I)\sum_x a_x\lvert x\rangle

=2ψxaxψxxaxx= 2\lvert\psi\rangle\sum_x a_x\langle\psi\rvert x\rangle - \sum_x a_x\lvert x\rangle

其中 ψx=1N\langle\psi\rvert x\rangle = \frac{1}{\sqrt{N}},所以 xaxψx=1Nxax=Naˉ\sum_x a_x\langle\psi\rvert x\rangle = \frac{1}{\sqrt{N}}\sum_x a_x = \sqrt{N}\cdot\bar{a},其中 aˉ=1Nxax\bar{a} = \frac{1}{N}\sum_x a_x 是平均振幅。

因此:

Udiffxaxx=2ψNaˉxaxxU_{\text{diff}}\sum_x a_x\lvert x\rangle = 2\lvert\psi\rangle\sqrt{N}\bar{a} - \sum_x a_x\lvert x\rangle

=2xaˉxxaxx= 2\sum_x \bar{a}\lvert x\rangle - \sum_x a_x\lvert x\rangle

=x(2aˉax)x= \sum_x (2\bar{a} - a_x)\lvert x\rangle

这正是”反演关于均值”:每个振幅 axa_x 被替换为 2aˉax2\bar{a} - a_x。如果 axa_x 低于均值,它会被”抬升”;如果 axa_x 高于均值,会被”压低”。这种操作的效果是:放大高于均值的振幅,缩小低于均值的振幅

电路实现Udiff=Hn(200I)HnU_{\text{diff}} = H^{\otimes n} (2\lvert0\rangle\langle0\rvert - I) H^{\otimes n},其中 200I2\lvert0\rangle\langle0\rvert - I 是”零态相位翻转”:对 0n\lvert0\rangle^{\otimes n} 施加 1-1 相位,其他不变。

完整扩散算子的电路

|x⟩ —H^⊗n—(·)—H^⊗n—
           |
           X—•—X
             |
            X—•—X
               |
              ...
             X—•—X    (n 个 X 门对)
               |
           H—⊗n—H

中间的 (200I)(2\lvert0\rangle\langle0\rvert - I) 可以进一步分解为:将所有比特翻转(施加 XnX^{\otimes n}),做多比特受控 ZZ 门(111\lvert11\ldots1\rangle 翻转相位),再翻转回来。


4.4 几何解释:二维旋转

Grover 算法最优雅的数学解释是将整个 2n2^n 维空间约化为一个二维子空间。这是理解 Grover 算法为何如此高效的关键。

定义两个正交态:

α=1N1xxx(非标记项的均匀叠加)\lvert\alpha\rangle = \frac{1}{\sqrt{N-1}}\sum_{x \neq x^*}\lvert x\rangle \quad \text{(非标记项的均匀叠加)}

β=x(标记项)\lvert\beta\rangle = \lvert x^*\rangle \quad \text{(标记项)}

注意 αβ=0\langle\alpha\rvert\beta\rangle = 0,且初始态 ψ=Hn0\lvert\psi\rangle = H^{\otimes n}\lvert0\rangle 可以写为:

ψ=1Nβ+N1Nα=sinθβ+cosθα\lvert\psi\rangle = \frac{1}{\sqrt{N}}\lvert\beta\rangle + \frac{\sqrt{N-1}}{\sqrt{N}}\lvert\alpha\rangle = \sin\theta\lvert\beta\rangle + \cos\theta\lvert\alpha\rangle

其中 θ=arcsin(1/N)\theta = \arcsin(1/\sqrt{N})。对于大 NNθ1/N\theta \approx 1/\sqrt{N}

Grover 迭代由两步组成:

  1. Oracle: Uoracle=I2ββU_{\text{oracle}} = I - 2\lvert\beta\rangle\langle\beta\rvert — 关于 β\lvert\beta\rangle 的反射
  2. Diffusion: Udiff=2ψψIU_{\text{diff}} = 2\lvert\psi\rangle\langle\psi\rvert - I — 关于 ψ\lvert\psi\rangle 的反射

Grover 迭代 G=UdiffUoracleG = U_{\text{diff}} \cdot U_{\text{oracle}} 是一个旋转: 在 {α,β}\{\lvert\alpha\rangle, \lvert\beta\rangle\} 平面中,GG 将态旋转 2θ2\theta 弧度:

Gkψ=sin((2k+1)θ)β+cos((2k+1)θ)αG^k\lvert\psi\rangle = \sin((2k+1)\theta)\lvert\beta\rangle + \cos((2k+1)\theta)\lvert\alpha\rangle

证明(归纳法或直接几何论证):

UoracleU_{\text{oracle}} 是关于 β\lvert\beta\rangle 轴的反射(保持 β\lvert\beta\rangle 不变,翻转其垂直分量)。 UdiffU_{\text{diff}} 是关于 ψ\lvert\psi\rangle 轴的反射。 两次反射的组合是一个旋转,旋转角为两反射轴夹角的 2 倍。

初始态 ψ\lvert\psi\rangleα\lvert\alpha\rangle 的夹角为 θ\theta(因为 ψα=cosθ\langle\psi\rvert\alpha\rangle = \cos\theta)。 UoracleU_{\text{oracle}}ψ\lvert\psi\rangle 反射为与 β\lvert\beta\rangle 夹角为 θ\theta 的对称点。 UdiffU_{\text{diff}} 再将结果反射,相对于 ψ\lvert\psi\rangle 轴。 两次反射的净效果:旋转 2θ2\theta 弧度。

因此 kk 次迭代后,态的角度为 (2k+1)θ(2k+1)\theta


4.5 最优迭代次数的推导

我们希望在经过 kk 次迭代后,β\lvert\beta\rangle(标记项)的振幅尽可能大。

βGkψ=sin((2k+1)θ)\langle\beta\rvert G^k\lvert\psi\rangle = \sin((2k+1)\theta)

最大化该振幅的条件:

sin((2k+1)θ)1    (2k+1)θπ2\sin((2k+1)\theta) \approx 1 \;\Longrightarrow\; (2k+1)\theta \approx \frac{\pi}{2}

因此最优迭代次数:

kopt=π4θ12k_{\text{opt}} = \left\lfloor\frac{\pi}{4\theta} - \frac{1}{2}\right\rfloor

由于 θ=arcsin(1/N)1/N\theta = \arcsin(1/\sqrt{N}) \approx 1/\sqrt{N}(对 N1N \gg 1):

koptπ4Nk_{\text{opt}} \approx \left\lfloor\frac{\pi}{4}\sqrt{N}\right\rfloor

最终振幅

βGkoptψ2=sin2((2kopt+1)arcsin1N)|\langle\beta\rvert G^{k_{\text{opt}}}\lvert\psi\rangle|^2 = \sin^2\left((2k_{\text{opt}}+1)\arcsin\frac{1}{\sqrt{N}}\right)

k=koptk = k_{\text{opt}} 时,这个值接近 1。更精确地,最小失败概率:

Pfail=cos2((2kopt+1)arcsin1N)1NP_{\text{fail}} = \cos^2\left((2k_{\text{opt}}+1)\arcsin\frac{1}{\sqrt{N}}\right) \leq \frac{1}{N}

示例N=4N=4θ=arcsin(1/2)=π/6\theta = \arcsin(1/2) = \pi/6kopt=π/(4π/6)1/2=1.50.5=1k_{\text{opt}} = \lfloor \pi/(4\cdot\pi/6) - 1/2 \rfloor = \lfloor 1.5 - 0.5 \rfloor = 1。经过 1 次迭代,sin(3θ)=sin(π/2)=1\sin(3\theta) = \sin(\pi/2) = 1,以概率 1 找到标记项。

示例N=8N=8θ=arcsin(1/8)0.3614\theta = \arcsin(1/\sqrt{8}) \approx 0.3614 rad。kopt=π/(40.3614)0.5=2.170.5=1k_{\text{opt}} = \lfloor \pi/(4\cdot 0.3614) - 0.5 \rfloor = \lfloor 2.17 - 0.5 \rfloor = 1。经过 1 次迭代,sin(3θ)=sin(1.0842)0.882\sin(3\theta) = \sin(1.0842) \approx 0.882,成功概率约 0.7780.778。如果做 2 次迭代,sin(5θ)=sin(1.807)0.971\sin(5\theta) = \sin(1.807) \approx 0.971,概率约 0.9430.943


4.6 完整工作示例:N=4N=4

N=4N=4n=2n=2 量子比特),标记项设为 x=10x^* = 10(二进制,即十进制 2)。

初始化ψ0=00\lvert\psi_0\rangle = \lvert00\rangle

第一次 H2H^{\otimes 2}ψ1=12(00+01+10+11)\lvert\psi_1\rangle = \frac{1}{2}(\lvert00\rangle + \lvert01\rangle + \lvert10\rangle + \lvert11\rangle)

Grover 迭代 1(也是唯一需要的迭代):

Step A — Oracle(标记 1010): Uoracle=I21010U_{\text{oracle}} = I - 2\lvert10\rangle\langle10\rvert

ψ2=12(00+0110+11)\lvert\psi_2\rangle = \frac{1}{2}(\lvert00\rangle + \lvert01\rangle - \lvert10\rangle + \lvert11\rangle)

Step B — Diffusion:

先计算均值 aˉ=14(1+11+1)=12\bar{a} = \frac{1}{4}(1+1-1+1) = \frac{1}{2}

反演关于均值:ax=2aˉax=1axa_x' = 2\bar{a} - a_x = 1 - a_x

xxaxa_x(Oracle后)axa_x'(Diffusion后)
001/21/211/2=1/21 - 1/2 = 1/2
011/21/211/2=1/21 - 1/2 = 1/2
101/2-1/21(1/2)=3/21 - (-1/2) = 3/2
111/21/211/2=1/21 - 1/2 = 1/2

ψ3=12(00+01+310+11)\lvert\psi_3\rangle = \frac{1}{2}(\lvert00\rangle + \lvert01\rangle + 3\lvert10\rangle + \lvert11\rangle)

归一化后(注意 14+14+94+14=3\frac{1}{4} + \frac{1}{4} + \frac{9}{4} + \frac{1}{4} = 3):

ψ3=112(00+01+310+11)\lvert\psi_3\rangle = \frac{1}{\sqrt{12}}(\lvert00\rangle + \lvert01\rangle + 3\lvert10\rangle + \lvert11\rangle)

测量P(x=10)=(312)2=912=34P(x^* = 10) = \left(\frac{3}{\sqrt{12}}\right)^2 = \frac{9}{12} = \frac{3}{4}

75%75\% 概率找到标记项。注意对于 N=4N=4,单次迭代的成功概率已经很高。如果再做第二次迭代(不必要),成功概率还会上升。

事实上,N=4N=4 是最特殊的例子,因为 θ=π/6\theta = \pi/6(21+1)π/6=π/2(2\cdot 1 + 1)\cdot\pi/6 = \pi/2,完美旋转到 β\lvert\beta\rangle,成功概率应该是 100%。上面的 75%75\% 差异来自我使用了 axa_x' 的表达式,但让我们更精确地计算 UdiffU_{\text{diff}}

Udiff=H2(20000I)H2U_{\text{diff}} = H^{\otimes 2}(2\lvert00\rangle\langle00\rvert - I)H^{\otimes 2}

让我们一步步算:

H2ψ2H^{\otimes 2}\lvert\psi_2\rangle:先计算 H2H^{\otimes 2}ψ2\lvert\psi_2\rangle 的作用。

ψ2=12(1111)\lvert\psi_2\rangle = \frac{1}{2}\begin{pmatrix}1\\1\\-1\\1\end{pmatrix}

H2=12(1111111111111111)H^{\otimes 2} = \frac{1}{2}\begin{pmatrix}1&1&1&1\\1&-1&1&-1\\1&1&-1&-1\\1&-1&-1&1\end{pmatrix}

H2ψ2=14(1+11+111111+1+1111+1+1)=14(2222)H^{\otimes 2}\lvert\psi_2\rangle = \frac{1}{4}\begin{pmatrix}1+1-1+1\\1-1-1-1\\1+1+1-1\\1-1+1+1\end{pmatrix} = \frac{1}{4}\begin{pmatrix}2\\-2\\2\\2\end{pmatrix}

应用 (20000I)(2\lvert00\rangle\langle00\rvert - I):将 00\lvert00\rangle 分量的相位翻转(21=12-1=1 变成 22424=2424=02\cdot\frac{2}{4} - \frac{2}{4} = \frac{2}{4} - \frac{2}{4} = 0,等等——实际上 (20000I)(2\lvert00\rangle\langle00\rvert - I)00\lvert00\rangle 的系数从 2/42/4 变成 2/4-2/4,不对!)

更直接的方法:20000I=diag(1,1,1,1)2\lvert00\rangle\langle00\rvert - I = \text{diag}(1, -1, -1, -1)

所以:14(2222)20000I14(2222)\frac{1}{4}\begin{pmatrix}2\\-2\\2\\2\end{pmatrix} \xrightarrow{2\lvert00\rangle\langle00\rvert - I} \frac{1}{4}\begin{pmatrix}2\\2\\-2\\-2\end{pmatrix}

最后施加 H2H^{\otimes 2}

ψ3=18(1111111111111111)(2222)=18(2+222222+22+2+2+222+22)=(0010)\lvert\psi_3\rangle = \frac{1}{8}\begin{pmatrix}1&1&1&1\\1&-1&1&-1\\1&1&-1&-1\\1&-1&-1&1\end{pmatrix}\begin{pmatrix}2\\2\\-2\\-2\end{pmatrix} = \frac{1}{8}\begin{pmatrix}2+2-2-2\\2-2-2+2\\2+2+2+2\\2-2+2-2\end{pmatrix} = \begin{pmatrix}0\\0\\1\\0\end{pmatrix}

完美!ψ3=10\lvert\psi_3\rangle = \lvert10\rangle,概率 1 找到标记项。N=4N=4 只需要 1 次 Grover 迭代就能以概率 1 找到标记项。


4.7 工作示例:N=8N=8n=3n=3

标记项设为 x=110x^* = 110(十进制 6)。

初始化ψ0=18x=07x\lvert\psi_0\rangle = \frac{1}{\sqrt{8}}\sum_{x=0}^{7}\lvert x\rangle

第 1 次迭代: Oracle 将 110\lvert110\rangle 振幅从 1/81/\sqrt{8} 翻转为 1/8-1/\sqrt{8}。 均值 aˉ=(71+(1))/88=6/(88)=3/(48)\bar{a} = (7\cdot 1 + (-1))/8\sqrt{8} = 6/(8\sqrt{8}) = 3/(4\sqrt{8})。 反演后,110\lvert110\rangle 振幅变为 23/(48)(1/8)=(6/4+1)/8=(5/2)/8=5/(28)2\cdot 3/(4\sqrt{8}) - (-1/\sqrt{8}) = (6/4 + 1)/\sqrt{8} = (5/2)/\sqrt{8} = 5/(2\sqrt{8})。 其他项振幅变为 23/(48)1/8=(3/21)/8=1/(28)2\cdot 3/(4\sqrt{8}) - 1/\sqrt{8} = (3/2 - 1)/\sqrt{8} = 1/(2\sqrt{8})。 成功概率 P=(5/(28))2=25/320.781P = (5/(2\sqrt{8}))^2 = 25/32 \approx 0.781

第 2 次迭代: 从 ψ=128(1,1,1,1,1,1,5,1)T\lvert\psi\rangle = \frac{1}{2\sqrt{8}}(1,1,1,1,1,1,5,1)^T 开始(归一化前的振幅向量)。 Oracle:将 110\lvert110\rangle 振幅翻转为 5-5(保持其他不变)。 新均值 aˉ=(71+(5))/(828)=2/(168)=1/(88)\bar{a}' = (7\cdot 1 + (-5))/(8\cdot 2\sqrt{8}) = 2/(16\sqrt{8}) = 1/(8\sqrt{8})。 反演后,110\lvert110\rangle 振幅变为 2aˉ(5/(28))=1/(48)+5/(28)=11/(48)2\bar{a}' - (-5/(2\sqrt{8})) = 1/(4\sqrt{8}) + 5/(2\sqrt{8}) = 11/(4\sqrt{8})。 其他项振幅变为 2aˉ1/(28)=1/(48)1/(28)=1/(48)2\bar{a}' - 1/(2\sqrt{8}) = 1/(4\sqrt{8}) - 1/(2\sqrt{8}) = -1/(4\sqrt{8})。 成功概率 P=(11/(48))2=121/1280.945P = (11/(4\sqrt{8}))^2 = 121/128 \approx 0.945

因此 N=8N=8 做 2 次迭代即可达到约 94.5% 的成功概率。与理论值 kopt=π8/4=2.22=2k_{\text{opt}} = \lfloor \pi\sqrt{8}/4 \rfloor = \lfloor 2.22 \rfloor = 2 一致。


4.8 最优性证明:BBBV 定理

问题:是否存在比 O(N)O(\sqrt{N}) 更快的量子搜索算法?

BBBV 定理(Bennett, Bernstein, Brassard & Vazirani, 1997):任何解决无结构搜索问题的量子算法必须调用 oracle 至少 Ω(N)\Omega(\sqrt{N}) 次。

证明思路(高维几何论证):

考虑一个量子算法的状态演化过程。初始态为 ψ0\lvert\psi_0\rangle,经过 TT 次 oracle 调用 UfU_fT+1T+1 次非 oracle 酉变换 U0,U1,,UTU_0, U_1, \ldots, U_T

ψT=UTUfUT1UfU1UfU0ψ0\lvert\psi_T\rangle = U_T U_f U_{T-1} U_f \cdots U_1 U_f U_0\lvert\psi_0\rangle

定义”查询状态”序列:ψt\lvert\psi^t\rangle 为第 tt 次 oracle 调用后的状态。关键想法是跟踪标记项的振幅如何随查询次数增长。

引入”无标记”oracle U0U_0(对所有 xx 都不翻转相位),定义 ϕt\lvert\phi^t\rangle 为使用 U0U_0 代替 UfU_f 时的状态序列。引理:ψtϕt12t2/N|\langle\psi^t\rvert\phi^t\rangle| \geq 1 - 2t^2/N 的某种形式——即很难区分搜索空间是否包含标记项。

更严格的论证如下:

定义 ψk\lvert\psi_k\ranglekk 次 oracle 查询后的态。定义 oracle 操作 O=I2xxO = I - 2\lvert x^*\rangle\langle x^*\rvert 和参考 oracle O0=IO_0 = I(无标记)。

考虑差异向量 Δk=ψkϕk\Delta_k = \lVert\lvert\psi_k\rangle - \lvert\phi_k\rangle\rVert,其中 ϕk\lvert\phi_k\rangle 使用 O0O_0

可以证明每次 oracle 调用最多能将 Δk\Delta_k 增加 2/N2/\sqrt{N}(因为每次 oracle 调用最多影响 O(1/N)O(1/\sqrt{N}) 的振幅)。经过 TT 次查询后:

ΔT2TN\Delta_T \leq \frac{2T}{\sqrt{N}}

另一方面,要成功区分有标记和无标记的情况(即找到标记项),需要 ΔT=Ω(1)\Delta_T = \Omega(1)。因此 T=Ω(N)T = \Omega(\sqrt{N})

直观理解:每次 oracle 查询只能”小幅”改变量子态,需要 N\sqrt{N} 次查询才能积累足够的改变以可靠地定位标记项。这就是 Grover 算法二次加速的最优性证明。


4.9 多标记项的推广

如果存在 MM 个标记项(而非 1 个),Grover 算法仍然有效。

定义 N=N/MN' = N/M,旋转角 θ=arcsin(M/N)\theta' = \arcsin(\sqrt{M/N})。最优迭代次数为:

koptπ4NM12k'_{\text{opt}} \approx \frac{\pi}{4}\sqrt{\frac{N}{M}} - \frac{1}{2}

成功概率接近 1。当 M=N/4M = N/4 时,kopt=1k'_{\text{opt}} = 1,一次迭代即可找到某个标记项。

特殊情形

  • M>N/2M > N/2,可以随机猜测的经典方法已经很快,量子加速减弱
  • MM 未知,需要量子计数 (Quantum Counting) 先估计 MM,再运行 Grover
  • 量子计数本身是 Grover 迭代和 QPE 的结合

5. Shor 因子分解算法

5.1 问题定义与经典复杂度

问题:给定一个 LL 比特的合数 N=pqN = pqp,qp,q 为素数),求 ppqq

经典复杂度(已知最优算法):

  • 一般数域筛法 (GNFS):O(exp((logN)1/3(loglogN)2/3))O\left(\exp\left((\log N)^{1/3}(\log\log N)^{2/3}\right)\right),亚指数但超多项式
  • 对大 NN(如 RSA-2048,L=2048L=2048),经典算法完全不可行

量子复杂度O((logN)3)O((\log N)^3),多项式时间!这是 RSA 加密安全性的根本威胁。


5.2 从因子分解到周期发现:详细归约

Shor 算法的核心洞察是将因子分解转化为周期发现问题。归约分为以下几步:

步骤 1:琐碎情况排除

如果 NN 是偶数,直接输出因子 2。如果 N=abN = a^b 对某些 a1,b2a \geq 1, b \geq 2,直接分解。这些情况可以在多项式时间内判定。

步骤 2:随机选择 aa

随机选择 a{2,3,,N1}a \in \{2, 3, \ldots, N-1\}。计算 gcd(a,N)\gcd(a, N)。若 gcd(a,N)>1\gcd(a, N) > 1,则已找到因子。否则 aaNN 互素。

步骤 3:周期发现问题

考虑模指数函数:

fa,N(x)=axmodNf_{a,N}(x) = a^x \bmod N

这个函数是周期性的,因为模运算的有限性保证存在最小的 r>0r > 0 使得 ar1(modN)a^r \equiv 1 \pmod{N}(费马小定理的推广,rraa 在乘法群 ZN\mathbb{Z}_N^* 中的阶)。周期 rr 就是函数 fa,Nf_{a,N} 的周期。

步骤 4:周期 rr 到因子的转化

定理:若 rr 是偶函数且 ar/2≢1(modN)a^{r/2} \not\equiv -1 \pmod{N},则 gcd(ar/21,N)\gcd(a^{r/2} - 1, N)gcd(ar/2+1,N)\gcd(a^{r/2} + 1, N) 都是 NN 的非平凡因子。

证明

ar1(modN)a^r \equiv 1 \pmod{N}ar10(modN)a^r - 1 \equiv 0 \pmod{N}

rr 是偶数,ar/21a^{r/2} - 1ar/2+1a^{r/2} + 1 的乘积被 NN 整除:

(ar/21)(ar/2+1)=ar10(modN)(a^{r/2} - 1)(a^{r/2} + 1) = a^r - 1 \equiv 0 \pmod{N}

如果 ar/2≢±1(modN)a^{r/2} \not\equiv \pm 1 \pmod{N},则 NN 有公共因子与 ar/2±1a^{r/2} \pm 1 共享——这些因子就是 ppqq\blacksquare

步骤 5:失败处理

如果 rr 是奇数,或 ar/21(modN)a^{r/2} \equiv -1 \pmod{N},选择另一个 aa 重复。

成功概率:对于随机选择的 aa,成功找到因子的概率至少 11/2k11 - 1/2^{k-1}(其中 kkNN 的不同素因子个数)。对 RSA 密钥(两个不同奇素数),k=2k=2,成功概率至少 1/21/2。因此期望在 O(1)O(1) 次尝试内成功。

完整归约示例

N=15N = 15:随机选 a=7a = 7gcd(7,15)=1\gcd(7,15)=1。计算 f(x)=7xmod15f(x) = 7^x \bmod 15

xx7x7^x7xmod157^x \bmod 15
011
177
2494
334313
424011

周期 r=4r = 4(偶数)。ar/2=72=494(mod15)a^{r/2} = 7^2 = 49 \equiv 4 \pmod{15}4≢±1(mod15)4 \not\equiv \pm 1 \pmod{15}。因此:

gcd(721,15)=gcd(48,15)=3\gcd(7^2 - 1, 15) = \gcd(48, 15) = 3 gcd(72+1,15)=gcd(50,15)=5\gcd(7^2 + 1, 15) = \gcd(50, 15) = 5

得到因子 3 和 5。✓


5.3 量子周期发现:电路设计

周期发现是 Shor 算法的”量子引擎”,它用量子相位估计 (QPE) 来找到 f(x)=axmodNf(x) = a^x \bmod N 的周期。

总电路图(两层寄存器):

|0⟩^⊗m —H^⊗m—•—•—•—•—•—•—•—•—H^⊗m—QFT†—[M]
               | | | | | | | |
|0⟩^⊗L ——————U U U U U U U U—————————[M]  (验证)
               a² a a¹ a⁸
               ⁸ ⁴ ²
  • 上层:m=2Lm = 2L 个量子比特的”计数”寄存器(L=log2NL = \lceil \log_2 N \rceil
  • 下层:LL 个量子比特的”工作”寄存器
  • UaU_a 是模乘算子:Uay=aymodNU_a\lvert y\rangle = \lvert ay \bmod N\rangle
  • 受控 Ua2jU_{a^{2^j}} 是模指数运算的构建块
模块化指数运算(Modular Exponentiation)

模指数运算 axmodNa^x \bmod N 是 Shor 算法中最昂贵的部分。它通过”平方-乘”法 (square-and-multiply) 分解为一系列受控模乘:

axmodN=aj=0m1xj2jmodN=j:xj=1(a2jmodN)modNa^x \bmod N = a^{\sum_{j=0}^{m-1} x_j 2^j} \bmod N = \prod_{j: x_j=1} \left(a^{2^j} \bmod N\right) \bmod N

电路实现需要:

  • 预计算 a2jmodNa^{2^j} \bmod N 对于 j=0,1,,m1j = 0, 1, \ldots, m-1(可以在经典计算机上预处理)
  • 使用受控模乘器 C-Ua2jC\text{-}U_{a^{2^j}}:当控制比特为 1\lvert1\rangle 时,施加 Ua2jU_{a^{2^j}}
  • 将这些受控操作级联起来

量子资源:模块化指数运算需要 O(L3)O(L^3) 个基本门(使用经典辅助的模乘算法),占用 O(L)O(L) 个辅助量子比特。

量子相位估计(QPE)用于周期发现

QPE 的详细电路实现在教程 Part 5.2 中。这里概述其在 Shor 算法中的应用:

输入寄存器:初始化为 +m=Hm0m\lvert+\rangle^{\otimes m} = H^{\otimes m}\lvert0\rangle^{\otimes m}(叠加态)。

工作寄存器:初始化为 1\lvert1\rangle

对每个比特 jj0j<m0 \leq j < m),施加受控 Ua2jU_{a^{2^j}} 操作:

Ua2jy=a2jymodNU_{a^{2^j}}\lvert y\rangle = \lvert a^{2^j} y \bmod N\rangle

QPE 的核心结果是:在计数寄存器上运行逆 QFT(QFTQFT^\dagger)后,测量得到周期 rr 的近似值 φ~\tilde{\varphi}

φ~kr\tilde{\varphi} \approx \frac{k}{r}

其中 kk00r1r-1 之间的随机整数。

关键数学推导UaU_a 的本征态包含: Uaus=e2πis/rusU_a\lvert u_s\rangle = e^{2\pi i s/r}\lvert u_s\rangle 其中 us=1rj=0r1e2πisj/rajmodN\lvert u_s\rangle = \frac{1}{\sqrt{r}}\sum_{j=0}^{r-1} e^{-2\pi i s j/r}\lvert a^j \bmod N\rangle

初始态 1\lvert1\rangle 可以展开为这些本征态的叠加: 1=1rs=0r1us\lvert1\rangle = \frac{1}{\sqrt{r}}\sum_{s=0}^{r-1}\lvert u_s\rangle

因此 QPE 测量得到相位 s/rs/r 的概率为 1/r1/r。从 φ~\tilde{\varphi} 提取 rr 需要连分数展开。


5.4 连分数算法(Continued Fractions Algorithm)

QPE 输出的是 φ=k/r\varphi = k/rmm 比特近似值 φ~\tilde{\varphi}(一个 mm 位二进制小数)。我们需要从 φ~\tilde{\varphi} 还原出分母 rr

连分数展开:任何有理数 φ\varphi 都可以唯一地表示为:

φ=a0+1a1+1a2+1+1an\varphi = a_0 + \cfrac{1}{a_1 + \cfrac{1}{a_2 + \cfrac{1}{\ddots + \cfrac{1}{a_n}}}}

简记为 [a0;a1,a2,,an][a_0; a_1, a_2, \ldots, a_n]

算法

φ0=φ~\varphi_0 = \tilde{\varphi} 开始,迭代:

  • aj=φja_j = \lfloor \varphi_j \rfloor(取整数部分)
  • φj+1=1/(φjaj)\varphi_{j+1} = 1/(\varphi_j - a_j)(取倒数的小数部分)

在第 jj 步,收敛分数 pj/qjp_j/q_jQFTQFT 测量值 φ~\tilde{\varphi} 的第 jj 个收敛)给出 k/rk/r 的近似:

[a0;a1,,aj]=pjqj[a_0; a_1, \ldots, a_j] = \frac{p_j}{q_j}

其中递推关系为: p2=0,  p1=1,  pj=ajpj1+pj2p_{-2} = 0,\; p_{-1} = 1,\; p_j = a_j p_{j-1} + p_{j-2} q2=1,  q1=0,  qj=ajqj1+qj2q_{-2} = 1,\; q_{-1} = 0,\; q_j = a_j q_{j-1} + q_{j-2}

收敛定理:若 pqφ<12q2\left|\frac{p}{q} - \varphi\right| < \frac{1}{2q^2},则 p/qp/qφ\varphi 的一个收敛分数。

应用于 Shor 算法:选择 m=2Lm = 2L(即 2log2N2\lceil\log_2 N\rceil),QPE 的精度保证:

krφ~<12m<12r2\left|\frac{k}{r} - \tilde{\varphi}\right| < \frac{1}{2^{m}} < \frac{1}{2r^2}

因此 φ~\tilde{\varphi} 的收敛分数之一就是 k/rk/r。检查分母 qjq_j:若 aqj1(modN)a^{q_j} \equiv 1 \pmod{N}(或找到了正确的因式分解),则 r=qjr = q_j

完整示例:从 φ~\tilde{\varphi} 提取周期

N=15N=15,实际周期 r=4r=4k=1k=1,所以 φ=1/4=0.25\varphi = 1/4 = 0.25

假定 QPE 输出 φ~=0.2501\tilde{\varphi} = 0.2501(微小误差)。连分数展开:

  • φ0=0.2501\varphi_0 = 0.2501a0=0a_0 = 0
  • φ1=1/0.25013.9984\varphi_1 = 1/0.2501 \approx 3.9984a1=3a_1 = 3
  • φ2=1/0.99841.0016\varphi_2 = 1/0.9984 \approx 1.0016a2=1a_2 = 1
  • φ3=1/0.0016=625\varphi_3 = 1/0.0016 = 625a3=625a_3 = 625

收敛分数:

  • [0;3]=0+1/3=1/3[0; 3] = 0 + 1/3 = 1/3r=3r=3(检查:73mod15=1317^3 \bmod 15 = 13 \neq 1,不是周期)
  • [0;3,1]=1/(3+1/1)=1/4[0; 3, 1] = 1/(3 + 1/1) = 1/4r=4r=4(检查:74mod15=17^4 \bmod 15 = 1,正确!)

因此 r=4r=4


5.5 完整工作示例:分解 15

这是 Shor 算法最经典的教学示例。

参数

  • N=15N = 15L=4L = 4(因为 15<24=1615 < 2^4 = 16
  • m=2L=8m = 2L = 8 个计数比特
  • 随机选择 a=7a = 7gcd(7,15)=1\gcd(7,15)=1

步骤 1:排除琐碎情况。15 是奇数,不是完全幂。

步骤 2:预计算 72jmod157^{2^j} \bmod 15

jj2j2^j72jmod157^{2^j} \bmod 15
0171mod15=77^1 \bmod 15 = 7
1272mod15=47^2 \bmod 15 = 4
2474mod15=17^4 \bmod 15 = 1
3878mod15=17^8 \bmod 15 = 1
416716mod15=17^{16} \bmod 15 = 1
\vdots\vdots\vdots

注意 j2j \geq 272j1(mod15)7^{2^j} \equiv 1 \pmod{15}

步骤 3:制备计数寄存器叠加态。

ψcount=1256x=0255x\lvert\psi_{\text{count}}\rangle = \frac{1}{\sqrt{256}}\sum_{x=0}^{255}\lvert x\rangle

工作寄存器为 1\lvert1\rangle

步骤 4:施加受控模乘。

由于 j2j\geq272jmod15=17^{2^j} \bmod 15 = 1,只有 j=0,1j=0,1 的受控模乘有实际效果。因此:

ψ=1256x=0255x7xmod15\lvert\psi\rangle = \frac{1}{\sqrt{256}}\sum_{x=0}^{255}\lvert x\rangle \lvert 7^x \bmod 15\rangle

步骤 5:施加逆 QFT 并测量。

假设测得 φ~=64/256=0.25\tilde{\varphi} = 64/256 = 0.25(即二进制 0.010.01)。

步骤 6:连分数展开。

0.25=[0;4]0.25 = [0; 4]r=4r = 4

步骤 7:验证与分解。

ar/2modN=72mod15=49mod15=4a^{r/2} \bmod N = 7^2 \bmod 15 = 49 \bmod 15 = 4gcd(721,15)=gcd(48,15)=3\gcd(7^2 - 1, 15) = \gcd(48, 15) = 3gcd(72+1,15)=gcd(50,15)=5\gcd(7^2 + 1, 15) = \gcd(50, 15) = 5

可能的历史差异:注意这里用 a=7a=7。在实际的量子实验中(如 Martin-López et al. 2012 的光量子实现),通常使用更简单的 a=11a=11(因为 112mod15=121mod15=111^2 \bmod 15 = 121 \bmod 15 = 1,周期 r=2r=2),或 a=4a=4(周期 r=2r=2)。周期 r=2r=2 的情况更简单,因为只需要受控模乘 j=0j=0 的一层操作。


5.6 复杂度分析与资源估算

量子门复杂度
子程序门复杂度说明
模指数运算O(L3)O(L^3)L 为比特数,使用经典模乘+量子受控操作
QFT / QFT†O(L2)O(L^2)m = 2L 比特的 QFT
总复杂度O(L3)O(L^3)= O((logN)3)O((\log N)^3)

总量子门数估算(以 RSA-2048 为例,L=2048L=2048):

  • 模指数运算:约 204838.6×1092048^3 \approx 8.6 \times 10^9 个基本门
  • 辅助比特:约 40004000
  • 所需物理量子比特(考虑纠错):数百万个
经典 vs 量子复杂度总表
操作经典复杂度量子复杂度加速
整数因子分解exp(O((logN)1/3(loglogN)2/3))\exp(O((\log N)^{1/3}(\log\log N)^{2/3}))O((logN)3)O((\log N)^3)指数
离散对数exp(O((logp)1/3(loglogp)2/3))\exp(O((\log p)^{1/3}(\log\log p)^{2/3}))O((logp)3)O((\log p)^3)指数
椭圆曲线离散对数O(exp(logp)1/3(loglogp)2/3)O(\exp(\log p)^{1/3}(\log\log p)^{2/3})O((logp)3)O((\log p)^3)指数(打破 ECC)

资源估算参考

要分解的数比特数经典 GNFS 时间量子门数量子时间估测*
154即时~10³微秒
215即时~10⁴微秒
1438微秒~10⁶毫秒
RSA-512512~10⁴ 年~10¹⁰小时
RSA-10241024~10⁶ 年~10¹¹
RSA-20482048~10¹¹ 年~10¹²

*量子时间估测基于 1 MHz 的物理门率,非常理想化。实际工程中还需考虑纠错开销(约 10³–10⁴ 倍资源增加)。


5.7 Shor 算法对 RSA 的威胁分析

短期(<10 年):不构成威胁。需要数百万高质量物理量子比特,远超当前水平。

中期(10-20 年):可能威胁 RSA-1024。需要约 100 万物理量子比特,纠错门保真度 > 99.9%。

长期(20-30 年):威胁所有 RSA 密钥长度。NIST 已经推进 PQC(后量子密码)标准化,FIPS 203/204/205 已于 2024 年发布。

量子比特资源:Gidney & Ekerå (2021) 的优化方案表明,分解 RSA-2048 需要约 2,000 个逻辑量子比特和约 2,000 万个物理量子比特门,运行时间约 8 小时。这比之前的估算(需要 10-100 倍资源)要乐观得多,但工程实现仍然面临巨大挑战。


6. 算法复杂度对比总表

算法问题经典复杂度量子复杂度加速类型实际应用
Deutsch-Jozsa常数/平衡判定Θ(2n)\Theta(2^n)O(n)O(n)指数教学
Bernstein-Vazirani隐藏比特串Θ(n)\Theta(n)O(1)O(1)线性教学
Simon隐藏周期(Z2n\mathbb{Z}_2^n)Θ(2n/2)\Theta(2^{n/2})O(n)O(n)指数教学,启发了 Shor
Grover无结构搜索Θ(N)\Theta(N)O(N)O(\sqrt{N})二次数据库搜索、优化
Shor整数因子分解亚指数O((logN)3)O((\log N)^3)指数破解 RSA

加速类型说明

  • 指数加速:量子比经典快指数级(问题规模翻倍时优势远超常数倍)
  • 二次加速:量子比经典快平方根级别
  • 线性加速:量子比经典快常数倍

7. 参考文献与延伸阅读

原始论文

  1. Deutsch, D. & Jozsa, R. (1992). “Rapid solution of problems by quantum computation”. Proceedings of the Royal Society A, 439(1907), 553–558.
  2. Bernstein, E. & Vazirani, U. (1997). “Quantum complexity theory”. SIAM Journal on Computing, 26(5), 1411–1473.
  3. Simon, D. R. (1997). “On the power of quantum computation”. SIAM Journal on Computing, 26(5), 1474–1483.
  4. Grover, L. K. (1996). “A fast quantum mechanical algorithm for database search”. Proceedings of STOC 1996, 212–219.
  5. Shor, P. W. (1997). “Polynomial-time algorithms for prime factorization and discrete logarithms on a quantum computer”. SIAM Journal on Computing, 26(5), 1484–1509.
  6. Bennett, C. H., Bernstein, E., Brassard, G., & Vazirani, U. (1997). “Strengths and weaknesses of quantum computing”. SIAM Journal on Computing, 26(5), 1510–1523. (BBBV 最优性证明)

教材与综述

  1. Nielsen, M. A. & Chuang, I. L. (2010). Quantum Computation and Quantum Information: 10th Anniversary Edition. Cambridge University Press. (第 6 章: 量子搜索算法; 第 5 章: QFT 与 Shor 算法)
  2. Kaye, P., Laflamme, R., & Mosca, M. (2007). An Introduction to Quantum Computing. Oxford University Press.
  3. Mermin, N. D. (2007). Quantum Computer Science: An Introduction. Cambridge University Press.
  4. Jozsa, R. (1997). “Quantum algorithms and the Fourier transform”. Proceedings of the Royal Society A, 454(1969), 323–337.

专题深入

  1. Grover, L. K. (1998). “Quantum computers can search arbitrarily large databases by a single query”. Physical Review Letters, 79(23), 4709.
  2. Boyer, M., Brassard, G., Høyer, P., & Tapp, A. (1998). “Tight bounds on quantum searching”. Fortschritte der Physik, 46(4-5), 493–506.
  3. Brassard, G., Høyer, P., Mosca, M., & Tapp, A. (2002). “Quantum amplitude amplification and estimation”. Contemporary Mathematics, 305, 53–74. (Grover 的推广)
  4. Kitaev, A. Y. (1995). “Quantum measurements and the Abelian stabilizer problem”. arXiv:quant-ph/9511026. (HSP 框架)
  5. Gidney, C. & Ekerå, M. (2021). “How to factor 2048 bit RSA integers in 8 hours using 20 million noisy qubits”. Quantum, 5, 433. (资源优化里程碑)
  6. Beauregard, S. (2003). “Circuit for Shor’s algorithm using 2n+3 qubits”. Quantum Information and Computation, 3(2), 175–185.
  7. Martin-López, E. et al. (2012). “Experimental realization of Shor’s quantum factoring algorithm using qubit recycling”. Nature Photonics, 6, 773–776. (首次完整演示 Shor 分解 15)

本教程内部参考

  1. Part 1.8: 离散傅里叶变换(DFT)的数学形式
  2. Part 5.1: QFT 电路的详细构造(受控相位门、级联结构)
  3. Part 5.2: 量子相位估计(QPE)算法——Shor 算法的量子引擎

文档版本:v1.0
最后更新:2026-06-03
关联教程:量子计算前置教程——从第一性原理出发 (Part 3.6, Part 5.1, Part 5.2)
编写原则:标准量子计算记法,与教程 Part 1–3 的数学物理基础一致