第二章 常微分方程的初值问题讲解

设解为指数形式 y=Arn，代入递推关系得 上述方程的解为
宽度为 h=1/N 2). 寻找一个递推关系，把 yn 同 {yn-1, yn-2, …} 联系起来。
欧拉法
在 xn 点，将微分方程左边的微分用向前差分公式替代
得到递推关系
局部误差为 O(h2) 全局误差为 NO(h2)≈O(h)
精度太低！
例子：牛顿动力学方程
一维运动的粒子方程为
用欧拉法写出具体的递推关系
练习：二阶 Runge-Kutta 法处理二阶常微分方程
2.4 算法的稳定性
对微分方程进行积分时，一个首要的考虑是所用算 法的数值稳定性，也就是说，舍入误差或数值计算 中的其它误差能被放大的程度。
以欧拉法的一个简单扩展为例
当把这一方法应用到下述问题
这个问题的解析解是 y=e-x 为了启动上面的递推关系，还需要 y1 的值，这可以由泰 勒展开获得
两式相比，有
可选解为 或
若取 得二阶Runge-Kutta法
二阶Runge-Kutta 法与二阶泰勒级数法比较
可以看出，相比二阶泰勒级数法而言，二阶Runge-Kutta法 适用性更广，使用也更为方便
三阶Runge-Kutta 法
最为常用的是四阶Runge-Kutta法
通常认为四阶Runge-Kutta法在效率和精度间达到了最好 的平衡
多步法
为了获得更高的精度，可以将 yn+1 不仅仅是同 yn，而且 还同更早的点，如 yn-1, yn-2 等点相联系 例如: 利用 f 在 xn-1 和 xn 处的值通过线型插值得到 f 在积分区 域的值
代入积分公式得Adams-Bashforth 二步积分法
类似的，Adams-Bashforth 四步积分法为
边界值问题
在自变量的两个端点 上对待求函数施加约 束
本征值问题
含有待定参数，只有 在参数取特定值时， 方程才有非零解
实际问题往往涉及不止一类问题，例如对偏微分方程分离变量 时，同时得到分离开来的初值问题与本征问题。 本章只讨论初值问题
2.1 简单方法
待求问题
求 x=1 处 y 的值
策略
1). 离散化：将区间 [0, 1分]为 N 个等间隔的子区间，每个区间
欧拉法的Matlab 实现
标量形式
for n = 1:N x(n+1) = x(n) + h * p(n) /m p(n+1) = p(n) + h * g(x(n), t(n))
end 矢量形式
f=@(x, t)[x(2)/m, g(x(1), t)]
for n = 1:N z(n+1,:) = z(n,:) + h*f(z(n,:), t(n))
end
泰勒级数法
将 yn+1 在 yn 附近做泰勒展开 而又已知
从而有 这个公式在已知 f 的解析形式时很好用。阶数增大时，变得复杂。
练习：推导三阶泰勒级数法
2.2 多步法和隐式法 在子区间 [xn, xn+1] 将微分方程写成积分形式
关键在于对积分号下的 f(x, y) 取合理的近似 欧拉法实际上将 f(x, y) 近似为 f(xn, yn)
练习：四阶 Runge-Kutta 法处理一阶常微分方程
f = @(x, y)(-x*y); y(1) = 1; for i=1:n
K1 =f(x(i), y(i)); K2 =f(x(i) + h/2, y(i) + h*K1/2 ); K3 =f(x(i) + h/2, y(i) + h*K2/2 ); K4(i)=f(x(i) + h, y(i) + h*K3 ); y(i+1)=y(i) + 1/6*h*(K1 + 2*K2 + 2*K3 + K4 ); end
一个常常使用的具有局域误差 O(h5)的预报校正算法：
显式的 Adams-Bashforth 四步法
Adams-Moulton三步法
练习：简谐振子方程
2.3 Runge-Kutta 法
思想：从欧拉法说起
在子区间 [xn, xn+1] 将微分方程写成积分形式
平均斜率 中值定理
欧拉法就是用 xn 点的斜率近似 [xn, xn+1] 区间的平均斜率 Kave
第二章 常微分方程的初值问题
本章要研究的物理问题：经典动力学中的有序和混沌
本章内容
1
简单方法
2
多步法和隐式法
3
龙格库塔法
4
稳定性问题
54
动力学中的有序和混沌
2.0 引子
常微分方程是物理学中经常碰到，以一维运动粒子为例
一般形式为M个藕合wk.baidu.com一阶方程
常微分方程定解问题的分类
初值问题
给定待求函数在某 个初始点上的值
将上式展开到 O(hm)，得到 m 个方程，而有 m+m(m?1)/2 个待定参数 αi , νij. 所以还有灵活选择的空间 以 m=2 为例 , Runge-Kutta 法公式为
其中
将 K2 做泰勒展开到 O(h2) 项，得 利用 将 yn+1在 xn 点附近做泰勒展开至 O(h2) 项，得
改进的欧拉法 用 xn xn+1 两点斜率的平均值来近似 [xn, xn+1] 平均斜率
K1、K2 分别是 xn xn+1 点处的斜率 但是由于 yn+1 待定，因此需要做“预报”
Runge-Kutta 法 思想：为了提高精度，多取几点的斜率值作为加权平均当作平 均斜率
其中 αi (i = 1,2,...,m) 和 νij (i = 2,3,...,m 且 j < i) 是待定参数
得到其显式解
Adams-Moulton 方法——多步隐式法 利用 fn-1, fn, fn+1 三个值在区间 [xn , xn+1] 对 f 进行二次多
项式插值，得到隐式递推关系
带入积分公式 得
用三次多项式插值，得到相应的三步法为
隐式法很少直接使用，它通常用在预报校正算法中
先通过显式法“预报”一个 yn+1的值，然后利用隐式法来校正 得到一个更精确的值。这样的算法有个优点，它可以持续监控 积分的精度。
多步法有个小麻烦：最初几个格点的启动值需要通过其它的 积分法如欧拉法、泰勒级数法来获得
隐式法 在区间 [xn , xn+1] 上过 fn， fn+1 两点对 f 插值，得
带入积分公式 得
隐式法意味着在每一个积分步都必须解一个方程，会非常耗时。 一个特别简单的情况是，如果 f 对于 y 是线性的，即 f(x,y)=g(x) y，则上面的方程可以为