Python中四阶龙格-库塔法(RK4)步长h与误差异常的技术问询

这是截断误差和舍入误差相互作用的结果，再加上初始条件的微小误差，共同导致了这个看似反常的现象，具体分析如下：

1. 先明确微分方程的解析解 #

你的微分方程是一阶线性非齐次方程：y' = 2x² -4x + y，结合初始条件x=1时y=0.7182818，对应的解析解为：

import math
def exact_solution(x):
    return -2 * x**2 + math.exp(x)

在目标点x=3处的真实值是exact_solution(3) ≈ 2.0855369232。

2. RK4的误差来源 #

RK4方法的总误差由两部分核心因素构成：

• 截断误差：来自用泰勒展开近似微分方程的固有误差，每步截断误差为O(h⁵)，累计后整体截断误差是O(h⁴)——理论上步长h越小，截断误差应该越小。
• 舍入误差：来自浮点数计算的精度限制（Python采用双精度浮点数，约15-17位有效数字），每一步计算都会引入微小的舍入误差，步数越多（h越小，总步数N=(b-x0)/h越大），累计的舍入误差就越显著。

3. 反常现象的具体原因 #

当你测试不同步长时，两种误差的主导地位在变化，叠加后出现了看似矛盾的结果：

• h=0.1：仅20步计算，舍入误差累计极小（≈20×1e-16=2e-15），此时截断误差主导（≈(0.1)^4=1e-4）。但巧合的是，你给出的初始值0.7182818比真实初始值0.718281828459小了约2.8e-8，这个初始误差的传播方向与截断误差相反，两者部分抵消，导致总误差看起来比更小的步长还小。
• h=0.01和h=0.001：步数分别为200和2000步，截断误差确实更小（分别为1e-8和1e-12），但此时舍入误差累计量上升，更关键的是，截断误差与初始误差的传播方向相同，误差叠加后，总误差反而超过了h=0.1的情况。
• h=0.0001：步数达到20000步，舍入误差累计为2e-12，但此时截断误差已经极小（(0.0001)^4=1e-16），几乎可以忽略，总误差主要由舍入误差主导，而这个量级远小于h=0.01时的截断误差，所以误差再次降低。

4. 验证建议 #

你可以通过以下操作验证这个结论：

• 使用更精确的初始值0.718281828459045，消除初始误差的干扰后，再对比不同步长的误差。
• 计算每个步长下数值解与解析解的绝对误差，绘制误差随h变化的曲线，你会看到误差先因截断误差减小而降低，随后因舍入误差累计而上升，当h继续减小到截断误差远小于舍入误差时，误差会趋近于舍入误差的量级（约1e-12）。

内容的提问来源于stack exchange，提问作者user161524