龙格-库塔(Runge-Kutta methods)是用于非线性常微分方程的解的重要的一类迭代法。 由数学家卡尔·龙格和马丁·威尔海姆·库塔于 1900 年左右发明。
# 原理
龙格-库塔方法实质上是简介的使用泰勒级数法的一种技术。
考虑插商 ,根据微分中值定理,存在 ,使得
于是,利用所给方程 得到
这里 称作区间 上的平均斜率。
如果设法在区间 内多预测几个点的斜率值,然后将他们加权平均作为平均斜率 ,则有可能构造出更高精度的计算公式,这就是龙格-库塔的基本思想。
# 显式龙格-库塔
# 四阶
# 经典龙格-库塔方法
在各种 Runge-Kutta 法当中有一个方法十分常用,以至于经常被称为 “RK4” (四阶四级)。
令初值问题表述如下。
则,对于该问题的 RK4 由如下方程给出:
其中
这样,下一个值()由现在的值()加上时间间隔()和一个估算的斜率的乘积所决定。 该斜率是以下斜率的加权平均:
- 是时间段开始时的斜率;
- 是时间段中点的斜率,通过欧拉法采用斜率 来决定 在点 的值;
- 也是中点的斜率,但是这次采用斜率 决定 值;
- 是时间段终点的斜率,其 值用 决定。
当四个斜率取平均时,中点的斜率有更大的权值:
RK4 法是四阶方法,也就是说每步的误差是 阶,而总积累误差为 阶。
注意上述公式对于标量或者向量函数( 可以是向量)都适用。
::: tips 值得注意的是龙格-库塔方法基于泰勒展开,因而他要求所具有的解具有较好的光滑性质。 反之,如果解的光滑性较差,那么使用龙格-库塔得到的解,其精度可能反而不如低阶方法。 因此应结合具体问题的特点选择合适的算法。 :::
# 基尔方法
四阶四级。
其中
# 其他
四阶四级。
其中
# 高阶
高阶相对于四阶而言,多了更多的计算量,却没有相对来说性价比较高的精度进步,故而一般还是采用四阶龙格-库塔,例如,6 级龙格-库塔拥有 5 阶精度,8 级龙格库塔拥有 6 级精度。
# 隐式
# 四阶龙格-库塔
显式 Runge-Kutta 法一般来讲不适用于求解刚性方程。 这是因为显式 Runge-Kutta 方法的稳定区域被局限在一个特定的区域里。 显式 Runge-Kutta 方法的这种缺陷使得人们开始研究隐式 Runge-Kutta 方法,一般而言,隐式 Runge-Kutta 方法具有以下形式:
其中,
在显式 Runge-Kutta 方法的框架里,定义参数 的矩阵是一个下三角矩阵,而隐式 Runge-Kutta 方法并没有这个性质,这是两个方法最直观的区别。
需要注意的是,与显式 Runge-Kutta 方法不同,隐式 Runge-Kutta 方法在每一步的计算里需要求解一个线性方程组,这相应的增加了计算的成本。
# 隐式中点法
该方法为一级二阶方法。
# Hammer-Hollingsworth 方法
这是一种二级四阶的方法,其系数表有
# 一种三级六阶的隐式龙格库塔方法
其系数表为
# 变步长龙格-库塔
步长越小,截断误差也越小,但随着步长的缩小,在一定范围内所要完成的步数就增加了,步数的增加不但引起了计算量的增大,而且可能导致舍入误差的严重积累,因此同积分的数值计算一样,微分方程的数值解法也有一个选者步长的问题。
这就引入了两个问题:
- 怎样衡量和检验计算结果的精度?
- 如何依据的精度处理步长?
从经典的四阶龙格-库塔方法出发,由于公式的局部误差为 ,故有
然后将步长折半, 即步长为 ,从 跨两步到 ,再求得一个近似值 ,每跨一步的截断误差为 ,因此有
比较 和 式,我们可以看到,步长折半以后误差大约只有原误差的 ,既有
于是得到事后估计式:
随后我们检查步长,折半前后的两次计算结果的偏差:
来判断所选的步长是否合适,即:
对于给定的精度 ,
- 如果 ,我们将反复将步长折半进行计算,直到 为止,这时取最终的 为结果。
- 如果 ,我们将反复将步长加倍进行计算,直到 为止,这时在将步长折半一次, 这时取最终的 为结果。
当然,必不可少的引入了更多的计算量,但总体考虑有可能是合算的。
# 单步法的绝对稳定性
可将方程局部线性化为
从 计算一步得
依赖于所选的方法,而且是 的一个近似值。 如果在 的计算有误差 则会引起 的误差 。 如果在 的计算有误差 则会引起 的误差 。
若某方法满足 ,称此方法为绝对稳定;在 复平面上复变量 满足 的区域称之为绝对稳定区域,它与实轴的交称之为绝对稳定区间。
例子:
- 将欧拉法用于 ,有
其中 。 是 的一阶泰勒近似,当 时,欧拉法是绝对稳定的,即欧拉法是一个以 为中心的单位圆,绝对稳定区间为 。
- 将四阶经典龙格-库塔用于 ,有
有 为 的四阶泰勒近似。
# 代码实现
TODO