流匹配 ODE 求解器 · 图解¶
采样一张图 = 数值积分一条 ODE。这篇笔记把五个经典求解器（Euler、Midpoint、Heun、RK4、Dopri5）拆开讲清楚——每个动画都是真实积分，每张数据图都来自可复现实验。
1 · 引言¶
你训练好了一个 flow matching 模型。采样的时候，代码里总有那么一行 x = solver.step(x, t)——但那一行里到底发生了什么？为什么 Stable Diffusion 3 默认走 28 步而不是 1000 步？为什么有人用 Euler、有人用 Heun、有人用听起来很高级的 DPM-Solver？换一个求解器，到底能省多少算力？
这篇笔记用一个完全可复现的 2D 玩具实验把这些问题讲透。我们训练一个真实的流匹配模型（两份 2D 数据、15000 步、固定种子），然后让五种经典 ODE 求解器——Euler、Midpoint、Heun、RK4、Dopri5——在同一个场上同台竞技。
阅读方式随意，但有两个约定：页面里的每个交互图都是浏览器内的真实数值积分（不是预渲染的动图），每张数据图的数字都来自配套代码 fm_solvers.py 的真实运行。读完你应该能回答：给定一个 NFE 预算，我该选哪个求解器，为什么。
2 · 动机：为什么生成模型工程师要懂求解器¶
三个理由，一个比一个实际：
- NFE 就是钱。 生成一张图的成本 ≈ 网络前向次数 × 单次前向的开销。模型定了之后，单次前向的开销就定了——求解器决定的是同样的质量要付多少次前向。
- 差距是 3 倍量级，而且是免费的。 本文实验里，同一个训练好的模型要达到接近采样地板的质量：Euler 要 NFE=40，RK4 只要 12（第 8 节）。不改模型、不重训练，只换几行求解器代码。
- 结论并不显然。 极低预算下高阶方法反而更差；场没训练好时 Euler 会反超 RK4；「自动挡」的自适应方法也不是全场最优。凭直觉选，很容易选错——所以本文所有结论都配数据。
3 · 预备知识¶
3.1 · 什么是 ODE 与初值问题¶
一条常微分方程（ODE）长这样：
它说的是一件很朴素的事：「一个点在位置 \(x\)、时刻 \(t\) 时，它的速度是 \(f(x,t)\)」。\(f\) 叫做向量场——空间里每一点都插着一支小箭头（页首题图里那些淡色短线就是）。再给定出发点 \(x(0)=x_0\)，就构成一个初值问题（IVP）：从 \(x_0\) 出发、每时每刻顺着箭头走，走出来的轨迹 \(x(t)\) 就是解。
只有极少数 \(f\) 能写出解析解（比如第 5 节用来做测试的 \(dx/dt=x^2\)）。当 \(f\) 是一个神经网络时，解析解无从谈起——只能数值求解，这正是本文的主题。
3.2 · 生成模型是怎么变成一条 ODE 的¶
把「生成」看成「求解 ODE」不是流匹配的发明，而是一条汇流的河：neural ODE / 连续正规化流直接把生成建模为积分一条学出来的场[1]；score-based 扩散模型[3]的随机采样过程被证明存在一条确定性的「probability flow ODE」，两者边际分布完全相同[2]；而 flow matching[4] / rectified flow[5] 干脆跳过 SDE，直接用回归目标训练一个速度场。本文用的就是最后这种（linear/OT 路径）：
训练时在噪声 \(x_0\) 与数据 \(x_1\) 的连线上随机取点，让网络 \(v_\theta(x,t)\) 回归「平均往哪儿流」。采样时反过来，从 \(x_0 \sim \mathcal{N}(0, I)\) 出发解初值问题：
\(t=1\) 时粒子落在哪里，样本就是什么。于是：
生成质量 = 场学得多准（模型误差） × 积分得多准（离散化误差）。 前者是训练的事，后者是求解器的事——本文只讲后者：同一个场，不同求解器的差距。
3.3 · 数值积分入门：离散化与「阶」¶
数值求解的思路只有一句话：把 \(t\in[0,1]\) 切成 \(N\) 小步，每一步用有限次场评估近似真轨迹。两个关键概念：
- 局部截断误差：单独一步犯的错。把方法的更新公式和真解的泰勒展开逐项对比，如果能匹配到 \(dt^p\) 项，这一步的误差就是 \(O(dt^{p+1})\)。
- 全局误差：走完全程累积的错。一共走 \(N=1/dt\) 步，误差大致累加：\(O(dt^{p+1})\cdot(1/dt) = O(dt^p)\)。这样的方法称为 \(p\) 阶方法。
「阶」是理解一切的钥匙：\(p\) 阶方法步长减半、误差降 \(2^p\) 倍，在 log-log 图上是一条斜率 \(-p\) 的直线——第 5 节我们会当场实测。多说一句边界：本文的流匹配场都不「刚性」（stiff），显式方法完全够用；刚性问题需要隐式方法，那是另一本书的事[7]。
3.4 · 记账单位：NFE¶
比较求解器时，横轴永远用 NFE（number of function evaluations，网络前向次数）而非步数：一步 RK4 要 4 次前向，一步 Euler 只要 1 次。「10 步 RK4 vs 10 步 Euler」不公平，「NFE=40 的 RK4 vs NFE=40 的 Euler」才公平——NFE 就是钱，本文所有图的横轴都是它。
4 · 固定步长求解器：一步怎么走¶
预备知识就绪，进入正题。所有「单步法」都在回答同一个问题：已知当前点的速度，往前跨一步 dt，落在哪？ 差别只在于「偷看」几次场，以及怎么组合这些偷看——偷看得越聪明，阶数越高。
与其看公式，不如亲手拆一步。图 1 的场是一个匀速旋转场（真解是圆弧，弯曲程度一目了然）：切换求解器、拖动 dt，看每种方法如何用 1–4 次「偷看」拼出这一步，以及单步误差怎样长大。
NFE 记账¶
| 求解器 | 阶 | 每步评估 | 一次评估的用途 | steps=10 时 NFE |
|---|---|---|---|---|
| Euler | 1 | 1 | 起点速度 k₁ 直接走满一步 | 10 |
| Midpoint | 2 | 2 | 先走半步探路，用中点速度 k₂ 走全程 | 20 |
| Heun | 2 | 2 | 先走满一步探路，用首尾平均 (k₁+k₂)/2 修正 | 20 |
| RK4 | 4 | 4 | 两次中点 + 一次终点，加权 (k₁+2k₂+2k₃+k₄)/6 | 40 |
5 · 收敛阶：为什么高阶「省钱」¶
3.3 节说过：\(p\) 阶方法步长减半、误差降 \(2^p\) 倍，log-log 图上是斜率 \(-p\) 的直线。这不只是理论——配套代码 test_solvers.py 的收敛阶单元测试断言的就是这个比例（Euler 2×、Heun/Midpoint 4×、RK4 16×），而图 2 让你亲眼看到它。
图 2 在你的浏览器里现场求解 \(dx/dt = x^2\), \(x(0)=\tfrac12\)（真解 \(x(1)=1\)），逐个 NFE 画出误差。注意横轴已经换算成 NFE：RK4 每步贵 4 倍，但斜率 −4 很快就把这笔成本挣回来了。
一个有趣的坑
为什么不用更顺手的 \(dx/dt = x\)？因为在线性 ODE 上 Midpoint 和 Heun 的一步更新代数恒等（都精确等于乘 \(1+h+h^2/2\)），两条线会逐位重合、看起来像图画错了——要看出「同阶、误差常数不同」，必须用非线性方程。这个 bug 我们真踩过。
6 · 自适应步长：Dopri5 如何自己配速¶
固定步长有个尴尬：场平缓的地方步子迈小了浪费，弯急的地方步子迈大了翻车。Dormand–Prince 5(4)[6]——SciPy RK45、MATLAB ode45、torchdiffeq 默认求解器的本体——解法是每一步同时算两个答案：7 级评估拼出一个 5 阶解 \(y_5\) 和一个 4 阶解 \(y_4\)，两者之差就是免费的误差估计：
- \(\mathrm{err} \le 1\)：接受这一步，且末级导数下步复用（FSAL，7 级只花 6 次新评估）；
- \(\mathrm{err} > 1\)：拒绝，缩小 dt 重来（这一步的评估白花）；
- 无论接受与否，按
dt ×= clip(0.9·err^(−1/5), 0.2, 5)更新步长（经典的步长控制器，安全系数 0.9，clip 防震荡[7]）。
图 3 在 8-Gaussians 的真实流场上跑一遍（512 个粒子共享步长，和 fm_solvers.py 的实现一致）。这个场前段近似直线、临近 \(t=1\) 时粒子要挤进 σ=0.08 的窄峰——看 dt 如何自动「前松后紧」。
Butcher 表（Dormand–Prince 5(4)）¶
这些系数就是 fm_solvers.py 里 _DP_A / _DP_B5 / _DP_B4 的来源。第 \(i\) 行表示第 \(i\) 级在 \(t + c_i\cdot dt\) 处、用前面各级的 \(a_{ij}\) 线性组合位置去评估场：
| 0 | ||||||
| 1/5 | 1/5 | |||||
| 3/10 | 3/40 | 9/40 | ||||
| 4/5 | 44/45 | −56/15 | 32/9 | |||
| 8/9 | 19372/6561 | −25360/2187 | 64448/6561 | −212/729 | ||
| 1 | 9017/3168 | −355/33 | 46732/5247 | 49/176 | −5103/18656 | |
| 1 | 35/384 | 0 | 500/1113 | 125/192 | −2187/6784 | 11/84 |
| b5 | 35/384 | 0 | 500/1113 | 125/192 | −2187/6784 | 11/84 |
| b4 | 5179/57600 | 0 | 7571/16695 | 393/640 | −92097/339200 | 187/2100 |
注意 \(b_5\) 恰好等于最后一级的 \(a\) 行——所以「本步的 5 阶解」正是「下一步的第 1 级评估点」，这就是 FSAL（First Same As Last）省一次评估的来历。b₄ 行末尾还有一个 1/40（第 7 级权重）。
rtol。代价是每步同时算 5 阶与 4 阶两个解来拿免费的误差估计；
FSAL 让 7 级评估只花 6 次新的网络前向，被拒绝的步则是白花的评估。7 · 亲手采样：同一个场，五种走法¶
现在把第 4–6 节合起来玩真的：1500 个粒子从 \(\mathcal{N}(0,I)\) 出发，在 8-Gaussians 的精确流场上积分到 \(t=1\)。换求解器、拖步数，看散点云如何从「一圈糊」凝聚成八个尖峰。推荐一条对比路线：Euler·2 步 → Euler·8 步 → Midpoint·4 步（NFE=8）→ RK4·2 步（NFE=8）——同样的 NFE，形状差很多，这就是「怎么花钱」的差别。
8 · 实验：真实训练模型上的对比¶
玩具场好懂，但结论要在真实训练出来的模型上成立才算数。实验设置：在 8-Gaussians 与 Two-Moons 两份 2D 数据上各训练一个流匹配 MLP（15000 步、EMA、固定种子），冻结后让五种求解器从同一批 2000 个评测噪声出发采样（配对比较），横轴 NFE。完整代码与数据见 fm_solvers.py / results.json。
8.1 · 分布质量 W2：很快就到「地板」¶
W2（2-Wasserstein 距离）衡量生成分布与真实数据的距离。图 5 里它先快速下降，然后饱和——剩下的差距来自模型误差与 2000 点采样噪声，再加 NFE 也买不动。求解器的差别只在「多快到地板」。
8.2 · 纯数值精度 endpoint_err：教科书斜率现身¶
把 RK4 走 500 步（NFE=2000）当作「真答案」，量每个求解器终点偏了多少。这个指标剥离了模型误差，只剩离散化误差——于是第 5 节的教科书斜率在训练出来的网络上原样复现（图 6）。
8.3 · 五条值得记住的结论¶
到地板的最小 NFE 差 3 倍
Two-Moons 上达到 ≤1.05× 地板：Euler 40 · Midpoint 16 · Heun 16 · RK4 12。高阶方法用约 1/3 的算力买到同样的分布质量。
极低 NFE 时高阶未必赢
NFE=4 处 8-Gaussians 的 W2：Euler 0.194 < Heun 0.280。步数太少时，高阶法「每步昂贵」的代价压过了「每步精确」的收益——在图 4 里亲手试。
W2 有地板，endpoint 没谎
分布指标会被采样噪声饱和（8G 地板 0.165），饱和后所有求解器看起来一样好。要比数值精度，得用配对的 endpoint 误差。
Dopri5 = 免调参的前沿
给一个 rtol，它自动配速：粗容差 NFE≈37 就有不错的样本，紧容差平滑逼近参考解，全程不需要人肉扫步数网格。
场坏了，谁也救不了
欠训练消融（1000 步训练）：W2 地板从 0.05 飙到 0.33+，加 NFE 无用；Two-Moons 上 Euler 反超 RK4——高阶法只是更精确地积分了一个错误的场。
9 · 实践指南¶
把全文压缩成一张决策表：
| 求解器 | 阶 | NFE/步 | 什么时候用 |
|---|---|---|---|
| Euler | 1 | 1 | 预算极低（NFE≤6）或场很糙（欠训练/蒸馏模型）时的保底选择 |
| Midpoint | 2 | 2 | 低 NFE 甜点区（8–16）；本实验里低预算下最稳的一档 |
| Heun | 2 | 2 | 同为 2 阶，误差常数略不同；扩散社区常用（EDM 的默认二阶法[9]） |
| RK4 | 4 | 4 | 中高预算（NFE≥12）要高保真终点；到参考精度最快 |
| Dopri5 | 5(4) | 自适应 | 不想调步数：设 rtol 即可，自动落在帕累托前沿 |
再往前一步就出了本文的范围：真实的图像/视频生成里，低 NFE 区间还有一批利用扩散 ODE 结构的专用求解器——DDIM[8]、DPM-Solver++[10]、UniPC[11]，以及 SD3/Flux 这类流匹配模型实际采用的「Euler + 移位时间表」[12]。它们与通用求解器的关系，见问答最后一条。
与代码对照：solver 实现在 fm_solvers.py 的 solve_euler / solve_midpoint / solve_heun / solve_rk4 / solve_dopri5；收敛阶与 NFE 记账由 test_solvers.py 把关；本页数据来自 summary.csv。
10 · 问答¶
来自学习过程中的真实提问，持续补充。点开每条查看解答。
图里的 y 轴是「精度」吗？它和生成质量是一回事吗？
图 2（第 5 节）与图 6（8.2 节）的 y 轴是终点误差：求解器算出的终点到真解（图 2 是解析解，图 6 是 RK4-500 参考解）的距离，log 刻度。它衡量的是纯数值积分精度——越低积分越准，log-log 斜率就是收敛阶。
但它不是生成质量。质量看图 5（8.1 节）的 W2（生成分布 vs 真实数据）。两者的分工正是结论三：W2 会被采样噪声「饱和」，饱和后所有求解器看起来一样好；endpoint 误差剥离了模型误差与采样噪声，永远能分出高下。一句话：W2 回答「样本像不像」，endpoint 误差回答「积分准不准」。
rtol / atol 到底是什么？为什么叫「相对」容差？
容差 = 你允许求解器每一步犯多大的错。有两种定法：atol（绝对）——「误差不许超过 0.0001，不管数值多大」，一把固定刻度的尺子；rtol（相对）——「误差不许超过数值本身的千分之一」，按比例要求精度。
为什么要「相对」？量珠穆朗玛峰差 1 米无所谓（占比 0.01%），量身高差 1 米就离谱——同样的误差严不严重，取决于被测的东西多大。但纯相对标准在数值趋近 0 时会失控（「0 的千分之一还是 0」，步长被逼到无限小），所以用 atol 兜底。实际公式：
数值大时 rtol 项主导（按比例），数值趋近 0 时 atol 项接管（绝对底线）。本实验取 atol = rtol/100。
用大白话解释 RK4？
把积分一步想成闭着眼开一段弯路，睁眼看一次路 = 一次网络前向（1 NFE）。
Euler：出发前看一眼方向，闭眼开完全程——路一弯就冲出去（图 4 里低步数的「弥散环」就是这么来的）。
RK4 聪明地看四次：出发一眼（k₁）→ 按 k₁ 开到半路再看（k₂）→ 退回起点按 k₂ 重开到半路又看（k₃）→ 按 k₃ 开到终点看最后一眼（k₄）。然后加权平均 (k₁+2k₂+2k₃+k₄)/6——半路两眼权重最大，因为半路的方向最能代表全程——用这个「平均方向」真正走完这一步。
本质：四次采样拼出整段路的平均速度（数学上是步内做辛普森积分），所以弯路也能贴着走。代价是每步 4 NFE，回报是步长减半、误差降 16 倍——很快回本。在图 1 里切到 RK4 拖大 dt，能亲眼看到这四眼的位置。
用大白话解释 Dopri5？
RK4 是个好司机，但步子多大得你提前定死。Dopri5 = 好司机 + 自动巡航，四个招数：
① 一步看七眼，拼出两个版本的答案——精细版（5 阶）与粗糙版（4 阶）。
② 两版之差 ≈ 这步错了多少——不知道真答案也能自我检查，而且免费。
③ 据此自动调速：两版几乎一样 → 路很直 → 下一步油门加大（dt 最多 ×5）；差太多 → 这步作废重开，步子收紧（最少 ×0.2）。像导航限速：高速 120，进村自动 30——图 3 的「前松后紧 + 红框重试」就是它在工作。
④ 尾巴接头（FSAL）：这步的最后一眼恰好落在下一步起点上，直接复用，七眼只花六眼的钱。
你只需要说一句「我要千分之一的精度」（rtol），哪里快哪里慢它全自己决定——这就是「免调参」的含义。
实际用 diffusion / flow-matching 模型时，首选是不是 Dopri5？
不是。 Dopri5 是「不想调参时的稳妥默认」（SciPy RK45、MATLAB ode45、torchdiffeq 的默认都是它），但生产级生成模型的低 NFE 区间被专用方法占据：
| 场景 | 主流选择 |
|---|---|
| Stable Diffusion 3 / Flux（流匹配） | Euler + 移位时间表，NFE 20–30[12] |
| 扩散模型低 NFE 采样 | DPM-Solver++ / UniPC（NFE 10–20）[10][11] |
| EDM / Karras 系 | Heun + 精心设计的 σ 网格[9] |
| 确定性快速采样的源头 | DDIM[8] |
| 蒸馏模型（LCM / Turbo / Lightning） | Euler 1–4 步 |
| 精确似然 / inversion / 参考解 | Dopri5 ✓[6] |
生产不用它的四个原因：① 低 NFE 打不过专用法——本页数据里 Dopri5 最粗容差也要 NFE≈37 起步，而 Midpoint 在 NFE=16 就到地板；拒绝步在大模型上白烧钱。② DPM-Solver 系利用扩散 ODE 的半线性结构（线性部分解析积分），通用 RK 把网络当黑盒，放弃了这个先验。③ 自适应步长意味着每张图 NFE 不可预测，固定步数延迟恒定、好部署。④ 时间表（σ/timestep schedule）往往比 solver 阶数更重要，固定步长方法可以直接调它。
Dopri5 真正的主场：面对陌生的场、没有调参预算、要高精度（算 likelihood、做 inversion、生成参考解）、NFE 预算 ≥30 且不在乎波动。一句话：免调参的帕累托前沿，不是全场最优——预算紧张的生产采样，前沿在专用低 NFE 方法手里。
参考文献¶
- R. T. Q. Chen, Y. Rubanova, J. Bettencourt, D. Duvenaud. Neural Ordinary Differential Equations. NeurIPS 2018.
- Y. Song, J. Sohl-Dickstein, D. P. Kingma, A. Kumar, S. Ermon, B. Poole. Score-Based Generative Modeling through Stochastic Differential Equations. ICLR 2021.（probability flow ODE）
- J. Ho, A. Jain, P. Abbeel. Denoising Diffusion Probabilistic Models. NeurIPS 2020.
- Y. Lipman, R. T. Q. Chen, H. Ben-Hamu, M. Nickel, M. Le. Flow Matching for Generative Modeling. ICLR 2023.
- X. Liu, C. Gong, Q. Liu. Flow Straight and Fast: Learning to Generate and Transfer Data with Rectified Flow. ICLR 2023.
- J. R. Dormand, P. J. Prince. A family of embedded Runge–Kutta formulae. Journal of Computational and Applied Mathematics, 6(1):19–26, 1980.
- E. Hairer, S. P. Nørsett, G. Wanner. Solving Ordinary Differential Equations I: Nonstiff Problems. Springer, 2nd ed., 1993.（步长控制与刚性问题的标准教材）
- J. Song, C. Meng, S. Ermon. Denoising Diffusion Implicit Models. ICLR 2021.
- T. Karras, M. Aittala, T. Aila, S. Laine. Elucidating the Design Space of Diffusion-Based Generative Models. NeurIPS 2022.（EDM，Heun 采样器与 σ 网格）
- C. Lu, Y. Zhou, F. Bao, J. Chen, C. Li, J. Zhu. DPM-Solver: A Fast ODE Solver for Diffusion Probabilistic Model Sampling in Around 10 Steps. NeurIPS 2022；及 DPM-Solver++, arXiv:2211.01095.
- W. Zhao, L. Bai, Y. Rao, J. Zhou, J. Lu. UniPC: A Unified Predictor-Corrector Framework for Fast Sampling of Diffusion Models. NeurIPS 2023.
- P. Esser, S. Kulal, A. Blattmann, et al. Scaling Rectified Flow Transformers for High-Resolution Image Synthesis. ICML 2024.（Stable Diffusion 3）
本文实验部分（图 5–6 数据、训练配置、指标定义）完全由配套代码 fm_solvers.py / plots.py / test_solvers.py 生成，可一键复现；本页交互图均为浏览器内真实数值积分。