因为生着脚,所以要走

# 量子搜索算法

# 黑箱(Oracle)

假设在 N 个元素的搜索空间中进行搜索。我们不关心元素是什么,只关心元素的编号 0~N-1。方便起见,规定N=2nN=2^n。假设搜索共有 M 个解,1MN1\leq M\leq N

  • 该搜索问题的一个特例可以方便地表示为一个输入 x 的函数 f,x 是从 0 到 N-1 的整数,f 的定义是,若 x 是搜索问题的一个解f(x)=1f(x)=1,而如果 x 不是搜索问题的解,则f(x)=0f(x)=0

  • 设有一个量子 oracle 可以识别搜索问题的解,Oracle 的内部工作方式将在后面讨论,识别的结果通过 oracle 的一个量子比特给出。更确切低,这个 oracle 是一个酉算子 O,满足:

    Oxqxqf(x)O|x\rang|q\rang\rightarrow|x\rang|q\oplus f(x)\rang

    其中\oplus 是模二加。因此我们可以通过制备x0|x\rang|0\rang 再应用 O 和测量第二寄存器来判断 x 是否是搜索问题的一个解。

  • 特别地,在状态012\frac{|0\rang-|1\rang}{\sqrt{2}} 上 oracle 的作用为

    Ox012(1)f(x)x012O|x\rang\frac{|0\rang-|1\rang}{\sqrt{2}}\rightarrow(-1)^{f(x)}|x\rang\frac{|0\rang-|1\rang}{\sqrt{2}}

    可以看作第二寄存器状态未发生改变,而第一寄存器状态乘了个(1)f(x)(-1)^{f(x)}。由于第二寄存器状态在量子搜索算法过程中保持为012\frac{|0\rang-|1\rang}{\sqrt{2}},因此我们可以简写 O 的作用为Ox(1)f(x)xO|x\rang\rightarrow (-1)^{f(x)}|x\rang

  • 我们说 oracle 通过改变解的相位标记了搜索问题的解。对有 M 个解的 N 元搜索问题,实际上只需要在量子计算机上应用搜索 oracleO(N/M)O(\sqrt{N/M}) 次。


  • 以上讨论是基于我们知道 oracle 是怎样工作的。那这些讨论有什么意义呢?事实上,有时候不需要解是哪些,我们就可以识别和判断解

# 过程

  • 搜索算法运行的框架如下图所示。

1

​ oracle 的工作可能需要附加的工作量子比特,但是算法的目的是用最少的 oracle 应用次数求出搜索问题的一个解。

  • 其中,量子搜索算法由反复应用记作GG 的被称为 Grover 迭代 Grover 算子的量子子程序组成。Grover 迭代的量子线路如下图:

    2

    可分为四步:

    • 应用 oracle O;
    • 应用 Hadamard 变换HnH^{\otimes n}
    • 在计算机上执行条件相移,使0|0\rang 以外每个计算基态获得 - 1 的相位:(200I)x(1)δx,0x(2|0\rang\lang 0|-I)|x\rang\rightarrow(-1)^{\delta_{x,0}}|x\rang
    • 应用 Hadamard 变换HnH^{\otimes n}
  • Grover 迭代中的每个运算都可在量子计算机上有效实现。第 2,4 步的应用 Hadamard 变换各需要n=logNn=logN 个门,第三步的条件相移利用之后介绍的技术也可以只用O(n)O(n) 个门。第一步 oracle 的调用消耗需要依赖特定的应用,这个例子中只需要一个单 oracle 调用。

  • 把 2,3,4 步的算子结合起来:

    Hn(200I)Hn=22nx,y=02n1xyIH^{\otimes n}(2|0\rang\lang 0|-I)H^{\otimes n}=\frac{2}{2^n}\sum_{x,y=0}^{2^n-1}|x\rang\lang y|-I

    ψ=1Nx=0N1x|\psi\rang=\frac{1}{\sqrt{N}}\sum_{x=0}^{N-1}|x\rang,则 Grover 迭代 G 可写作

    G=(2ψψI)OG=(2|\psi\rang\lang\psi|-I)O

    应用运算2ψψI2|\psi\rang\lang\psi|-I 到一般状态k=0N1akk\sum_{k=0}^{N-1}a_k|k\rang 会得到状态

    (2ψψI)k=0N1akk=k=0N1[ak+2a]k(2|\psi\rang\lang\psi|-I)\sum_{k=0}^{N-1}a_k|k\rang=\sum_{k=0}^{N-1}[-a_k+2\lang a\rang]|k\rang

    其中a=1Nk=0N1ak\lang a\rang=\frac{1}{N}\sum_{k=0}^{N-1}a_k 是均值,因此2ψψI2|\psi\rang\lang\psi|-I 有时称为关于均值的反演运算。

# 几何可视化

  • 事实上,Grover 迭代可视为在由开始向量ψ|\psi\rang 和搜索问题解组成的均匀叠加态张成的二维空间中的一个旋转。

    首先,我们定义x\sum_{x}' 求和符号中的 x 枚举的是搜索问题所有的解,x\sum_{x}'' 枚举的是非搜索问题的解。定义归一化状态:

    α=1NMxxβ=1Mxx|\alpha\rang=\frac{1}{\sqrt{N-M}}\sum_x''|x\rang\\ |\beta\rang=\frac{1}{\sqrt{M}}\sum_x'|x\rang

    则有x=0N1=x+x\sum_{x=0}^{N-1}=\sum_x'+\sum_{x}'',因此均匀叠加态:

    ψ=1Nx=0N1x=NMNα+MNβ|\psi\rang=\frac{1}{\sqrt{N}}\sum_{x=0}^{N-1}|x\rang=\sqrt{\frac{N-M}{N}}|\alpha\rang+\sqrt{\frac{M}{N}}|\beta\rang

    故量子计算机初态属于α|\alpha\rangβ|\beta\rang 张成的向量空间。

  • 可以认识到α|\alpha\rangβ|\beta\rang 是正交的,对于任何一个在α,β|\alpha\rang,|\beta\rang 所在平面上的向量aα+bβa|\alpha\rang+b|\beta\rang,O 作用在上面就是:

    O(aα+bβ)=aαbβO(a|\alpha\rang+b|\beta\rang)=a|\alpha\rang-b|\beta\rang

    几何上上就是把aα+bβa|\alpha\rang+b|\beta\rang 关于α|\alpha\rang 对称一下。不妨令cosθ2=NMNcos\frac{\theta}{2}=\sqrt{\frac{N-M}{N}},则ψ=cosθ2α+sinθ2β|\psi\rang=cos\frac{\theta}{2}|\alpha\rang+sin\frac{\theta}{2}|\beta\rangOψ=cosθ2αsinθ2βO|\psi\rang=cos\frac{\theta}{2}|\alpha\rang-sin\frac{\theta}{2}|\beta\rang

    (2ψψI)(cosθ2αsinθ2β)=(2cosθ1)cosθ2α+(2cosθ+1)sinθ2b=cos3θ2α+sin3θ2β(2|\psi\rang\lang\psi|-I)(cos\frac{\theta}{2}|\alpha\rang-sin\frac{\theta}{2}|\beta\rang)\\=(2cos\theta-1)cos\frac{\theta}{2}|\alpha\rang+(2cos\theta+1)sin\frac{\theta}{2}|b\rang\\ =cos\frac{3\theta}{2}|\alpha\rang+sin\frac{3\theta}{2}|\beta\rang

    就是一个三倍角θ2\frac{\theta}{2}。画到图上如下:

    3

    Gkψ=cos(2k+12θ)α+sin(2k+12θ)βG^k|\psi\rang=cos(\frac{2k+1}{2}\theta)|\alpha\rang+sin(\frac{2k+1}{2}\theta)|\beta\rang 不难发现每次迭代就是把状态向量逆时针转θ\theta 角,而θ\theta 是二倍的ψ|\psi\rangα|\alpha\rang 的夹角。经过O(N/M)O(\sqrt{N/M}) 次后,状态向量将接近β|\beta\rang(至多偏离βπ4|\beta\rang\frac{\pi}{4} 角)。然后其实对β|\beta\rang 做计算基下的测量得出的就是搜索问题的解,而状态向量离β|\beta\rang 越近,对其测量越大概率得到与β|\beta\rang 重合的解,即搜索问题的解。


    例子:下面是在 N=4 的搜索空间上的量子搜索算法中的 Grover 迭代一次对应的线路。oracle 对除x0x_0 外所有的 x 有f(x)=0f(x)=0,而f(x0)=1f(x_0)=1。通过不同的x0x_0oracle 线路也不同。下图是从左到右对应x0=0,1,2,3x_0=0,1,2,3 的 oracle:

    4

    起初,上面两个量子比特处于状态00|0\rang|0\rang,下面的量子比特为1|1\rang。中间的用点框起来的部分实现的是条件相移20000I2|00\rang\lang 00|-I

    为得到x0x_0,我们需要应用上图线路多少次?

    其实可以看出θ=2arctanMNM=2arctan13=π3\theta=2arctan\sqrt{\frac{M}{N-M}}=2arctan\frac{1}{\sqrt{3}}=\frac{\pi}{3},所以迭代一次就可以获得β|\beta\rang。迭代一次后对上面两个量子比特进行测量就可以得出x0x_0。与之对照,经典计算机四选一的 oracle 平均需要 2.25 次调用。

    注:其实这里可能产生疑问,最后一个量子比特全程没有参与感,为啥要用它?实际上,它参与了最后的测量。最后对整个系统测量算子为001001,011011,101101,111111|001\rang\lang001|,|011\rang\lang011|,|101\rang\lang101|,|111\rang\lang 111|,即我们需要保证下面一个量子比特状态为1|1\rang(和输入值一样)的情况下去测量前两个量子比特。


# 性能

  • 系统初态为ψ=NMNα+MNβ|\psi\rang=\sqrt{\frac{N-M}{N}}|\alpha\rang+\sqrt{\frac{M}{N}}|\beta\rang,与β|\beta\rang 夹角成arccosMNarccos\sqrt{\frac{M}{N}},所以实际上需要迭代arccosM/Nθ\lfloor\frac{arccos\sqrt{M/N}}{\theta}\rceil 次,把ψ|\psi\rang 转到离β|\beta\rangθ/2π/4\theta/2\leq\pi/4 范围内。于是对状态在计算基中的观测,将至少以1/21/2 的概率是在β|\beta\rang 部分,即得到搜索的一个解。

  • 特别地,当M<<NM<<N 时,我们有θsinθ=2M/N\theta\approx sin\theta=2\sqrt{M/N},故最终误差至多是θ/2M/N\theta/2\approx\sqrt{M/N}。事实上,需要迭代的次数和精度分析并不取决于解具体是什么只取决于解的个数,因此知道 M 就可以搜索算法分析误差。

  • R=arccosM/Nθπ2θθ2sinθ2=MNRπ4NMR=\lfloor\frac{arccos\sqrt{M/N}}{\theta}\rceil\leq\lceil\frac{\pi}{2\theta}\rceil\\ \because \frac{\theta}{2}\geq sin\frac{\theta}{2}=\sqrt{\frac{M}{N}}\\ \therefore R\leq\lceil\frac{\pi}{4}\sqrt{\frac{N}{M}}\rceil

    给出了需要迭代次数一个非常漂亮的上界,即R=O(N/M)R=O(\sqrt{N/M}) 次 Grover 迭代,oracle 调用,才能以最高的概率得到搜索问题的一个解。

  • 特别地,当MN/2M\geq N/2 时,即θ/2π/4\theta/2\geq\pi/4 时,我们可以一开始就随机从搜索空间中选择一个元检测它是否为一个解。此时成功概率已经至少为1/21/2。同样,我们也可以增加 N 个均不为解的元到搜索空间,使得解的个数一定不到新搜索空间大小的一半。这可以通过增加一个单量子比特q|q\rang 都搜索指标,把要搜索的元数目加倍到2N2N 实现。构造一个增广 oracle OO',因为引进q|q\rang 而使得OO' 增加的元素全部置为 0 即可。

  • 然而在实际问题中,我们必须懂得如何实现 oracle。因此在我们考虑的每个实际问题中,我们都给出 oracle 实现的具体描述。

# 作为量子仿真的量子搜索

  • 我们考虑构造一个 Hamilton 量,使得在自然选择的演化下,在一定时间内就可以实现初态ψ|\psi\rang 到解x|x\rang 的演化。

    例如H=xx+ψψH=|x\rang\lang x|+|\psi\rang\lang\psi|,根据薛定谔方程的解ψ=exp(iHt)ψ|\psi'\rang=exp(-iHt)|\psi\rang,在极短时间tt 内,将exp(iHt)exp(-iHt) 泰勒展开保留前两项:

    (IiHt)ψ=(1it)ψitxψx(I-iHt)|\psi\rang=(1-it)|\psi\rang-it\lang x|\psi\rang|x\rang

    ψ|\psi\rang 稍微转了一部分去x|x\rang

  • 但当 t 不再是小量时,我们做一下全程分析看是否存在 t,使得exp(iHt)ψ=xexp(-iHt)|\psi\rang=|x\rang

    考虑ψ|\psi\rangx|x\rang 张成的二维空间,y|y\rang 是二维空间上和x|x\rang 正交的单位向量。不妨设ψ=αx+βy|\psi\rang=\alpha|x\rang+\beta|y\rang,假定选择合适的x,y|x\rang,|y\rang 相位,使得α,β\alpha,\beta 是非负实数(即把α=reiθ\alpha=re^{i\theta}eiθe^{i\theta} 部分作为相位乘到x|x\rang 上),且α2+β2=1\alpha^2+\beta^2=1。此时,我们考虑在x,y|x\rang,|y\rang 基中,(即令x=(10),y=(01)|x\rang=\begin{pmatrix}1\\0\end{pmatrix},|y\rang=\begin{pmatrix}0\\1\end{pmatrix},有:

    H=xx+ψψ=(1000)+(α2αβαββ2)=(1+α2αβαβ1α2)=I+α(βX+αZ)H=|x\rang\lang x|+|\psi\rang\lang\psi|=\begin{pmatrix}1&0\\0&0\end{pmatrix}+\begin{pmatrix}\alpha^2&\alpha\beta\\\alpha\beta&\beta^2\end{pmatrix}=\begin{pmatrix}1+\alpha^2&\alpha\beta\\ \alpha\beta&1-\alpha^2\end{pmatrix}=I+\alpha(\beta X+\alpha Z)

    exp(iHt)=eiIteiαt(βX+αZ)Rn^(θ)=exp(iθn^σ/2)=cos(θ2)Iisinθ2(n^σ)eiαt(βX+αZ)=cos(αt)Iisin(αt)(βX+αZ)exp(iHt)=eit[cos(αt)Iisin(αt)(βX+αZ)]exp(-iHt)=e^{-iIt}e^{-i\alpha t(\beta X+\alpha Z)}\\ \because R_{\hat{n}}(\theta)=exp(-i\theta\hat{n}\cdot\overrightarrow{\sigma}/2)=cos(\frac{\theta}{2})I-isin\frac{\theta}{2}(\hat{n}\cdot\overrightarrow{\sigma})\\ \therefore e^{-i\alpha t(\beta X+\alpha Z)}=cos(\alpha t)I-isin(\alpha t)(\beta X+\alpha Z)\\ \therefore exp(-iHt)=e^{-it}[cos(\alpha t)I-isin(\alpha t)(\beta X+\alpha Z)]

    所以

    exp(iHt)ψ=eit[cos(αt)ψisin(αt)(βX+αZ)ψ]exp(-iHt)|\psi\rang=e^{-it}[cos(\alpha t)|\psi\rang-isin(\alpha t)(\beta X+\alpha Z)|\psi\rang]

    x,y|x\rang,|y\rang 基下,ψ=(αβ)|\psi\rang=\begin{pmatrix}\alpha\\ \beta\end{pmatrix},所以有(βX+αZ)ψ=(10)=x(\beta X+\alpha Z)|\psi\rang=\begin{pmatrix}1\\0\end{pmatrix}=|x\rang。忽略全局相位eite^{-it},于是过了时间 t,系统进入

    cos(αt)ψisin(αt)xcos(\alpha t)|\psi\rang-isin(\alpha t)|x\rang

    于是在t=π2αt=\frac{\pi}{2\alpha} 时刻,对系统的观测将以概率 1 产生结果x|x\rang,即我们找到了解。


    不幸的是,观测时间tt 依赖ψ|\psi\rangx|x\rang 方向上的分量α\alpha,从而依赖x|x\rang,而x|x\rang 又是我们要求的。

    实际上,考虑到搜索解xx 一定是一个长度为 n 的二进制串,因此x|x\rang 不可能是叠加态,而是nn 量子比特标准正交计算基下的一个(例如x|x\rang 可以为010|010\rang,但不能为12110+12100\frac{1}{\sqrt{2}}|110\rang+\frac{1}{\sqrt{2}}|100\rang)。所以我们可以制备均匀叠加态ψ=1Nx=0N1x|\psi\rang=\frac{1}{\sqrt{N}}\sum_{x=0}^{N-1}|x\rang,这样对于任意一个x|x\rangψ|\psi\rangx|x\rang 上的投影αψx=1N\alpha\equiv\lang\psi|x\rang=\frac{1}{\sqrt{N}}t=πN/2t=\pi\sqrt{N/2} 不依赖于xx 的值。而且均匀叠加态是很好用HH 门制造出来。


  • 现在我们用一个量子线路来模拟H=xx+ψψH=|x\rang\lang x|+|\psi\rang\lang\psi|。下面两个图便分别实现了IiΔtxxexp(ixxΔt)I-i\Delta t|x\rang\lang x|\approx exp(-i|x\rang\lang x|\Delta t)IiΔtψψexp(iψψΔt)I-i\Delta t|\psi\rang\lang\psi |\approx exp(-i|\psi\rang\lang\psi|\Delta t),其中Δt\Delta t 是个小量。在作业中我验证了其正确性。注意图中eiΔte^{i\Delta t} 应换为eiΔte^{-i\Delta t}

    5

    又因为eiHt=eixxteiψψte^{-iHt}=e^{-i|x\rang\lang x|t}e^{-i|\psi\rang\lang\psi| t},所以我们可以以Δt\Delta t 的步长交替使用eixxΔte^{-i|x\rang\lang x|\Delta t}eiψψΔte^{-i|\psi\rang\lang\psi |\Delta t},每次使用上面线路时实际上和真正的eixxΔt,eiψψΔte^{-i|x\rang\lang x|\Delta t},e^{-i|\psi\rang\lang\psi |\Delta t}O(Δt2)O(\Delta t^2) 的误差(因为我只把eixxΔt=IiΔtxxe^{-i|x\rang\lang x|\Delta t}=I-i\Delta t|x\rang\lang x 展开到一阶小量)

  • 而需要的总步数为t/Δt=πN/2/Δt=O(N/Δt)t/\Delta t=\pi\sqrt{N/2}/\Delta t=O(\sqrt{N}/\Delta t),误差在O(Δt2×N/Δt)=O(ΔtN)O(\Delta t^2\times \sqrt{N}/\Delta t)=O(\Delta t\sqrt{N})

    因此如果想让误差为O(1)O(1),可以选择Δt\Delta t~O(1/N)O(1/\sqrt{N})

  • 进一步地,如果我们实现的线路更为出色,可以以O(Δt3)O(\Delta t^3) 的精度去实现eixxΔte^{-i|x\rang\lang x|\Delta t}eiψψΔte^{-i|\psi\rang\lang\psi |\Delta t},那么积累误差会更小。


  • 那么考虑算子U(Δt)=exp(ixxΔt)exp(iψψΔt)U(\Delta t)=exp(-i|x\rang\lang x|\Delta t)exp(-i|\psi\rang\lang\psi|\Delta t),显然它只有作用在x,ψ|x\rang,|\psi\rang 张成的二维平面上时是不平凡的。

    在基x,y|x\rang,|y\rang 讨论下,xx=(I+Z)/2,ψψ=(I+ψσ)/2|x\rang\lang x|=(I+Z)/2,|\psi\rang\lang\psi|=(I+\overrightarrow{\psi}\cdot\overrightarrow{\sigma})/2,其中ψ=(2αβ,0,α2β2)\overrightarrow{\psi}=(2\alpha\beta,0,\alpha^2-\beta^2)。所以在 Bloch 球面上,x|x\rang 对应(0,0,1)(0,0,1)y|y\rang 对应(0,0,1)(0,0,-1)ψ|\psi\rang 对应(2αβ,0,α2β2)(2\alpha\beta,0,\alpha^2-\beta^2)。`

    不难验证(好吧其实挺烦),U(Δt)=exp(ixxΔt)exp(iψψΔt)U(\Delta t)=exp(-i|x\rang\lang x|\Delta t)exp(-i|\psi\rang\lang\psi|\Delta t) 实际上是先绕z^\hat{z} 旋转Δt\Delta t,再绕ψ\overrightarrow{\psi}Δt\Delta t。合成后:

    U(Δt)=(cos2Δt2sin2Δt2ψz^)I2isinΔt2(cosΔt2ψ+z^2+sinΔt2ψ×z^2)σU(\Delta t)=(cos^2\frac{\Delta t}{2}-sin^2\frac{\Delta t}{2}\overrightarrow{\psi}\cdot\hat{z})I-2isin\frac{\Delta t}{2}(cos\frac{\Delta t}{2}\frac{\overrightarrow{\psi}+\hat{z}}{2}+sin\frac{\Delta t}{2}\frac{\overrightarrow{\psi}\times\hat{z}}{2})\cdot\overrightarrow{\sigma}

    所以就是绕r=cosΔt2ψ+z^2+sinΔt2ψ×z^2\overrightarrow{r}=cos\frac{\Delta t}{2}\frac{\overrightarrow{\psi}+\hat{z}}{2}+sin\frac{\Delta t}{2}\frac{\overrightarrow{\psi}\times\hat{z}}{2}(提个 2 是为了归一化)

    2arccos(cos2Δt2sin2Δt2ψz^)=2arccos(12Nsin2Δt2)2arccos(cos^2\frac{\Delta t}{2}-sin^2\frac{\Delta t}{2}\overrightarrow{\psi}\cdot\hat{z})=2arccos(1-\frac{2}{N}sin^2\frac{\Delta t}{2})。下图为一个形象描述:

    6

  • 为了使每次转的角度尽量大,我们可以选择Δt=π\Delta t=\pi,这这样每次就转2arccos(12/N)4/N2arccos(1-2/N)\approx 4/\sqrt{N}。次数需要的总步数为πN/2/π=O(N)\pi\sqrt{N/2}/\pi=O(\sqrt{N})。实际上,取Δt=π\Delta t=\pi 时已经等价于原始的量子搜索算法。因为eiπψψ=I2ψψe^{-i\pi|\psi\rang\lang\psi|}=I-2|\psi\rang\lang\psi |。并且除了一个不重要的相位因子,这些算子就等同于构成 Grover 迭代的步骤。


  • 给定要解的问题,包含量子算法输入和输出的描述

    猜测一个 Hamilton 量来解决问题,并验证其可用

    寻找一个模拟该 Hamilton 量的过程

    分析消耗的资源

  • 量子算法就是量子仿真。

# 量子计数

在事先未知的情况下,多快可以确定 N 元搜索问题的解的个数 M?

  • 显然经典计算机上需要Θ(N)\varTheta(N) 次对 oracle 的调用来确定 M。

    在量子计算机上,有可能通过把 Grover 迭代和基于量子 Fourier 变换的相位估计结合,以比经典计算机快得多的方式估计解的数目。

  • a,b|a\rang,|b\rangα,β|\alpha\rang,|\beta\rang 张成的平面中 Grover 迭代算子 G 对应的两个特征向量,令θ\theta 表示 G 作用效果的旋转角,因为G=(cosθsinθsinθcosθ)G=\begin{pmatrix}cos\theta&-sin\theta\\sin\theta&cos\theta\end{pmatrix},所以其两个特征值为eiθ,ei(2πθ)e^{i\theta},e^{i(2\pi-\theta)}。假设 O 已经被扩充(但 G 和 O 还是只作用于α,β|\alpha\rang,|\beta\rang 张成的平面上的向量是非平凡的),即搜索空间规模被扩大到2N2N,并保证sinθ2=M2Nsin\frac{\theta}{2}=\sqrt{\frac{M}{2N}}

    事实上,通过相位估计在 G 上的应用,可以让我们估计θ\theta,从而估计解的个数。

    下图是一个至少以1ϵ1-\epsilon 的成功概率估计θ\theta 到 m 比特精度。

    7

    其中,第一寄存器共有t=m+log(2+12ϵ)t=m+\lceil log(2+\frac{1}{2\epsilon})\rceil 个量子比特,第二个寄存器包含n+1n+1 个量子比特。注意到第二寄存器上来也制备了均匀叠加态ψ|\psi\rang,因为保证了ψ|\psi\rang 一定在α,β|\alpha\rang,|\beta\rang 平面上,因此它一定是a,b|a\rang,|b\rang 的叠加,因此最后测出来以一定概率是θ\theta 的近似,一定概率是2πθ2\pi-\theta 的近似。而且精度高达Δθ2m|\Delta\theta|\leq 2^{-m},且sinθ2=sin2πθ2sin\frac{\theta}{2}=sin\frac{2\pi-\theta}{2}


    又因为sin2θ2=M2N\sin^2\frac{\theta}{2}=\frac{M}{2N},所以ΔM2N=sin2θ+Δθ2sin2θ2=(sinθ+Δθ2+sinθ2)sinθ+Δθ2sinθ2\frac{|\Delta M|}{2N}=|sin^2\frac{\theta+\Delta\theta}{2}-sin^2\frac{\theta}{2}|=(\sin\frac{\theta+\Delta\theta}{2}+sin\frac{\theta}{2})|sin\frac{\theta+\Delta\theta}{2}-sin\frac{\theta}{2}|

    又因为sinθ+Δθ2sinθ2Δθ2,sinθ+Δθ2<sin(θ/2)+Δθ/2|sin\frac{\theta+\Delta\theta}{2}-sin\frac{\theta}{2}|\leq\frac{|\Delta\theta|}{2},|sin\frac{\theta+\Delta\theta}{2}|<sin(\theta/2)+|\Delta\theta|/2,故\frac{|\Delta M|}{2N}<(2sin\frac{\theta}{2}+\frac{|\Delta\theta|}{2})\frac{|\Delta\theta|}< p>

    再把sin2θ2=M2N,Δθ2msin^2\frac{\theta}{2}=\frac{M}{2N},|\Delta\theta|\leq 2^{-m} 代入,有|\Delta M|<(\sqrt{2MN}+\frac{N}{2^{m+1}})2^< p>

    • 例如我们取m=n2+1m=\lceil\frac{n}{2}\rceil+1,且ϵ=1/6\epsilon=1/6,则有t=n2+3t=\lceil\frac{n}{2}\rceil+3。注意此时线路图中虽然好像有快速幂的样子,但实则和前面不同(前面的算子 U 满足Ux=xyU|x\rang=|xy\rang 是个乘法)。实际上G2iG^{2^i} 就是很朴素地连续作用2i2^i 个 G 门。因此总共需要O(2t)O(2^{t}) 个 G 门,即O(N)O(\sqrt{N}) 个 Grover 迭代、oracle 调用。此时ΔM<M/2+1/4=O(M)|\Delta M|<\sqrt{M/2}+1/4=O(\sqrt{M})

    • 特别地,若M=0M=0,则ΔM<1/4|\Delta M|<1/4,于是算法至少以1ϵ=5/61-\epsilon=5/6 的概率给出估计 0。反之,若M0M\neq 0,则以5/65/6 的概率给出 M 不是 0 的估计。这样在判断 NP 问题是否有解时非常有用。


  • 量子计数的另一个应用是当解数目 M 未知时,找出问题的一个解。之前朴素的黑箱量子搜索的 Grover 迭代次数需要依赖 M,但是我们现在可以先用相位估计估计 M,再代入量子搜索。

    在之前的量子搜索中,我们分析了通过不断地 Grover 迭代和旋转,初状态ψ|\psi\rang 可以旋转至离末状态x|x\rang 角度θ/2π/4\theta/2\leq\pi/4 范围里。此时因为量子技术过程中的相位估计会对θ\theta 估计不准确,因此ψ|\psi\rang 旋转后最多最多偏离x|x\rangπ4(1+Δθθ)\frac{\pi}{4}(1+\frac{\Delta \theta}{\theta})

# 判断 Hamilton 圈问题的加速

HC 问题为确定给定的 n 个顶点的图是否具有一个 Hamilton 圈。这个问题属于 NP-complete 问题,一般认为在经典计算机上不可解。

  • 朴素算法是遍历顶点的全排列进行分析。但这显然很差劲。事实上,NP 类中的任何问题都可以以类似的方式来解:如果有ω(n)\omega(n) 比特给定的证据,其中ω(n)\omega(n) 为 n 的多项式,如果解存在,则可以遍历2ω(n)2^{\omega(n)} 则一定可以得出解。

  • 然而使用量子算法,我们可以令mlognm\equiv\lceil logn\rceil,然后将搜索空间用mnmn 个量子比特状态来表示:v1...vn|v_1...v_n\rang,其中vi|v_i\rang 是一个 m 量子比特状态用于表示第 i 个顶点的编号。那么黑箱 oracle 就要实现:

    Ov1...vn={v1...vnv1,...,vn不是Hamiltonv1...vnv1,...,vnHamiltonO|v_1...v_n\rang=\begin{cases}|v_1...v_n\rang&v_1,...,v_n不是Hamilton圈\\ -|v_1...v_n\rang&v_1,...,v_n是Hamilton圈\end{cases}

    当我们有图的描述时,这样的 oracle 是容易实现的。可以取一个能判断图中 Hamilton 圈的多项式规模的电路,然后使用量子搜索算法。


    • 事实上,此时我们需要调用 oracle O(N)=O(2t/2)=O(2mn/2)O(\sqrt{N})=O(2^{t/2})=O(2^{mn/2}) 次来判断是否有 Hamilton 圈。

      量子搜索算法渐进于经典算法平方根的运算次数。

# 非结构化数据库的量子搜索

给你一个包含一千个名字的清单,问你 “Perth Rose” 出现在名单上第几个?

  • 如果遍历搜索,平均需要检查 500 个名字才能找到。但用量子数据库搜索算法可以加速这个过程。

  • 设我们有一个包含N2nN\equiv 2^n 个元的数据库,每个元的长度为 l 比特,这些元记为d1,...,dNd_1,...,d_N,我们希望确定一个特定的 l 比特串 s 是否在数据库中。

  • 假设我们的量子计算机由 CPU 和内存这两个单元组成,其中 CPU 有四个寄存器:

    • 一个初始化为0|0\rang 的 n 量子比特索引寄存器
    • 一个初始化为s|s\rang,并在整个过程中保持该状态的 l 量子比特寄存器
    • 一个初始化为0|0\rang 的 l 量子比特数据寄存器
    • 一个初始化为(01)/2(|0\rang-|1\rang)/\sqrt{2} 的 1 量子比特寄存器

    而内存则为2n2^n 个单元的经典内存,每个单元有 l 比特。然而与经典内存不同的是,每个单元的索引 x 可以处于多重值的叠加,即他不再是离散的地址数字,而且这个索引x|x\rang 对应的也是单元值的叠加dx|d_x\rang

  • 例如现在 CPU 的索引寄存器在状态x|x\rang,且数据寄存器处在状态d|d\rang。则x|x\rang 索引的内容dxd_x 加载到数据寄存器则是dddx|d\rang\rightarrow|d\oplus d_x\rang,其中\oplus 表示异或。实现量子搜索算法的关键在于 oracle 的实现,它必须翻转在内存中定位 s 索引的相位。

  • 设 CPU 处于状态xs0012|x\rang|s\rang|0\rang\frac{|0\rang-|1\rang}{\sqrt{2}}

    应用 LOAD 操作,载入数据使进入状态xsdx012\rightarrow |x\rang|s\rang|d_x\rang\frac{|0\rang-|1\rang}{\sqrt{2}}

    这时应用 oracle 比较第二,第三寄存器,如果相同则应用比特翻转到第四寄存器,否则不变:

    (1)[s==dx]xsdx012\rightarrow(-1)^{[s==d_x]}|x\rang|s\rang|d_x\rang\frac{|0\rang-|1\rang}{\sqrt{2}}

    然后再应用 LOAD 操作,把数据寄存器还原xs0012\rightarrow|x\rang|s\rang|0\rang\frac{|0\rang-|1\rang}{\sqrt{2}}


    其实整个过程我们发现寄存器之间都是独立的没有纠缠在一起,而且二三四寄存器始终没有变化就一直是在改变第一寄存器的相位。此时我们就有了一个很好的 oracle,然后我们可以把第一寄存器制备初态ψ|\psi\rang,然后利用 Grover 迭代量子搜索,在O(2n)O(\sqrt{2^n}) 次 LOAD 操作后第一寄存器旋转至 s 对应的索引xs|x_s\rang 的附近,从而得出解。

  • 不难发现上述过程中其实内存部分是可以经典实现的(索引也可以是经典的,只不过在 LOAD 使 LOAD 一个叠加态的索引)。说明只要量子地实现 CPU 和 LOAD 过程量子搜索算法也可以用于搜索经典数据库。如果维持单调字典序的话,一个经典的数据库可以在O(logN)O(logN) 时间内定位一个元。然后我们给出的方案需要O(NlogN)O(NlogN) 个量子开关,实际上和经典数据结构没啥区别,所以实际上,量子搜索算法的主要用途还是用于困难问题的解。数据库的话经典数据结构已经很成熟了。

# 搜索算法的最优性

  • 我们要好像吃饱了撑的证明不可能调用 oracle 少于O(N)O(\sqrt{N}) 次就在 N 元数据库里搜索,因此上面说的算法是最优的。

  • 设算法从状态ψ|\psi\rang 出发,考虑搜索问题只有唯一解的情况。我们被允许使用一个对解产生 (-1) 相移而其他状态不变的 oracle Ox=I2xxO_x=I-2|x\rang\lang x|。我们假设算法从ψ|\psi\rang 出发并恰好应用 k 次OxO_x,并在 oracle 运算之间穿插任意酉运算U1,U2,...,UkU_1,U_2,...,U_k,定义:

    ψkxUkOxUk1Ox....U1Oxψψk=UkUk1...U1ψ|\psi_k^x\rang\equiv U_kO_xU_{k-1}O_x....U_1O_x|\psi\rang\\ |\psi_k\rang=U_kU_{k-1}...U_1|\psi\rang

    我们的目的是估计量Dk=xψkxψk2D_k=\sum_x|||\psi_k^x\rang-|\psi_k\rang||^2,即同一个状态在酉运算的操作下,加不加 oracle 的差别是否很大。如果DkD_k 很小,则说明我们的 oracle 没起到什么效果,那么其实在搜索过程中,ψ|\psi\rang 也就以很小的概率最后转到离x|x\rang 很近。

    下面的证明过程省略狄拉克符号的 ket,用ψ\psi 表示向量ψ|\psi\rang


  • 首先我们证明Dk4k2D_k\leq 4k^2。显然这对 k=0 成立,其中Dk=0D_k=0,然后:

    Dk+1=xUk+1(Oxψkxψk)2=xOxψkxψk2=x(Ox(ψkxψk)+(OxI)ψk2b+c2b2+2bc+c2letb=Ox(ψkxψk),c=(OxI)ψk=2xψkxDk+1x(ψkxψk2+4ψkxψkxψk+4ψkx2)xxψk2=1Dk+1Dk+4+4xψkxψkxψkD_{k+1}=\sum_x||U_{k+1}(O_x\psi_k^x-\psi_k)||^2=\sum_x||O_x\psi_k^x-\psi_k||^2\\ =\sum_x||(O_x(\psi_k^x-\psi_k)+(O_x-I)\psi_k||^2\\ \because ||b+c||^2\leq||b||^2+2||b||||c||+||c||^2\\ let\quad b=O_x(\psi_k^x-\psi_k),c=(O_x-I)\psi_k=-2\lang x|\psi_k\rang|x\rang\\ \therefore D_{k+1}\leq\sum_x(||\psi_k^x-\psi_k||^2+4||\psi_k^x-\psi_k|||\lang x|\psi_k\rang|+4|\lang\psi_k|x\rang|^2)\\ \because \sum_x|\lang x|\psi_k\rang|^2=1\\ \therefore D_{k+1}\leq D_k+4+4\sum_x||\psi_k^x-\psi_k|||\lang x|\psi_k\rang|

    用柯西不等式D_{k+1}\leq D_k+4+4\sqrt{\sum_x{||\psi_k^x-\psi_k||^2\sum_x|\lang x|\psi_k\rang|^2}}=D_k+4+4\sqrt

    根据归纳假设Dk+14k2+4+8k=4(k+1)2D_{k+1}\leq 4k^2+4+8k=4(k+1)^2


  • 然后我们需要证明,事实上只有当DkD_kΩ(N)\Omega(N) 时,成功概率才是高的。

    假设要求对所有 x 都有xψkx21/2|\lang x|\psi_k^x\rang|^2\geq1/2,即要求无论解是什么,在搜索后我们的ψkx|\psi_k^x\rang 必须至少以1/21/2 的概率测量得到x|x\rang

    ψkxx2=22xψkx22||\psi_k^x-x||^2=2-2|\lang x|\psi_k^x\rang|\leq 2-\sqrt{2}

    定义Ekxψkxx2E_k\equiv\sum_x|||\psi_k^x-x||^2,则Ek(22)NE_k\leq(2-\sqrt{2})N,下面证明EkE_kΩ(N)\Omega(N) 的。

    定义Fkxxψk2F_k\equiv\sum_x||x-\psi_k||^2,有

    Dk=x(ψkxx)+(xψk)2xψkxx22xψkxxxψk+xxψk2=Ek+Fk2xψkxxxψkD_k=\sum_x||(\psi_k^x-x)+(x-\psi_k)||^2\\ \geq\sum_x||\psi_k^x-x||^2-2\sum_x||\psi_k^x-x||||x-\psi_k||+\sum_x||x-\psi_k||^2\\ =E_k+F_k-2\sum_x||\psi_k^x-x||||x-\psi_k||

    由柯西不等式:xψkxxxψkEkFk\sum_x||\psi_k^x-x||||x-\psi_k||\leq\sqrt{E_kF_k},于是

    DkEk+Fk2EkFk=(FkEk)2D_k\geq E_k+F_k-2\sqrt{E_kF_k}=(\sqrt{F_k}-\sqrt{E_k})^2

    而其实可以证明,Ek(22)N<2N,Fk2N2NE_k\leq(2-\sqrt{2})N<2N,F_k\geq 2N-2\sqrt{N},于是Dk>(222)2ND_k> (\sqrt{2}-\sqrt{2-\sqrt{2}})^2N

    所以对于任意常数c<(222)2c<(\sqrt{2}-\sqrt{2-\sqrt{2}})^2,都有Dk>cND_k>cN。又因为Dk4k2D_k\leq 4k^2,所以

    kcN/4k\geq\sqrt{cN/4}

    所以 oracle 调用次数 k 是Ω(N)\Omega(\sqrt{N}) 的我们才可以保证DkD_k 足够大,ψkx|\psi_k^x\rang1/2\geq 1/2 的概率得到解。


  • 然而这也证明了简单量子搜索的极限,即它的调用 oracle 次数仍然是O(2ω(n))O(2^{\omega(n)}) 的,其中ω(n)\omega(n) 是关于 n 的一个多项式。

  • 许多研究者相信,NP 完全问题苦难的本质在于它们搜索空间本质上是无结构的,并且解决这类问题的(精确到一个多项式因子)的最好方法是采用搜索的办法。然而这也只是片面之言,或许 NP 完全问题仍可能具有某些我们还未发现的结构。

  • 被普遍视为属于难度介于 P 和 NP 完全之间的 NPI 类的因子问题,量子计算机可以很好地求解其的关键在于利用了该问题的一个隐含结构,即可归约为求阶问题。然而经典算法至今也没有很好的求阶算法。还有许多被怀疑属于 NPI 类的问题(如图同构)甚至 NP 完全问题都被怀疑有类似结构。

# 黑箱算法的极限

给定一个函数 f,判断是否存在 x,使得f(x)=1f(x)=1

  • 实际上,这个问题的难度和找出一个 x 使得f(x)=1f(x)=1 是一样的。其实我们就是判断 Boolean 函数F(f)=f(0)...f(N1)F(f)=f(0)\vee...\vee f(N-1) 是否为 1。推广下,我们希望F(f)F(f) 是个任意 Boolean 函数,即里面的元为f(0),...,f(N1)f(0),...,f(N-1),运算符可以是与,或,非,异或等等。给定实现 f 的黑箱,如何快速确定 F 的值?复杂度分析基于 oracle 的调用次数。

  • 确定性查询复杂性D(F)D(F) 定义为经典计算机完成 F 的确定性计算,需要调用 oracle 的最少次数。

    相应的量子复杂性QE(F)Q_E(F) 是一台量子计算机为确定性计算 F 需要的最少 oracle 调用次数。因为量子计算机的概率输出特性,一个更有意义的量是有界误差性Q2(F)Q_2(F),即量子计算机以至少2/32/3 的概率输出 F 的值,需要最少调用 oracle 的次数。(2/32/3 只是任取的一个大于 0.5 的数,这也就可以通过重复来逼近 1)我们可以介绍一种多项式的方法计算 F。

  • 多项式方法基于表示 Boolean 函数的最小阶实多线性多项式的性质。下面要考虑的所有多项式是f(k){0,1}f(k)\in\{0,1\} 的函数,且因为f(k)2=f(k)f(k)^2=f(k),所以多项式是多线性的(即多项式中每个变量都是一次的,因为幂等),我们构造一个多项式:

    p(X)=Y{0,1}NF(Y)k=0N1[1(YkXk)2]p(X)=\sum_{Y\in\{0,1\}^N}F(Y)\prod_{k=0}^{N-1}[1-(Y_k-X_k)^2]

    其中,Xk=f(k),XX_k=f(k),X 表示函数的取值集合,Y 是一个长度为 N 的二进制串,即也可以表示为函数 f 的取值集合。上面右式表示当且仅当YkXkY_k\equiv X_k 才等于F(Y)=F(X)F(Y)=F(X)。于是我们构造出了个p:RNRp:R^N\rightarrow R 来表示 F,且X,p(X)=F(X)\forall X,p(X)=F(X)。实际上我们就是把 Boolean 函数用一般多项式表示。例如F(X)=X0X1,p(X)=1(X01)(X11),p(X)=F(X)F(X)=X_0\vee X_1,p(X)=1-(X_0-1)(X_1-1),p(X)=F(X)

  • 其实,Boolean 函数的多项式表示的最小次数是唯一的(上面的 p 不一定是最小次数)。譬如F(X)=(X0X1)X2F(X)=(X_0\vee X_1)\oplus X_2 的话,就有p(X)=(1(1X0)(1X1))+X2p(X)=(1-(1-X_0)(1-X_1))+X_2 为二次。这个可以证明但我不会(。可以感性理解一下,事实上位运算和加法,乘法是等价的。而且一个位运算对应次数翻倍。即有:

    ab=1(1a)(1b),ab=ab,ab=(a+b)(2ab)a\vee b=1-(1-a)(1-b),a\wedge b=ab,a\oplus b=(a+b)(2-a-b)

    想一想应该加法乘法域和位运算域是同构的。所以我们可以把加法乘法和位运算相互转换,例如就把与,或,异或转换为加法和乘法后的多项式次数是很固定的,一个与,或,非就代表次数需要左右相加。

  • 但由于这里 F 我们不知道里面有什么位运算和结构,所以我记其对应的多项式最小次数为deg(F)deg(F)。事实上,若 F 里面全是或,即把 N 个数或起来,则deg(OR)=Ndeg(OR)=N,同理deg(AND)=N,deg(XOR)=Ndeg(AND)=N,deg(XOR)=N

    然后已经证明了D(F)2deg(F)4D(F)\leq 2deg(F)^4(这个是事实,但相当难证明)。思考下,发现其实可以拓展p(X)p(X)。因为我们知道F(X)F(X) 的值一定是 0 或 1,所以如果X,p(X)F(X)1/3\forall X,|p(X)-F(X)|\leq 1/3,我们也可以通过p(X)p(X) 来得出F(X)F(X)。此时我们说pp 近似FF。用deg~(F)\widetilde{deg}(F) 表示这种近似多项式的最小次数,目前已经得到证明的有:

    deg~(OR)=Θ(N),deg~(AND)=Θ(N),deg~(XOR)=ND(F)216deg~(F)6\widetilde{deg}(OR)=\varTheta(\sqrt{N}),\widetilde{deg}(AND)=\varTheta(\sqrt{N}),\widetilde{deg}(XOR)=N\\ D(F)\leq 216\widetilde{deg}(F)^6

    这是目前已知最好结果,然后证明过程都略略略了。


  • 回到黑箱。我们本来想判断黑箱实现的函数 f 是否有解,然后转化为求F(X)F(X)。其实我们可以证明,F(X)F(X) 的近似多项式p(X)p(X) 可以用量子算法得到。

    记执行 T 次 oracle 调用的量子算法 Q 的输出为k=0N1ckk\sum_{k=0}^{N-1}c_k|k\rang,而 Q 实际上是任意一个算法,也就是任意作用了一些门和 oracle。其一般形式为:Q=U0OU1O...UT1OUTQ=U_0OU_1O...U_{T-1}OU_T,如下图:

    8

    考虑算法 Q 作用在前n+1n+1 个量子比特上,不失一般性还可以允许其引入一些工作量子比特(虽然不影响什么)。

    事实上,我们考虑任意一个状态ψin=i,jai0ji0j+ai1ji1j|\psi_{in}\rang=\sum_{i,j}a_{i0j}|i\rang|0\rang|j\rang+a_{i1j}|i\rang|1\rang|j\rang,它在经过一个 oracle 后变为:

    i,j[((1Xi)ai0j+Xiai1j)i0+((1Xi)ai1j+Xiai0j)i1]j\rightarrow \sum_{i,j}[((1-X_i)a_{i0j}+X_ia_{i1j})|i\rang|0\rang+((1-X_i)a_{i1j}+X_ia_{i0j})|i\rang|1\rang]|j\rang

    我们发现状态i0j,i1j|i0j\rang,|i1j\rang 前面的概率幅度关于XiX_i 的次数从 0 便变成了 1,实际上就是涨了个次数。而作用任何酉运算都不会改变概率幅度关于XiX_i 的次数(所以说引入什么量子比特都不重要,最后面那个j|j\rang 属实没啥贡献)

    于是经过了 T 次调用,概率幅度的阶已经为 T 了,而且正好为 T。我们考虑对算法输出k=0N1ckk\sum_{k=0}^{N-1}c_k|k\rang 进行测量,则概率Pk(X)=ck2P_k(X)=|c_k|^2 就是关于 X 的 2T 次实多项式。


    假设我们有个量子算法 Q 可以使得对算法输出k=0N1ckk\sum_{k=0}^{N-1}c_k|k\rang 进行测量可以正确地得出F(X)F(X) 的值(即输出不是0|0\rang 就是1|1\rang)此时F(X)F(X) 就等于得到状态0...01|0...01\rang 的概率,要么是 0 要么是 1。故其也可以表示为c12|c_1|^2,关于 X 的次数2T\leq 2T。所以有

    QE(F)=Tdeg(F)2Q_E(F)=T\geq\frac{deg(F)}{2}

    上式表示调用 oracle 次数。同样,若QQ 的输出结果是0,1|0\rang,|1\rang 的叠加态,即会以很小的误差概率测出错误的F(X)F(X),此时c12|c_1|^2F(X)F(X) 的近似,即有:

    Q2(F)deg~(F)2Q_2(F)\geq\frac{\widetilde{deg}(F)}{2}

    再结合D(F)2deg(F)4,D(F)216deg~(F)6D(F)\leq 2deg(F)^4,D(F)\leq 216\widetilde{deg}(F)^6,得到

    QE(F)D(F)321/4,Q2(F)D(F)138241/6Q_E(F)\geq\frac{D(F)}{32}^{1/4},Q_2(F)\geq\frac{D(F)}{13824}^{1/6}

    这意味着黑箱算法对于经典算法最多只能进行多项式加速,因为如果希望得出结果的概率>1/2>1/2 根据之前的推导D(F)(222)2ND(F)\geq(\sqrt{2}-\sqrt{2-\sqrt{2}})^2N 所以QE(F),Q2(F)Q_E(F),Q_2(F) 也至少是关于 N 多项式级别的,原因也在于常见的位运算函数(与或异或)的deg(F),deg~(F)deg(F),\widetilde{deg}(F)Ω(N)\Omega(N) 的。只有对实现 f 的 oracle 黑箱有更进一步的结构了解才有可能对经典算法进行指数加速。