Runge–Kutta 方法

形式陈述 ​

作为常微分方程时间步进方法公理库常微分方程时间步进框架Time-stepping method for ODE · ODE integrator · Numerical time integration将常微分方程初值问题离散为时间网格上的数值状态与一步或多步更新，并分离流映射、误差层和退出状态。，一个 $s$-stage Runge–Kutta 格式对初值问题 $y^{\prime }=f(t,y)$ 在步长 $h_n$ 内计算

$$k_i=f(t_n+c_i{h}_n,y_n+h_n\sum _{j=1}^s{a}_{ij}{k}_j),i=1,\dots ,s,$$

再更新

$$y_{n+1}=y_n+h_n\sum _{i=1}^s{b}_i{k}_i.$$

系数用 Butcher tableau 记录：

$$\begin{array}{cc}c& A\\ & b^{\mathsf{T}}\end{array}.$$

通常要求内部一致性

$$c_i=\sum _{j=1}^s{a}_{ij}.$$

若 $a_{ij}=0$ 对所有 $j\ge i$ 成立，stage 可按 $i=1,\dots ,s$ 依次直接计算，方法是显式 RK。若 $A$ 为下三角且对角可非零，得到 diagonally implicit RK，每个 stage 依次解一个状态维数为 $d$ 的隐式方程；一般稠密 $A$ 则把全部 $sd$ 个 stage 耦合成一个非线性系统。stages 因而不是任意并行、彼此独立的函数采样，其依赖由 $A$ 精确规定。

阶数来自把数值一步映射与精确流的多元Taylor 展开公理库多元 Taylor 定理Multivariable Taylor theorem用各阶导数给出多元函数的局部多项式展开及余项控制。逐项匹配。在 $c=A\mathbf{1}$ 约定下，一至三阶的必要阶条件依次包括

$$\sum _i{b}_i=1,$$$$\sum _i{b}_i{c}_i=\frac{1}{2},$$

以及

$$\sum _i{b}_i{c}_i^2=\frac{1}{3},\sum _{i,j}{b}_i{a}_{ij}{c}_j=\frac{1}{6}.$$

第一式给一致性；满足到第 $p$ 阶的全部条件，才可在相应光滑性与误差传播条件下得到局部一步缺陷 $O(h^{p+1})$ 和固定终点上的全局误差 $O(h^p)$。高阶条件按不同导数组合分叉，可由 rooted trees 系统组织；基础算法不能凭 stage 数或几条低阶等式猜测更高阶。

对测试方程 $y^{\prime }=\lambda y$，令 $z=h\lambda$。消去 stages 得稳定函数

$$R(z)=1+z{b}^{\mathsf{T}}(I-zA{)}^{-1}\mathbf{1},y_{n+1}=R(z)y_n.$$

显式 RK 的 $A$ 严格下三角，$R$ 是有限次多项式；隐式 RK 通常得到有理函数。$R$ 的具体衰减区域与方法阶数不同，高阶显式方法仍只拥有有限稳定区域。

算法输入应包括 IVP、时间网格和 tableau；自适应误差控制与隐式 stage 求解再按需提供各自容差和工作预算。输出包括状态、函数与 Jacobian 求值数、非线性迭代数及退出状态。显式 $s$-stage 方法直接实现时每步至多需要 $s$ 次新的函数求值；若方法具有 FSAL 性质且前一步成功接受，可以复用端点导数并少算一次。向量组合成本通常为 $O(s^{2}d)$，稀疏 $A$ 可相应减少；隐式方法还需反复求解 $sd$ 维耦合系统或 $s$ 个顺序 stage 系统，成本取决于 Jacobian 和线性代数结构。

固定网格显式方法在所有 stages 返回有限值并完成更新后结束该步。隐式方法还必须让 stage 方程的尺度化残差达到容差；若非线性求解失败、函数值非有限、时间加法不再前进或预算耗尽，必须返回失败。嵌入式误差估计、步长接受与控制器属于自适应 RK 的后续层，不应隐含在某个 tableau 名称中。

直觉

Euler 只在一步起点读取一个斜率，RK 方法则在同一步中安排多个试探位置。早期 stage 预测轨道会走向哪里，后续 stage 在这些预测位置重新读取斜率，最终按 $b_i$ 汇总。tableau 是这张依赖图的紧凑坐标，而不是一串可随意互换的系数。

阶条件要求这些试探斜率在所有低阶导数结构上共同复制精确流。两个方法可以使用相同 stage 数并具有不同阶，也可以阶数相同却拥有不同误差常数、存储需求和稳定函数；“更多 stage”不是单调的质量标签。

例子与边界

显式 midpoint 的 tableau 为

$$\begin{array}{ccc}0& 0& 0\\ \frac{1}{2}& \frac{1}{2}& 0\\ & 0& 1\end{array},$$

即

$$k_1=f(t_n,y_n),k_2=f(t_n+\frac{h}{2},y_n+\frac{h}{2}{k}_1),y_{n+1}=y_n+h{k}_2.$$

Heun 二阶法则取

$$\begin{array}{ccc}0& 0& 0\\ 1& 1& 0\\ & \frac{1}{2}& \frac{1}{2}\end{array}.$$

它先用 $k_1$ 预测终点斜率，再平均起终斜率。两者都有两个 stages，且满足二阶条件，但采样位置、权重和局部误差常数不同；“同阶”不表示逐步结果相同。

classical RK4 使用

$$\begin{array}{ccccc}0& 0& 0& 0& 0\\ \frac{1}{2}& \frac{1}{2}& 0& 0& 0\\ \frac{1}{2}& 0& \frac{1}{2}& 0& 0\\ 1& 0& 0& 1& 0\\ & \frac{1}{6}& \frac{1}{3}& \frac{1}{3}& \frac{1}{6}\end{array}$$

并对测试方程给出

$$R(z)=1+z+\frac{z^2}{2}+\frac{z^3}{6}+\frac{z^4}{24}.$$

这正是 $e^z$ 到四次项的截断，验证了 RK4 在线性测试方程上的四阶匹配。一般光滑非线性问题的四阶准确性仍来自前述完整的 rooted-tree 阶条件，不能由这一条稳定多项式单独推出。它仍是多项式，负实轴上的稳定区间有限；遇到很快衰减的刚性模式时，为保持稳定而被迫选择的步长可能远小于精度本身所需，增加阶数并未消除显式稳定限制。

隐式 RK 的 stages 也不能简单各自做一次固定点更新后接受。若 stage 残差没有达到与时间离散误差相容的尺度，代数求解误差会污染声明的阶数；全隐式方法还需处理 stage 之间的耦合 Jacobian。

推论与应用

显式 midpoint、Heun 和 classical RK4 展示同一 stage 模型如何形成不同阶与稳定函数；Gauss、Radau 和 DIRK 则通过隐式系数追求不同的稳定或结构性质。这里列出它们的共同计算骨架，不把方法名称堆成推荐清单。

生产求解器常用共享 stages 的嵌入式公式同时给出高低阶近似，再根据差值控制步长。那一层需要误差尺度、接受/拒绝逻辑和控制器保护；本页只把 $A,b,c$、stage 方程、成本和失败接口固定下来。

参考资料

• John C. Butcher, Numerical Methods for Ordinary Differential Equations, 3rd ed., Wiley, 2016.
• Ernst Hairer, Syvert P. Nørsett, and Gerhard Wanner, Solving Ordinary Differential Equations I: Nonstiff Problems, 2nd rev. ed., Springer, 1993.
• NIST Digital Library of Mathematical Functions, §3.7: Ordinary Differential Equations, accessed 2026.