中国矿业大学：《数值分析》课程教学课件（讲稿，研究生）第九章 常微分方程数值解法

第九章常微分方程数值解法S1Euler方法 2 Runge-Kutta法S3 单步法的绝对稳定性S4线性多步法85一阶方程组与高阶方程的初值问题

-1- 第九章 常微分方程数值解法 §1 Euler方法 §2 Runge-Kutta法 §3 单步法的绝对稳定性 §4 线性多步法 §5 一阶方程组与高阶方程的初值问题

常微分方程数值解法口必要性在工程和科学技术的实际问题中，常需要求解微分方程。只有简单的和典型的微分方程可以求出解析解，而在实际问题中的微分方程往往无法求出解析解。y'=1-2xy如微分方程初值问题其解析解（精确解）为：y(0) = 0(x)=e-" ["e" dt但y(1)、y(1.5)等值却无法直接计算。2-

-2-  必要性 在工程和科学技术的实际问题中，常需要求解微分方程。 只有简单的和典型的微分方程可以求出解析解，而在实际问题 中的微分方程往往无法求出解析解。 1 2 (0) 0 y xy y        如微分方程初值问题 ，其解析解(精确解)为： 但 、 等值却无法直接计算。 y y ( ) ( ) 1 1.5 2 2 0 ( ) x x t y x e e dt    常微分方程数值解法

国什么叫微分方程数值解就是求微分方程解函数y(x)在区间[a,b]上的一系列离散点x,：a=x<x, <x <...<x,=b上函数值y(x,)的近似值y,(k=1,,n)，称y,为问题的数值解。心哪些微分方程的数值解？y'=f(x,y) a≤x≤b一阶方程初值问题[y(a)= yoy"= f(x,y,y') a≤x≤b一高阶方程初值问题(a) = yo , y'(a) = y'[=(x,Ji,J2) (xo)=--方程组初值问题y2 = f(x,y1,y2) yz(xo) = y2P186-定理1-3-微分方程"解析解”存在的条件国

- 3 -  什么叫微分方程数值解 . ( ) [ , : ] k 就是求微分方程 在区间 上的一系列离散点 解函数y x a b x 0 1 2 n a b   x x x x     ( ) ( 1, , ) k k 上函数值 的近似值 ， y x y k n  k 称 为问题的 y 数值解。 哪些微分方程的数值解？ . 0 ( , ) ( ) y f x y a x b y a y           一阶方程初值问题 0 0 ( , , ) ( ) , ( ) y f x y y a x b y a y y a y               高阶方程初值问题 1 1 1 2 1 0 1 2 2 1 2 2 0 2 ( , , ) ( ) ( , , ) ( ) y f x y y y x y y f x y y y x y            方程组初值问题  微分方程 存在的条件 . "解析解" P186 1 定理

$1欧拉方法一 问题已知初值问题y'=f(x,y) a≤x≤bb-a, h=Ny(a)= yo求其解函数y= y(x)在等距节点x,=a+nh(n=0,,N)上的近似值 y,?-4-

-4- §1 欧拉方法 0 ( ) ( 0 ( , ) ( ) , , ) ? n k y f x y a x y y x x a nh b b a h y N N y a y n                 已知初值问题 ， ， 求其解函数 在等距节点 上的近似值 一 问题

[y'=f(x,)a≤x≤b二方法(y(a)= yo显式公式1.Euler方法将初值问题的解函数y=y(x)在x.点Tarloy展开，有：y"()y(x)= y(x,)+y'(x,)(x-x,)+c-x.)-2!而y'=f(x,y),所以 y(x,)= f(x,y(x,), 代入上式:y"()J(x)= y(x,)+ f(xn,J(x,)(x-x,)+ -x2!令x= xn+1 :-y(xn+1) = y(x,)+ f(x,,y(x,)(xnt1 -x,)2!y"(5n)h, 其中5,E(xn,Xn+1)= y(x.)+ f(x,,y(x,)h+2!(5)n，得y(x1)的近似值y/)满足：截去T：2!--Euler公式yn+ = y, +hf(x,,y,)-5-

-5- 二 方法1. Euler方法 ( ) Tarloy n 将初值问题的解函数 在 点 展开，有： y y  x x 2 ( ) ( ) ( ) ( )( ) ( ) 2 ! n n n n y y x y x y x x x x x         而 ， y f x y   ( , ) ( ) ( , ( )) n n n 所以 ，代入上式： y x f x y x   2 ( ) ( ) ( ) ( , ( ))( ) ( ) 2 ! n n n n n y y x y x f x y x x x x x        1 : n x x 令   2 1 1 1 ( ) ( ) ( ) ( , ( ))( ) ( ) 2 ! n n n n n n n n n y y x y x f x y x x x x x           2 1 ( ) ( ) ( , ( )) , 2 ! ( , ) n n n n n n n y y x f x y x h h x x        其中  2 1 ( ) , 2! n y T h   截去  1 1 ( ) n n y x y 得 的近似值 满足：   1 ( , ) n n n n y y hf x y    0 ( , ) ( ) y f x y a x b y a y          Euler公式 显式公式

2x0≤x≤1y例1用Euler公式求解初值问题(h = 0.1)y(0)=1解2由题意知：2xx。= a = 0,n= 10,b=1,yo=1f(x,y)= yy根据Euler公式：2x(n = 0,1,...,9)n+1 = yn + hf(x,,y,)= y, + h(y,yn代入数据：2x2×0= 1.1yi = yo + h(yo1+ 0.1(1yo2 ×0.1y2 = yi +h(y)- 1.19181.1 + 0.1(1.11.1y1依次类推注方程的精确解：y= /1+2x6

-6- 解 2 ( , ) , x f x y y y   0 0 x a n b y      0, 10, 1, 1 由题意知： 根据 公式： Euler 1 ( , ) n n n n y y hf x y    2 ( ) n n n n x y h y y    ( 0,1, ,9) n  0 1 0 0 0 2 ( ) x y y h y y    2 0 1 0.1(1 ) 1      1.1 1 2 1 1 1 2 ( ) x y y h y y    2 0.1 1.1 0.1(1.1 ) 1.1      1.1918 代入数据： 依次类推 . 注 方程的精确解：y x  1 2 例1 用 公式求解初值问题 Euler 2 0 1 x y y x y      y(0) 1     （ ） h  0.1

2x0≤x≤1y'=yy例1(续）(h= 0.1)(0)=11.8xn数值解Jn精确解y(xn）n1.700.01.00001.0000Euler公式0.111.10001.09541.620.21.19181.18321.530.31.27741.264940.41.35821.34161.450.51.43511.41421.360.61.50901.483270.71.58031.54921.2精确解：y=V1+2x80.81.64981.61251.190.91.71781.6733101.01.78481.73210.1020.30.40.50.60.70.80.90-7-

-7- 例1 续( ) 2 0 1 x y y x y      y(0) 1     （ ） h  0.1

Euler方法的几何意义y不y(xt)y(xn)Yn+1y(x2)yopyny(x)1y21y1-10xoX2xnXn+1Xix-8-

-8-  n 1 y  . Euler方法的几何意义 n y 2 y 1 y y x n1  y x n  y x 2  y x 1  . 0 y n 1 x n  x 2 x 1 x0 x

隐式公式2.后退的Euler方法将y=y(x)在x,点Tarloy展开：y"()y(x) = y(xn+1)+ y'(xu+)(x - xn+1)r2!y"(5)= y(xn+1)+ f(xn+1, y(xn+)(x - xn2!令x=x,:2y(x,) = y(xn+1)+ f(xn+1, y(xn+1)(x, - Xn+12!"(5)), 其中引 e(x,X+1)= y(xu+1)-hf(xu+1,y(xn+)+2!y"(n)即:y(xn+) = y(x,)+ hf(xu+1, y(x2!y"(5')n,截去 =得y(x,)的近似值y,满足：2!一后退Euler公式Yu+1 = yn + hf(xn+1,yn+1)-9-

-9- 2. Euler 后退的 方法 ( ) 1 Tarl yo n y y x x 将 在 点 展开：   2 1 1 1 1 ( ) ( ) ( ) ( ) ) 2 ( ) ( ! n n n n y y x y x x y x x x x             1 1 1 2 1 1 ( , ( ) ( ) ( ) ( ) 2 ( )) ! n n n n n y y x x f x y x x x x             : n 令x x  2 1 1 1 1 1 ( ) ( ) ( ) ( , ( ))( ) ( ) 2 ! n n n n n n n n n y y y x f x y x x x x x x              2 1 1 1 1 ( ) ( ) ( , ( )) , ( , ) 2 ! n n n n n n n y y x h h f x y x x x              其中 2 2 ( ) , 2! n y T h    截去   ( ) n n 得 的近似值 满足： y x y 1 1 1 ( , ) n n n n y y hf x y       后退 公式 Euler 即: 2 1 1 1 ( ) ( ) ( ) ( , ( )) 2 ! n n n n n h y y x y x hf x y x          隐式公式

3.梯形公式隐式公式y"(5)n: y(xut)=y(x,)+hf(xn,y(x,))+2!"(E)ny(xn+1)= y(x,)+ hf(xn+1,J(xn+1))2!y"(5.)y"(5')h和T, =注意到：T2的"符号"相反12!2!所以，两式相加并截去"T+T"得：一梯形公式+[f(x,yn)+ f(xn+1,+1)Yn+1 = ynX-4.改进的Euler公式梯形公式为隐式公式，求解时往往需要求解非线性方程，实际计算中通常由Euler公式对y进行"预测"，利用梯形公式进行校正”Jn+i = yn +hf(x,,yn)一改进的Euler公式hyn+1 = y, +=If(x,yn)+ f(xn+1,Jn+1))-10-

-10- 1 2 ( ) ( ) ( ) ( , ( )) 2 ! n n n n n y y x y x hf x y x h        1 2 2 2 ( ) 2 ) ! ( 2 ! n n y T y T h h       注意到：  和  的"符号"相反 1 1 2 1 ( ) ( ) ( , ( ( )) ) 2 ! n n n n n y x y x hf x y x y h          所以，两式相加并截去" + "得： T1 T2 1 1 1 ( , ) ( , ) 2 n n n n n n h y y f x y f x y           3.梯形公式 隐式公式 Euler n 1 y  梯形公式为 ，求解时往往需要求解非线性方程，实际计算 中通常由 公式对 进行"预测"，利用梯形公式进行 隐式公式 "校正" 1 ( , ) n n n n y y hf x y    1 1 1 [ ( , ) ( , )] 2 n n n n n n h y y f x y f x y            4. Euler 改进的 公式  梯形公式  改进的 公式 Euler