已知微分方程
y′(t)=f(t,y),y(t0)=y0
当这个方程很难直接求出解析解时,RK4 可以一步一步近似算出 y(t) 的值。
1. 它在做什么
最简单的欧拉法是:
yn+1=yn+hf(tn,yn)
也就是用当前点的斜率直接往前走一步。
但欧拉法精度比较低,因为它只看了起点处一个斜率。
RK4 更聪明:它在一步长度 h 内,取 4 次斜率样本,然后加权平均,所以精度高很多。
2. 四阶龙格–库塔公式
对
y′=f(t,y)
从 (tn,yn)走到 (tn+1,yn+1),步长为 h,计算:
k1=f(tn,yn)
k2=f(tn+2h,;yn+2hk1)
k3=f(tn+2h,;yn+2hk2)
k4=f(tn+h,;yn+hk3)
然后更新:
yn+1=yn+6h(k1+2k2+2k3+k4)
同时
tn+1=tn+h
可以把 k1,k2,k3,k4 理解为这一步区间内对斜率的四次估计:
- k1:起点斜率
- k2:用 k1 预测到中点后,在中点算一次斜率
- k3:再用 k2 预测中点,再算一次斜率
- k4:用 k3 预测终点,在终点算斜率
最后用
61(1,2,2,1)
这个权重平均,得到更好的步进结果。
. 一个简单例子
考虑初值问题:
y′=y,y(0)=1
它的解析解是
y=et
现在用 RK4 从 t=0 算到 t=0.1,取 h=0.1。
因为 f(t,y)=y,所以:
k1=f(0,1)=1
k2=f(0.05,1+0.05×1)=1.05
k3=f(0.05,1+0.05×1.05)=1.0525
k4=f(0.1,1+0.1×1.0525)=1.10525
于是
y1=1+60.1(1+2×1.05+2×1.0525+1.10525)
所以
y1=1+60.1×6.31025=1.105170833⋯
事实上,这种方法是对欧拉法(RK1)的一种改进,在适当增加复杂度后获得了更高的精度,并且这种方法也能进一步推广到更高阶,但是RK4是一个复杂度-精确度平衡点。