§6 常微分方程的数值解法
一、 一阶微分方程初值问题的数值解
要求微分方程
在初始条件
下的数值解,就是要在解存在的区间的一系列点
上从初值y0出发,逐个求出的近似值
.
[改进的欧拉方法(预报校正法)] 计算公式为
式中yn表示y(xn)的近似值,表示步长.这里截断误差为
这个方法中,(1)式用折线法提供初值,称为预报公式.(2)式用梯形法给出较精确的值,称为校正公式.合称预报校正公式.
[龙格-库塔方法] 计算公式为
式中
截断误差为
R=O(h5)
手算时按下表自上而下进行.
|
|
|
k |
x0 |
y0 |
|
|
[阿达姆斯方法]
1° 内插公式
这是关于的隐式方程,只要h比较小,可用迭代法求解.
2° 外推公式
这是关于yn+1的显式方程,只要知道前几点的值,就可从公式中直接算出yn+1.
3° 预报校正公式
截断误差都为
R=O(h5)
阿达姆斯方法可以单独采用外推公式计算,每算一个yn+1,只要计算一次f(x,y)的值,计算量比龙格-库塔法小,而截断误差同阶,所以计算量小是外推法优点之一.缺点是前几个yi值不能用此法计算,计算y1 ,y2 ,y3需采用其他方法(一般可用龙格-库塔法).
补充说明 1° 阿达姆斯方法中,计算yn+1的值,要用到所以称为多步法.而改进的欧拉法、龙格-库塔法,只要用到yn,所以称为一步法.
2° 一步法中途改变步长方便,多步法中途改变步长麻烦,因要用一步法重新算出开头几项.
3° 多步法还有一个优点,即可顺便估计出截断误差.
用yn+1(E)和yn+1(I)分别表示用外推和内插公式算得的值,有
式中.若计算时规定的允许误差不超过ε,就将ε与
比较,如果,可继续计算,否则说明误差太大,应当缩小步长.若比ε小得多,不妨将步长放大,提高计算速度.
二、 一阶微分方程组初值问题的数值解
这里为书写简便,只讨论含两个未知函数的微分方程组,含多个未知函数的微分方程组,计算公式类同.
和初始条件
[改进的欧拉方法(预报校正法)]
预报公式
校正公式
[龙格-库塔方法] 计算公式为
其中
先计算再计算yn+1,zn+1.
[阿达姆斯方法的预报校正公式]
预报值(外推公式)
校正值(内插公式)
前几个yj和zj的计算与本节一中的龙格-库塔方法相同.
三、 边值问题
这里只讨论二阶线性常微分方程的边值问题.
(1)
[差分方法]
将区间[a,b]分成n等份,步长,分点x0=a, x1=a+h,
…,
xk=a+kh,
…,xn=b 称为节点.把微商用差商代替,边值问题化为下面差分方程组的求解问题.
式中yk=y(xk),pk=p(xk),qk=q(xk),fk=f(xk)
().上面的差分方程组可看成n+1个未知量y0,y1,y2,…,yn 的线性代数方程组,方程个数也是n+1.整理合并同类项,差分方程组可改写成
()
式中
这是一种特殊形式的线性代数方程组,除了可用消元法,迭代法等求解外,还可用更简便有效的方法—追赶法(见第四章,
§3).
在应用上,还可能遇到下面形式的边界条件:
式中α1,α2,β1,β2是已知常数.这时,上面差分方程组中相应于边界条件的那两个方程要换成:
[化为初值问题的数值解] 先求出初值问题
的数值解y0(x),再求出满足初始条件y(a)=0,y'(a)=1的相应齐次方程
的数值解y1(x).原边值问题(1)的解可表示为
四、 小参数法
当微分方程含有绝对值很小的参数ε时,将解写成ε的幂级数以求得近似解析解的方法,称为小参数法.下面就常微分方程初值问题作简略介绍.
给定微分方程组
(1)
假设所有作为t的函数是充分光滑的,作为的函数是解析的,并且在 ε=0可以展开成ε的幂级数
系数又都是x1,…,xn的解析函数.
又给定初始条件
t=t0时,xi=xi(ε)=xi(0) +εxi(1) +…
(2)
于是方程组(1)的满足初始条件(2)的解在t与ε的某区域上存在,并且可以展开成ε 的幂级数
(3)
而(3)的任何部分和便是(1)的满足(2)的近似解析解.在具体计算时,只需把级数(3)代入(1),把两边化为ε的幂级数形状,再比较系数.
例 求包含小参数ε的黎卡提方程
满足初始条件
t=0, x=0
的解.
解 设
代入上述微分方程,得
即
令两边ε的同次幂系数相等,得到
另外初始条件可改写为
于是又可列出
求出
而
是所求解的近似表达式(即近似解析解).