量子线性求解器详解:从 HHL 到现代改进¶
求解线性方程组 \(A\vec{x} = \vec{b}\) 是科学计算中最基础、最高频的问题。HHL 算法(2009)首次证明量子计算机可在多项式时间内完成这一任务,但原始 HHL 有多项重大局限。本文系统梳理从 HHL 原始算法到近年来基于 Qubitization、QSVT、变分方法等现代改进的完整演进,并补全各方法中关键的推导环节。
问题定义¶
给定 \(A \in \mathbb{C}^{N \times N}\)(\(N = 2^n\)),\(\vec{b} \in \mathbb{C}^N\),我们要求 \(\vec{x} = A^{-1}\vec{b}\)。记 \(s\) 为稀疏度(每行非零元个数的上界),\(\kappa\) 为条件数(condition number)。
在经典计算一侧,主要有以下两类方法:
高斯消元:\(O(N^3)\)
共轭梯度法(稀疏 \(A\)):\(O(Ns\kappa \log(1/\varepsilon))\)
我们的量子目标是:以 \(O(\text{poly}(n) \cdot \text{poly}(\kappa) \cdot \log(1/\varepsilon))\) 的代价制备 \(|x\rangle \propto A^{-1}|b\rangle\)(原始 HHL 对精度的依赖是多项式级的,见下文)。
第一代:原始 HHL(2009)¶
算法回顾¶
Harrow, Hassidim, Lloyd (2009) 提出的流程为:
QPE:对 \(|b\rangle\) 执行量子相位估计,将 \(A\) 的特征值估计 \(\tilde\lambda_j\) 读入时钟寄存器;
受控旋转:以寄存器中的 \(\tilde\lambda_j\) 为控制条件,对辅助比特作用 \(R_y(\theta_j)\),其中 \(\theta_j = 2\arcsin(C/\tilde\lambda_j)\)、\(C \le \lambda_{\min}\),使辅助比特 \(|1\rangle\) 分量的振幅为 \(C/\tilde\lambda_j\),从而把因子 \(1/\lambda_j\) 写入振幅(注意旋转的"角度"不等于 \(C/\lambda_j\),是角度的正弦值等于 \(C/\lambda_j\);推导见本系列的 HHL 一篇);
逆 QPE:消去时钟寄存器中的特征值信息;
后选择(post-selection):测量辅助比特,得 \(|1\rangle\) 时系统态坍缩为 \(|x\rangle \propto A^{-1}|b\rangle\)。
局限性¶
局限 |
描述 |
|---|---|
\(\kappa^2\) 依赖 |
成功概率 \(P \ge 1/\kappa^2\),需 \(O(\kappa^2)\) 次重复 |
稀疏性要求 |
\(A\) 必须是 \(s\)-稀疏的,QPE 需要 \(e^{iA}\) 的高效实现 |
输出问题 |
只能提取 $\langle x |
输入问题 |
$ |
条件数限制 |
\(\kappa\) 不能指数级大 |
精度依赖 |
相位估计与旋转的精度开销使总代价带 \(O(1/\varepsilon)\) 因子 |
第二代:最优精度改进(2017)¶
Childs-Kothari-Somma (SIAM J. Comp. 2017)¶
关键改进:将总代价中对 \(\kappa\) 的依赖从 \(O(\kappa^2)\) 降至 \(O(\kappa\log(\kappa/\varepsilon))\),同时把对精度的依赖从 \(O(1/\varepsilon)\) 降至 \(O(\log(1/\varepsilon))\)——这正是论文标题"对精度指数改进"的含义。
方法:用多项式滤波(eigenvalue filtering)替代简单的受控旋转。将 \(\lambda_{\max}\) 归一化为 \(1\),则全部特征值落在 \([1/\kappa, 1]\) 内。我们不再逐特征值做受控旋转,而是构造多项式 \(P\) 使 \(P(x) \approx 1/x\)(在 \([1/\kappa, 1]\) 上,配合适当的整体缩放),再通过"QPE + 受控多项式函数 + 逆 QPE"把每个特征分量的振幅乘以 \(P(\lambda_j) \approx 1/\lambda_j\)。
多项式次数为什么是 \(d = O(\kappa\log(\kappa/\varepsilon))\):我们可以用一个初等推导说明所需的量级(CKS 的具体构造用切比雪夫多项式与线性组合 of 酉算子实现,达到同一量级)。对 \(x\in[1/\kappa,1]\),令 \(y = \frac{1-x}{1+x} \in [0, \frac{\kappa-1}{\kappa+1}]\);反解得 \(x = \frac{1-y}{1+y}\),代入可验证恒等式
(最后一个等号用了等比级数 \(1/(1-y^2) = \sum_j y^{2j}\)。)截断到 \(J\) 项的误差不超过 \(\frac{2\,y^{2J+2}}{1-y^2}\);在 \(y\) 的上界 \(y_{\max} = \frac{\kappa-1}{\kappa+1}\) 处估计
合并得截断误差 \(\lesssim \frac{\kappa}{2}\,e^{-4J/(\kappa+1)}\),故取 \(J = O(\kappa\log(\kappa/\varepsilon))\) 即可把它压到 \(\varepsilon\)。这正是次数中 \(\kappa\log(\kappa/\varepsilon)\) 因子的来源。
总复杂度:滤波过程需要 \(d\) 次对特征值预言机(等价地,\(e^{iA}\) 的实现)的调用,每次调用的代价记为 \(C_{U_A}\),故
其中在稀疏存取模型下 \(C_{U_A} = O(s\cdot n)\) 是 \(e^{iA}\) 的一次实现代价。
优势:多项式滤波直接缩放振幅,输出态的范数损失通过把多项式整体缩放并配合不依赖初态的振幅放大(oblivious amplitude amplification)补回,因此不需要 \(O(\kappa^2)\) 次后选择重复;对 \(\kappa\) 的依赖完全来自多项式次数 \(d\),即线性依赖。
第三代:Qubitization + QSVT 框架(2019)¶
Gilyén, Su, Low, Wiebe (STOC 2019)¶
QSVT(量子奇异值变换,Quantum Singular Value Transformation)为量子线性求解提供了最通用、最优雅的框架。
核心思想:将 \(A^{-1}\) 的实现归约为对 \(A\) 的奇异值施加多项式变换 \(P(\sigma_j) \approx c/\sigma_j\)(\(c\) 为缩放常数)。
步骤:
块编码(block encoding):构造 \(A\) 的 \((\alpha, m, \varepsilon)\)-块编码 \(U_A\),即满足 \(\big\|A - \alpha\,(\langle 0^m|\otimes I)\,U_A\,(|0^m\rangle\otimes I)\big\| \le \varepsilon\) 的酉算符——矩阵 \(A/\alpha\) 被嵌入在 \(U_A\) 的某个块中,归一化因子 \(\alpha\) 衡量编码的"代价"。
Qubitization:由 \(U_A\) 构造量子行走(quantum walk)酉算符 \(W(A)\)。标准结论(此处只陈述不证明)是:在由 \(|0^m\rangle\otimes|u_j\rangle\) 与一个正交态张成的二维不变子空间中,\(W\) 的作用等价于平面旋转,旋转角 \(\theta_j\) 满足 \(\cos\theta_j = \sigma_j/\alpha\)(\(\sigma_j\) 为 \(A\) 的奇异值)。注意归一化因子 \(\alpha\) 必须出现在此处:块编码的左上块是 \(A/\alpha\),因此是 \(\sigma_j/\alpha\) 而非 \(\sigma_j\) 本身决定旋转角。
QSP 多项式变换:量子信号处理(Quantum Signal Processing, QSP)告诉我们,选取相位序列 \(\vec\Phi = (\phi_0,\dots,\phi_d)\) 并交替作用 \(W\) 与单比特旋转后,所得电路的相应块是 \(\cos\theta\) 的 \(d\) 次多项式 \(P_{\vec\Phi}(\cos\theta)\),可实现的多项式有清晰的代数刻画(给定宇称与 \(|P|\le 1\) 型约束)。取 \(x = \sigma/\alpha = \cos\theta \in [1/\kappa', 1]\)(\(\kappa' = \alpha/\sigma_{\min}\) 为有效条件数),用与第二代相同的分析选 \(P(x)\approx c/x\),则 \(P(\cos\theta_j)\propto 1/\sigma_j\)。
执行:\(d\) 次 \(W(A)\) 调用,电路深度 \(O(d)\)。
多项式次数:\(d = O(\kappa'\log(\kappa'/\varepsilon))\);当块编码的归一化 \(\alpha\) 与 \(\sigma_{\max}\) 同阶时,\(\kappa' = \alpha/\sigma_{\min}\) 与通常的条件数 \(\kappa = \sigma_{\max}/\sigma_{\min}\) 同阶,即 \(d = O(\kappa\log(\kappa/\varepsilon))\)。
总复杂度:
对比 HHL:
特性 |
原始 HHL |
Childs (2017) |
QSVT (2019) |
|---|---|---|---|
\(\kappa\) 依赖 |
\(O(\kappa^2)\) |
\(O(\kappa)\) |
\(O(\kappa)\) |
精度依赖 |
\(O(1/\varepsilon)\) |
\(O(\log(1/\varepsilon))\) |
\(O(\log(1/\varepsilon))\) |
后选择 |
有 |
无 |
无 |
块编码 |
不需要 |
隐式使用(LCU) |
需要(显式框架) |
通用性 |
QPE-based |
QPE-based |
统一框架 |
第四代:变分量子线性求解器(VQLS,2020)¶
Bravo-Prieto et al. (Quantum 2023)¶
VQLS 是面向 NISQ 设备的量子-经典混合算法,不使用 QPE,无需深层电路。
核心思想:将 \(A\vec{x} = \vec{b}\) 转化为变分优化问题。直接的优化目标是残差范数:
由于 \(|x(\vec{\theta})\rangle\) 已被归一化,无法靠整体放大缩小残差,实际采用的是与之同解的归一化损失函数:
两者的关系可以严格说清。由 Cauchy–Schwarz 不等式与 \(\langle b|b\rangle = 1\),
因此 \(\mathcal{L} \ge 0\);且等号成立当且仅当 \(A|x\rangle\) 与 \(|b\rangle\) 平行(\(A\) 可逆时即 \(|x\rangle \propto A^{-1}|b\rangle\),整体相位可吸收进拟设)。于是 \(\mathcal{L}\) 的全局最小点恰是方程的解,且 \(\mathcal{L} = 0\) 与残差为零同时发生。两个目标函数在一般参数处取值不同,但零点集一致,因此作为优化目标互相等价。
算法步骤:
参数化电路:\(|x(\vec{\theta})\rangle = U(\vec{\theta})|0\rangle^{\otimes n}\)
损失函数估计:用 Hadamard 测试估计 \(\langle x|A^\dagger A|x\rangle\)、\(\langle b|A|x\rangle\) 等
经典优化:Adam/COBYLA 更新 \(\vec{\theta}\)
收敛后测量:测量 \(|x(\vec{\theta})\rangle\) 得到解的统计信息
Hadamard Test 电路:
|0⟩ ─ H ─ ctrl-A|x(θ)⟩ ─ H ─ 测量
|b⟩ ──────────────────────
测量期望值 \(\text{Re}\langle b|A|x(\vec{\theta})\rangle\)。该输出可由三行推导得到:控制比特经 \(H\)、受控 \(U\)、再经 \(H\) 后,末态为 \(\tfrac12|0\rangle(1+U)|\psi\rangle + \tfrac12|1\rangle(1-U)|\psi\rangle\);两个分支的模方差为
即 \(Z\) 期望值正是交叠的实部。需要注意 \(A\) 一般不是酉算子,不能直接受控作用;VQLS 先把 \(A\) 分解为酉算子的线性组合(Linear Combination of Unitaries, LCU)\(A = \sum_l c_l U_l\)(例如 Pauli 分解),对每个 \(U_l\) 分别做 Hadamard 测试,再按系数 \(c_l\) 线性组合即得所需交叠。
复杂度:
量子电路深度:\(O(C_A \cdot D)\)(\(C_A\) 为 \(A\) 的电路深度,\(D\) 为拟设深度)
测量次数:\(O(\text{poly}(n) / \varepsilon^2)\)(振幅估计到精度 \(\varepsilon\) 的采样代价)
优化迭代:\(O(T_{\text{opt}})\)(经验)
优势:
无 QPE,电路浅
适用于 NISQ 设备
可处理非稀疏矩阵
局限:
无收敛保证(可能陷入局部最小值)
Barren Plateaus(贫瘠高原)问题
精度有限
LDSE(Linear Differential Equation Solver)¶
Xu, Sun, Yuan (PRL 2023) 将 VQLS 推广到微分方程求解,用变分方法实现 \(A(t)\vec{x}(t) = \vec{b}(t)\) 的时间演化,思路是把各时刻的线性系统沿时间轴逐点变分求解并拼接轨迹。
第五代:基于 Schrödingerization 的方法(2022-2024)¶
Jin, Liu, Yu (2022-2024)¶
Schrödingerization 框架将线性方程组 \(A\vec{x} = \vec{b}\) 的求解转化为薛定谔方程的模拟。
思想:考虑动力学
其解可以显式写出:令误差 \(\vec{e}(t) = \vec{x}(t)-\vec{x}_{\text{eq}}\),则 \(\dot{\vec{e}} = -A\vec{e}\),在 \(A\) 的特征基下逐分量积分得 \(\vec{e}(t) = e^{-At}\vec{e}(0)\),即
取 \(\vec{x}(0) = \vec{b}\)。若 \(A\) 的特征值均有正实部,则 \(e^{-At}\to 0\)、\(\vec{x}(t)\to A^{-1}\vec{b}\)。收敛速度由最慢模 \(e^{-\lambda_{\min}t}\) 决定:由 \(\|\vec{x}(0)-\vec{x}_{\text{eq}}\| \le O(\kappa)\|\vec{b}\|\)、\(\|\vec{x}_{\text{eq}}\|\ge \|\vec{b}\|/\lambda_{\max}\),要把相对误差压到 \(\varepsilon\) 需要 \(t = \Omega(\kappa\log(\kappa/\varepsilon))\)(\(\lambda_{\min} = \lambda_{\max}/\kappa\))。
困难与方法:演化 \(e^{-At}\) 是耗散收缩而非酉的,不能直接作为量子门实现。Schrödingerization 引入一个附加维度并对其实施傅里叶变换("变形变换"),把非酉演化嵌入更高维空间中的薛定谔方程,从而可用标准的哈密顿量模拟技术模拟;这也使方法不要求 \(A\) 厄米,非对称性由附加维度处理。
代价:非酉嵌入的归一化因子约为 \(e^{\|A\|t}\),因此门数为 \(O(sn\cdot e^{\|A\|t})\),随演化时间指数增长;而收敛又要求 \(t = \Omega(\kappa\log(\kappa/\varepsilon))\),两者相互制约。这是此类"直接演化到稳态"方案的已知瓶颈。
其他现代方法¶
量子梯度下降法(2025)¶
arXiv:2502.13630 提出基于梯度下降的量子线性求解器,迭代格式为最小二乘目标 \(f(\vec{x}) = \tfrac12\|A\vec{x}-\vec{b}\|^2\) 的梯度步:
收敛性可以完整推导。记精确解 \(\vec{x}^* = A^{-1}\vec{b}\)(满足 \(A^TA\vec{x}^* = A^T\vec{b}\))与误差 \(\vec{e}_k = \vec{x}_k - \vec{x}^*\),代入迭代式得
在 \(A^TA\) 的特征基下,每个模态每步乘以 \(|1-\eta\sigma_j^2|\),其中 \(\sigma_j^2\in[\sigma_{\min}^2,\sigma_{\max}^2]\)(\(\sigma_j\) 为奇异值)。\(|1-\eta\sigma^2|\) 是 \(\sigma^2\) 的线性函数,最大值在区间端点取得;取最优步长 \(\eta = 2/(\sigma_{\min}^2+\sigma_{\max}^2)\),两端取值相同,故最坏收缩因子为
因此需要 \(O(\kappa^2\log(1/\varepsilon))\) 步收敛。每步在量子电路上实现,避免了 QPE 的深度要求。
量子随机 Kaczmarz 方法¶
将经典 Kaczmarz 迭代(逐行投影)量子化。每步以概率 \(p_i = \|a_i\|^2/\|A\|_F^2\) 随机选取一行 \(i\)(\(\|\cdot\|_F\) 为 Frobenius 范数),把当前解投影到超平面 \(\{\vec{x}:\langle a_i,\vec{x}\rangle = b_i\}\) 上:
投影的正确性可以直接验证:\(\langle a_i,\vec{x}_{k+1}\rangle = \langle a_i,\vec{x}_k\rangle + (b_i-\langle a_i,\vec{x}_k\rangle) = b_i\)。收敛速度同样可以推导:记 \(\vec{e}_k = \vec{x}_k-\vec{x}^*\),由 \(\langle a_i,\vec{x}^*\rangle = b_i\),投影公式化为 \(\vec{e}_{k+1} = \vec{e}_k - \frac{\langle a_i,\vec{e}_k\rangle}{\|a_i\|^2}a_i\)(投影后的误差与 \(a_i\) 正交);取模方并对 \(i\) 按上述概率平均,得
其中用了 \(\|A\vec{e}\| \ge \|\vec{e}\|/\|A^{-1}\|\)(\(A\) 可逆),并记 \(\kappa_F = \|A^{-1}\|\,\|A\|_F\);由 \(\|A\|_F \ge \|A\|\) 与 \(\|A\|_F^2\le Ns\|A\|^2\)(每个元素满足 \(|a_{ij}|\le\|A\|\),非零元至多 \(Ns\) 个)知 \(\kappa \le \kappa_F \le \sqrt{Ns}\,\kappa\)。故步数为 \(O(\kappa_F^2\log(1/\varepsilon))\),在 \(\kappa_F\) 与 \(\kappa\) 同阶的常见情形即 \(O(\kappa^2\log(1/\varepsilon))\) 步。在相应的量子化方案中,每步可用 \(O(n)\) 个门实现。
绝热量子线性求解¶
将 \(A\vec{x} = \vec{b}\) 编码为哈密顿量基态问题,用绝热演化(adiabatic evolution)求解。一个正确的构造是取
(\(|b\rangle\) 已归一化)。验证其基态恰为解:\(I-|b\rangle\langle b|\) 是到 \(|b\rangle\) 正交补空间的投影算子,半正定,且 \((I-|b\rangle\langle b|)|y\rangle = 0\) 当且仅当 \(|y\rangle \propto |b\rangle\)。于是 \(H_P \succeq 0\),且 \(H_P|x\rangle = 0\) 当且仅当 \(A|x\rangle \propto |b\rangle\);由 \(A\) 可逆即得 \(|x\rangle \propto A^{-1}|b\rangle\),且该零能基态非简并。从基态易于制备的初值哈密顿量 \(H_0\) 出发,绝热演化 \(H(s) = (1-s)H_0 + sH_P\) 到 \(H_P\);由绝热定理,只要沿途谱隙不小于 \(1/\mathrm{poly}(\kappa)\),运行时间就是 \(\mathrm{poly}(n,\kappa)\)。
复杂度对比¶
方法 |
门数 |
\(\kappa\) 依赖 |
精度依赖 |
后选择 |
适用设备 |
|---|---|---|---|---|---|
HHL (2009) |
\(O(s^2 \kappa^2 \text{poly}(n) / \varepsilon)\) |
\(O(\kappa^2)\) |
\(O(1/\varepsilon)\) |
有 |
容错 |
Childs (2017) |
\(O(s n \kappa \log(\kappa/\varepsilon))\) |
\(O(\kappa)\) |
\(O(\log(1/\varepsilon))\) |
无 |
容错 |
QSVT (2019) |
\(O(s n \kappa \log(\kappa/\varepsilon))\) |
\(O(\kappa)\) |
\(O(\log(1/\varepsilon))\) |
无 |
容错 |
VQLS (2020) |
\(O(C_A D T_{\text{opt}})\) |
隐式 |
\(O(1/\varepsilon^2)\) |
无 |
NISQ |
Schrödingerization |
\(O(s n \cdot e^{|A|t})\) |
隐式 |
\(O(\log(1/\varepsilon))\) |
有 |
容错 |
梯度下降 (2025) |
\(O(s n \kappa^2 \log(1/\varepsilon))\) |
\(O(\kappa^2)\) |
\(O(\log(1/\varepsilon))\) |
无 |
NISQ |
各 suppressed 因子(如 \(\mathrm{poly}\log\) 项)随访问模型与实现细节而有所不同;HHL 一行采用振幅放大后的形式,门数与重复次数的详细分解见本系列的 HHL 一篇。
实验进展¶
年份 |
平台 |
方法 |
规模 |
|---|---|---|---|
2019 |
IBM |
HHL |
4 量子比特 |
2020 |
光量子 |
HHL |
4 量子比特 |
2021 |
超导 |
VQLS |
6 量子比特 |
2023 |
超导 |
变分方法 |
12 量子比特 |
2024 |
GPU 模拟 |
QSVT |
20+ 量子比特 |
关键限制与开放问题¶
输入问题(Input Problem)¶
将经典向量 \(\vec{b}\) 加载为量子态 \(|b\rangle\) 需要 \(O(N)\) 门——抵消量子加速。对结构化 \(\vec{b}\)(稀疏、低秩),可用稀疏态制备降低代价。
输出问题(Output Problem)¶
\(|x\rangle\) 编码了全部解分量,但读出单个 \(x_j\) 需要以常数概率观察到基态 \(|j\rangle\),即 \(O(1/|x_j|^2)\) 次测量;把振幅估计到精度 \(\delta\) 需要 \(O(1/\delta^2)\) 次测量。完整读出所有分量需要 \(O(N)\) 次——与经典相同。
量子优势场景:
计算 \(\langle x|M|x\rangle\)(期望值):\(O(1)\) 次测量
计算 \(\|x\|^2\):\(O(1)\) 次测量
从 \(|x\rangle\) 采样:\(O(\text{poly}(n))\) 次
QRAM 的角色¶
许多"量子线性代数"算法假设数据可通过 QRAM 高效访问。若没有 QRAM,输入/输出的开销可能抵消理论加速。这是当前争议的焦点。
Barren Plateaus(对 VQLS)¶
对深拟设,VQLS 的损失景观存在贫瘠高原(barren plateau)——梯度随规模指数衰减。缓解策略包括局域损失函数、结构化拟设等。
应用场景¶
1. 有限元分析¶
将偏微分方程离散化为 \(A\vec{x} = \vec{b}\),用 QSVT 求解。对稀疏、良态的 \(A\),有望实现指数加速。
2. 机器学习¶
最小二乘回归:\(A = X^T X\),\(\vec{b} = X^T \vec{y}\)。用 VQLS 在 NISQ 设备上求解。
3. 金融¶
Black-Scholes 方程的有限差分离散化产生大型线性系统,量子线性求解器有望加速风险计算。
4. 电力系统¶
潮流计算的核心是求解非线性方程组的线性化版本,每步迭代涉及线性方程组。
总结¶
量子线性求解器从 HHL (2009) 的开创性工作,经过 Childs-Kothari-Somma (2017) 的精度优化,到 QSVT (2019) 的统一框架,再到 VQLS (2020) 的 NISQ 适用方案,经历了五代演进。现代方法将 \(\kappa\) 依赖从 \(\kappa^2\) 改进到 \(\kappa\)(理论最优),并消除了后选择需求。在实际应用中,方法的选择取决于问题规模、精度要求和可用量子硬件。
参考文献:
Harrow, A. W., Hassidim, A., & Lloyd, S. (2009). Quantum algorithm for linear systems of equations. Physical Review Letters, 103(15), 150502.
Childs, A. M., Kothari, R., & Somma, R. D. (2017). Quantum algorithm for systems of linear equations with exponentially improved dependence on precision. SIAM Journal on Computing, 46(6), 1920-1950.
Gilyén, A., Su, Y., Low, G. H., & Wiebe, N. (2019). Quantum singular value transformation and beyond. STOC 2019.
Bravo-Prieto, C., et al. (2023). Variational Quantum Linear Solver. Quantum, 7, 1188.
Jin, S., Liu, N., & Yu, Y. (2022). Quantum simulation of partial differential equations via Schrödingerization. arXiv:2212.13969.
返回目录:量子计算算法教程系列