量子常微分方程（ODE）求解是量子科学计算的核心方向之一。从线性 ODE 的哈密顿量模拟，到非线性 ODE 的 Carleman 线性化和 Schrödingerization，近年来取得了一系列理论突破。本文系统介绍量子 ODE 求解的各方法及其复杂度演进。

问题定义

一般 ODE

$$\frac{d\mathbf{u}}{dt} = f(\mathbf{u}, t), \quad \mathbf{u}(0) = \mathbf{u}_0 \in \mathbb{R}^d$$

分类：

• 线性 ODE：$\dot{\mathbf{u}} = A(t)\mathbf{u} + \mathbf{b}(t)$
• 非线性 ODE：$\dot{\mathbf{u}} = f(\mathbf{u})$（$f$ 含 $\mathbf{u}$ 的非线性项）

经典方法

• 显式/隐式 Euler：$O(N/\varepsilon)$ 步（$N$ 为空间维度，$\varepsilon$ 为目标精度）
• Runge-Kutta（RK4）：$O(N/\varepsilon^{1/4})$ 步
• 多步法（Adams-Bashforth）：$O(N/\varepsilon^{1/p})$ 步

经典下界：对 $d$ 维 ODE，达到精度 $\varepsilon$ 至少需要 $\Omega(\log(1/\varepsilon))$ 次函数评估（信息论下界）。

线性 ODE 的量子求解

情况一：常系数线性 ODE

$\dot{\mathbf{u}} = A\mathbf{u} + \mathbf{b}$，$A$ 为常数矩阵。

解：$\mathbf{u}(t) = e^{At}\mathbf{u}_0 + \int_0^t e^{A(t-s)}\mathbf{b}\,ds$

量子算法：对 $e^{At}$ 施加哈密顿量模拟（Trotter、LCHS、QSVT 等），见哈密顿量模拟相关教程。

复杂度（QSVT）：

$$O\left((\|A\|t + \log(1/\varepsilon)) \cdot C_U\right)$$

其中 $C_U$ 为 $A$ 的块编码代价。

情况二：时变线性 ODE

$\dot{\mathbf{u}} = A(t)\mathbf{u} + \mathbf{b}(t)$。

困难：$A(t)$ 随时间变化，$e^{A(t)}$ 不再简单地可分解。

方法一：Magnus 展开

$$\mathbf{u}(t) = \exp\left(\Omega(t)\right) \mathbf{u}_0$$

其中 $\Omega(t) = \int_0^t A(s)ds + \frac{1}{2}\int_0^t\int_0^s [A(s), A(\tau)]d\tau ds + \cdots$

截断 Magnus 展开到 $p$ 阶，每项可由 Trotter 分解实现。

方法二：Schrödingerization

将 $\dot{\mathbf{u}} = A(t)\mathbf{u}$ 直接转化为薛定谔方程（见 Schrödingerization 教程）。对 $A(t)$ 非 Hermitian 的情况，通过辅助维度实现 Hermitian 化。

复杂度：

$$O\left(\|A\|_{\max} t \cdot \log(1/\varepsilon) \cdot C_U\right)$$

情况三：非齐次线性 ODE

$\dot{\mathbf{u}} = A\mathbf{u} + \mathbf{b}(t)$，$\mathbf{b}(t) \neq 0$。

Duhamel 公式：

$$\mathbf{u}(t) = e^{At}\mathbf{u}_0 + \int_0^t e^{A(t-s)}\mathbf{b}(s)ds$$

量子实现：用量子积分（quantum quadrature）实现卷积积分。将积分转化为受控旋转，用 QPE 提取相位。

非线性 ODE 的量子求解

非线性 ODE 的量子求解是量子科学计算的核心挑战。量子力学是线性的，不能直接模拟非线性动力学。目前有三大类方法。

方法一：Carleman 线性化

核心思想（Liu, Kolda, Krovi, Loureiro, Trivisa, Childs, PNAS 2021; arXiv:2011.03185；Krovi, Quantum 7, 913, 2023；arXiv:2202.01054）：

将非线性 ODE $\dot{u} = f(u)$ 通过引入新变量 $v_k = u^k$（$k = 1, 2, \ldots, K$）线性化为无限维线性 ODE，然后截断到有限维。

对二次非线性 ODE：$\dot{u} = au + bu^2 + c$

引入 $v_1 = u$，$v_2 = u^2$，$v_3 = u^3$，…，则：

$$\dot{v}_1 = av_1 + bv_2 + c$$ $$\dot{v}_2 = 2v_1(av_1 + bv_2 + c) = 2av_2 + 2bv_3 + 2cv_1$$ $$\dot{v}_k = k v_{k-1}(av_1 + bv_2 + c) = kav_k + kbv_{k+1} + kcv_{k-1}$$

这是一个有限带宽的无限维线性 ODE（$v_k$ 只耦合到 $v_{k-1}, v_k, v_{k+1}$）。

截断到 $K$ 阶：忽略 $v_{K+1}$ 以上的项，得到 $K$ 维线性 ODE：

$$\dot{\mathbf{V}} = M_K \mathbf{V} + \mathbf{c}$$

其中 $M_K$ 是 $K \times K$ 三对角矩阵。

复杂度分析（Krovi, Quantum 7, 913, 2023；arXiv:2202.01054）：

Krovi 将矩阵指数的范数特征化引入量子 ODE 算法的运行时分析，替代了先前工作对矩阵可对角化性的要求。对 Carleman 线性化后的非线性 ODE，其误差对精度 $\varepsilon$ 的依赖从多项式改进为对数级：

$$O\left(\text{poly}(\kappa, \log(1/\varepsilon)) \cdot T^{1+o(1)}\right)$$

其中 $T$ 为模拟时间，$\kappa$ 为 Carleman 矩阵的条件数。该算法可处理更广泛类型的稀疏矩阵（包括不可对角化的奇异矩阵），只要矩阵具有负对数范数（即 $\mu(A) < 0$，对应耗散系统）。

局限性：

• Carleman 截断阶数 $K$ 随模拟时间和非线性强度增长，可能指数爆炸
• 条件数 $\kappa$ 可能指数增长
• 仅适用于多项式非线性（$f(u) = \sum_k a_k u^k$）
• 需要耗散条件：非线性强度与耗散之比 $R < 1$（对高维系统，$R$ 可能随空间离散点数增长）

2025 年进一步改进（Costa, Schleich, Morales, Berry, npj Quantum Information 11, 141, 2025；arXiv:2312.09518）：

• 高阶 Carleman 截断：改进 Carleman 线性化的误差界
• 重缩放技术（rescaling）：显著降低代价，避免 Carleman 阶数上的指数增长（否则会阻止 PDE 的量子加速）
• 对线性化方程的高精度求解技术：实现复杂度对误差的对数依赖和对时间的近线性依赖
• 注意：重缩放（max-norm）与离散化稳定性（2-norm）之间可能存在冲突，需在高阶离散化下平衡

方法二：非齐次线性 ODE 的量子求解（LCHS）

核心思想（An, Liu, Lin, PRL 131, 150603, 2023；arXiv:2303.01029）：

对非齐次线性 ODE $\dot{\mathbf{u}} = A\mathbf{u} + \mathbf{b}(t)$，用线性组合哈密顿量模拟（LCHS）框架求解。

方法：Duhamel 公式给出

$$\mathbf{u}(T) = e^{AT}\mathbf{u}_0 + \int_0^T e^{A(T-s)}\mathbf{b}(s)\,ds$$

LCHS 将非酉演化（$e^{At}$，当 $A$ 非反 Hermitian 时）表示为多个酉哈密顿量模拟的线性组合，从而在态制备代价上达到最优。

注意：LCHS 本身处理的是线性非酉动力学，不能直接求解非线性 ODE。要求解非线性 ODE，需先通过 Carleman 线性化将非线性问题转化为高维线性问题，再对线性化后的系统应用 LCHS。这两个步骤分别来自不同论文。

后续改进（An, Childs, Lin, Comm. Math. Phys. 2026；arXiv:2312.03916）：

对线性 ODE 的 LCHS 框架实现指数级精度提升，在矩阵查询的所有参数上达到近最优缩放。

方法三：Schrödingerization + 线性化

核心思想（Jin, Liu, Yu, PRL 133, 230602, 2024；arXiv:2212.13969；及后续多篇论文）：

1. 将非线性 ODE 线性化（Carleman 或 Taylor 线性化）
2. 对线性化后的系统施加 Schrödingerization
3. 用量子电路实现 Hermitian 化后的薛定谔方程

完整流程（以 $\dot{u} = -u + u^2$ 为例）：

1. Carleman 化：引入 $v_k = u^k$，得到 $\dot{\mathbf{V}} = M_K \mathbf{V}$
2. Schrödingerization：$\mathcal{H} = (M_K + M_K^\dagger)/2 \otimes I + (M_K - M_K^\dagger)/(2i) \otimes P$
3. 量子模拟：Trotter/QSVT 实现 $e^{-i\mathcal{H}t}$
4. 后选择 + 测量

方法四：量子 Fourier ODE 求解器

Xiao, Yang, Shu, Du, Song（arXiv:2504.10218, 2025）：

将驱动函数 $f(x)$ 展开为傅里叶级数，简化量子电路中的 oracle 构造。传统量子 ODE 求解器使用振幅估计并要求 $f(x) \in [0,1]$，而 Fourier 方法移除了这一限制。

该方法已在 Navier-Stokes 方程等线性和非线性 PDE 上进行了测试，数值结果与解析解及参考解吻合。

注意：该论文是 2025 年的预印本，复杂度分析和理论保证尚需同行评审验证。

各方法复杂度对比

说明：$\|A\|$ 为矩阵范数，$C_U$ 为块编码单次查询代价，$K$ 为 Carleman 截断阶数，$\kappa$ 为条件数，$n_p$ 为 Schrödingerization 辅助量子比特数。以上复杂度均为理论分析结果，实际门数取决于具体实现和问题规模。

具体例子

例子一：Lotka-Volterra 方程（捕食者-被捕食者）

$$\dot{x} = \alpha x - \beta xy, \quad \dot{y} = \delta xy - \gamma y$$

二次非线性。Carleman 化到 $K$ 阶，得到 $K$ 维线性系统。

对 $\alpha = 1, \beta = 0.1, \delta = 0.075, \gamma = 1.5$，截断 $K = 10$ 已足够精确（对 $T = 10$）。

例子二：Van der Pol 振子

$$\dot{x} = y, \quad \dot{y} = \mu(1 - x^2)y - x$$

三次非线性（$x^2 y$ 项）。Carleman 化需要更高阶截断。

例子三：Lorenz 系统（混沌）

$$\dot{x} = \sigma(y - x), \quad \dot{y} = x(\rho - z) - y, \quad \dot{z} = xy - \beta z$$

二次非线性，但对初值极度敏感（混沌）。量子算法只能在有限时间内有效预测，长期行为因蝴蝶效应不可计算。

代码示例（uniqc 实现）

类型：核心子程序演示 —— 情况一（常系数线性 ODE）的量子算法本质就是哈密顿量模拟。用 uniqc 构建一阶 Trotter 电路 $e^{-i H \Delta t}$, 重复 $N$ 步逼近 $e^{-i H T}$；以 $H = X_0 + X_1 + Z_0Z_1$ 为例与 scipy.expm 精确解比较。

import numpy as np
from scipy.linalg import expm
from uniqc import Circuit

def trotter_circuit(t: float, n_steps: int) -> Circuit:
    """一阶 Trotter: e^{-iHt} ≈ (e^{-iX0·dt} e^{-iX1·dt} e^{-iZ0Z1·dt})^N"""
    dt = t / n_steps
    c = Circuit()
    for _ in range(n_steps):
        c.rx(0, 2 * dt)      # rx(θ) = exp(-i θ/2 X), 故 angle = 2·dt
        c.rx(1, 2 * dt)
        c.zz(0, 1, 2 * dt)   # zz(θ) = exp(-i θ/2 Z·Z)
    return c

I, X, Z = np.eye(2), np.array([[0,1],[1,0]]), np.diag([1,-1])
H = np.kron(X,I) + np.kron(I,X) + np.kron(Z,Z)
T = 1.0
psi0 = np.zeros(4, complex); psi0[0] = 1.0
psi_exact = expm(-1j * H * T) @ psi0

for N in [1, 2, 4, 8, 16, 32, 64, 128]:
    psi_tr = trotter_circuit(T, N).get_matrix() @ psi0
    err = np.linalg.norm(psi_tr - psi_exact)
    F = abs(np.vdot(psi_exact, psi_tr)) ** 2
    print(f'N={N:>3}err={err:.4e}F={F:.6f}')

运行结果：

■ 哈密顿量 H 维度: (4, 4),  演化时间 T = 1.0
  精确解 ψ(T) (前 4 分量):
    [-0.0385-0.5967j -0.    -0.3518j  0.    -0.3518j -0.5788+0.2448j]

■ 用 uniqc Trotter 电路逼近 e^{-i H T} (不同步数 N):
  N (Trotter 步) |   Trotter 误差   | 保真度 F
  ---------------+-------------------+---------
         1       |    7.9919e-01    |  0.478068
         2       |    2.8472e-01    |  0.933483
         4       |    1.2867e-01    |  0.986998
         8       |    6.2725e-02    |  0.996976
        16       |    3.1165e-02    |  0.999261
        32       |    1.5558e-02    |  0.999817
        64       |    7.7758e-03    |  0.999954
       128       |    3.8875e-03    |  0.999989

观察：

1. Trotter 误差近似呈 $\mathcal O(\|H\|^2 T^2 / N)$ 一阶下降：$N$ 翻倍, 误差减半 — 与表格中 8e-1, 2.8e-1, 1.3e-1, 6.3e-2… 的下降一致
2. $N=128$ 步即可达到 $F > 0.9999$, 演示了情况一（常系数线性 ODE）通过 Hamilton 量模拟可在多项式门数下高精度求解
3. 对非线性 ODE（如本节例子里的 Lotka-Volterra）需要先做 Carleman 线性化得到形如 $\dot V = M V$ 的高维线性系统, 然后再对 $M$ 做（可能要先 Hermitian 化的）量子模拟; 本演示故意聚焦最简单的”线性 + Hermitian”情形, 复杂度链路最清

局限性与开放问题

1. Carleman 爆炸

对强非线性（大系数或高次非线性），Carleman 截断阶数 $K$ 可能需要指数大才能保持精度。

2. 混沌系统

对混沌 ODE（如 Lorenz），量子算法不能超越经典算法的长期预测能力——蝴蝶效应是信息论的，而非计算复杂度的。

3. 守恒量

量子模拟天然保持范数（幺正演化），但 ODE 的解不一定守恒范数。需要后选择或归一化处理。

4. 边界条件

ODE 的边界条件在量子编码中可能难以实现（特别是边值问题，而非初值问题）。

总结

量子 ODE 求解已从线性 ODE 的标准哈密顿量模拟，发展到非线性 ODE 的 Carleman 线性化和 Schrödingerization 等多种方法。2023—2025 年的改进（矩阵指数范数方法、重缩放技术、高阶 Carleman 截断）显著改善了复杂度对精度和时间的依赖。然而，非线性 ODE 的量子求解仍是开放问题：Carleman 截断阶数随模拟时间增长、耗散条件 $R < 1$ 在高维系统中可能难以满足、混沌系统的长期不可预测性是主要障碍。

参考文献：

1. Berry, D. W. (2014). High-order quantum algorithm for solving linear differential equations. Journal of Physics A: Mathematical and Theoretical, 47, 105301. [arXiv:1010.2745]
2. Liu, J.-P., Kolda, H. Ø., Krovi, H. K., Loureiro, N. F., Trivisa, K., & Childs, A. M. (2021). Efficient quantum algorithm for dissipative nonlinear differential equations. Proceedings of the National Academy of Sciences (PNAS). [arXiv:2011.03185]
3. Krovi, H. (2023). Improved quantum algorithms for linear and nonlinear differential equations. Quantum, 7, 913. [arXiv:2202.01054]
4. An, D., Liu, J.-P., & Lin, L. (2023). Linear combination of Hamiltonian simulation for nonunitary dynamics with optimal state preparation cost. Physical Review Letters, 131, 150603. [arXiv:2303.01029]
5. An, D., Childs, A. M., & Lin, L. (2026). Quantum algorithm for linear ODEs with near-optimal dependence on all parameters. Communications in Mathematical Physics. [arXiv:2312.03916]
6. Jin, S., Liu, N., & Yu, Y. (2024). Quantum simulation of partial differential equations via Schrödingerization. Physical Review Letters, 133, 230602. [arXiv:2212.13969]
7. Costa, P. C. S., Schleich, P., Morales, M. E. S., & Berry, D. W. (2025). Further improving quantum algorithms for nonlinear differential equations via higher-order methods and rescaling. npj Quantum Information, 11, 141. [arXiv:2312.09518]
8. Xiao, Y., Yang, L., Shu, C., Du, Y., & Song, Y. (2025). A novel quantum Fourier ordinary differential equation solver for solving linear and nonlinear partial differential equations. [arXiv:2504.10218]

返回目录：量子计算算法教程系列