第 16 章 · 常微分方程积分器选型

解 $\dot{\mathbf{y}} = \mathbf{f}(t, \mathbf{y})$：数值方法一大家族，选型的核心问题只有两个——刚性吗、要精度还是要速度。这一章是 ODE 求解器的导购指南。

非刚性首选：RK 与多步

• 四阶 Runge-Kutta：教科书经典，每步 4 次函数求值，固定步长的简单场景够用。
• 自适应 RK（RKCash-Karp/Dormand-Prince）：内嵌误差估计自动调步——默认选择。精度要求 $10^{-8}$ 量级的光滑问题，它几乎总是对的。
• Bulirsch-Stoer：外推法，光滑问题超高精度（$10^{-14}$）的冠军，但遇到不可导/尖峰就翻车。

刚性：换一整套逻辑

刚性 = 解里含快慢悬殊的时间尺度（化学反应、电路、化学-扩散耦合）。显式法被迫用最慢尺度的步长，算一万年。诊断信号：步长被稳定性（不是精度）卡死，减小误差容差反而更慢。专药是隐式方法：

每步解（非线性）方程组，换无条件（或极大）稳定域——backward Euler、Gear/BDF 家族（变阶变步隐式多步）、隐式 RK（Radau）。隐式法的代价是每步要解雅可比方程组，刚性场景下这笔账划算。

工程决策表

• 光滑非刚性 + 常规精度：自适应 RK（DOPRI）。
• 超高精度 + 非常光滑：Bulirsch-Stoer。
• 刚性（或疑似）：BDF/Radau，提供雅可比更好。
• 二阶方程（天体力学）：辛积分器（蛙跳/速度 Verlet）保长期能量行为——普通 RK 长期能量漂移，辛格式不漂。

老师的提醒

• 长时间积分（轨道、分子动力学）精度指标不是单步误差而是守恒量漂移——这时辛格式或 Symplectic 处理优先于经典精度阶。
• 事件定位（如"何时相撞"）用求解器的 dense output + 求根，不要靠加密步长扫描。

自检

• 刚性的诊断信号是什么？
• 显式与隐式方法的代价各花在哪？
• 什么场景必须辛积分器？

更深一层：ODE 求解器的"体检指标"

选 ODE 求解器前先给问题做体检，四个指标定终身：刚性比（最快/最慢时间尺度之比，> 1000 即刚性）、光滑性（右端函数 $f$ 可微到几阶——外推法的前提）、守恒结构（能量/动量是否该守恒——辛算法的适用判据）、事件（是否需要在特定状态处停下——dense output 的需求）。scipy.solve_ivp 的方法参数（RK45/DOP853/Radau/BDF/LSODA）正是这四指标的映射：LSODA 自动在非刚/刚之间切换（用刚性检测），是"不确定时"的懒人选择，代价是切换开销与黑盒性。

刚性的识别还有一个实操信号值得记：步长被稳定性卡死时，减小容差反而更慢——精度要求提高 → 容差收紧 → 显式法被迫用更小步长。出现这个悖论就该换隐式（BDF/Radau）并提供雅可比（解析或稀疏模式），速度往往快两个量级。

推导细读：辛积分器的"长期保真"

经典 RK 积分单步精度高，但长期积分保守系统（行星轨道、分子动力学）时能量缓慢漂移——每步的小误差有系统性偏置，累积不可逆。辛积分器（蛙跳、速度 Verlet）不追求单步最准，而是精确保真离散版的哈密顿结构：能量在真值附近有界震荡而不漂移。数学身份：它是对"修正哈密顿量"（back error analysis 的产物）精确积分——轨迹严格守恒一个略偏离的能量泛函。工程含义：$10^9$ 步的分子动力学、$10^6$ 年的天体演化，只有辛格式（或其对称变体）可信。选型口诀：短时高精度用 RK 类，长时保结构用辛类——两句话划分 ODE 应用的半壁江山。

在真实世界里：事件定位与参数扫描的工程技巧

两个高频工程需求的标准解法：事件定位（"轨迹何时到达某阈值"）——用求解器的 dense output + 布伦特求根精确夹逼时刻，切忌靠加密步长扫描（既慢又不准）；参数扫描（"解对参数的依赖"）——暖启动（上一参数的解作为下一参数初值）+ 雅可比连续性，比逐参数独立求解快数倍。这两个技巧的共同哲学：利用解的连续结构，不要把每次求解当独立问题——这也是任何仿真工程师与"调包侠"的分水岭。

专题深入：刚性的"诊断信号"与隐式方法家族

刚性的诊断信号值得做成体检表：步长被稳定性卡死（减小容差反而更慢——精度要求提高 → 容差收紧 → 显式法被迫用更小步长）、快慢时间尺度悬殊（化学反应的纳秒与小时、电路的快充慢放）、$\lambda_{max}/\lambda_{min}$ 巨大（雅可比特征值比）。刚性 = 显式方法的死刑——隐式方法（backward Euler、Gear/BDF 家族、隐式 RK Radau）的无条件（或极大）稳定域是唯一出路。隐式法的代价（每步解非线性方程组）在刚性场景下划算。

数学细节展开：辛积分器的"修正哈密顿量"思想

辛积分器（蛙跳/速度 Verlet）不追求单步最准，而是精确保真离散版的哈密顿结构：能量在真值附近有界震荡而不漂移。数学身份：它是对"修正哈密顿量"（back error analysis 的产物）精确积分——轨迹严格守恒一个略偏离的能量泛函。$10^9$ 步的分子动力学、$10^6$ 年的天体演化，只有辛格式可信——普通 RK 长期能量漂移。长时间积分的精度指标不是单步误差而是守恒量漂移——这是辛几何数值分析的纲领。

工程案例细节：ODE 求解器的"体检指标"决策表

ODE 求解器的体检指标决策表：刚性比 > 1000 → 隐式（BDF/Radau）；光滑非刚性 + 常规精度 → 自适应 RK（DOPRI 默认）；超高精度 + 非常光滑 → Bulirsch-Stoer；二阶保守系统（天体力学、分子动力学）→ 辛积分器；不确定时 → LSODA（自动切换非刚/刚，代价是切换开销与黑盒性）。事件定位（轨迹何时到达某阈值）用 dense output + 求根，不靠加密步长扫描。