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

常微分方程初值问题的数值解法在自然科学、工程技术、经济和医学等领域中，常常会遇到一阶常微分方程初值问题：(,),,(),y f x y a x b y a y '=≤≤⎧⎨=⎩ (1) 此处f 为,x y 的已知函数，0y 是给定的初始值。

本章讨论该问题的数值解法，要求f 在区域{(,)|,}G x y a x b y =≤≤<∞内连续，并对y 满足Lipschitz 条件，从而初值问题(1)有唯一的连续可微解()y y x =，且它是适定的。

1 几个简单的数值积分法1.1 Euler 方法（1）向前Euler 公式（显式Euler 公式）10(,),0,1,2,,(),n n n n y y hf x y n y y a +=+=⎧⎨=⎩(2) 其中h 为步长。

由此便可由初值0y 逐步算出一阶常微分方程初值问题(1)的解()y y x =在节点12,,x x 处的近似值12,,y y 。

该公式的局部截断误差为2()O h ，是一阶方法。

（2）向后Euler 公式（隐式Euler 公式）1110(,),0,1,2,,(),n n n n y y hf x y n y y a +++=+=⎧⎨=⎩(3) 这是一个隐格式，也是一阶方法。

这类隐格式的计算比显格式困难，一般采用迭代法求解。

首先用向前Euler 公式提供迭代初值，然后迭代计算：(0)1(1)()111(,),(,),0,1,2,n n n n k k n n n n y y hf x y y y hf x y k +++++⎧=+⎨=+=⎩ (4)1.2 梯形方法1110[(,)(,)],2(),(0,1,2,)n n n n n n h y y f x y f x y y y a n +++⎧=++⎪⎨⎪=⎩= (5) 这也是一个隐格式，是二阶方法。

一般也采用迭代法求解。

迭代公式如下：(0)1(1)()111(,),[(,)(,)],0,1,2,2n n n n k k n n n n n n y y hf x y h y y f x y f x y k +++++⎧=+⎪⎨=++=⎪⎩ (6)1.3 改进的Euler 方法11110(,),[(,)(,)],0,1,2,,2(),n n n n n n n n n n y y hf x y h y y f x y f x y n y y a ++++⎧=+⎪⎪=++=⎨⎪=⎪⎩(7) 为了便于上机编程计算，(7)可改写为110(,),(,),0,1,2,,1(),2(),p n n n cn n p n p c y y hf x y y y hf x y n y y y y y a ++=+⎧⎪=+⎪⎪=⎨=+⎪⎪=⎪⎩(8) 该格式是显式，也是二阶方法。

例1 函数f ( t , y ) = t y 在区域D0 = {( t , y ) | 1 ≤ t ≤ 2, −3 ≤ y ≤ 4}
关于y满足Lipschitz条件，相应的Lipschitz常数可取为L = 2
3 存在性定理 定理1 设函数f ( t , y )在凸集D ⊂ R 2中有定义，若存在常数
(7.2.7)
称为显式Runge-Kutta(龙格-库塔 )方法，简称R-K方法，
其中正整数N 称为R-K方法的级，所有ci , ai , bij 都是待定 常数。
根据定义(7.2.7)，N 级R-K方法(7.9)的局部截断误差为
Rn+1 = y( t n+1 ) − y( t n ) − h∑ ci ki
dy 其斜率为 = f ( t0 , y0 ) dt ( t0 , y0 ) 由 点 斜 式 写 出 切线 方 程 dy y = y0 + ( t − t0 ) = y0 + ( t − t0 ) f ( t0 , y0 ) dt ( x0 , y )
0
等步长为h,则t1 - t0 = h, 可由切线算出 y1 ： 则 y1 = y0 + hf ( t0 , y0 ) 按此逐步计算y( tn ), 在tn +1处的值 ： yn+1 = yn + hf ( tn , yn ) y 注意： 这是“ 注意 ： 这是 “ 折 yN 线法” 而非“ 线法 ” 而非 “ 切 线法” 线法 ” 除第一个 点是曲线切线外， 点是曲线切线外 ， 其他点不是切线 y2 而是折线(如右 y1 y0 图所示)。 图所示 。
பைடு நூலகம்
则称数值解法(7.5)为显式方法。否则，称数值解法(7.3) 为隐式方法。

根据微分方程的性质和初始条件，常 微分方程初值问题可以分为多种类型， 如一阶、高阶、线性、非线性等。
数值解法的必要性
实际应用需求
许多实际问题需要求解常微分方程初值问题，如物理、 化学、生物、工程等领域。
解析解的局限性
对于复杂问题，解析解难以求得或不存在，因此需要 采用数值方法近似求解。
数值解法的优势
未来发展的方向与挑战
高精度算法
研究和发展更高精度的算法，以提高数值解的准确性和稳定性。
并行计算
利用并行计算技术，提高计算效率，处理大规模问题。
自适应方法
研究自适应算法，根据问题特性自动调整计算精度和步长。
计算机技术的发展对数值解法的影响
1 2
硬件升级
计算机硬件的升级为数值解法提供了更强大的计 算能力。
它首先使用预估方法（如欧拉方法）得到一个 初步解，然后使用校正方法（如龙格-库塔方法） 对初步解进行修正，以提高精度。
预估校正方法的优点是精度较高，且计算量相 对较小，适用于各种复杂问题。
步长与误差控制
01
在离散化过程中，步长是一个重要的参数，它决定 了离散化的精度和计算量。
02
误差控制是数值逼近的一个重要环节，它通过设定 误差阈值来控制计算的精度和稳定性。
能够给出近似解的近似值，方便快捷，适用范围广。
数值解法的历史与发展
早期发展
早在17世纪，科学家就开始尝 试用数值方法求解常微分方程。
重要进展
随着计算机技术的发展，数值 解法在20世纪取得了重要进展， 如欧拉法、龙格-库塔法等。
当前研究热点
目前，常微分方程初值问题的 数值解法仍有许多研究热点和 挑战，如高精度算法、并行计
软件优化
软件技术的发展为数值解法提供了更多的优化手 段和工具。

则有
yn 1 yn hf ( xn , yn )
( 5.2 ) Euler格式
例5.1 用Euler格式解初值问题
2x y y y y (0) 1
取步长h=0.1.
(0 x 1)
Euler格式的具体形式为
y n 1 y n hf ( x n , y n ) 2 xn yn 0.1( yn ) yn 0.2 xn 1.1 yn yn
计算公式的精度 常以Taylor展开为工具来分析计算公式的精度. 为简化分析，假定yn是准确的，即在 yn y ( xn ) 的前提下估计误差 y ( xn 1 ) yn 1 Euler格式的局部截断误差 由 从而 局部截断误差
f ( xn , yn ) f ( xn , y ( xn )) y '( xn ) y ( xn 1 ) yn 1 y ( xn 1 ) ( yn hf ( xn , yn )) y ( xn 1 ) y ( xn ) hy '( xn )
y ( xn ), y ( xn 1 ), 的近似值 y1 , y2 , , yn , yn 1 ,
相邻两个节点的间距 h xi 1 xi 称为步长，步 长可以相等，也可以不等.本章总是假定h为定数， 称为定步长，这时节点可表示为
xn x0 nh , n 0,1, 2,
由f ( xn 1 , yn 1 ) f ( xn 1 , y ( xn 1 )) f y ( xn 1 , )( yn 1 y ( xn 1 )) f ( xn 1 , y ( xn 1 )) y '( xn 1 )（在xn点Taylor展开） h2 y '( xn ) hy ''( xn ) y '''( xn ) ... 2 3 2 h h 因此yn 1 y ( xn ) hy '( xn ) y ''( xn ) y '''( xn ) 2 4 h f y ( xn 1 , )( yn 1 y ( xn 1 )) 2 h2 h3 y ( xn 1 ) y ( xn ) hy '( xn ) y ''( xn ) y '''( xn ) 2 3!