四阶龙格-库塔法——以2022年高教杯数学建模A题为例

提示：文章写完后，目录可以自动生成，如何生成可参考右边的帮助文档

---

前言

在微分方程的数值求解当中，四阶龙格－库塔法（Runge-Kutta methods）是用于非线性常微分方程的解的重要的一类隐式或显式迭代法。本文以2022年高教杯数学建模A题为例，详细讲解四阶龙格－库塔法的原理及其python实现。

---

一、原理

假设我们求解以下微分方程
$$\{\begin{array}{rl}& \frac{dy}{dt}=f(y,t)\\ & y(0)=y_0\end{array}$$
由中值定理可得
$$y_{n+1}=y_n+f(y_{\epsilon },t_{\epsilon })\Delta t$$
其中，$f(y_{\epsilon },t_{\epsilon })$ 为区间 $(t_n,t_{n+1})$ 上的平均斜率。四阶龙格－库塔法的基本思想，就是找一个容易计算且与 $f(y_{\epsilon },t_{\epsilon })$ 相差极小的值来代替 $f(y_{\epsilon },t_{\epsilon })$ 进行运算。

四阶龙格－库塔法的表达式为
$$y_{n+1}=y_n+\frac{\Delta t}{6}(k_1+2k_2+2k_3+k_4)$$
其中
k 1 = f ( y n , t n ) k 2 = f ( y n + k 1 Δ t 2 , t n + Δ t 2 ) k 3 = f ( y n + k 2 Δ t 2 , t n + Δ t 2 ) k 4 = f ( y n + k 3 Δ t , t n + Δ t ) \begin{aligned} &k_1=f(y_n,t_n) \\ &k_2=f\left(y_n+k_1\frac {\Delta t}2,t_n+\frac {\Delta t}2\right) \\ &k_3=f\left(y_n+k_2\frac {\Delta t}2,t_n+\frac {\Delta t}2\right) \\ &k_4=f\left(y_n+k_3\Delta t,t_n+\Delta t\right) \end{aligned} ​k1​=f(yn​,tn​)k2​=f(yn​+k1​2Δt​,tn​+