第8章 常微分方程.ppt

《第8章 常微分方程.ppt》由会员分享，可在线阅读，更多相关《第8章 常微分方程.ppt（41页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS第8章常微分方程 实际中，很多问题的数学模型都是微分方程。我们可以研究它们的一些性质。但是，只有极少数特殊的方程有解析解。对于绝大部分的微分方程是没有解析解的。 常微分方程作为微分方程的基本类型之一，在自然界与工程界有很广泛的应用。很多问题的数学表述都可以归结为常微分方程的定解问题。很多偏微分方程问题，也可以化为常微分方程问题来近似求解。本章讨论常微分方程的数值解法数 学 系University of Science and Technology

2、 of ChinaDEPARTMENT OF MATHEMATICS对于一个常微分方程：, , ),(baxyxfdxdyy通常会有无穷个解。如：Raaxyxdxdy,)sin( )cos(因此，我们要加入一个限定条件。通常会在端点出给出，如下面的初值问题：0)(, , ),(yaybaxyxfdxdy为了使解存在唯一，一般，要加限制条件在f上，要求f对y满足Lipschitz条件：2121),(),(yyLyxfyxf数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS 常微分方程的解是一个函数

3、，但是，计算机没有办法对函数进行运算。因此，常微分方程的数值解并不是求函数的近似，而是求解函数在某些节点的近似值。例：我们对区间做等距分割： , ()/ixai hhbam设解函数在节点的近似为iy由数值微分公式，我们有iixxxxyxfdxdy),(，则：),(1iiiiyxfhyy向前差商公式),( 1iiiiyxfhyy可以看到，给出初值，就可以用上式求出所有的iy数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS基本步骤如下： 解差分方程，求出格点函数 对区间作分割：bxxxamI10:

4、 求y(x)在xi上的近似值yi。iy称为分割I上的格点函数 由微分方程出发，建立求格点函数的差分方程。这个方程满足：A、解存在唯一；B、稳定，收敛；C、相容数值方法，主要研究步骤，即如何建立差分方程，并研究差分方程的性质。 这种方法 ，称为数值离散方法。求的是在一系列离散点列上，求未知函数y在这些点上的值的近似。我们的目的，就是求这个格点函数数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS 为了考察数值方法提供的数值解，是否有实用价值，需要知道如下几个结论： 收敛性问题 误差估计 稳定性问题

5、步长充分小时，所得到的数值解能否逼近问题的真解；舍入误差，在以后各步的计算中，是否会无限制扩大；数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS8.1 Euler公式公式做等距分割，利用数值微分代替导数项，建立差分方程。imabaxiI :1、向前差商公式)( 2)( )()(1iiiiyhxyhxyxy)( 2)(,()()(1iiiiiyhxyxfhxyxy)( 2)(,()()(21iiiiiyhxyxhfxyxy所以，可以构造差分方程),(1iiiiyxhfyy称为局部截断误差。显然，

6、这个误差在逐步计算过程中会传播，积累。因此还要估计这种积累数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS定义定义在假设在假设 yi = y(xi)，即第，即第 i 步计算是精确的前提下，考虑步计算是精确的前提下，考虑的截断误差的截断误差 Ri = y(xi+1) yi+1 称为称为局部截断误差局部截断误差 /* local truncation error */。记为2、收敛性)( 2)(,()()(21iiiiiyhxyxhfxyxy),(1iiiiyxhfyy1ihT考察局部误差的传播和

7、积累111)(iiiyxye1),()(,()(iiiiiiiThyxfxyxfhyxy数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS1)(iiiiThyxyhLejjiTThTehLmax , )1 (hThTehLhLi1)1 ()1 (hThLehLi1)1 ()1 (12hThLhTehLhLi1)1 ()1 ()1 (22hThLhLehLi1)1 ()1 ()1 (223hThLhLehLii1)1 ()1 ()1 (01hThLhLehLii)1 (1)1 (1)1 (101数

8、 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS称为整体截断误差hThLhLehLii101)1 ()1 (LTehLi01)1 (LTeLTeLabhLi)()1()( hOei)()1 (00hOTexenxn是1阶方法定义定义若某算法的局部截断误差为若某算法的局部截断误差为O(hp+1)，则称该算法有，则称该算法有p 阶精度。阶精度。数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS3、稳定

9、性误差在以后各步的计算中不会无限制扩大。我们考虑简单情况：仅初值有误差，而其他计算步骤无误差。设iz是初值有误差后的计算值，则),(),(11iiiiiiiizxhfzzyxhfyy所以，我们有：),(),(111iiiiiiiizxfyxfhezye)1 (hLezyhLeiiiihLiieehLe)1(010)1 ( 可以看出，向前差商公式关于初值是稳定的。当初始误差充分小，以后各步的误差也充分小数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS4、向后差商公式11()()()()2nnnn

10、y xy xhy xyh)( 2)(,()()(111nnnnnyhxyxfhxyxy)( 2)(,()()(2111nnnnnyhxyxhfxyxy),(111nnnnyxhfyy是隐格式，要迭代求解)0(1)(11)1(1),(nknnnknyyxhfyy可以由向前差商公式求出数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS5、中心差商公式)( 6)( 2)()(211nnnnyhxyhxyxy),(211nnnnyxhfyy是多步，2阶格式，该格式不稳定6、梯形法基于数值积分的公式对微分

11、方程, , ),(baxyxfdxdyy做积分，则：1),(nnxxyxfdxdy数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS)(21hOen类似，可以算出其误差估计式：2阶的方法所以，有格式为：),(),(2111nnnnnnyxfyxfhyy是个隐式的方法，要用迭代法求解11 ()()( , ( )nnxnnxy xy xf x y x dx)( 12)(,()(,(2)()(3111fhxyxfxyxfhxyxynnnnnn局部截断误差数 学 系University of Scien

12、ce and Technology of ChinaDEPARTMENT OF MATHEMATICS 基于数值积分的构造法基于数值积分的构造法将将 在在 上积分，得到上积分，得到),(yxfy ,1npnxx1)(,()()(1npnxxpnndxxyxfxyxy只要只要近似地算出右边的积分近似地算出右边的积分 ，则可通，则可通过过 近似近似y(xn+1) 。而。而选用不同近似式选用不同近似式 Ik，可得到不，可得到不同的计算公式同的计算公式。1)(,(npnxxkdxxyxfIkpnnIyy1数 学 系University of Science and Technology of Chin

13、aDEPARTMENT OF MATHEMATICS1)()!1()()2(npnxxqqdxxqy1)(,(npnxxdxxyxf若积分用节点qnnnxxx,1作为积分点，则1110)()()()(,(1nqnqnnxxhTxfaxfaxfahdxxyxfnpn)(,()(nnnxyxfxf101)(,()()(nqjjnjnjnpnnhTxyxfahxyxy积分系数1)(npnxxjjndxxlha这是显格式，q+1阶r+1步格式。r=maxp,q若以xn+1, xn+1, xn-q+1 为积分节点，可以构造r+1步q+1阶隐格式局部截断误差数 学 系University of Scien

14、ce and Technology of ChinaDEPARTMENT OF MATHEMATICS例：建立p=1,q=2的显格式p=1，q=2，显格式，11)(,(nnxxdxxyxf积分区间为积分节点为21,nnnxxxhdxxxxxxxxxhannxxnnnnnn37)()(1121210所以hdxxxxxxxxxhannxxnnnnnn32)()(1121121hdxxxxxxxxxhannxxnnnnnn31)()(1112212数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS),

15、(31),(32),(3722111nnnnnnnnyxfyxfyxfhyy)(,(31)(,(32)(,(37)()(22111nnnnnnnnxyxfxyxfxyxfhxyxy误差分析)( )(,(),(nnnnnxyxyxfyxf数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS例：建立p=2,q=2的隐格式p=2，q=2，隐格式，12)(,(nnxxdxxyxf积分区间为积分节点为11,nnnxxxhdxxxxxxxxxhannxxnnnnnn43)()(1211110所以0)()(12

16、11111nnxxnnnnnndxxxxxxxxxhahdxxxxxxxxxhannxxnnnnnn49)()(1211112数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS它的截断误差较 显格式 小，通常也具有更好的稳定性。 Adams公式公式 p=0 时候的多步法时候的多步法参见书数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS线性多步法线性多步法用用若干若干节点处的节点处的 y 及及 y

17、 值的值的线性组合线性组合来近似来近似y(xn+1)。).(.110111101knknnnknknnnffffhyyyy b bb bb bb ba aa aa a其通式可写为：其通式可写为：),(nnnyxff 当当 b b 1 0 时，为时，为隐式公隐式公式式; b b 1=0 则为则为显式公式显式公式。数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS由Taylor展开)()!1()(!)( )()()1(1)(1nkknkknnnykhxykhxhyxyxy记为1nhT),()( yx

18、fxy),(),()( yyxfyxfxyyx)( xy所以，可以构造格式),(),(),(! 2),(21nnnnynnxnnnnyxfyxfyxfhyxhfyy这种格式使用到了各阶偏导数，使用不便。如何起步？如何得到高精度、单步的格式？数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS从另一个角度看，11),(,()()(nnnnnhTfxyxhhFxyxy取(xn,y(xn)及其附近的点做线性组合，表示F，问题就好办了。当然，要求此时的展开精度相同。这种方法称为RungeKutta法8.2

19、 RungeKutta法法数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS在(xn,y(xn)处展开，)()(,()(,()(,()(,()(,(),(,(22122221hOxyxfxyxhfbcxyxhfacxyxfcxyxfcfxyxhFnnynnnnxnnnnnn比较11),(),(),(! 2),()()(nnnnnynnxnnnnhTyxfyxfyxfhyxfhxyxy以2阶为例，设),(,(),(),(21221yxhfbyhaxfcyxfcfyxhF数 学 系Universit

20、y of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS有：2/12/112122221bcaccc112/121221bacc1、改进的Euler公式),(),( 2/121211hKyhxfKyxfKKKhyynnnnnn3/23/24/3, 4/121221bacc2、Heun公式数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS一般的RungeKutta法构造),(),(),(),(111212212211miimimmm

21、mKbhyhaxfKhKbyhaxfKyxfKKcKcKcfyxhF常见的为3阶，4阶公式数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS),(),(),(),()22(34222312221432161hKyhxfKKyxfKKyxfKyxfKKKKKyyiihihihihiiihii 最常用为四级最常用为四级4阶阶经典龙格经典龙格-库塔法库塔法 ：数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMAT

22、ICSLab07 常微分方程3.用如上程序求常微分方程22 (01.5)(0)3yx yxy 3( )3/(1)y xx分别取步长h=0.1,0.1/2,0.1/4,0.1/8计算y(1.5)，并与精确解比较1.经典4阶Runge-Kutta方法解常微分方程的通用程序2.Adams隐式3阶方法解常微分方程的通用程序(由1提供初值)4.给出误差和误差阶。简单分析数据数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICSSample Output ( represents a space)Runge-K

23、utta法，误差和误差阶为k=0,0.244934066848e00k=1,0.534607244904e-01, 3.90.Adams隐格式，误差和误差阶为k=0,0.244934066848e00k=1,0.534607244904e-01 , 3.01.数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS8.4 方程组和高阶方程的数值解法方程组和高阶方程的数值解法bxaayayyyxfdxdyyyxfdxdymmmmmm , )()(),(),(111111写成向量的形式：)(),(aYYx

24、FdxdY数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS各种方法都可以直接运用过来。bxazazyayzyxgdxdzzyxfdxdy , )()(),(),(00Euler公式),(),(11nnnnnnnnnnzyxhgzzzyxhfyy以两个方程的方程组为例数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICSRunge-Kutta公式112341(22)6nnnnyyhKKKKzz)2,2

25、(12KhYhxFK1(1)(2)112(1)(2)11(,)(,)(,)222(,)222nnnnnnnnnnnnf xyzKg xyzhhhf xyKzKKhhhg xyKzK 数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS(1)(1)223(1)(1)22(1)(1)334(1)(1)33(,)222(,)222(,)(,)nnnnnnnnnnnnhhhf xyKzKKhhhg xyKzKf xh yhKzhKKg xh yhKzhK数 学 系University of Science

26、 and Technology of ChinaDEPARTMENT OF MATHEMATICS0.05 (1)0.002200.09 (1)0.1515(0)0.193(0)0.083duuuuvdtdvvvuvdtuv1、( , , )0.05 (1/20) 0.002( , , )0.09 (1/15) 0.15f u v tuuuvFg u v tvvuv2、确定方法，然后求解（0.20276 0.0881157）（0.213007 0.0934037）（0.223763 0.0988499）（0.235052 0.104437）（0.246902 0.110146）4阶Runge-

27、Kutta法，h=1数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS高阶方程bxaayayayyyyxfdxydmmmmm , )()( )(), ,()1(21)1(1数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS则有：bxaayayyyxfdxdyydxdymmmm , )()(),(11121令mmydxdyydxdyyy1211数 学 系University of Science a

28、nd Technology of ChinaDEPARTMENT OF MATHEMATICS例：例：考察初值问题考察初值问题 在区间在区间0, 0.5上的解。上的解。分别用欧拉显、隐式格式和改进的欧拉格式计算数值解。分别用欧拉显、隐式格式和改进的欧拉格式计算数值解。 1)0()(30)(yxyxy0.00.10.20.30.40.5精确解精确解改进欧拉法改进欧拉法 欧拉隐式欧拉隐式欧拉显式欧拉显式 节点节点 xixey30 1.0000 2.0000 4.0000 8.0000 1.6000 101 3.2000 101 1.00002.5000 10 1 6.2500 10 21.5625

29、 10 23.9063 10 39.7656 10 41.00002.50006.25001.5626 1013.9063 1019.7656 1011.00004.9787 10 22.4788 10 31.2341 10 46.1442 10 63.0590 10 7出了什么问题出了什么问题 ?!数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS8.5 差分方程的绝对稳定性差分方程的绝对稳定性 显然，这个误差不仅于差分格式有关，而且与微分方程本身有关。如果微分方程本身是不稳定，那就没理由要求

30、这2组解充分接近。因此，差分方程的稳定性概念是建立在微分方程稳定的基础上的。考虑最简单的模型：只有初值产生误差，看看这个误差的传播。0)()0(Re yayydxdy把这个典型微分方程规定为：数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS差分方程运用到如上的微分方程后，可以得到0)Re( , 00kjjnjkjjnjybhya对于给定的初始误差110,keee，误差方程具有一样的形式0)Re( , 00kjjnjkjjnjebhea对于一般的差分方程kjjnjnjkjjnjyxfbhya00

31、),(数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS定义：差分方程称为绝对稳定的，若差分方程作用到微分方程)0(Re ydxdy时，对任意的初值，总存在左半复平面上的一个区域，当 在这个区域时，差分方程的解趋于0。这个区域称为稳定区域h例：Euler公式的稳定性nnnhyyy1误差方程：nnnheee1heenn1110-1-2ReImg数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS考察隐

32、式欧拉法考察隐式欧拉法11 iiiyhyy 111iiyyh11011iih可见绝对稳定区域为：可见绝对稳定区域为：|1|1h210ReImg注：注：一般来说，隐式欧拉法的绝对稳定性比同阶的显式一般来说，隐式欧拉法的绝对稳定性比同阶的显式法的好。法的好。数 学 系University of Science and Technology of ChinaDEPARTMENT OF MATHEMATICS3阶Runge-Kutta232132111)211 (46hhyKhyKyKKKKhyynnnnn32161211hhhyynn显式显式 1 4 阶方法的绝对稳定区域为阶方法的绝对稳定区域为k=1k=2k=3k=4-1-2-3-123ReImg