稀疏矩阵幂对角元

给一个维数 $N$ 可能指数大的稀疏实对称矩阵 $A$，我们只关心它的幂的一个对角元

$$(A^m{)}_{jj}=⟨j|A^m|j⟩.$$

这看起来像是矩阵计算里最小的局部问题——不读整张矩阵、不输出整行整列，只要一个数。但本词条将看到：在适当的归一化口径下，这个”只算一个数”的问题精确刻画了量子计算的全部能力——它是 PromiseBQP-complete 的。一方面，量子计算机可以从基态 $|j⟩$ 的谱测度 (spectral measure) 中采样，把该对角元解释为这个测度的第 $m$ 阶矩来估计；另一方面，Feynman clock 构造可以把任意量子电路的接受概率编码进某个稀疏矩阵的谱矩，因此任何能解该问题的经典算法都能模拟通用量子计算。

这个结果由 Janzing 与 Wocjan 于 2006 年证明（见文末参考文献，Zoo 编号 60）。它与字符串重写一节的结果互为表里：两者都把”局部转移规则生成的隐式稀疏矩阵”作为载体，把通用量子计算编码进高阶谱矩。

1. 问题陈述与输入模型

1.1 稀疏矩阵预言机

设 $N=2^n$，即矩阵的行（列）用 $n$ 个比特索引。输入是一个实对称矩阵 $A\in {\mathbb{R}}^{N\times N}$，满足稀疏性假设：每行至多有

$$s=\mathrm{polylog}N=\mathrm{poly}(n)$$

个非零矩阵元。由于 $N$ 指数大，$A$ 不可能以显式形式给出；我们通过 sparse oracle 访问它。预言机回答两类查询，每次查询的代价是 $\mathrm{poly}(n)$：

1. 位置查询：给定行号 $r$ 与序号 $\ell \in \{1,\dots ,s\}$，返回第 $r$ 行第 $\ell$ 个非零元的列号 $c$；
2. 取值查询：给定 $(r,c)$，返回实数值 $A_{r,c}$。

这就是 Hamiltonian simulation 文献中的标准稀疏访问模型（参见量子模拟一章同款假设—— sparse Hamiltonian simulation 用的正是这样的预言机）。此外已知谱范数的上界

$$\Vert A\Vert \le b,$$

其中 $\Vert A\Vert ={\max }_r|{\lambda }_r(A)|$ 是算子范数，$b=\mathrm{poly}(n)$。

1.2 要估计什么：归一化是问题定义的核心

问题（稀疏矩阵幂对角元估计）. 给定上述稀疏预言机、行号 $j$、幂次

$$m=\mathrm{polylog}N,$$

要求以加性误差

$$ϵb^m,1/ϵ=\mathrm{polylog}N$$

估计 $(A^m{)}_{jj}$，即输出一个数 $\tilde{\mu }$，使得 $|\tilde{\mu }-(A^m{)}_{jj}|\le ϵb^m$。

为什么误差按 $b^m$ 缩放？因为 $b^m$ 是这个量的自然尺度。事实上，由算子范数的次可乘性（$\Vert XY\Vert \le \Vert X\Vert \Vert Y\Vert$，对 $m$ 归纳得 $\Vert A^m\Vert \le \Vert A{\Vert }^m$）与算子范数的定义（$|⟨j|M|j⟩|\le \Vert M\Vert$，因为单位矢量的像的范数不超过 $\Vert M\Vert$），

$$\left|(A^m{)}_{jj}\right|=\left|⟨j|A^m|j⟩\right|\le \Vert A^m\Vert \le \Vert A{\Vert }^m\le b^m.$$

也就是说，$(A^m{)}_{jj}$ 本身可以达到 $b^m$ 量级。除以 $b^m$ 之后，问题等价于估计

$${\mu }_m=⟨j|H^m|j⟩,H:=A/b,\Vert H\Vert \le 1,$$

到加性误差 $ϵ$。下文一律在归一化矩阵 $H$ 上工作，最后再把因子 $b^m$ 乘回去。

这个归一化是问题定义的核心，而不是技术细节。 若改为要求”原值的绝对误差 $1$“，而 $b^m$ 随 $n$ 指数大，那么我们是在要求相对误差 $b^{-m}$——复杂度的口径完全不同，量子算法不再高效。本词条（以及其余各节）反复强调的原则在这里同样适用：一个量的可估计性永远要连同它的归一化一起陈述。

1.3 动机：谱矩在哪里出现

矩阵幂的对角元绝不是人为造出来的问题。下面三个例子说明它的来源：

• 路径计数。若 $A$ 是某个稀疏图的邻接矩阵，$(A^m{)}_{jj}$ 是从顶点 $j$ 出发、恰好 $m$ 步回到 $j$ 的闭合路径条数（我们将在 §2.1 证明这一点）。更一般地，带权矩阵给出带权路径和。统计物理中的配分函数、组合学中的生成函数都涉及这类计数。
• 谱测度与矩问题。数值线性代数中，估计大规模对称矩阵的谱密度（例如凝聚态物理中的态密度）的标准做法正是通过 $⟨v|f(A)|v⟩$ 型的量提取谱信息；取 $f(x)=x^m$ 就是谱矩。
• 量子计算自身。§7 将看到，任意量子电路的演化历史可以被写成一个稀疏 Hermitian 矩阵，其幂的对角元编码了电路的接受概率——这是 BQP-hardness 的来源。

2. 经典算法的瓶颈

在给出量子算法之前，先弄清楚经典算法卡在哪里。这能帮我们理解量子加速究竟来自哪一步。

2.1 路径展开：矩阵幂的组合含义

引理（路径展开恒等式）. 对任意矩阵 $A$ 与任意正整数 $m$，

$$(A^m{)}_{jj}=\sum _{i_1,\dots ,i_{m-1}}{A}_{j,i_1}{A}_{i_1,i_2}\cdots A_{i_{m-1},j}.$$

若 $A$ 是图 $G$ 的邻接矩阵（$A_{ik}=1$ 当且仅当 $(i,k)$ 是一条边），则 $(A^m{)}_{jj}$ 等于从 $j$ 出发长度为 $m$ 的闭合路径（closed walk）条数。

证明。$m=1$ 时右端没有求和指标，就是 $A_{jj}$ 本身。设公式对 $m$ 成立，则由矩阵乘法的定义 $(XY{)}_{ik}={\sum }_l{X}_{il}{Y}_{lk}$，

$$(A^{m+1}{)}_{jj}=\sum _{i_1}{A}_{j,i_1}(A^m{)}_{i_1,j}=\sum _{i_1}{A}_{j,i_1}\sum _{i_2,\dots ,i_m}{A}_{i_1,i_2}\cdots A_{i_m,j},$$

把两个求和合并即得 $m+1$ 的情形，归纳完成。对邻接矩阵，乘积 $A_{j,i_1}{A}_{i_1,i_2}\cdots A_{i_{m-1},j}$ 等于 $1$ 当且仅当 $(j,i_1,\dots ,i_{m-1},j)$ 是一条合法的闭合路径（每条边都存在），否则为 $0$；求和就是数路径。Q.E.D.

这个恒等式给我们两个重要信息：

• 非零项的数量。从 $j$ 出发，每一步至多有 $s$ 个非零元可选（稀疏性），所以长度为 $m$ 的路径至多有 $s^m$ 条，即求和中至多有 $s^m$ 个非零项。
• 与干涉的对应。量子力学中”振幅 = 路径振幅之和”的 Feynman 图像在这里完全精确：$(A^m{)}_{jj}$ 是所有闭合路径贡献的代数和。§7 的 Feynman clock 构造正是把量子电路的演化历史嵌入这些路径。

2.2 直接计算与快速幂为什么都不行

策略一：逐条枚举路径。 由 §2.1，用深度优先搜索枚举所有从 $j$ 出发的长度 $m$ 路径并累加贡献，代价是 $O(s^m)$ 次预言机查询。代入问题的参数尺度：

$$s^m=(\mathrm{polylog}N{)}^{\mathrm{polylog}N}=2^{\mathrm{polylog}N\cdot \log \mathrm{polylog}N},$$

这是关于 $n=\log N$ 的超多项式（quasi-polynomial 以上）增长。$m$ 哪怕只是 $n^2$，$s^m$ 也已经是 $2^{\Theta (n^2\log n)}$。路径枚举不可用。

策略二：把 $A^m$ 整个算出来。 这更糟：$A^m$ 通常是稠密 (dense) 的，有 $\Theta (N^2)$ 个非零矩阵元，仅仅把输出写下来就需要 $\Omega (N^2)$ 时间，而 $N=2^n$。

策略三：快速幂（repeated squaring）。 经典计算中，算数的 $m$ 次幂只需 $O(\log m)$ 次乘法。对矩阵同样可以做 $A\to A^2\to A^4\to \cdots$，只需 $O(\log m)$ 次矩阵乘法。但这里有两个致命伤：

1. 稀疏性在乘方后不保持（fill-in）。若 $A$ 每行至多 $s$ 个非零元，则 $(A^2{)}_{ik}={\sum }_l{A}_{il}{A}_{lk}$ 中非零项出现在”存在共同邻居 $l$“的位置，每行至多 $s^2$ 个；归纳地看，$A^{2^k}$ 每行的非零元至多 $s^{2^k}$ 个。平方 $\log m$ 次之后上界已是 $s^m$——矩阵早已稠密，稀疏预言机的优势荡然无存。
2. 维数不可承受。快速幂把乘法次数降到 $O(\log m)$，但每次乘法的对象是 $N\times N$ 矩阵，$N=2^n$。即使最朴素的稠密乘法也是 $\Theta (N^3)$，任何已知乘法都不可能低于输出规模的 $\Omega (N^2)$。

结论：经典快速幂优化的是乘法次数，而我们的瓶颈是矩阵维数。两者正交。

策略四：Monte Carlo 采样路径。 若矩阵元全部非负（例如随机游走的转移矩阵），可以从 $j$ 出发随机游走 $m$ 步，用”回到 $j$ 的频率”估计对角元——这正是 Markov chain Monte Carlo 的思想，经典上是可行的。但一般的实对称矩阵有正有负，路径贡献带符号：此时估计量是大量正负数相消后的小差值，方差与所求量的比值可以指数大——这就是著名的符号问题 (sign problem)。量子算法的相干干涉恰好擅长处理这种带符号的求和；这也是 BQP-hardness 结果（§7）告诉我们”一般情形经典难”并不矛盾于”非负特例经典易”的原因。

3. 谱测度：把对角元看成矩

量子算法的第一步不是线路，而是一个重新表述：把 ${\mu }_m$ 解释为某个概率分布的矩。这一步是纯粹的经典数学。

3.1 从谱分解到矩公式

$H$ 是实对称矩阵，由谱定理它有完备的正交归一特征矢量基：

$$H=\sum _r{\lambda }_r|u_r⟩⟨u_r|,{\lambda }_r\in [-1,1],$$

其中特征值的范围来自 $\Vert H\Vert \le 1$（§1.2 的归一化）。

第一步：$H$ 的幂的谱分解。 我们证明

$$H^m=\sum _r{\lambda }_r^m|u_r⟩⟨u_r|.$$

$m=1$ 时就是谱分解本身。设上式对 $m$ 成立，则

$$H^{m+1}=H\cdot H^m=\left(\sum _s{\lambda }_s|u_s⟩⟨u_s|\right)\left(\sum _r{\lambda }_r^m|u_r⟩⟨u_r|\right)=\sum _{s,r}{\lambda }_s{\lambda }_r^m|u_s⟩⟨u_s|u_r⟩⟨u_r|.$$

由正交归一性 $⟨u_s|u_r⟩={\delta }_{sr}$，双重求和退化为 $s=r$ 的单重求和，系数为 ${\lambda }_r\cdot {\lambda }_r^m={\lambda }_r^{m+1}$，归纳完成。这一步只用到了”投影算子 $|u_r⟩⟨u_r|$ 互相正交”——同一个特征矢量基同时对角化 $H$ 的所有幂。

第二步：夹在 $|j⟩$ 之间。 两边取对角元：

$${\mu }_m=⟨j|H^m|j⟩=\sum _r{\lambda }_r^m⟨j|u_r⟩⟨u_r|j⟩=\sum _r|⟨u_r|j⟩{|}^2{\lambda }_r^m,$$

最后一步用了 $⟨j|u_r⟩=\overline{⟨u_r|j⟩}$（内积的共轭对称性），故两者的乘积是模方。

3.2 谱测度是概率分布

定义

$$p_r:=|⟨u_r|j⟩{|}^2,{\nu }_j:=\sum _r{p}_r{\delta }_{{\lambda }_r},$$

其中 ${\delta }_{\lambda }$ 是在点 $\lambda$ 处的 Dirac 测度。称 ${\nu }_j$ 为 $H$ 在基态 $|j⟩$ 上看到的谱测度。它确实是概率测度：

• 非负性：$p_r=|⟨u_r|j⟩{|}^2\ge 0$；
• 归一性：由完备性关系 ${\sum }_r|u_r⟩⟨u_r|=I$，

$$\sum _r{p}_r=\sum _r⟨j|u_r⟩⟨u_r|j⟩=⟨j|\left(\sum _r|u_r⟩⟨u_r|\right)|j⟩=⟨j|j⟩=1.$$

（这就是 Parseval 恒等式：$|j⟩$ 在任何正交归一基下展开的系数模方和等于它的范数平方。）

于是 §3.1 的公式可以写成期望的形式：

$$\boxed{{\mu }_m={\mathbb{E}}_{\lambda \sim {\nu }_j}\left[{\lambda }^m\right]}$$

即 ${\mu }_m$ 是谱测度 ${\nu }_j$ 的第 $m$ 阶矩。记号 $\lambda \sim {\nu }_j$ 表示”以概率 $p_r$ 取值 ${\lambda }_r$“的随机变量。

这个表述立即解释了一些基本事实。例如 $|{\mu }_m|\le 1$：由三角不等式与 $|{\lambda }_r|\le 1$，

$$|{\mu }_m|\le \sum _r{p}_r|{\lambda }_r{|}^m\le \sum _r{p}_r=1,$$

与 §1.2 中 $|(A^m{)}_{jj}|\le b^m$ 的归一化版本一致（换回原矩阵就是乘以 $b^m$）。

3.3 直觉：局部谱信息

谱测度 ${\nu }_j$ 回答的问题是：“基态 $|j⟩$ 坐在 $H$ 的哪些特征方向上，各坐多重？” 把 $|j⟩$ 在特征基下展开，$|j⟩={\sum }_r⟨u_r|j⟩|u_r⟩$，谱测度记录的正是各分量的权重与对应的特征值——但不记录相位与方向信息。

这个观察解释了为什么这个问题天然地适合量子计算机。相位估计（）的工作方式恰恰是：对 $e^{iHt}$、输入态 $|j⟩$ 做测量，输出以概率 $p_r=|⟨u_r|j⟩{|}^2$ 抽到特征值 ${\lambda }_r$ 的估计——也就是说，相位估计天然地从谱测度 ${\nu }_j$ 中采样。经典算法若不读取指数多的矩阵结构，很难获得同一个局部谱测度的样本：算出 $p_r$ 需要知道特征矢量，而对角化一个 $2^n$ 维矩阵无从谈起。

值得强调：单个矩 ${\mu }_m$ 只是谱测度的一个线性泛函，远不能确定整个测度。我们也不需要整个测度——这正是问题”只问一个数”的省力之处，也是它最终落入 BQP 而非更难复杂类的原因。

4. 量子算法一：谱采样

4.1 算法流程

第一个实现把 §3 的概率解释直接变成采样算法：

1. 制备 $|j⟩$（就是一个计算基态，$n$ 个比特的制备代价 $O(n)$）；
2. 相位估计：对酉算子 $e^{iHt}$ 与输入 $|j⟩$ 运行相位估计（），得到某个特征值的估计 $\tilde{\lambda }$。由相位估计的分析，输出 ${\tilde{\lambda }}_r$ 的概率是 $p_r=|⟨u_r|j⟩{|}^2$，即恰好从 ${\nu }_j$ 采样；演化算子 $e^{iHt}$ 本身由 sparse Hamiltonian simulation 实现（稀疏预言机 + Trotter 化或更优方法），其代价关于 $s$、$t$、$\log N$ 是多项式；
3. 经典后处理：计算 ${\tilde{\lambda }}^m\in [-1,1]$；
4. 重复平均：独立重复 $R$ 次，输出样本均值 $\frac{1}{R}{\sum }_i{\tilde{\lambda }}_i^m$。

由期望的定义，单次输出的期望是 ${\mathbb{E}}_{{\nu }_j}[{\lambda }^m]={\mu }_m$（相差一个由估计精度造成的系统误差，见下），所以样本均值确实收敛到目标。剩下的问题是：特征值需要多精确？样本需要多少个？ 这是两个独立的误差来源，需要分别控制再做预算平衡。

4.2 关键引理：幂函数的 Lipschitz 常数

系统误差来自相位估计的分辨率：我们拿到的是 $\tilde{\lambda }\approx \lambda$，但代入的是 $m$ 次幂。误差如何放大？

引理. 对任意 $|x|,|y|\le 1$ 与正整数 $m$，

$$|x^m-y^m|\le m|x-y|.$$

证明。用有限几何级数的因式分解

$$x^m-y^m=(x-y)\sum _{k=0}^{m-1}{x}^{m-1-k}{y}^k,$$

这可以通过把右端展开验证：求和式乘以 $(x-y)$ 后是一个错位相消（telescoping）的和，${\sum }_k(x^{m-k}{y}^k-x^{m-1-k}{y}^{k+1})=x^m-y^m$。右端的和共有 $m$ 项，每项满足 $|x^{m-1-k}{y}^k|=|x{|}^{m-1-k}|y{|}^k\le 1$（因为 $|x|,|y|\le 1$）。由三角不等式，

$$|x^m-y^m|\le |x-y|\cdot \sum _{k=0}^{m-1}1=m|x-y|.$$

Q.E.D.

注意条件 $|x|,|y|\le 1$ 必不可少：若 $|x|>1$，因子 $|x{|}^{m-1-k}$ 会随 $m$ 指数放大。这正是 §1.2 把矩阵归一化到 $\Vert H\Vert \le 1$ 的又一个原因——归一化把幂函数变成了 $m$-Lipschitz 函数，特征值误差最多被放大 $m$ 倍，而 $m=\mathrm{polylog}N$ 是可承受的。

4.3 误差预算与参数平衡

总误差有两个来源，我们要求它们各自不超过 $ϵ/2$，从而总和不超过 $ϵ$。

系统误差（分辨率）。设相位估计给出的特征值误差为 $|\tilde{\lambda }-\lambda |\le \delta$。由引理，每次采样的偏差满足 $|{\tilde{\lambda }}^m-{\lambda }^m|\le m\delta$。这个偏差是系统的：它不会在平均中抵消，但因为每次采样都受同一个界控制，均值也受它控制：

$$\left|\mathbb{E}[{\tilde{\lambda }}^m]-\mathbb{E}[{\lambda }^m]\right|\le \mathbb{E}\left|{\tilde{\lambda }}^m-{\lambda }^m\right|\le m\delta ,$$

其中第一步是期望的三角不等式（Jensen 不等式的特例）。要求 $m\delta \le ϵ/2$，解出所需分辨率

$$\delta =O\left(\frac{ϵ}{m}\right).$$

相位估计要把特征值分辨到 $\delta$，需要的演化时间是 $t=O(1/\delta )=O(m/ϵ)$（相位估计的标准标度：精度与时间成反比）。

统计误差（采样）。每个样本 $X_i={\tilde{\lambda }}_i^m$ 满足 $|X_i|\le 1$，故方差 $\mathrm{Var}(X_i)\le \mathbb{E}[X_i^2]\le 1$。$R$ 个独立样本的均值方差为 $\mathrm{Var}/R\le 1/R$。由 Chebyshev 不等式，

$$\mathrm{Pr}\left[\left|\frac{1}{R}{\sum }_i{X}_i-\mathbb{E}[X]\right|\ge \frac{ϵ}{2}\right]\le \frac{1/R}{(ϵ/2{)}^2}=\frac{4}{R{ϵ}^2}.$$

要让失败概率为常数（例如 $1/3$），取

$$R=O\left(\frac{1}{{ϵ}^2}\right)$$

个样本即可。这就是普通的 Monte Carlo 标度：$[-1,1]$ 上有界随机变量的均值估计，代价 $O(1/{ϵ}^2)$。

平衡的结果。把两部分合起来：$R=O(1/{ϵ}^2)$ 次独立重复，每次的 Hamiltonian simulation 时间为 $t=O(m/ϵ)$，每次模拟的线路代价关于 $s$、$\log N$ 与 $t$ 为多项式。总代价

$$\underset{\text{样本数}}{\underbrace{O\left(\frac{1}{{ϵ}^2}\right)}}\times \underset{\text{单次模拟}}{\underbrace{\mathrm{poly}\left(s,\log N,\frac{m}{ϵ}\right)}}=\mathrm{poly}\left(s,m,\frac{1}{ϵ},\log N\right).$$

代入问题的参数尺度 $s,m,1/ϵ=\mathrm{polylog}N=\mathrm{poly}(n)$，总代价是 $\mathrm{poly}(n)$——该问题属于 PromiseBQP。注意这里的参数平衡是标准的”对分误差预算”：系统误差与统计误差各分到 $ϵ/2$，再分别反解出分辨率 $\delta$ 与样本数 $R$；不存在需要联立求解的隐式平衡。

4.4 复杂度表达式中每个因子的来源

把总代价拆开，逐项对照来源：

• 因子 $m$：来自引理中幂函数的 Lipschitz 常数——幂次越高，特征值需要分辨得越细，演化时间越长；
• 因子 $1/{ϵ}^2$：Monte Carlo 平均的统计涨落，有界随机变量均值估计的通用代价；
• 因子 $s$ 与 $\log N$：sparse Hamiltonian simulation 每一步访问预言机的代价，以及 $n=\log N$ 个系统比特上的基本门操作；
• 隐含的 $\log (1/{\delta }_{\mathrm{fail}})$：若要求成功概率 $1-{\delta }_{\mathrm{fail}}$，Chebyshev 换成 Hoeffding 不等式（样本有界时可用）会引入对数因子，标准的成功概率放大技术。

5. 量子算法二：Block Encoding 与 QSVT

第二个实现更贴近块编码的语言，也更直接地暴露”归一化”与”次数”这两个关键参数。

5.1 为什么 $p(x)=x^m$ 是合法多项式

由稀疏预言机可以构造 $H$ 的 block encoding（标准构造，见 block-encoding 教程）：一个作用在系统加少量辅助比特上的酉算子 $U_H$，其左上角块就是 $H$（相差一个已知的归一化常数，构造时吸收进误差预算）。

回忆 QSVT 对多项式 $p$ 的三个可实施性条件（）：(i) $p$ 的次数为 $d$；(ii) $p$ 有确定的奇偶性（奇或偶）；(iii) 在 $[-1,1]$ 上 $|p(x)|\le 1$。逐条核对 $p(x)=x^m$：

• 次数为 $m$；
• 奇偶性：$(-x{)}^m=(-1{)}^m{x}^m$，奇偶与 $m$ 一致，确定；
• 有界性：$|x|\le 1$ 时 $|x^m|=|x{|}^m\le 1$。

三条都满足——而且条件 (iii) 恰好又是 $\Vert H\Vert \le 1$ 这个归一化在起作用。于是 QSVT 给出 $H^m$ 的 block encoding，其代价是 $m$ 次对 $U_H$ 及其逆的调用（QSVT 的查询数等于多项式次数）。在 QSVT 出现之前，原论文的表述是”测量 $A$ 并估计 spectral moment”（即 §4 的谱采样路线）；两种表述等价，但 block-encoding 语言把 normalization（块编码的比例因子）与 degree（多项式次数 = 查询次数）显示在账面上。

5.2 从块编码到 Hadamard test

有了 $H^m$ 的 block encoding $U_{H^m}$，估计 $⟨j|H^m|j⟩$ 是标准的 Hadamard test：用控制比特叠加地施加 $U_{H^m}$，测量控制比特得到 $\mathrm{Re}⟨j|H^m|j⟩$（本问题中矩阵为实，虚部为零，一次 Hadamard test 即可）。块编码的比例因子是已知常数，除以它即得 ${\mu }_m$ 的估计。

Hadamard test 的单次输出是 $\pm 1$ 随机比特，均值估计同样需要 $O(1/{ϵ}^2)$ 次重复；但这里可以做得更好：把”控制比特为 $+1$“视为好事件做振幅估计 (amplitude estimation)（相位估计的又一应用），把精度代价从 $O(1/{ϵ}^2)$ 降到 $O(1/ϵ)$。

5.3 复杂度对比

$$\underset{\text{QSVT 查询}}{\underbrace{O(m)}}\times \underset{\text{块编码代价}}{\underbrace{\mathrm{poly}(s,\log N)}}\times \underset{\text{振幅估计}}{\underbrace{O\left(\frac{1}{ϵ}\right)}}=\mathrm{poly}\left(s,m,\frac{1}{ϵ},\log N\right).$$

与 §4 的谱采样相比，$ϵ$ 的依赖从 $1/{ϵ}^3$（样本数 $1/{ϵ}^2$ 乘以单次模拟 $m/ϵ$）改进到接近 $1/ϵ$，代价是动用了更重的机制（块编码、QSVT、振幅估计）。两条路线在参数尺度 $\mathrm{polylog}N$ 下都是 $\mathrm{poly}(n)$，所以 BQP 成员资格的结论相同。

6. 为什么不是”快速输出 $A^m$“

6.1 输出本身就太大

值得把 §2 的结论再强调一遍，因为它划定了这类算法的边界：量子算法不输出 $A^m$，也不输出它的某一行。$A^m$ 通常稠密，有 $\Theta (N^2)$ 个矩阵元；输出它的任何经典表示都至少需要 $\Omega (N^2)$ 时间。量子算法回答的只是”一个矩阵元”或更一般的少量线性泛函 $⟨\psi |A^m|\phi ⟩$。把问题理解成”快速计算矩阵幂”会得到错误的复杂度预期；正确的理解是”从局部稀疏预言机直接访问谱函数的值，避免物化任何中间矩阵”。

同样的边界条件也适用于读者熟悉的其他算法：HHL 不输出解向量本身，只输出关于解的测量统计；量子相位估计不输出整个谱，只输出被采样到的那个特征值。

6.2 非对角元：polarization 恒等式

对角元算法其实够用于所有矩阵元。对非对角元，用 polarization 恒等式把问题归约到对角元。给定 $i\ne j$，定义叠加态（这里的 $|+⟩$ 依赖于 $i,j$，是临时记号，不要与 Hadamard 基混淆）

$$|+⟩:=\frac{|i⟩+|j⟩}{\sqrt{2}}.$$

对实对称矩阵 $B$（例如 $B=H^m$），有恒等式

$$⟨i|B|j⟩=⟨+|B|+⟩-\frac{1}{2}⟨i|B|i⟩-\frac{1}{2}⟨j|B|j⟩.$$

验证。把 $|+⟩$ 的定义代入并展开（四个交叉项）：

$$⟨+|B|+⟩=\frac{1}{2}\left(⟨i|B|i⟩+⟨i|B|j⟩+⟨j|B|i⟩+⟨j|B|j⟩\right).$$

因为 $B$ 实对称，$⟨j|B|i⟩=⟨i|B|j⟩$（实对称的内积 $⟨j|B|i⟩=(B^T{)}_{ij}=B_{ij}$，就是同一个数），两个交叉项合并：

$$⟨+|B|+⟩=\frac{1}{2}⟨i|B|i⟩+\frac{1}{2}⟨j|B|j⟩+⟨i|B|j⟩.$$

把对角项移到左端即得恒等式。Q.E.D.

右端三项都是”某个易制备态上的对角元”：$|i⟩,|j⟩$ 是计算基态，$|+⟩$ 只需对 $i,j$ 的地址比特做一次受控叠加即可制备。因此对角元算法的任何实现都自动覆盖全部矩阵元，代价至多是三倍。注意这个归约依赖 $B$ 实对称；一般 Hermitian 情形需要额外的 $i$ 相位版本区分实部与虚部（留作习题）。

7. PromiseBQP-hardness：Feynman clock

现在转向反方向：为什么这个问题难。我们将看到任意 BQP 计算的接受概率可以被编码成某个稀疏 Hermitian 矩阵的谱矩，因此该估计问题是 PromiseBQP-hard 的；结合 §4–§5 的成员资格，它是 PromiseBQP-complete 的。

构造的思想可以追溯到 Feynman 的 clock 寄存器与 Kitaev 的 history state：给量子电路配一个”时钟”，让矩阵的幂把演化历史当作路径来数。

7.1 propagation matrix 的构造与基本性质

给定作用在 $n$ 个比特上的量子电路

$$U=U_T\cdots U_1,T=\mathrm{poly}(n),$$

引入时钟寄存器（维数 $T+1$，基态 $|0⟩,\dots ,|T⟩$），在”时钟 $\otimes$ 数据”空间上定义 propagation matrix

$$K=\sum _{t=0}^{T-1}\left(|t+1⟩⟨t|\otimes U_{t+1}+|t⟩⟨t+1|\otimes U_{t+1}^{†}\right).$$

逐项理解：第一项把时钟从 $t$ 拨到 $t+1$，同时在数据寄存器上施加第 $t+1$ 个门；第二项是它的逆操作（时钟回拨、施加逆门）。$K$ 有三个直接可验证的性质。

Hermitian 性。对和式中的每一项取伴随：

$${\left(|t+1⟩⟨t|\otimes U_{t+1}\right)}^{†}=|t⟩⟨t+1|\otimes U_{t+1}^{†},$$

恰好是和式中的另一项。所以 $K^{†}=K$，它有实谱，谱矩语言适用。

稀疏性。固定基态 $|t⟩|x⟩$（$x$ 是数据的计算基态），

$$K|t⟩|x⟩=|t+1⟩\otimes U_{t+1}|x⟩+|t-1⟩\otimes U_t^{†}|x⟩$$

（边界处只有一项）。每个时钟位形只连接前后两步；而一个局部门（一至两个比特）把计算基态映到至多常数个基态（例如两比特门至多 $4$ 个）。所以 $K$ 的每列只有 $O(1)$ 个非零元，且非零元就是门的矩阵元——位置与取值都可以从电路描述在 $\mathrm{poly}(n)$ 时间内算出，即稀疏预言机可以有效实现。

范数有界。记”前进部分”

$$W:=\sum _{t=0}^{T-1}|t+1⟩⟨t|\otimes U_{t+1},$$

则 $K=W+W^{†}$。计算

$$W^{†}W=\sum _{t,t^{\prime }}\left(|t^{\prime }⟩⟨t^{\prime }+1|\otimes U_{t^{\prime }+1}^{†}\right)\left(|t+1⟩⟨t|\otimes U_{t+1}\right)=\sum _t|t⟩⟨t|\otimes U_{t+1}^{†}{U}_{t+1},$$

其中第二个等号用了时钟基的正交性 $⟨t^{\prime }+1|t+1⟩={\delta }_{t^{\prime }t}$。由 $U_{t+1}^{†}{U}_{t+1}=I$，$W^{†}W$ 是一个投影算子（向时钟 $t\le T-1$ 子空间的投影），故 $\Vert W\Vert =1$（$W$ 是部分等距）。于是

$$\Vert K\Vert \le \Vert W\Vert +\Vert W^{†}\Vert =2,$$

即归一化常数 $b\le 2$，与 $n$ 无关。

7.2 $K^m$ 的闭合路径展开

取时钟为 $0$、数据为某个计算基态 $x$ 的联合基态 $|j⟩=|0⟩|x⟩$，考察 $⟨j|K^m|j⟩$。由 §2.1 的路径展开恒等式（把基态写成 $(t_i,x_i)$ 对），

$$⟨0,x|K^m|0,x⟩=\sum _{\begin{array}{c}{t}_0=\cdots =t_m=0\\ x_0=x_m=x\end{array}}\prod _{i=0}^{m-1}⟨t_{i+1},x_{i+1}|K|t_i,x_i⟩.$$

由 §7.1，因子 $⟨t_{i+1},x_{i+1}|K|t_i,x_i⟩$ 非零仅当 $|t_{i+1}-t_i|=1$：时钟序列必须是一条在线段 $\{0,1,\dots ,T\}$ 上、每步走一格的路径。这条路径从 $0$ 出发又回到 $0$，所以前进与后退的步数必须相等，$m$ 必为偶数；奇数阶矩恒为零——这与对称谱 ${\nu }_j(\lambda )={\nu }_j(-\lambda )$ 的图像一致。

每条时钟路径对应一个数据寄存器上的算子乘积：前进经过 $t\to t+1$ 贡献 $U_{t+1}$，后退经过 $t+1\to t$ 贡献 $U_{t+1}^{†}$。对中间数据基态 $x_1,\dots ,x_{m-1}$ 求和，由路径展开恒等式，恰好拼出算子乘积的对角元：

$$⟨0,x|K^m|0,x⟩=\sum _{\text{闭合时钟路径 }(t_0,\dots ,t_m)}⟨x|V_{t_{m-1}\to t_m}\cdots V_{t_0\to t_1}|x⟩,V_{t\to t+1}=U_{t+1},V_{t+1\to t}=U_{t+1}^{†}.$$

例如，“单调走到 $T$ 再走回来”的 excursion（要求 $m\ge 2T$）贡献

$$⟨x|U_1^{†}{U}_2^{†}\cdots U_T^{†}\cdot U_T\cdots U_2{U}_1|x⟩=1,$$

即这类路径只贡献平凡的计数。这说明裸的 $K$ 还不够：有用的信号必须来自被门振幅非平凡加权的闭合路径，其干涉权重包含电路各门的矩阵元。

7.3 接受概率如何进入谱矩

原论文通过三个补充构造让非平凡信号分离出来：

• 输入 flag：在时钟起点附近加入标志位与投影项，强制参与求和的路径从合法初态（例如数据寄存器为 $|0^n⟩$）出发，排除垃圾历史；
• 输出 gadget：在时钟终点处放置与接受输出位耦合的结构，使穿过终点的路径带上正比于电路接受振幅的权重；
• padding（恒等步）：在电路末端拼接若干恒等”门”，延长时钟线段，使得可以选择一个 $\mathrm{poly}(n)$ 大小的 $m$，让含信号的路径族恰好被第 $m$ 阶矩分离出来。

最终的结论是：对任意 BQP 语言与任意输入，可以构造（$\mathrm{poly}(n)$ 时间内构造）一个稀疏 Hermitian 矩阵与一个基态 $|j⟩$，使得电路接受时 $⟨j|K^m|j⟩$ 落在一个区间、拒绝时落在与之不相交的另一区间，两区间之间隔着逆多项式（在归一化口径下逆多对数）的 gap。于是，任何能以 $ϵb^m$ 加性误差估计稀疏矩阵幂对角元的算法，都能判定任意 BQP 语言——PromiseBQP-hardness 得证。

精确的 gadget 构造技术性较强，超出本词条范围；这里的要点是机制：矩阵幂把演化历史展开成闭合路径的和，路径的代数干涉权重里藏着门振幅，而 flag 与 padding 负责把”完整穿过计算”的路径从平凡计数中挑出来。 这与 §2.1 路径展开恒等式给出的 Feynman 图像是同一枚硬币的两面。

还有一个耐人寻味的细节：通过把每个门替换为只含 $\pm 1$ 矩阵元的小 gadget（门的”离散化”），构造甚至可以把矩阵元限制为 $-1,0,1$。这说明该问题的困难不来自高精度输入数值——即使矩阵元粗糙到只有三个取值，问题依然 BQP-hard；困难完全来自隐式图的组合结构与相位干涉。

7.4 保留条款

• 结论针对的是 promise problem：gap（接受/拒绝两区间的分离）是输入承诺的一部分。没有这个 gap，“矩落在哪个区间”本身不是一个良定义的判定问题。
• 复杂度口径是 §1.2 的归一化口径：误差 $ϵb^m$、$1/ϵ=\mathrm{polylog}N$。改成绝对误差或相对误差，结论不再成立。
• “稀疏”与”实对称”两个条件都在构造中起到实质作用：稀疏性保证预言机高效、保证量子算法一端的多项式代价；对称（Hermitian）性保证谱矩语言与谱采样成立。

8. 与字符串重写的统一

字符串重写一节考虑的问题看似完全不同：按局部替换规则把字符串 $s$ 改成 $t$，问”恰好 $m$ 步的重写序列数”之差的符号。但那里的重写图邻接矩阵同样是稀疏 Hermitian 的，重写路径数同样是其幂的矩阵元（那里用的是 §2.1 的路径展开恒等式）。两个 BQP-complete 结果共享同一条流水线：

$$\text{局部 transition rules}⟶\text{隐式 sparse matrix}⟶\text{高阶谱矩/路径振幅}.$$

差别只在符号的来源：矩阵幂问题允许一般的正负矩阵元，符号直接来自权重；字符串重写中只有非负的路径计数，符号靠”两个目标路径数之差”人为制造。对应的 Feynman clock / propagation Hamiltonian 构造在两边几乎逐字相同。把两节对照阅读，可以看到”局部规则 $\Rightarrow$ 隐式稀疏矩阵 $\Rightarrow$ 谱矩”是反复出现的范式。

9. 两个小例子

9.1 Pauli $X$：三种算法互相验证

取

$$H=\left(\begin{array}{cc}0& 1\\ 1& 0\end{array}\right)=X,|j⟩=|0⟩.$$

这里 $N=2$、$s=1$、$\Vert H\Vert =1$。我们用三种方法计算 ${\mu }_m=⟨0|H^m|0⟩$。

方法一：直接乘方。计算 $X^2=\left(\begin{array}{cc}0& 1\\ 1& 0\end{array}\right)\left(\begin{array}{cc}0& 1\\ 1& 0\end{array}\right)=\left(\begin{array}{cc}1& 0\\ 0& 1\end{array}\right)=I$。于是

$$H^{2r}=(X^2{)}^r=I,H^{2r+1}=X\cdot X^{2r}=X,$$

故

$$⟨0|H^m|0⟩=\{\begin{array}{ll}⟨0|I|0⟩=1,& m\text{ 偶},\\ ⟨0|X|0⟩=0,& m\text{ 奇}.\end{array}$$

方法二：谱测度。$X$ 的特征分解为 $X=|+⟩⟨+|-|-⟩⟨-|$，特征值 ${\lambda }_{\pm }=\pm 1$，特征矢量 $|\pm ⟩=(|0⟩\pm |1⟩)/\sqrt{2}$。权重

$$p_{+}=|⟨+|0⟩{|}^2=\frac{1}{2},p_{-}=|⟨-|0⟩{|}^2=\frac{1}{2},$$

即谱测度在 $\lambda =\pm 1$ 上各有质量 $1/2$。于是

$${\mu }_m=\frac{1}{2}\cdot 1^m+\frac{1}{2}\cdot (-1{)}^m=\frac{1+(-1{)}^m}{2}=\{\begin{array}{ll}1,& m\text{ 偶},\\ 0,& m\text{ 奇}.\end{array}$$

方法三：路径计数。把 $X$ 看成两个顶点、两条平行边（或一条无向边）的图的邻接矩阵：从 $0$ 出发每步只能走到 $1$、再走回来，路径被唯一决定。$m$ 为奇数时不可能回到 $0$（零条闭合路径）；$m$ 为偶数时恰有一条（来回交替）。与方法一、二一致。

三种视角——算子代数、谱测度、组合路径——给出同一个答案。这个例子也验证了一般理论：谱测度关于 $0$ 对称（$\pm 1$ 等权），所以奇数阶矩为零，与 §7.2 中”奇数阶矩恒零”的时钟奇偶性论证同源。

9.2 三顶点路径图：归一化的作用

取三个顶点的路径图的邻接矩阵，$j$ 为中间顶点：

$$A=\left(\begin{array}{ccc}0& 1& 0\\ 1& 0& 1\\ 0& 1& 0\end{array}\right),|j⟩=(0,1,0{)}^T.$$

路径计数。从中间顶点出发走两步：先到左端再回来，或先到右端再回来，共 $2$ 条闭合路径，故 $(A^2{)}_{jj}=2$。走四步时，每条两步路径可以继续分裂：逐步数下去得 $(A^4{)}_{jj}=4$（先左/右、回中、再左/右、回中，共 $2\times 2=4$ 条）。

谱方法验证。$A$ 的特征值为 $\sqrt{2},0,-\sqrt{2}$（特征方程 $\det (\lambda I-A)=\lambda ({\lambda }^2-2)=0$：展开行列式，${\lambda }^3-\lambda -\lambda ={\lambda }^3-2\lambda$）。对应 $\lambda =\pm \sqrt{2}$ 的归一化特征矢量为 $\frac{1}{2}(1,\pm \sqrt{2},1{)}^T$，与 $|j⟩$ 的内积为 $\pm \sqrt{2}/2$，权重各为 $1/2$；$\lambda =0$ 的特征矢量权重为 $0$。于是

$${\mu }_m=\frac{1}{2}(\sqrt{2}{)}^m+\frac{1}{2}(-\sqrt{2}{)}^m=\{\begin{array}{ll}{2}^{m/2},& m\text{ 偶},\\ 0,& m\text{ 奇},\end{array}$$

与路径计数一致（$m=2$ 得 $2$，$m=4$ 得 $4$）。

归一化。这里 $\Vert A\Vert =\sqrt{2}>1$，所以矩随 $m$ 增长：${\mu }_m=\Vert A{\Vert }^m$（$m$ 偶）。算法运行时必须先归一化 $H=A/\sqrt{2}$，此时 $⟨j|H^2|j⟩=1$，落在 $[-1,1]$ 内，谱采样与 QSVT 的有界性条件才成立；最后把答案乘回 $(\sqrt{2}{)}^m=b^m$。这个小例子把 §1.2 的抽象口径落到了具体数字上：$b^m$ 不是人为的缩放，而是这个量真实的增长尺度。

10. 小结与习题

• 稀疏矩阵幂的对角元 $⟨j|H^m|j⟩$ 是基态 $|j⟩$ 的谱测度 ${\nu }_j={\sum }_r|⟨u_r|j⟩{|}^2{\delta }_{{\lambda }_r}$ 的第 $m$ 阶矩；它同时是从 $j$ 出发长度 $m$ 的带权闭合路径之和。
• 量子算法有两条等价路线：谱采样（sparse Hamiltonian simulation + 相位估计 + Monte Carlo，代价 $\mathrm{poly}(s,m,1/ϵ,\log N)$）与 block encoding / QSVT（次数 $m$ 的多项式变换 + Hadamard test / 振幅估计）。归一化 $\Vert H\Vert \le 1$ 同时保证幂函数的 $m$-Lipschitz 性与 QSVT 多项式的有界性。
• 算法不输出稠密的 $A^m$，只回答单个矩阵元或少量线性泛函；误差口径是 $ϵb^m$，归一化是问题定义的核心。
• Feynman clock 把任意量子电路编码为稀疏 Hermitian propagation matrix，其幂的对角元把演化历史展开成闭合路径的干涉求和；配合 flag 与 padding，接受/拒绝电路给出落在不相交区间的矩，从而该 promise problem 是 PromiseBQP-complete 的。
• 与字符串重写问题共享”局部规则 $\to$ 隐式稀疏矩阵 $\to$ 高阶谱矩”的范式。

习题

1. （谱矩公式）设 $H={\sum }_r{\lambda }_r|u_r⟩⟨u_r|$ 是 $\Vert H\Vert \le 1$ 的 Hermitian 矩阵的谱分解。证明 $H^m={\sum }_r{\lambda }_r^m|u_r⟩⟨u_r|$，并由此推导 ${\mu }_m={\sum }_r|⟨u_r|j⟩{|}^2{\lambda }_r^m$。说明证明中哪一步用到了特征矢量的正交归一性。
2. （Lipschitz 界）证明对 $|x|,|y|\le 1$ 有 $|x^m-y^m|\le m|x-y|$（提示：因式分解 $x^m-y^m$）。给出一个 $|x|>1$ 时该界失效的例子，并解释这与”先把 $A$ 归一化到 $\Vert H\Vert \le 1$“之间的关系。
3. （polarization）展开 $⟨+|B|+⟩$ 验证 §6.2 的恒等式，并指出”实对称”用在了哪一步。若 $B$ 只是 Hermitian（允许复矩阵元），写出用 $|+⟩=(|i⟩+|j⟩)/\sqrt{2}$ 与 $|+i⟩=(|i⟩+i|j⟩)/\sqrt{2}$ 分别恢复 $\mathrm{Re}⟨i|B|j⟩$ 与 $\mathrm{Im}⟨i|B|j⟩$ 的恒等式。
4. （手算谱矩）对 §9.2 的三顶点路径图，用路径计数直接数出 $(A^6{)}_{jj}$，再用谱测度公式 ${\mu }_m=\frac{1}{2}(\sqrt{2}{)}^m+\frac{1}{2}(-\sqrt{2}{)}^m$ 核对。中间顶点的谱测度为什么不含 $\lambda =0$ 的成分？
5. （Feynman clock 的奇偶性）证明 §7.2 的断言：对时钟基态 $|j⟩=|0⟩|x⟩$，当 $m$ 为奇数时 $⟨j|K^m|j⟩=0$。进一步证明 $\Vert K\Vert \le 2$（提示：把 $K$ 写成 $W+W^{†}$ 并计算 $W^{†}W$）。
6. （复杂度口径）本题要求把参数对应关系讲清楚。(a) 设 $N=2^n$、$m=\mathrm{poly}(n)$、$1/ϵ=\mathrm{poly}(n)$、$s=\mathrm{poly}(n)$，说明谱采样算法的总代价是 $\mathrm{poly}(n)$，逐因子对应。(b) 设 $b$ 为大于 $1$ 的常数（例如 $b=2$），则 $b^m=2^{\mathrm{poly}(n)}$ 指数大；此时若把误差要求从 $ϵb^m$ 改为原值的绝对误差 $1$，相当于要求归一化量 ${\mu }_m$ 的误差达到 $b^{-m}$。说明谱采样算法的代价如何变化，为什么 §1.2 说”复杂度口径完全不同”。

参考文献

• Zoo 编号 60：Dominik Janzing 与 Pawel Wocjan, Estimating Diagonal Entries of Powers of Sparse Symmetric Matrices Is BQP-Complete.

---

参考资料

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