内容发布更新时间 : 2024/10/21 13:26:26星期一 下面是文章的全部内容请认真阅读。

i.常微分方程初值问题数值解法

常微分方程初值问题的真解可以看成是从给定初始点出发的一条连续曲线。差分法是常微分方程初值问题的主要数值解法，其目的是得到若干个离散点来逼近这条解曲线。有两个基本途径。一个是用离散点上的差商近似替代微商。另一个是先对微分方程积分得到积分方程，再利用离散点作数值积分。 i.1 常微分方程差分法

考虑常微分方程初值问题：求函数u(t)满足

du?f(t,u), 0?t?T （i.1a） dtu(0)?u0 (i.1b)

其中f(t,u)是定义在区域G: 0?t?T, u??上的连续函数，u0和T是给定的常数。我们假设f(t,u)对u满足Lipschitz条件，即存在常数L使得

f(t,u1)?f(t,u2)?Lu1?u2, ?t?[0,T]; u1,u2?(??,?) （i.2） 这一条件保证了（i.1）的解是适定的，即存在，唯一，而且连续依赖于初值u0。

通常情况下，（i.1）的精确解不可能用简单的解析表达式给出，只能求近似解。本章讨论常微分方程最常用的近似数值解法-差分方法。先来讨论最简单的Euler法。为此，首先将求解区域[0,T]离散化为若干个离散点： 0?t0?t1?其中tn?hn，h?0称为步长。

在微积分课程中我们熟知，微商（即导数）是差商的极限。反过来，差商就是微商的近

?tN?1?tN?T （i.3）

似。在t?t0处，在（i.1a）中用向前差商

u(t1)?u(t0)du代替微商，便得

dth u(t1)?u(t0)?hf(t0,u(t0))??0

如果忽略误差项?0，再换个记号，用ui代替u(ti)便得到 u1?u0?hf(t0,u0) 一般地，我们有

Euler 方法：un?1?un?hf(tn,un), n?0,1,

从(i.1b) 给出的初始值u0出发，由上式可以依次算出t1,,N?1 （i.4）

,tN上的差分解u1,,uN。

下面我们用数值积分法重新导出 Euler法以及其它几种方法。为此，在区间[tn,tn?1]上

积分常微分方程（i.1a），得

u(tn?1)?u(tn)??tn?1tn f(t,u(t))dt （i.5）

用各种数值积分公式计算（i.5）中的积分，便导致各种不同的差分法。例如，若用左矩形公式就得到 Euler 法（i.4）。如果用右矩形公式，便得到下面的： 隐式Euler方法：un?1?un?hf(tn?1,un?1), n?0,1,类似地，如果用梯形公式，就得到 改进的Eule方法r un?1?un?,N?1 （i.6）

h2f[t,?)f?1t(n(unn?1nu,)?]n, (i.7) 0,1,?N,1当f(t,u)关于u是非线性函数的时候，不能由（i.6）或 (i.7) 从un直接算出un?1，称这一类方法为隐式，通常采用某种迭代法求解。例如，将一般的隐式方法写成

un?1?F(tn,un,un?1) (i.8) 则可以利用如下的迭代法由un算出un?1：

k?1k?un?F(t,u,unnn?1), k?0,1, ??10un?1?un? (i.9)

关于k的迭代通常只需进行很少几步就可以满足精度要求了。

为了避免对隐式方法进行迭代的麻烦，比如说对于改进的Euler方法（i .7），可以采用某种预估法近似算出f(tn?1,un?1)，然后再用（i .7）作校正，这就导致所谓预估校正法。下面给出一个例子：

??f(tn,un)??un???预估： ?un?1?un?1?2hun???u??f(t,u) （i.10） ? n?1n?1?n?1??校正：u?u?h(u??u?)n?1nn?1n??2这是一个多步法，即计算节点tn?1上的近似值un?1时，除了用到前一点的近似值un之外，还要用到un?1，甚至可能用到un?2,。而用前面的各种Euler法计算节点tn?1上的近似值un?1时，只用到un，这样的方法称之为单步法。

下面给出另一个多步法的例子。在区间[tn,tn?2]上积分（i.1a），得

u(tn?2)?u(tn)??tn?2tnf(t,u(t))dt

用Simpson公式(即把被积函数看作二次函数)近似计算积分，便得到 Miln法：e un?2?un?h3fn(?f4?n1?f?n2),n? （i.11） 0,1,N? , 2用多步法（i.10）或（i.11）计算时，必须先用某种单步法由u0计算出u1，称为造表头。然后再逐次算出u2,u3,,uN。

一般说来，多步法比Euler法等简单的单步法精度要高一些。下面我们讨论一类所谓

Runge-Kutta法。他们是单步法，但是其精度可以与多步法比美。最常用的是下面的标准Runge-Kutta法和Gear法：

标准Lunge-Kutta法

?K1?f(tn,un)?hK1h?K2?f(tn?,un?)22??hK2hK3?f(tn?,un?)?（i.12） 22?K4?f(tn?h,un?hK3)??h?un?1?un?(K1?2K2?2K3?K4)6?

?K1?hf(tn,un)??K2?hf(tn?h,un?K1)22??h?111?)K1?(1?)K)?K3?hf(tn?,un?(?2 Gear方法： ? (i.13) 2222??K4?hf(tn?h,un?(?1)K2?(1?1)K)322???un?1?un?1(K1?2(1?1)K2?2(1?1)K?3K)4?622?从几何上， Runge-Kutta法可以粗略地解释为：在区间[tn,tn?1]中选取若干个点(可以重复)tn??1??2???k?1??k?tn?1，仅仅利用在区间[tn,tn?1]内可以得到的所有信息，依

然后把它们组合起,Kk，

次给出函数f(t,u(t))在这些点上尽可能精确的的近似值K1,K2,来，尽可能精确地近似计算（i.5）中的积分。下面介绍Runge-Kutta法的一般构造方式。