量子概率采样

量子计算机运行结束时做的事只有一件：测量。因此它的天然输出不是某个数，而是一个随机样本——每次运行从输出分布中独立抽取一个比特串。本词条讨论的两类工作正是从这个事实出发，但方向恰好相反：

• 第一类（BosonSampling、IQP）利用采样这件”免费”的事：构造一个量子器件能直接产生、而任何经典算法都难以高效模仿的输出分布。这类结果是”量子计算优势”论证的主要形式之一。
• 第二类（log-concave sampling）加速采样这件”有用”的事：统计物理、贝叶斯推断、机器学习里处处需要从某个目标分布 $\pi$ 抽样，量子 walk 可以把某些经典 Markov chain 的混合时间开方加速。

两类结果不矛盾：前者故意挑选由 permanent、指数和这类 p-hard 量控制的干涉分布，后者则挑选有凸性、快速混合等良好结构的目标分布。“量子采样容易还是难”完全取决于分布本身长什么样。

预备知识方面，我们假设读者熟悉本词条前面的词条的内容：量子力学与量子线路基础、Grover 与振幅放大、相位估计，以及 Markov 链和凸分析的最基本概念（会在用到时简要复习）。

历史背景。Shor 算法证明了量子计算机能在判定/函数问题上指数加速，但因子分解的困难性本身只是一个（广为相信的）假设。2000 年代末，研究者开始寻找一条更直接的路线来证明”量子设备能做经典设备做不到的事”：不去解某个有用的问题，而是让量子设备做它最自然的事——从一个分布采样——然后论证任何经典算法都无法高效模仿这个输出分布。这条路线的优点是采样任务直接对应实验上”运行设备、记录测量结果”这一动作，门槛远低于通用容错量子计算。Shepherd 与 Bremner 2009 年的 IQP（Zoo 474）和 Aaronson 与 Arkhipov 2011 年的 BosonSampling（Zoo 473）是这条路线的两个奠基模型；此后多个光量子与超导实验沿着这类思路展示了采样优势。反过来，“采样”在计算科学中本来就是一个经典老大难问题（MCMC 是其标准工具），2022 年 Childs 等人（Zoo 475）的工作表明量子 walk 对其中结构良好的一类——log-concave 分布——能给出来路清楚的 polynomial 加速与近最优的精度标度。本词条把这两条线放在一起讲。

1. Sampling 与 probability estimation 是两件事

先给出贯穿全文的定义。

定义 1（采样问题）. 给定目标分布 $p$（通过某种隐式描述，例如一个量子线路、一个未归一化的密度），一个 sampler（采样器） 是一个随机算法，每次运行输出一个样本 $x\sim \tilde{p}$，要求输出分布满足

$$\Vert p-\tilde{p}{\Vert }_{TV}=\frac{1}{2}\sum _x|p(x)-\tilde{p}(x)|\le ϵ.$$

左边的量称为 total variation distance（全变差距离）：它等于”用 $\tilde{p}$ 冒充 $p$ 时，任何统计检验能察觉的最大优势”。$ϵ$ 越小，两个分布在操作意义上越不可区分。一个最小例子帮助建立量感：两枚硬币，$p$ 是均匀分布（正反面各 $\frac{1}{2}$），$\tilde{p}$ 正面概率 $\frac{1}{2}+\delta$，则 $\Vert p-\tilde{p}{\Vert }_{TV}=\delta$——这意味着单次抽样时，再聪明的观测者也只有至多 $\delta$ 的优势分辨用的是哪枚硬币；要可靠地区分，需要 $\Theta (1/{\delta }^2)$ 个独立样本。采样定理用 TV 距离而不是逐点误差，正是因为它刻画”实验者用有限统计能观察到的一切”。

注意这个定义里没有要求算法会计算 $p(x)$ 的数值。sampler 只需”掷出一枚偏向正确的骰子”，不需要知道每个面朝上的概率是多少。这个区别是本质的：

• 一个方向：目标分布的概率本身可能是 p-hard 的（下文会看到 BosonSampling 的概率正比于 permanent 的模方），不存在高效算法能逐个算出 $p(x)$；但物理量子过程照样每次运行”吐”出一个样本。大自然不计算概率，它只演化。
• 另一个方向：能对某些边缘概率做加性估计（例如用量子线路加采样去估 $|⟨z|U|0^n⟩{|}^2$ 到加性精度 $2^{-n}/\mathrm{poly}$），也不等于能从联合分布高效采样到 TV 精度 $ϵ$——后者要求对所有输出整体的统计保真。

Stockmeyer 计数：困难性论证的核心工具

“经典难采样”的定理几乎都是反证法，骨架如下。假设存在一个高效经典 sampler，输出分布与某量子器件的输出分布在 TV 距离 $ϵ$ 内一致。注意经典 sampler 本身是一台随机多项式时间机器：给它随机种子 $r$，它确定性地输出 $x=C(r)$。于是”输出 $x$ 的概率”就是”使 $C(r)=x$ 的种子所占比例”——这是一个计数问题。

Stockmeyer 近似计数定理（定性陈述）：给定一个多项式规模的布尔电路 $f$，统计满足 $f(r)=1$ 的输入个数可以在 ${\mathrm{BPP}}^{NP}$（带 NP oracle 的随机多项式时间）内做到乘法近似，即输出落在真值的 $(1\pm 1/\mathrm{poly})$ 倍以内。直观地说，NP oracle 允许我们用”哈希到随机子集再询问是否为空”的办法，把计数化归为一串 NP 判定问题。

于是链条闭合：经典 sampler 存在 $\Rightarrow$ 它的输出概率可在 ${\mathrm{BPP}}^{NP}$ 内乘法近似 $\Rightarrow$ 若这些概率编码了 p-hard 量（如 permanent），则 p-hard 问题落入多项式层级的低层 $\Rightarrow$ polynomial hierarchy 发生非预期的塌缩。最后一步值得说透：Toda 定理给出 $\mathrm{PH}\subseteq P^{\#P}$，即多项式层级的每一层都可以用”一次 p 计数”来解决；如果 p-hard 的量能反过来在 ${\mathrm{BPP}}^{NP}$（约当于 PH 的第三层以内）里被计算，那么 PH 的每一层都被压回这个低层，整个无限层级坍缩到有限层。这与"$\mathrm{P}=\mathrm{NP}$"一样属于被普遍否定的情形，因此反设的经典 sampler 不存在。注意论证并不要求经典 sampler 自己会”算概率”——它只要存在，其输出概率就自动成为 Stockmeyer 机制可以处理的计数对象。

这条链的强度取决于每一环的假设：exact sampler（TV 距离 $0$）的论证最干净；实验上现实的 approximate sampler（TV 距离 $ϵ>0$）需要额外的平均情形困难性与 anti-concentration 猜想（见第 3 节）。读者应始终把”采样困难”读作条件性命题，而不是无条件证明。

2. BosonSampling：从线性光学到 permanent

物理模型

BosonSampling 由 Aaronson 与 Arkhipov 于 2011 年提出（Zoo 473），问题是：一个只做线性变换的光学网络，加上单光子源和光子计数探测器，能不能做出经典计算机做不到的事？

装置有三个要素：

• $m$ 个光学模式（spatial mode），第 $j$ 个模式的产生算符记为 $a_j^{†}$。Fock 态 $|s_1,\dots ,s_m⟩$ 表示模式 $j$ 中有 $s_j$ 个光子。
• 一个无源线性光学网络（分束器与相移器的任意组合）。它对单光子 Hilbert 空间的作用是一个 $m\times m$ 酉矩阵 $U$，产生算符按

$$a_i^{†}\mapsto \sum _{j=1}^m{U}_{ji}{b}_j^{†}$$

变换，其中 $b_j^{†}$ 是输出模式的产生算符。注意：这只是对光子的单粒子态做了 $U$；当 $n>1$ 个光子同时在网络里时，多体演化由 $U$ 的所有”路径分配”相干叠加而成，这正是复杂性的来源。

• 输入 $n$ 个不可分辨的单光子（例如前 $n$ 个模式各一个），输出端用光子计数探测器读出占据数图样

$$S=(s_1,\dots ,s_m),\sum _{j=1}^m{s}_j=n.$$

输出概率的推导

我们来一步步算出输出图样 $S$ 的概率。为简单起见，设输入是前 $n$ 个模式各一个光子，即初态

$$|{\psi }_{\mathrm{in}}⟩=a_1^{†}{a}_2^{†}\cdots a_n^{†}|0⟩.$$

第一步，把线性变换逐因子代入：

$$|{\psi }_{\mathrm{out}}⟩=\prod _{i=1}^n\left(\sum _{j=1}^m{U}_{ji}{b}_j^{†}\right)|0⟩.$$

第二步，把乘积展开。展开式的每一项对应一个分配函数 $\sigma :\{1,\dots ,n\}\to \{1,\dots ,m\}$，含义是”第 $i$ 个输入光子走去了模式 $\sigma (i)$”：

$$|{\psi }_{\mathrm{out}}⟩=\sum _{\sigma }\left(\prod _{i=1}^n{U}_{\sigma (i),i}\right)b_{\sigma (1)}^{†}\cdots b_{\sigma (n)}^{†}|0⟩.$$

第三步，按输出图样归并。一个固定的输出图样 $S$ 对应许多不同的 $\sigma$：只要每个模式 $j$ 恰好被击中 $s_j$ 次。由于光子不可分辨，这些 $\sigma$ 给出的是同一个 Fock 态，振幅必须相干相加。先处理”无重复占据”（$s_j\in \{0,1\}$）的情形：此时击中 $S$ 中 $n$ 个模式的 $\sigma$ 恰好是 $n$ 个输入光子到 $n$ 个输出模式的一一对应，即置换 $\pi \in S_n$，于是

$$⟨S|{\psi }_{\mathrm{out}}⟩=\sum _{\pi \in S_n}\prod _{i=1}^n{U}_{S\text{ 第 }i\text{ 个被占模式},\text{第 }i\text{ 个输入}}=\mathrm{Per}(U_S),$$

其中 $U_S$ 是把 $U$ 按”输入取前 $n$ 列、输出取 $S$ 中被占据的行”抽出的 $n\times n$ 子矩阵，$\mathrm{Per}$ 是 permanent（定义见下）。

第四步，处理重复占据。若模式 $j$ 中有 $s_j>1$ 个光子，Fock 态的归一化带来因子 $\sqrt{{\prod }_j{s}_j!}$（因为 $(b_j^{†}{)}^{s_j}|0⟩=\sqrt{s_j!}|\dots ,s_j,\dots ⟩$）。把 $U_S$ 的定义推广为”第 $j$ 行重复 $s_j$ 次”的 $n\times n$ 矩阵后，同样的归并给出（输入若也有重复占据，分母再乘上输入的 factorial）：

$$\mathrm{Pr}(S)=\frac{|\mathrm{Per}(U_S){|}^2}{{\prod }_j{s}_j!}.$$

至此采样问题的规模也清楚了：满足 ${\sum }_j{s}_j=n$ 的图样共有 $\left(\genfrac{}{}{0ex}{}{n+m-1}{n}\right)$ 个，随 $n$ 指数增长。经典模拟的困境因此是双重的——既要面对指数大的输出空间，每个输出的振幅本身又是一个 p-hard 的 permanent。量子设备则绕开两者：它不枚举图样，只是让干涉发生，然后读出结果。

一个实现层面的注记：理论分析与实验设计通常取模式数远大于光子数（粗略地说 $m$ 随 $n$ 的超线性增长）。这样在典型输出中每个模式至多一个光子（“无碰撞”图样占多数），分母的 factorial 通常为 $1$，输出概率直接就是某个 $n\times n$ 子矩阵 permanent 的模方。这让困难性论证可以完全围绕 permanent 展开，而不必处理占据数修正。

Permanent：没有符号的行列式

矩阵 $A=(A_{ij}{)}_{n\times n}$ 的 permanent 定义为

$$\mathrm{Per}(A)=\sum _{\pi \in S_n}\prod _{i=1}^n{A}_{i,\pi (i)},$$

它与行列式逐项相同，唯一的差别是没有置换符号 $\mathrm{sgn}(\pi )$。这一点点差别后果巨大：

• 行列式里不同路径的贡献带符号、彼此相消，高斯消元利用这种结构把计算压到 $O(n^3)$；
• permanent 没有相消，也没有对应的乘法结构（$\mathrm{Per}(AB)\ne \mathrm{Per}(A)\mathrm{Per}(B)$），Valiant 的经典定理表明精确计算 $0$–$1$ 矩阵的 permanent 是 #P-complete 的，即与”数 SAT 的解的个数”同样难。

物理图像因此非常干净：$n$ 个不可分辨光子从输入到输出的每条”分配路径”贡献一个振幅，玻色统计让所有路径同号相加——光学网络是一台”permanent 的物理求和器”。费米子（自由费米子线路）则给出带符号的行列式，所以反而经典易模拟；采样困难的根源正是这种 bosonic 多路径相长。

小例子：两个光子的 Hong–Ou–Mandel 效应

取 $m=2$ 个模式、$50/50$ 分束器

$$U=\frac{1}{\sqrt{2}}\left(\begin{array}{cc}1& 1\\ 1& -1\end{array}\right),$$

输入 $|1,1⟩$（每个模式一个光子）。先用产生算符直接算：

$$|{\psi }_{\mathrm{out}}⟩=\frac{1}{2}(b_1^{†}+b_2^{†})(b_1^{†}-b_2^{†})|0⟩=\frac{1}{2}\left[(b_1^{†}{)}^2-b_1^{†}{b}_2^{†}+b_2^{†}{b}_1^{†}-(b_2^{†}{)}^2\right]|0⟩.$$

玻色产生算符彼此对易，中间两项 $-b_1^{†}{b}_2^{†}+b_2^{†}{b}_1^{†}$ 恰好相消；再用 $(b_j^{†}{)}^2|0⟩=\sqrt{2}|2{⟩}_j$ 得

$$|{\psi }_{\mathrm{out}}⟩=\frac{1}{\sqrt{2}}\left(|2,0⟩-|0,2⟩\right).$$

即两光子总是从同一个输出端口成对出来， coincidence（两端各一个）概率为零。这就是 Hong–Ou–Mandel 效应。

再用上面的 permanent 公式验证一遍。三个可能的输出图样：

• $S=(1,1)$：$U_S=U$，$\mathrm{Per}(U)=U_{11}{U}_{22}+U_{12}{U}_{21}=\frac{1}{\sqrt{2}}\cdot \frac{-1}{\sqrt{2}}+\frac{1}{\sqrt{2}}\cdot \frac{1}{\sqrt{2}}=-\frac{1}{2}+\frac{1}{2}=0$，故 $\mathrm{Pr}(1,1)=0$。与直接计算的相消一致。
• $S=(2,0)$：$U_S$ 的两行都是 $U$ 的第一行，四个元素全为 $\frac{1}{\sqrt{2}}$，$\mathrm{Per}(U_S)=2\cdot \frac{1}{2}=1$，分母 $s_1!=2$，故 $\mathrm{Pr}(2,0)=\frac{|1{|}^2}{2}=\frac{1}{2}$。
• $S=(0,2)$：同理得 $\frac{1}{2}$。

三者相加 $0+\frac{1}{2}+\frac{1}{2}=1$，概率归一。这个例子也显示：即使 $m,n$ 很小，输出分布也已经由路径间的相消干涉决定——这正是大规模时经典模拟跟不上的原因。

3. BosonSampling 困难性结论的条件与保留条款

上一节给出了物理模型；本节把”它很难被经典模拟”这句话的确切含义与前提写清楚。这是全文最容易被误读的部分，请特别留意每条保留条款。

Exact sampling 的情形

若存在一个经典多项式时间算法，能从某个 $n$ 光子、$m$ 模式线性光学网络的精确输出分布采样，则由第 1 节的 Stockmeyer 论证：该 sampler 的输出概率（即 $|\mathrm{Per}(U_S){|}^2$ 型量）可在 ${\mathrm{BPP}}^{NP}$ 中被近似处理，进而可用来近似 permanent 相关的 p-hard 量，导致 polynomial hierarchy 塌缩。因此：除非 PH 塌缩，否则不存在高效经典 exact sampler。注意这个结论甚至不需要任何额外猜想——它只用到 exact sampler 的概率能被 Stockmeyer 机制”抓住”这一点。

Approximate sampling 的情形

实验装置永远有噪声，真正的输出分布与理想分布总有正的 TV 距离。对”存在 TV 距离 $ϵ$ 内的经典近似 sampler”这一更强的假设，Aaronson–Arkhipov 的论证还要额外依赖两条至今未证的猜想：

• Gaussian permanent 平均情形困难性（Permanent-of-Gaussians Conjecture）：对矩阵元素独立取自（适当归一化的）复高斯分布的矩阵 $X$，把 $|\mathrm{Per}(X){|}^2$ 近似到乘法（或适当的加性）精度是 p-hard 的。注意已证的 p-hardness 是最坏情形结果，而采样反证需要平均情形版本——这就是必须引入猜想的地方。
• Permanent anti-concentration（PACC）：高斯矩阵的 permanent 不会以不可忽略的概率指数小。为什么需要它？sampler 配合 Stockmeyer 只能给出输出概率的加性小误差估计；如果典型概率本身就指数小，加性估计就毫无信息量，反证链在这里断裂。anti-concentration 保证典型概率足够大，使加性估计能反推出乘法信息。

噪声条款

定理的正确陈述形式大致是：“若光子损耗率、光子间部分可区分性、器件误差都足够小，使真实输出与理想输出的 TV 距离保持在允许范围内，则经典模拟仍然困难。“这些不是技术细节而是定理的一部分：若损耗过高（绝大多数光子丢失）或光子高度可区分，已知存在高效经典模拟算法。实验中”是否进入了经典难模拟的区间”，正是围绕这些条款的定量核对。

最后强调两点定位：

• 以上全部是条件性结论，没有人无条件证明过”任何经典算法永远不可能模拟 BosonSampling”；其价值在于把采样困难性化归到少数几条表述清晰、可被独立检验的猜想上。
• BosonSampling 不是通用量子计算模型：没有证据表明它能做通用线路模拟。它的意义恰恰是”受限到只剩线性光学”的物理系统仍可能拥有采样优势——这大幅降低了展示量子优势的实验门槛。

4. IQP：时间上无序的 commuting circuits

模型定义

IQP（Instantaneous Quantum Polytime）由 Shepherd 与 Bremner 于 2009 年提出（Zoo 474）。原始表述是”X-program”：只允许对 $|+⟩$ 态作用彼此对易的、在 $X$ 基下对角的门，最后在 $X$ 基测量。经一层 Hadamard 共轭换到计算基，等价的典型形式是

$$U_{IQP}=H^{\otimes n}D{H}^{\otimes n},$$

其中 $D$ 由彼此交换的 $Z$-diagonal 门（如 $Z$ 旋转、$ZZ$ 耦合、$CCZ$ 等）组成。因为 $D$ 中所有门两两对易，它们可以以任意顺序执行、原则上甚至同时执行——这就是”temporally unstructured / instantaneous”的含义：电路没有有意义的时间结构，复杂性全部来自全局干涉。

输出振幅的逐步推导

从 $|0^n⟩$ 出发逐层计算。第一层 Hadamard 产生均匀叠加：

$$H^{\otimes n}|0^n⟩=\frac{1}{2^{n/2}}\sum _{x\in \{0,1{\}}^n}|x⟩.$$

第二层是对角电路 $D$。设 $D|x⟩=e^{i\varphi (x)}|x⟩$（对角门的作用就是给每个基矢一个依赖 $x$ 的相位，$\varphi (x)$ 是各门贡献的相位之和，例如对 $ZZ$ 门有 $\varphi$ 含 $\theta x_i{x}_j$ 这样的项）：

$$D{H}^{\otimes n}|0^n⟩=\frac{1}{2^{n/2}}\sum _x{e}^{i\varphi (x)}|x⟩.$$

第三层 Hadamard 用标准恒等式 $H^{\otimes n}|x⟩=2^{-n/2}{\sum }_z(-1{)}^{x\cdot z}|z⟩$（逐比特做 $H|x_i⟩=\frac{1}{\sqrt{2}}{\sum }_{z_i}(-1{)}^{x_i{z}_i}|z_i⟩$ 再张量积即得，$x\cdot z={\sum }_i{x}_i{z}_i\mathrm{mod}2$）：

$$U_{IQP}|0^n⟩=\frac{1}{2^n}\sum _z\left(\sum _{x\in \{0,1{\}}^n}(-1{)}^{z\cdot x}{e}^{i\varphi (x)}\right)|z⟩.$$

于是输出振幅

$$⟨z|U_{IQP}|0^n⟩=2^{-n}\sum _{x\in \{0,1{\}}^n}(-1{)}^{z\cdot x}{e}^{i\varphi (x)}.$$

这个和是什么

上式是函数 $e^{i\varphi (x)}$ 的 Boolean Fourier 变换（在 $z$ 处的取值）。当 $\varphi (x)$ 是低次布尔多项式时，它也就是文献中常说的 Boolean polynomial exponential sum——与数论中指数和同构的组合对象，其精确计算通常是 p-hard 的。

另一个有用视角：若 $D$ 只含单比特 $Z$ 旋转与两比特 $ZZ$ 耦合，把 $x_i\in \{0,1\}$ 换成自旋变量 $s_i=1-2x_i\in \{\pm 1\}$，相位 $e^{i\varphi (x)}$ 就写成 $e^{-\beta H(s)}$ 的形状，其中 $H(s)={\sum }_i{h}_i{s}_i+{\sum }_{ij}{J}_{ij}{s}_i{s}_j$ 是经典 Ising 哈密顿量，只是耦合与”温度”成了复数。于是输出振幅正比于复参数 Ising 模型的配分函数。这把 IQP 与本词条其他教程（配分函数的量子算法）联系起来：同一个数学对象，在”实温度估计”与”复参数采样”两种提问方式下难度完全不同。

小例子：受控 S 门的输出分布

取 $n=2$，$D$ 只含一个受控 S 门（controlled-$S$）：$D|x_1{x}_2⟩=i^{x_1{x}_2}|x_1{x}_2⟩$，即 $\varphi (x)=\frac{\pi }{2}{x}_1{x}_2$。它是 $Z$-diagonal 的（$D=\mathrm{diag}(1,1,1,i)$），与自身当然对易，是一个合法的（极小的）IQP 电路。用上面推出的振幅公式逐个算四个输出：

$$⟨z|U_{IQP}|00⟩=\frac{1}{4}\left[1+(-1{)}^{z_1}+(-1{)}^{z_2}+i(-1{)}^{z_1+z_2}\right],$$

其中方括号内四项依次来自 $x=00,10,01,11$。代入：

• $z=00$：振幅 $\frac{1}{4}(3+i)$，概率 $\frac{|3+i{|}^2}{16}=\frac{10}{16}=\frac{5}{8}$；
• $z=01$：振幅 $\frac{1}{4}(1-i)$，概率 $\frac{2}{16}=\frac{1}{8}$；
• $z=10$：振幅 $\frac{1}{4}(1-i)$，概率 $\frac{1}{8}$；
• $z=11$：振幅 $\frac{1}{4}(-1+i)$，概率 $\frac{1}{8}$。

总概率 $\frac{5}{8}+\frac{1}{8}+\frac{1}{8}+\frac{1}{8}=1$，归一正确。注意输出不是均匀分布：$00$ 的概率是其他串的五倍。偏向完全来自 $x=11$ 那一项的相位 $i$ 与其他三项的干涉——改变 $D$ 中的门角，干涉图样随之改变。这就是”单一 commuting 层也能产生非平凡干涉分布”的最小实例；把它换成含非 Clifford 门角的大规模相互作用图，就得到困难性定理所讨论的线路族。

困难性结论及其保留条款

与 BosonSampling 平行的结论：若经典算法能以足够强的乘法或 TV 精度，对某一族 IQP 线路 weakly sample（只输出样本，不计算概率），则 Stockmeyer 机制加上相应的平均情形困难与 anti-concentration 假设会导致 polynomial hierarchy 塌缩。这里的术语值得停下来区分一下：复杂度文献把”算出每个输出概率”称为 strong simulation，把”仅从输出分布抽样”称为 weak simulation。strong simulation 高效蕴含 weak simulation 高效，反之不然——IQP（和 BosonSampling）的输出概率本身是 p-hard 的，所以 strong simulation 困难几乎没有悬念；真正不平凡的结论恰恰是更弱的 weak simulation 也困难，这才对应”实验设备跑出来的结果无法被经典复现”。保留条款同样重要：

• 结论依赖于具体的门族：门角取值、相互作用图的结构、TV 误差 $ϵ$ 的定义方式都影响论证是否成立，不能笼统地说”一切 commuting 电路都难模拟”。
• 反面也成立：若 $D$ 中所有门都是 Clifford（例如只含 $S$、$CZ$），则整个 $H^{\otimes n}D{H}^{\otimes n}$ 仍是 Clifford 电路，由 Gottesman–Knill 定理可经典高效模拟。“commuting”与”难模拟”之间没有等号，非 Clifford 相位是必需的。

5. Log-concave sampling 是一个算法问题

现在转换方向：不再问”量子采样能否难倒经典计算机”，而是问”量子力学能否帮我们更快地从有用的分布采样”。

问题设定

给定凸函数 $f:{\mathbb{R}}^d\to \mathbb{R}$（通过 oracle 访问，见下），目标是从密度

$$\pi (x)=\frac{e^{-f(x)}}{Z},Z={\int }_{{\mathbb{R}}^d}{e}^{-f(x)}dx$$

采样。这类密度称为 log-concave（对数凹），因为 $-\log \pi =f$ 是凸函数。它是连续世界里最核心的一类”良性”分布：高斯、限制在凸体上的均匀分布、logistic 后验、许多贝叶斯模型的后验都属于此类。归一化常数 $Z$ 通常没有闭式解，本身也是计算目标（统计物理的配分函数、贝叶斯的证据）；对一般的非凸 $f$，即使近似 $Z$ 也是 p-hard 的（与本词条配分函数教程的主题一致），凸性假设正是把问题拉回多项式时间的关键。

为什么这个问题值得专门研究？因为它是连续优化与积分的”原语”：体积估计、贝叶斯模型选择、统计物理 expectation 的计算，最终都化归为”从某个 log-concave 分布抽样”或”估计其归一化常数”。经典上这几乎是 MCMC 的同义词，而 MCMC 的混合时间正是整个管线的瓶颈——所以混合时间的任何 polynomial 改善都有广泛的下游影响。

输入模型：算法拿不到 $f$ 的解析式，只能通过 oracle 查询。最弱的是 zeroth-order oracle：给 $x$ 返回 $f(x)$；更强的是 first-order（返回 $\nabla f(x)$）。量子算法的一个看点正是：只用 zeroth-order 量子 oracle，也能达到使用梯度 oracle 的经典算法的查询复杂度（见第 6 节）。

光滑性与条件数

两个标准正则性参数：

• $f$ 是 $\mu$-strongly convex：对所有 $x,y$， $f(y)\ge f(x)+\nabla f(x)\cdot (y-x)+\frac{\mu }{2}\Vert y-x{\Vert }^2,$ 即函数处处不劣于一个曲率为 $\mu$ 的抛物面（$\mu >0$ 保证 $\pi$ 的尾部至少像高斯一样衰减）。
• $f$ 是 $L$-smooth：$\Vert \nabla f(x)-\nabla f(y)\Vert \le L\Vert x-y\Vert$，即梯度不会变化太快，等价地（对二次可微 $f$）Hessian 的最大特征值不超过 $L$。

定义条件数

$$\kappa =\frac{L}{\mu }.$$

它衡量目标分布的”各向异性程度”：$\kappa =1$ 时等高线是球面（各方向同样容易采样），$\kappa$ 大时分布在某个方向被拉得很扁，Markov chain 需要更多步才能探索到窄方向。由定义可见 $\kappa \ge 1$：在最小值点 $x^{\text{*}}$ 处（$\nabla f(x^{\text{*}})=0$），strong convexity 给出 $f(x)\ge f(x^{\text{*}})+\frac{\mu }{2}\Vert x-x^{\text{*}}{\Vert }^2$，而 $L$-smoothness 的标准推论（沿梯度流的二次上界）给出 $f(x)\le f(x^{\text{*}})+\frac{L}{2}\Vert x-x^{\text{*}}{\Vert }^2$，两式比较即得 $\mu \le L$。所有经典与量子算法的复杂度都以 $\kappa$ 和 $d$ 的多项式出现，问题只在于次数。

经典基准：Langevin diffusion

从 log-concave 分布采样的经典工作马是 Langevin 动力学——模拟粒子在势场 $f$ 中的过阻尼布朗运动：

$$d{X}_t=-\nabla f(X_t)dt+\sqrt{2}d{B}_t,$$

其中 $B_t$ 是 $d$ 维标准布朗运动。为什么它以 $\pi$ 为稳态？把密度演化写成 Fokker–Planck 方程

$${\partial }_t\rho =\nabla \cdot (\rho \nabla f)+\Delta \rho .$$

代入 $\rho =e^{-f}$ 验证右边为零：第一项保持不动；第二项用链式法则，$\nabla e^{-f}=-e^{-f}\nabla f$，所以

$$\Delta e^{-f}=\nabla \cdot (\nabla e^{-f})=\nabla \cdot (-e^{-f}\nabla f)=-\nabla \cdot (e^{-f}\nabla f),$$

与第一项恰好相消，${\partial }_t\rho =0$。归一化不影响这个计算，故 $\pi$ 确为平稳密度。物理直觉：漂移项 $-\nabla f$ 把粒子推向低势能区，扩散项 $\sqrt{2}d{B}_t$ 把粒子往外推，二者在密度 $e^{-f}$ 处达到细致平衡。

实际算法是对该扩散做离散化（如 unadjusted Langevin algorithm，$X_{k+1}=X_k-h\nabla f(X_k)+\sqrt{2h}{\xi }_k$），这引入离散偏差；常用 Metropolis–Hastings 接受/拒绝步修正（即 MALA），再配合 warm start（从一个与 $\pi$ 不太远的初始分布出发）控制混合时间。对 $\mu$-strongly log-concave、$L$-smooth 的目标，经典 MALA 类算法的查询复杂度是 $\kappa$ 与 $d$ 的低次多项式（典型依赖形如 $\kappa \cdot d$ 乘以对数因子）——已经是多项式，但 $\kappa$ 和 $ϵ$ 上的次数正是量子算法要改善的对象。

离散化的误差预算。把连续扩散换成有限步长 $h$ 的迭代会引入偏差：步长太大则离散轨迹偏离连续扩散太远（偏差随 $h$ 增长），步长太小则需要更多步才能覆盖同样的物理时间（步数随 $1/h$ 增长）。MALA 的 Metropolis 修正以每次一步额外的 $f$（或 $\nabla f$）求值为代价，把这一步偏差精确吸收进接受概率，使离散链的平稳分布严格等于 $\pi$ 而不是某个近似——这是”corrected”与”unadjusted”版本的本质区别，也是量子算法选择相干化 MALA 核的原因：量子 walk 的开方加速作用于谱隙，前提是核的平稳分布必须是对的。

6. Quantum Langevin / MALA 的加速来源

Childs、Li、Liu、Wang 与 Zhang（Zoo 475）给出了 log-concave 采样与归一化常数估计的量子算法。加速由三个部件叠加而成，逐个解释。

部件一：相干化的 Markov 链与 Szegedy walk

经典 MCMC 的混合速度由转移核 $P$ 的谱隙 $\delta$ 控制：$P$ 是（关于 $\pi$）可逆的，特征值为 $1={\lambda }_0>{\lambda }_1\ge \cdots$，$\delta =1-{\lambda }_1$；混合时间大致为 $\tau \approx \frac{1}{\delta }\log \frac{1}{ϵ{\pi }_{\min }}$，即与谱隙成反比。谱隙小意味着链的态空间里存在”瓶颈”，随机游走要很久才能穿过。

Szegedy 量子 walk 的构造把任意可逆 Markov 核 $P$ 提升为一个酉算子 $W(P)$（作用在两份状态寄存器上，可看作 $P$ 的”量子化”），其关键性质是：$P$ 的谱隙 $\delta$ 变成 $W(P)$ 的相位隙 $\Theta (\sqrt{\delta })$。相位估计（读者在学过）分辨相位 $\theta$ 的代价是 $O(1/\theta )$，于是”制备平稳分布样本”的代价从 $O(1/\delta )$ 降为 $O(1/\sqrt{\delta })$——混合时间的平方根。这与 Grover 把 $\frac{1}{p}$ 次尝试开方为 $\frac{1}{\sqrt{p}}$ 是同一个”振幅比概率开方”的现象；事实上 Grover 正是把”均匀抽样直到命中”这条平凡 Markov 链量子化的特例，Szegedy walk 是它在任意可逆链上的推广。

量子算法的做法是：把 Langevin/MALA 的 proposal 与 accept/reject 核相干实现（每一步是一次对 $f$ 的量子查询加 d 维算术），再对其 Szegedy walk 做相位估计式的高效实现。对 strongly log-concave 目标，MALA 核的谱隙含因子 $1/\kappa$（扁分布的瓶颈），开方后 $\kappa$ 的依赖从 $\kappa$ 降为 ${\kappa }^{1/2}$。

部件二：从 zeroth-order oracle 相干提取梯度

漂移项需要 $\nabla f$，但最弱的输入模型只给 $f(x)$ 的量子 oracle ${\mathcal{O}}_f|x⟩|0⟩=|x⟩|f(x)⟩$。经典世界里，用函数值估计 $d$ 维梯度需要 $\Theta (d)$ 次查询（逐坐标有限差分）。量子世界里有一个叠加技巧（Jordan 梯度估计的思想）：在网格点的均匀叠加上调用相位化的 oracle，

$$\sum _x|x⟩⟼\sum _x{e}^{i\lambda f(x)}|x⟩,$$

其中 $\lambda$ 是适当选取的缩放因子。在光滑性保证的局部区域里，$f(x)\approx f(x_0)+\nabla f(x_0)\cdot (x-x_0)$（一阶泰勒展开，误差由 $L$-smoothness 控制），于是

$$e^{i\lambda f(x)}\approx e^{i\lambda f(x_0)}{e}^{i\lambda \nabla f(x_0)\cdot (x-x_0)}:$$

相位是 $x$ 的线性函数，其”频率”正是 $\lambda \nabla f(x_0)$。对该叠加做逆 QFT，测量即得梯度各分量的数字编码——$O(1)$ 次（而非 $d$ 次）oracle 查询拿到整个梯度向量（到相应精度）。关于这种”相位 kickback + Fourier 读出”的机制，读者可对照本词条梯度估计一章。结论：zeroth-order 量子查询在构造 drift 这一项上不吃亏，量子算法的查询复杂度得以匹配使用梯度 oracle 的经典算法。

部件三：振幅估计用于均值

估计归一化常数时需要估计许多形如 ${\mathbb{E}}_{x\sim {\pi }_k}[g_k(x)]$ 的期望值（见下文）。经典 Monte Carlo 用 $N$ 个样本估均值，误差 $\sim \sigma /\sqrt{N}$，达到误差 $ϵ$ 需 $O(1/{ϵ}^2)$ 个样本；量子 mean/amplitude estimation（振幅估计，相位估计对 Grover 型算子的应用）把同样的精度做到 $O(1/ϵ)$ 次调用——又是平方根级的 $ϵ$ 改善。

复杂度账目

把三个部件合起来（省略 polylog 因子与 warmness 参数，记为 $\tilde{O}$）：

采样：生成一个与 $\pi$ 的 TV 距离为 $ϵ$ 的样本需要

$$\tilde{O}({\kappa }^{1/2}d)$$

次对 $f$ 的查询。逐项读这个表达式：

• ${\kappa }^{1/2}$：来自 Szegedy walk 对谱隙的开方——经典链的谱隙含 $1/\kappa$ 因子，相位隙为其平方根；
• $d$：每步 walk 需要更新 $d$ 维坐标（对 $f$ 的查询、d 维向量算术），维度依赖以近线性的形式进入；
• $\tilde{O}$ 吸收了离散化步长、warm start 质量、目标 TV 精度 $ϵ$ 带来的 polylog 因子。

归一化常数：估计 $Z$ 到乘法精度 $1\pm ϵ$ 需要

$$\tilde{O}\left(\frac{{\kappa }^{1/2}{d}^{3/2}}{ϵ}\right)$$

次查询。框架是经典的退火/telescoping：取一列容易到难分布 $f_0,f_1,\dots ,f_K=f$（$f_0$ 取 $Z_0$ 有闭式解的良态二次型，如适当缩放的高斯），把 $Z$ 写成比值连乘

$$\frac{Z}{Z_0}=\prod _{k=0}^{K-1}\frac{Z_{k+1}}{Z_k},\frac{Z_{k+1}}{Z_k}=\int \frac{e^{-f_{k+1}}}{Z_k}dx={\mathbb{E}}_{x\sim {\pi }_k}\left[e^{-(f_{k+1}(x)-f_k(x))}\right].$$

关键观察：每个比值是 ${\pi }_k$ 下某个有界随机变量的期望——正好是部件三能处理的对象。粗略账目：schedule 长度 $K$ 与每个比值所需样本数共同贡献约 $d^{1/2}$ 量级的因子（并把误差预算 $1/ϵ$ 摊到 $K$ 个比值上只会再产生 polylog 因子）；每个样本的制备成本是上面的 $\tilde{O}({\kappa }^{1/2}d)$；振幅估计给出 $1/ϵ$ 而非 Monte Carlo 的 $1/{ϵ}^2$。三者相乘即为 ${\kappa }^{1/2}{d}^{3/2}/ϵ$ 的形状（polylog 吸收进 $\tilde{O}$）。

两点关于定理适用范围的说明。其一，上述复杂度都隐藏了 warmness 参数：算法需要一个与目标分布不太差的初始分布（warm start），其质量以对数或低次多项式因子进入；对 strongly log-concave 目标，以最小值点为中心的高斯通常是现成的 warm start。其二，$\mu$、$L$ 需已知或有可靠估计——它们决定离散步长、schedule 长度与误差预算的分配；把这些参数当成黑箱自适应估计会带来额外开销，不在上述表达式内。

小例子：一维高斯的 telescoping。取 $f_{\beta }(x)=\frac{\beta }{2}{x}^2$，则 ${\pi }_{\beta }$ 是方差 $1/\beta$ 的高斯，归一化常数有闭式解 $Z(\beta )=\int e^{-\beta x^2/2}dx=\sqrt{2\pi /\beta }$。取 schedule ${\beta }_0=1>{\beta }_1>\cdots >{\beta }_K$（每步乘以一个小于 $1$ 的常数），第 $k$ 个比值是

$$\frac{Z({\beta }_{k+1})}{Z({\beta }_k)}={\mathbb{E}}_{x\sim {\pi }_{{\beta }_k}}\left[e^{-({\beta }_{k+1}-{\beta }_k)x^2/2}\right].$$

这个期望可以手算验证：代入高斯密度，

$${\mathbb{E}}_{{\pi }_{{\beta }_k}}\left[e^{-({\beta }_{k+1}-{\beta }_k)x^2/2}\right]=\frac{1}{Z({\beta }_k)}\int e^{-{\beta }_{k+1}{x}^2/2}dx=\frac{Z({\beta }_{k+1})}{Z({\beta }_k)}=\sqrt{\frac{{\beta }_k}{{\beta }_{k+1}}},$$

即被估计量是一个量级为 $O(1)$ 的有界比值（被积函数取值在 $(0,1]$）。这正是振幅估计适用的场景：每个比值用 $\tilde{O}(1/ϵ)$ 次对 ${\pi }_{{\beta }_k}$ 的采样制备调用估到相对精度，$K$ 个比值相乘的误差只需把单步误差预算除以 $K$（代价被吸收进 polylog）。经典 Monte Carlo 做同一件事，每个比值要 $O(1/{ϵ}^2)$ 个样本——$1/ϵ$ 与 $1/{ϵ}^2$ 的差距在最终的复杂度表达式里原样保留下来。

近最优性：同一工作还证明了 $ϵ$ 维度的查询下界 $1/{ϵ}^{1-o(1)}$，即任何量子算法至少要以近 $1/ϵ$ 的标度随精度增长。因此上述上界在 $ϵ$ 依赖上近最优——平方根加速是这个问题的信息论极限，不可能再有”另一个平方根”。

7. 输出与验证

本节的方法论与本词条另一篇教程 quantum-simulated-annealing.md 互补：那里讨论的是用退火 schedule 寻找基态（优化），这里讨论的是从固定（或缓变）的平衡分布中抽样（积分/计数）。二者共享”量子 walk 开方谱隙”这一核心机制，复杂度分析的套路也相同。

两个实践层面的提醒。

输出是经典样本。Log-concave sampler 的最后一步是测量位置寄存器，得到一个经典向量 $x$。若应用需要 $M$ 个独立样本，就要把状态制备重复 $M$ 次——量子叠加不能”一次打印指数多个样本”，测量会摧毁叠加。量子算法的优势体现在每个样本的制备成本上，而不是样本数量的并行化。

验证与保证不对称。理论上的 TV 距离保证依赖一串前提：oracle 实现正确、$\mu$ 与 $L$ 的估计准确、离散化与相位估计的误差预算分配正确。反过来，给定一个声称的 sampler 和有限多样本，要认证其输出确实接近某个高维分布本身是统计上困难的任务（高维 goodness-of-fit 的样本复杂度随维数爆炸）。实验中可行的是检查低阶矩、能量分布等低维统计量是否匹配——这是必要而非充分条件。这也是为什么第 3、4 节的困难性定理要把”TV 距离 $ϵ$ 以内”写进假设，而不是声称认证了某个具体装置。

最后再强调一次两类结果的关系：BosonSampling/IQP 证明”采样”这种量子力学白送的能力，在刻意构造的干涉分布上经典计算机追不上；log-concave 结果证明量子 walk 能把有良好凸结构的经典采样问题做得更快。前者靠 permanent/指数和的 p-hard 性，后者靠凸性带来的大谱隙——看似相反，实则统一于”分布的结构决定采样的难度”。

8. 小结与习题

• Sampling 只要求输出样本的分布接近目标（TV 距离），不要求计算任何 $p(x)$ 的数值；这使”概率 p-hard”与”物理可采样”可以并存。
• 采样困难性定理是反证链：经典 sampler + Stockmeyer 计数 $\Rightarrow$ p-hard 量落入 PH 低层 $\Rightarrow$ PH 塌缩。exact 情形无条件；approximate 情形依赖平均情形困难与 anti-concentration 猜想。
• BosonSampling 的输出概率是 permanent 模方除以占据数 factorial；bosonic 多路径相长是物理本质，两个光子的 Hong–Ou–Mandel 效应是最小实例。
• IQP 线路 $H^{\otimes n}D{H}^{\otimes n}$ 的输出振幅是 Boolean exponential sum，也可写成复参数 Ising 配分函数；困难性依赖具体门族与误差模型，Clifford 子族可经典模拟。
• Log-concave 量子采样 = 相干 MALA 核 + Szegedy walk 开方谱隙 + 相位 kickback 提取梯度 + 振幅估计加速均值；复杂度显式依赖 $\kappa ,d,ϵ$，其中 $ϵ$ 依赖近最优。

习题（按难度大致递增；前四题有正文可对照，后两题需要综合多节内容）：

1. （推广 HOM 例子）把第 2 节的 $50/50$ 分束器换成一般分束器 $U=\left(\begin{array}{cc}\cos \theta & \sin \theta \\ \sin \theta & -\cos \theta \end{array}\right)$，输入仍为 $|1,1⟩$。用 permanent 公式证明 coincidence 概率 $\mathrm{Pr}(1,1)={\cos }^{2}2\theta$，并检查 $\theta =\pi /4$ 时回到正文结果。这个公式说明”干涉相消”如何随分束比连续变化。
2. （IQP 推导）设 $D={\prod }_j\exp (i{\theta }_j{Z}_j){\prod }_{(j,k)}\exp (i{\varphi }_{jk}{Z}_j{Z}_k)$，写出 $\varphi (x)$ 的显式形式，并完整重现第 4 节从 $|0^n⟩$ 到振幅公式 $2^{-n}{\sum }_x(-1{)}^{z\cdot x}{e}^{i\varphi (x)}$ 的三步推导，指出每一步用到的恒等式。
3. （Stockmeyer 的角色）第 1 节的反证链中，为什么必须假设经典 sampler 是”多项式时间随机算法”而不能是任意黑箱？如果把 Stockmeyer 计数换成精确的 $\#P$ oracle，论证会在哪一步变得平凡、从而失去意义？
4. （条件数）证明 $\kappa =L/\mu \ge 1$（提示：在 $f$ 的最小值点比较 strong convexity 下界与 smoothness 二次上界），并举一个 $\kappa$ 大的二维高斯例子，直观解释为什么条件数大时 Langevin 链混合慢。
5. （Fokker–Planck 验证）对一维情形 $f(x)=\frac{1}{2}{x}^2$，直接计算 $\nabla \cdot (\rho \nabla f)+\Delta \rho$ 在 $\rho =e^{-x^2/2}$ 处为零，确认正文第 5 节稳态验证的每一步；再说明为什么把 $\rho$ 换成 $e^{-x^2/4}$ 就不再平稳。
6. （精度标度）经典 Monte Carlo 估均值达到加性误差 $ϵ$ 需要 $O(1/{ϵ}^2)$ 样本，量子振幅估计需要 $O(1/ϵ)$ 次调用。用量子线路输出伯努利变量的例子说明：相位估计把”统计误差 $1/\sqrt{N}$“换成”旋转分辨误差 $1/N$“，从而解释归一化常数复杂度中 $1/ϵ$ 因子的来源；并结合 $1/{ϵ}^{1-o(1)}$ 下界说明为什么不可能再改善为 $1/\sqrt{ϵ}$。

参考文献

• Zoo 编号 473：Aaronson 与 Arkhipov, The Computational Complexity of Linear Optics.
• Zoo 编号 474：Shepherd 与 Bremner, Temporally Unstructured Quantum Computation.
• Zoo 编号 475：Childs、Li、Liu、Wang 与 Zhang, Quantum Algorithms for Sampling Log-Concave Distributions and Estimating Normalizing Constants.

---

参考资料

• 本词条整理自《量子计算算法教程》原文：ch13-topology-statistical-physics/probabilistic-sampling
• 内容遵循 CC BY-NC-SA 4.0 许可协议