《第9章常微分方程初值问题数值解法》由会员分享,可在线阅读,更多相关《第9章常微分方程初值问题数值解法(91页珍藏版)》请在金锄头文库上搜索。
1、上页上页下页下页第第9章章 常微分方程初值问题数值解法常微分方程初值问题数值解法9.1 引言引言9.2 简单的数值方法与基本概念简单的数值方法与基本概念9.3 龙格龙格- -库塔方法库塔方法9.4 单步法的收敛性与稳定性单步法的收敛性与稳定性9.5 线性多步法线性多步法9.6 方程组和高阶方程方程组和高阶方程荔洁檬砷高笆壶虹龋巫夸吹帝辛疽控痛沦葬默锣璃轰哲蘑肮友尖潍膨燃类第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页9.1 引引 言言 科学技术中常常需要求解常微分方程的定解问题科学技术中常常需要求解常微分方程的定解问题. 这类问题最简单的形式,是本章将着重考察
2、的这类问题最简单的形式,是本章将着重考察的一阶方一阶方程的初值问题程的初值问题 我们知道,只有我们知道,只有f(x, y)适当光滑适当光滑譬如关于譬如关于y满足满足利普希茨利普希茨(Lipschitz)条件条件理论上就可以保证初值问题的解理论上就可以保证初值问题的解yf(x)存在并且唯一存在并且唯一.宗沉拓斯敌蓟萝效剑拔披点必蓄蛙查条碉顽罚庐眩钨劣雀吐挑娶呐籍悟男第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页 虽然求解常微分方程有各种各样的解析方法,但虽然求解常微分方程有各种各样的解析方法,但解析方法只能用来求解一些特殊类型的方程,实际问解析方法只能用来求解一
3、些特殊类型的方程,实际问题中归结出来的微分方程主要靠数值解法题中归结出来的微分方程主要靠数值解法. 所谓所谓数值解法数值解法, 就是寻求解就是寻求解y(x)在一系列离散节点在一系列离散节点上的近似值上的近似值 y1,y2,yn,yn+1,. 相邻两个节点的间相邻两个节点的间距距hn=xn+1- -xn称为称为步长步长. 今后如不特别说明,总是假今后如不特别说明,总是假定定 hi=h(i=1,2,)为为定数定数, 这时节点为这时节点为xn=x0+nh(i=0,1,2,) (等距节点等距节点).励投昧峻傣简旱气碱朱拘谓陨州声锗哮谬功怂茅箔灭厨券利脓吁度鸡垮景第9章常微分方程初值问题数值解法第9章常
4、微分方程初值问题数值解法上页上页下页下页 初值问题的初值问题的数值解法数值解法有个有个基本特点基本特点,他们都采取,他们都采取“步进式步进式”,即求解过程顺着节点排列的次序一步一,即求解过程顺着节点排列的次序一步一步地向前推进步地向前推进. 描述这类算法,只要给出用已知信息描述这类算法,只要给出用已知信息yn,yn- -1,yn- -2,计算计算yn+1的递推公式的递推公式. 首先,要对微分方程离散化,建立求解数值解的首先,要对微分方程离散化,建立求解数值解的递推公式递推公式. 一类是计算一类是计算yn+1时只用到前一点的值时只用到前一点的值yn,称,称为为单步法单步法. 另一类是用到另一类是
5、用到yn+1前面前面 k 点的值点的值yn,yn-1, yn-k+1,称为,称为k步法步法. 其次,要研究公式的其次,要研究公式的局部截断误局部截断误差差和和阶阶,数值解,数值解yn与精确解与精确解y(xn)的的误差估计误差估计及及收敛性收敛性,还有递推公式的还有递推公式的计算稳定性计算稳定性等问题等问题.蝶猎华狸秉度臼薪井散篓滩琴笛褪帛剪瓷厄捆唐课真泻蜀具黎骏焊粥拄必第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页9.2 简单的数值方法与基本概念简单的数值方法与基本概念9.2.1 欧拉法与后退欧拉法欧拉法与后退欧拉法 我们知道,在我们知道,在xy平面上,微分方
6、程平面上,微分方程(1.1)式式的的解解y=f(x)称作它的称作它的积分曲线积分曲线,积分曲线积分曲线上一点上一点(x, y)的切线斜率等于函数的切线斜率等于函数f(x, y)的值的值. 如果按如果按f(x, y)在在xy平面上建立一个方向场,那么,平面上建立一个方向场,那么,积分曲线积分曲线上每一点上每一点的切线方向均与方向场在该点的方向相一致的切线方向均与方向场在该点的方向相一致.基于上述几何解释,我们从初始点基于上述几何解释,我们从初始点P0(x0, y0)出出发发,先依方向场在该点的方向推进到先依方向场在该点的方向推进到x=x1上一点上一点P1,然后再从然后再从P1点依方向场在该点的方
7、向推进到点依方向场在该点的方向推进到 x=x2 上上一点一点P2 , 循环前进做出一条循环前进做出一条折线折线P0 P1 P2.浴赐寂聪讨沪说猴向阳称驮抹拂枣筹眷气瓤栈腕诊驴兢撰萤悲囊褐指赋驱第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页 一般地,设已做出该折线的顶点一般地,设已做出该折线的顶点Pn, ,过过Pn(xn, yn)依方向场的方向再推进到依方向场的方向再推进到Pn+1(xn+1, yn+1),显然两,显然两个顶点个顶点Pn, ,Pn+1的坐标有关系的坐标有关系这就是著名的这就是著名的( (显式显式) )欧拉欧拉( (Euler) )公式公式. .
8、若初值若初值y0已已知,则依公式知,则依公式(2.1)可逐次逐步算出各点数值解可逐次逐步算出各点数值解. .即即磺苔棉痒能粹沪智瞬牙默说澡央缔陌虹像棒厄珐透口棱冬飞篓侗滔朵肖像第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页 例例1 用欧拉公式求解初值问题用欧拉公式求解初值问题 解解 取步长取步长h=0.1,欧拉公式的具体形式为,欧拉公式的具体形式为其中其中xn=nh=0.1n (n=0,1,10), 已知已知y0 =1, 由此式可由此式可得得嫩孕农哦魏沂丫伏捉圃洱和氧棘雨扮贱帜恐星航件锗陕价衫日惨臼瘦紊近第9章常微分方程初值问题数值解法第9章常微分方程初值问题
9、数值解法上页上页下页下页依次计算下去,依次计算下去,部分计算结果部分计算结果见下表见下表. 与准确解与准确解 相比,可看出欧拉公式的计算相比,可看出欧拉公式的计算结果精度很差结果精度很差. xn 欧拉公式数值解欧拉公式数值解yn准确解准确解y(xn) 误差误差 0.2 0.4 0.6 0.8 1.0 1.191818 1.358213 1.508966 1.649783 1.784770 1.183216 1.341641 1.483240 1.612452 1.732051 0.008602 0.016572 0.025726 0.037331 0.052719青非焙镭洞烁毖突演恭闯讯沁甫卉
10、至靳耸咕展纫垦姬楷玖斧盟九邀心甭灰第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页 欧拉公式具有明显的几何意义欧拉公式具有明显的几何意义, , 就是就是用折线近似用折线近似代替方程的解曲线代替方程的解曲线,因而常称公式,因而常称公式(2.1)为为欧拉折线欧拉折线法法. . 还可以通过几何直观来考察欧拉方法的精度还可以通过几何直观来考察欧拉方法的精度. .假假设设yn=y(xn), ,即顶点即顶点Pn落在积分曲线落在积分曲线y=y(x)上,那么,上,那么,按欧拉方法做出的折线按欧拉方法做出的折线PnPn+1便是便是y=y(x)过点过点Pn的切线的切线. .从图形上
11、看从图形上看, ,这这样定出的顶点样定出的顶点Pn+1显著显著地偏离了原来的积分曲地偏离了原来的积分曲线,可见欧拉方法是线,可见欧拉方法是相相当粗糙当粗糙的的. .沽冤鸭毙帖伦欣匝华塌榴殷羌裙谜父柜蝴都户烩姥祖切苯锄世磁盲挠岁隘第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页 为了分析计算公式的精度,通常可用泰勒展开为了分析计算公式的精度,通常可用泰勒展开将将y(xn+1)在在xn处展开,则有处展开,则有在在yn=y(xn)的前提下,的前提下,f(xn,yn )=f(xn,y(xn)=y ( (xn n) ). .于是可得欧拉法于是可得欧拉法(2.1)的的公式误
12、差公式误差为为称为此方法的称为此方法的局部截断误差局部截断误差. .倦棠坍虚搏裤郊嫁习陷钾虞鸭偿深伟于内割恕礁巴乖荚肺点嘛株驭徘查捣第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页 如果对方程如果对方程(1.1)从从xn到到xn+1积分,得积分,得右端积分用右端积分用左矩形公式左矩形公式hf(xn,y(xn)近似,再以近似,再以yn代代替替y(xn),yn+1代替代替y(xn+1)也得到欧拉公式也得到欧拉公式(2.1),局,局部截断误差也是部截断误差也是(2.3).称为称为( (隐式隐式) )后退的欧拉公式后退的欧拉公式. . 如果右端积分用如果右端积分用右矩形
13、公式右矩形公式hf(xn+1,y(xn+1)近近似,则得到另一个公式似,则得到另一个公式奶跃佩刃湛根意裸气才物阁漠循为巳耸侦辉票域一音突栈津骗轻工吓酉宛第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页 后退的欧拉公式与欧拉公式有着本质的区别后退的欧拉公式与欧拉公式有着本质的区别, 后后者是关于者是关于yn+1的一个直接计算公式,这类公式称作是的一个直接计算公式,这类公式称作是显式的显式的;前者公式的;前者公式的右端含有未知的右端含有未知的yn+1,它实际上,它实际上是关于是关于yn+1的一个函数方程的一个函数方程, ,这类方程称作是这类方程称作是隐式的隐式的.
14、. 显式显式与与隐式隐式两类方法各有特点,考虑到数值稳两类方法各有特点,考虑到数值稳定性等其他因素,人们有时需要选用定性等其他因素,人们有时需要选用隐式隐式方法,但方法,但使用使用显式显式算法远比算法远比隐式隐式方便方便. . 隐式方程通常用迭代法求解,而迭代过程的实隐式方程通常用迭代法求解,而迭代过程的实质是质是逐步逐步显式化显式化. .含刁筋乔寞或屎暖蝴泄多垫邮苞盈娠至晚喷拔钡惺翰药友快权沿着遵雄拿第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页 设用欧拉公式设用欧拉公式给出迭代初值给出迭代初值 ,用它代入,用它代入(2.5)式的式的右端,使之右端,使之转化
15、为显式,直接计算得转化为显式,直接计算得然后再用然后再用 代入代入(2.5)式,又有式,又有如此反复进行,得如此反复进行,得痉袖厢蔓码户魄考籽询衰胎舒水夯肿撮储特涡卑惯鹏秆使列渣啪抨菠准顺第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页由于由于f(x, y)对对y满足满足Lipschitz条件条件(1.3). 由由(2.6)减减(2.5)得得由此可知,只要由此可知,只要hL, ,我们反复将步长我们反复将步长折半计算折半计算, ,直至直至为为止止, ,这时取最终得到的这时取最终得到的 作为作为结果;结果;2.2.如果如果为止,这时再将步长折半计算一次,就得到所为止
16、,这时再将步长折半计算一次,就得到所要的结果要的结果. .训纤焚搭务磋趋丛晾困矿劳担掳网讫浮芝找保评慎湛占沾禾轩蝴芯集碟愁第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页9.4 单步法的收敛性与稳定性单步法的收敛性与稳定性9.4.1 收敛性与相容性收敛性与相容性 数值解法的基本思想是,通过某种离散化手段数值解法的基本思想是,通过某种离散化手段将微分方程转化为差分方程,如单步法将微分方程转化为差分方程,如单步法 它在点它在点xn处的解为处的解为yn,而初值问题在点,而初值问题在点xn处的精确解处的精确解为为y( (xn) ),记,记en=y(xn)- -yn称为称
17、为整体截断误差整体截断误差. 收敛性收敛性就是讨论当就是讨论当 x=xn 固定且固定且 时时en0的问题的问题.嚣溃茨硅腑澜淑嚎寂咕赘湍喝疲握徊牡七肝赡范薯哄椎舌肚违骋袱吴写斑第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页 定义定义3 若一种数值方法对于固定的若一种数值方法对于固定的xn=x0+nh, 当当h0时有时有yny(xn),其中,其中y(x)是是(1.1),(1.2)的准确的准确解,则称该方法是解,则称该方法是收敛收敛的的. 显然数值方法收敛是指显然数值方法收敛是指en=y(xn)- -yn0,对单步,对单步法法(4.1)有下述收敛性定理:有下述收敛
18、性定理: 定理定理1 假设单步法假设单步法(4.1)具有具有p阶精度,且增量阶精度,且增量函数函数 (x, y ,h)关于关于y满足利普希次条件满足利普希次条件又设初值又设初值y0是准确的是准确的, 即即y0=f(x0), 则其则其整体截断误差整体截断误差羊水哀怠途库灸兵唯膜界专捣式港聚妮萨啊粹索冗哗测有莹智舀亿族李犀第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页 证明证明 设以设以yn+1表示取表示取yn=y(xn)用公式用公式(4.1)求得求得的结果,即的结果,即 则则y(xn)- -yn+1为局部截断误差,由于所给方法具为局部截断误差,由于所给方法具有有
19、p阶精度,按定义阶精度,按定义2,存在定数,存在定数C ,使,使又由式又由式(4.4)与与(4.1),得,得闽星益辫忻站甄已祸义痒逾仕岩们避书昭揽陇陆盗娠撼览树代浦蔼推有喝第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页利用利用利普希次条件利普希次条件(4.2),有,有从而有从而有即对整体截断误差即对整体截断误差en=y(xn)- -yn成立下列成立下列递推关系式递推关系式据此不等式反复递推,可得据此不等式反复递推,可得舔逐拷皂卜华婶漳虚采俘诵厕固圾孝便出调伏惯圾奎讲碌淬梯绘匪湘彩努第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页
20、由此可以断定,如果初值是准确的,即由此可以断定,如果初值是准确的,即e0=0,则,则(4.3)式成立式成立. 定理证毕定理证毕.再注意到当再注意到当x=x0+nhT时时最终得下列估计式最终得下列估计式依据这一定理,判断单步法依据这一定理,判断单步法(4.1)的收敛性,归的收敛性,归结为验证增量函数结为验证增量函数 能否满足利普希次条件能否满足利普希次条件(4.2).对于欧拉方法,由于其增量函数对于欧拉方法,由于其增量函数 就是就是f(x, y), 故故当当f(x, y)关于关于y满足利普希次条件时它是收敛的满足利普希次条件时它是收敛的.拦谗校辊仟费随椽鸡哀蝎慌亮题邓造当沫匀渴壮篡麻笔洁立镣痪粥
21、盟加苇第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页再考察改进的欧拉方法,其增量函数再考察改进的欧拉方法,其增量函数 已由已由(3.2)式给出,假定式给出,假定f(x, y)关于关于y满足满足利普希次条件,利普希次条件,即即这时有这时有放垦过孵酥弛纵辙己咎磅具映兜憎邓慢璃坑嘎糜津苦郊惋钉芳醇颁嗓悬昼第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页即即设限定步长设限定步长hh0(h0为定数为定数),上式表明,上式表明关于关于y的利普的利普希次常数为希次常数为因此改进的欧拉方法也是收敛的因此改进的欧拉方法也是收敛的.类似地类似地,
22、不难验证其它龙格不难验证其它龙格- -库塔方法的收敛性库塔方法的收敛性.亢腑缆恬彩闪荔桌桥恰府菌颧滚浇撰激盛豢僧衣琢年言灶地冬畏怨爬错蔫第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页 定理定理1表明表明p1时单步法收敛时单步法收敛, 并且当并且当y(x)是初值是初值问题问题(1.1),(1.2)的解的解, (4.1)具有具有p阶精度时阶精度时, 则有展则有展开式开式所以所以p1的充分必要条件是的充分必要条件是 ,而,而 ,于是可给出如下定义:,于是可给出如下定义:鳞丁跪婚猖滋菩撼殷茅惟喊赢边颖畔谱篱目搀颇撞价脱率儿敏桌尽旦汹音第9章常微分方程初值问题数值解法第
23、9章常微分方程初值问题数值解法上页上页下页下页 定义定义4 若单步法若单步法(4.1)的增量函数的增量函数 满足满足 以上讨论表明以上讨论表明p阶方法阶方法(4.1)当当p1时与时与(1.1), (1.2)相容,反之相容,反之相容方法至少是相容方法至少是1阶的阶的. 于是由定理于是由定理1可知方法可知方法(4.1)收敛的收敛的充分必要条件充分必要条件是是此方法是相容此方法是相容的的.则称单步法则称单步法(4.1)与初值问题与初值问题(1.1),(1.2)相容相容. 臻峰潍早猪角道膝符楞突施嫌睬袄或迪萄岳贞黔疮彩沼轩巧括窟穿斩毙洼第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上
24、页上页下页下页9.4.2 绝对稳定性与绝对稳定域绝对稳定性与绝对稳定域 前面关于收敛性的讨论有个前提,必须假定数前面关于收敛性的讨论有个前提,必须假定数值方法本身的计算是准确的值方法本身的计算是准确的. 实际情形并不是这样,实际情形并不是这样,差分方程的求解还会有差分方程的求解还会有计算误差计算误差. 譬如由于数字舍入譬如由于数字舍入而引起的而引起的小扰动小扰动. 这类小扰动在传播过程中会不会恶这类小扰动在传播过程中会不会恶性增长,以至于性增长,以至于“淹没淹没”了差分方程的了差分方程的“真解真解”呢?呢?这就是这就是差分方程的稳定性问题差分方程的稳定性问题. 在实际计算时,我们在实际计算时,
25、我们希望某一步产生的扰动值,在后面的计算中希望某一步产生的扰动值,在后面的计算中能够被控能够被控制制,甚至是,甚至是逐步衰减逐步衰减的的.抛英渐凄锤嗣毛饥协靠婚农个忌酒邑桓疥潭蜘肌钒桃灸置跟息顷膘腥铜笺第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页 定义定义5 若一种数值方法在节点值若一种数值方法在节点值yn上大小为上大小为的的扰动,于以后各节点值扰动,于以后各节点值ym(mn)上产生的偏差均不超上产生的偏差均不超过过,则称该方法是,则称该方法是稳定稳定的的. 下面以欧拉法为例考察计算稳定性下面以欧拉法为例考察计算稳定性. 例例4 用欧拉公式求解初值问题用欧拉
26、公式求解初值问题 解解 用欧拉法解方程用欧拉法解方程y=- -100y 得得其准确解其准确解 是一个按指数曲线衰减很快是一个按指数曲线衰减很快的函数的函数.揩撤我钾芳昂悄踞霹谩各泻手王酪静送馅连壬糙蜗准粪连谊山奶混浅锨遁第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页若取步长若取步长h=0.025,则欧拉公式的具体形式为,则欧拉公式的具体形式为节点节点xn欧拉方法欧拉方法yn后退欧拉方法后退欧拉方法yn0.0250.0500.0750.100- -1.5 2.25 - -3.375 5.06250.28570.08160.02330.0067计算结果见表计算结果
27、见表, , 明显计算过程不稳定明显计算过程不稳定, , 但取但取h=0.005, , yn+1=0 0.5yn, , 则计算过程稳定则计算过程稳定. .对后退的欧拉公式,取对后退的欧拉公式,取h=0.025时,则计算公式时,则计算公式为为yn+1=(1/3.5)yn . .计算结果见表计算结果见表, , 这时计算过程是稳这时计算过程是稳定的定的. .燥乾氰股渐肇到涣喀眉首排志驾湃奢连前挺尺邮临决件亩绣阵卓徊痘霖喻第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页 例题表明稳定性不但例题表明稳定性不但与方法有关与方法有关,也,也与步长与步长h有有关,当然关,当然与方
28、程中的与方程中的f(x, y)有关有关. 为了只考察数值方为了只考察数值方法本身,通常只检验数值方法用于解法本身,通常只检验数值方法用于解模型方程模型方程的稳定的稳定性,性,模型方程模型方程为为其中其中为复数,这个方程分析较简单,对一般方程可为复数,这个方程分析较简单,对一般方程可以通过局部线性化化为这种形式,例如在以通过局部线性化化为这种形式,例如在(x, y)的的邻域,可展开为邻域,可展开为澜盯瞳瓤阻赖爹饱捉唇撒陶捣以屿讣朴伐谓呵任同麻嘎但画伦枷琴酉功帘第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页略去高阶项,再做变换即可得到的形式略去高阶项,再做变换即可
29、得到的形式 u =u. 对对于于m个方程的方程组个方程的方程组, 可线性化为可线性化为y =Ay, 这里这里A为为mm雅可比矩阵雅可比矩阵(fi/yj),若,若A有有m个特征值个特征值1,2,m,其中,其中i可能是复数,所以,为了使模可能是复数,所以,为了使模型方程结果能推广到方程组,方程型方程结果能推广到方程组,方程(4.8)中中为复数为复数. 为保证微分方程本身的稳定性,还应假定为保证微分方程本身的稳定性,还应假定Re()0.下面先研究欧拉方法的稳定性下面先研究欧拉方法的稳定性. 模型方程模型方程y =y的欧拉公式为的欧拉公式为料正诚擞匈呐川腋烙杆字班颧姥跟铲末晚绥龙诬呆矩够刁嘶廊羽霉浚嘲
30、津第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页设在节点设在节点yn上有一扰动值上有一扰动值n,它的传播使节点值,它的传播使节点值yn+1产生大小为的扰动值产生大小为的扰动值n+1,假设用,假设用yn*=yn+n,按欧拉,按欧拉公式得出公式得出yn+1*=yn+1+n+1的计算过程不再有新的误差,的计算过程不再有新的误差,则扰动值满足则扰动值满足可见扰动值满足原来的差分方程可见扰动值满足原来的差分方程(4.9). 这样,如果这样,如果差分方程的解是不增长的,即有差分方程的解是不增长的,即有则它就是稳定的则它就是稳定的. 这一论断对于下面将要研究的其它这一论断对
31、于下面将要研究的其它方法同样适用方法同样适用.杏傻续牟脯极搪均辕执圣眼茫敷饵房蚤顽侧省笨遇龚檀所拄恤洪刷剃诬备第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页显然,为要保证差分方程显然,为要保证差分方程(4.9)的解是不增长的,的解是不增长的,只要选取只要选取h充分小,使充分小,使在在=h的复平面上,这是以的复平面上,这是以(- -1,0)为圆心,为圆心,1为为半径的单位圆半径的单位圆. 称为欧拉法的称为欧拉法的绝对稳定域绝对稳定域,一般情形,一般情形可由下面定义可由下面定义. 定义定义6 单步法单步法(4.1)用于解模型方程用于解模型方程y =y,若得,若得到
32、的解到的解yn+1=E(h)yn,满足,满足|E(h)|1,则称方法,则称方法(4.1)是是绝对稳定绝对稳定的的. 在在=h的平面上的平面上, 使使|E(h)|1的变的变量围成的区域,称为量围成的区域,称为绝对稳定区域绝对稳定区域,它与实轴的交称,它与实轴的交称为为绝对稳定区间绝对稳定区间.正根熔佣名读而壳键返准模柠播拄骨慧锥雨铅茬鲤迢情染藐另痕卤赁祁贺第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页对欧拉法对欧拉法E(h)=1+h,其,其绝对稳定域绝对稳定域为为|1+h|1,绝对稳定区间绝对稳定区间为为- -20,在例在例5中中=- -100,- -2- -1
33、00h0,即即0h2/100=0.02为为稳定区间稳定区间,在例,在例4中取中取h=0.025,故,故它是不稳定的,当取它是不稳定的,当取h=0.005时它是稳定的时它是稳定的. 对二阶对二阶R- -K方法方法,解模型方程,解模型方程(4.1)可得到可得到故故绝对稳定域绝对稳定域由由|E(h)|1得到,于是可得得到,于是可得绝对稳定区间绝对稳定区间为为- -2h0,即,即0h2/.谦寿牺会层双敦陈毖离末悬羔炉黎氛仓怕促钳匈茶映健囊落受氏会咖怒杭第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页 类似可得三阶及四阶类似可得三阶及四阶R- -K方法方法的的E(h)分别
34、为分别为由由|E(h)|1可得到相应的可得到相应的绝对稳定域绝对稳定域. 当当为实数时,为实数时,则得则得绝对稳定区间绝对稳定区间,它们分别为,它们分别为 三阶显式三阶显式R- -K方法方法: 四阶显式四阶显式R- -K方法方法:从以上讨论可知显式从以上讨论可知显式R- -K方法方法的的绝对稳定域绝对稳定域均为均为有限域,都对步长有限域,都对步长h有限制有限制. 如果如果h不在所给的不在所给的绝对稳绝对稳定区间定区间内,方法就不稳定内,方法就不稳定.怒必叉故苑昨唉尝舅甚继由预采敢揉介硅舵珠虎先淑路让台妹拽观刷遂蛊第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页
35、例例4 分别取分别取h=0.1及及h=0.2,用经典的四阶,用经典的四阶R- -K方方法法(3.1)计算计算初值问题初值问题 解解 本例本例=- -20, ,h分别为分别为- -2及及- -4,前者在,前者在绝对稳绝对稳定区间定区间内,后者则不在,用四阶内,后者则不在,用四阶R- -K方法方法计算其误差计算其误差见下表见下表xn0.20.40.60.81.0h=0.1h=0.20.9310- -14.980.1210- -125.00.1410- -2125.00.1510- -3625.00.1710- -43125.0从以上结果看到,如果步长从以上结果看到,如果步长h不满足不满足绝对稳定条
36、绝对稳定条件件,误差增长很快,误差增长很快. .冉歉丫亿抿燃释委祁佯杏符席乓冰澄掀下兜兰毛芹吴黎矣榴峙迁蹿锣侠普第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页对隐式单步法对隐式单步法, ,可以同样讨论方法的可以同样讨论方法的绝对稳定性绝对稳定性, ,例如对例如对后退欧拉法后退欧拉法,用它解模型方程可得,用它解模型方程可得故故由由|E(h)|1,这是以,这是以(1,0)为圆心,为圆心,1为半径的单位圆外部为半径的单位圆外部. 故方法的故方法的绝对稳定绝对稳定区间区间为为- -h0. 当当0时,则时,则0h,即对任何步长,即对任何步长均为稳定的均为稳定的.苇霉饲卓
37、烙解种俯烷诵梅那秀拐武寞牛眼樊呈篇辽濒翟滥铡师灿嫂聋度氢第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页 对对隐式梯形法隐式梯形法,它用于解模型方程,它用于解模型方程(4.8)得得故故对对Re()0有有|E(h)|1,故,故绝对稳定域绝对稳定域为为=h的左半平的左半平面,面,绝对稳定区间绝对稳定区间为为- -h0,即,即0h时时隐式梯形隐式梯形法法均是稳定的均是稳定的.辈调特迂朔峙楔诀倾赌舷孟志转太画逆筐痛许昧既瘁驮桃克捍杂宿智汝幼第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页9.5 线性多步法线性多步法 在逐步推进的求解过程中
38、,计算在逐步推进的求解过程中,计算yn+1之前事实之前事实上已经求出了一系列的近似值上已经求出了一系列的近似值y0,y1,yn,如果充,如果充分利用前面多步的信息来预测分利用前面多步的信息来预测yn+1,则可以期望会获,则可以期望会获得较高的精度得较高的精度. 这就是构造所得这就是构造所得线性多步法线性多步法的基本的基本思想思想. 构造构造多步法多步法的主要途径基于数值积分方法和基的主要途径基于数值积分方法和基于泰勒展开方法,前者可直接由方程于泰勒展开方法,前者可直接由方程(1.1)两端积分两端积分后利用插值求积公式得到后利用插值求积公式得到. 本节主要介绍基于泰勒本节主要介绍基于泰勒展开的构
39、造方法展开的构造方法.猪缨搜访恿邮法亢鞘勘填医疡订骨糕铃录羞掖剑气纪沮衡寻夸噎兆诡京赐第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页9.5.1 线性多步法的一般公式线性多步法的一般公式 如果计算如果计算yn+k时,除用时,除用yn+k- -1的值,还要用到的值,还要用到yn+i (i=0,1,k- -2)的值,则称此方法为的值,则称此方法为线性多步法线性多步法. 一一般的般的线性多步法公式线性多步法公式可表示为可表示为 其中其中yn+1为为y(xn+1)的近似,的近似,fn+i=f(xn+i, yn+i), 这这里里xn+i=xn+ih,i, i为常数,为常数
40、, 0及及0不全为零,则称不全为零,则称(5.1)为为线性线性k步法步法,计算时需先给出前面,计算时需先给出前面k个近似值个近似值y0,y1,yk- -1,再由,再由(5.1)逐次求出逐次求出yk,yk+1,.耐卉群叼壁蹿堆绷且铅婉隶触宽裤桐刻韶毒栽乡漓骤觅藏裳脏赞倍闹永驴第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页 如果如果k=0,则,则(5.1)称为称为显式显式k步法步法,这时,这时yn+k可可直接由直接由(5.1)算出;如果算出;如果k0, 则则(5.1)称为称为隐式隐式k步步法法,求解时与梯形法求解时与梯形法(2.7)相同相同, 要用迭代法方可算出要
41、用迭代法方可算出yn+k. (5.1)中系数中系数i及及i可根据方法的局部截断误差可根据方法的局部截断误差及阶确定,其定义为及阶确定,其定义为 定义定义7 设设y(x)是初值问题是初值问题(1.1), (1.2)的准确解,的准确解,线性多步法线性多步法(5.1)在在xn+k上上局部截断误差局部截断误差为为若若Tn+k=O(hp+1),则称方法,则称方法(5.1)是是p阶的阶的,p1则称方则称方法法(5.1)与方程与方程(1.1)是是相容的相容的.党兹泉瞒撇砚葫桌脑言钎煎铱蕉隶迅元汐闽鞍个堑影印囚柬褐萤庸腔畏靡第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页由定义
42、由定义7,对,对Tn+k在在xn处泰勒展开,由于处泰勒展开,由于代入代入(5.2)得得其中其中阶袭掉淤胚洁漏砚晾稍亭拙搞复戒奎捆咏挎意颁嘿钩身荔焉蛆薪撬允枢活第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页若在公式若在公式(5.1)中选择系数中选择系数i及及i,使它满足,使它满足由定义可知此时所构造的多步法是由定义可知此时所构造的多步法是p阶的,且阶的,且称右端第一项为称右端第一项为局部截断误差主项局部截断误差主项, cp+1称为称为误差常数误差常数.郝浩赣笺满拧闻尝财纤尹客褐跳腥袄辽紧辗若飘州离糯决净宇旨蛛弟衅鹰第9章常微分方程初值问题数值解法第9章常微分方程
43、初值问题数值解法上页上页下页下页根据相容性定义根据相容性定义, p1,即,即c0=c1=0,由,由(5.4)得得故线性多步法故线性多步法(5.1)与微分方程与微分方程(1.1)相容的充分必要相容的充分必要条件是条件是(5.6)成立成立.显然,当显然,当k=1时,若时,若1=0,则由,则由(5.6)可求得可求得0=1,0=1.此时公式此时公式(5.1)为为沾评贝筹硝渔板蔡毙砌叉邻氓辕袋舶滚坛赋动岩谬盏环锑强妖妆晰隐圾皖第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页即为欧拉法即为欧拉法. 从从(5.4)可求得可求得c2=1/20,故方法为,故方法为1阶阶精度,且局
44、部截断误差为精度,且局部截断误差为这和第这和第2节给出的定义及结果是一致的节给出的定义及结果是一致的.对对k=1,若,若10,此时方法为隐式公式,为了确,此时方法为隐式公式,为了确定系数定系数0,0,1,可由,可由c0=c1=c2=0解得解得0=1, 0=1=1/2.于是得到公式于是得到公式即为梯形公式即为梯形公式.彦毛膀癸饺骋禾斯梧辽炬胶杂梁桨阴线党捅柜凌着犹占贺观媒俩康吠咙帽第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页由由(5.4)可求得可求得c2=- -1/12,故,故p=2,所以梯形法是,所以梯形法是二阶方法,其局部截断误差主项是二阶方法,其局部截断
45、误差主项是这与第这与第2节中讨论是一致的节中讨论是一致的.对对k2的多步法公式都可利用的多步法公式都可利用(5.4) 确定系数确定系数i,i,并由并由(5.5)给出局部截断误差,下面只就若干给出局部截断误差,下面只就若干常用的多步法导出具体公式常用的多步法导出具体公式.遥搞邮泛趣湖萝陈齐挽捕蜂泣气侠侍花争抨蛙烂芋丽菲弹皑禽狠研钓痈没第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页9.5.2 阿当姆斯显式与隐式公式阿当姆斯显式与隐式公式p362自学自学.9.5.3 米尔尼方法与辛普森方法米尔尼方法与辛普森方法p366自学自学.9.5.4 汉明方法汉明方法p367自
46、学自学.9.5.5 预测预测- -校正方法校正方法p368自学自学.9.5.6 构造多步法公式的注记和例构造多步法公式的注记和例p371自学自学.撇妥嘉刨磷鬃硒舅犬卒斩牢弗行足捞喉雁土蛛曹拂碱酵益壹皱沽细净耍柔第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页9.6 方程组和高阶方程方程组和高阶方程9.6.1 一阶方程组一阶方程组 前面我们研究了单个方程前面我们研究了单个方程y = =f 的数值解,只的数值解,只要把要把y 和和f 理解为向量,那么,所提供的各种计算公理解为向量,那么,所提供的各种计算公式即可应用到一阶方程组的情形式即可应用到一阶方程组的情形. 考
47、察考察一阶方程组一阶方程组的初值问题,初始条件给为的初值问题,初始条件给为队菌侈凸反探蔽菲炭史迎靡窍胀邹搐跋伍屡砍佰磋兵冀阐坡灸肖鳖矩透何第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页若采用向量的记号,记若采用向量的记号,记(向量向量)则上述方程组的初值问题可表示为则上述方程组的初值问题可表示为催挥颓钥仗构板妮奋溃电受虽娃办妈桔事设七狈栗翟熬间开迅韭稗摇撕辑第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页求解这一初值问题的求解这一初值问题的四阶龙格四阶龙格- -库塔公式库塔公式为为(向量向量)式中式中(向量向量)资云伯吧叉罪借铆
48、糖郊烂元航君哉祖栈脏邵贡铁鲍和胺跺还闸誊旗甸希判第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页或表示为或表示为(分量分量)其中其中(分量分量)这里这里yin是第是第i个因变量个因变量yi(x)在节点在节点xn=x0+nh的近似值的近似值.痛芝翌茅动撬届央滑斜哄仰苹达猛酷烈睛柠坷勾抽梁昏至眠噶梗茎米剔茹第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页为了帮助理解这一公式的计算过程,我们再考察为了帮助理解这一公式的计算过程,我们再考察两个方程的特殊情形两个方程的特殊情形这时这时四阶龙格四阶龙格- -库塔公式库塔公式具有形式具有形式嚏
49、苍债囊凯陡余假净淋劳姿久香靳梅妙蝶谣淹销拷省剃茵换垒斗献分卷匪第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页其中其中熙膘灵豢臀榨骚录桌耶揍猫涣烫码茸扰坝瓜汲帛凿莱考浅锹咸踏历末岔寅第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页这是一步法,利用节点这是一步法,利用节点xn上的值上的值yn, zn,由,由(6.3)式顺序计算式顺序计算K1,L1,K2,L2,K3,L3,K4,L4,然后代入,然后代入(6.2)式即可求得节点式即可求得节点xn+1上的上的yn+1, zn+1.衣粮缠蚕编筹需施智雍泛挟它椿侄锤稽谣彤罐氮吉矿高蝗缔肇钨凭
50、屁嘘撵第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页9.6.2 化高阶方程为一阶方程组化高阶方程为一阶方程组关于高阶微分方程关于高阶微分方程(或方程组或方程组)的初值问题,原则的初值问题,原则上总可以归结为一阶方程组来求解上总可以归结为一阶方程组来求解. 例如,考察下列例如,考察下列m阶微分方程阶微分方程初始条件为初始条件为只要引进新的变量只要引进新的变量即可将即可将m阶方程阶方程(6.4)化为如下的一阶方程组:化为如下的一阶方程组:祁棚咱巾范穿这戊嗡籽刺笺父蒙钳衙茄圃脚芜晨勉罗神崎矫宵拴娘绒钓工第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法
51、上页上页下页下页初始条件初始条件(6.5)则相应地化为则相应地化为不难证明初始条件不难证明初始条件(6.4),(6.5)和和(6.6),(6.7)是是彼此等价的彼此等价的.胯纳荡旅丝玛鬃拐熔屿型勤刑毙潍陌呈笺档塘畸侍明歹铅毖鄙各孵柔尔蓉第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页特别地,对于下列二阶方程的问题特别地,对于下列二阶方程的问题引进新变量引进新变量z= =y ,即可化为下列一阶方程组的初值问题即可化为下列一阶方程组的初值问题:芯擞加琅杂鸟鞋歹井糟病缔古鸥龋防鲸什刚钎赊认月处香剑濒撵筏荐鞍象第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值
52、解法上页上页下页下页针对这个问题应用针对这个问题应用四阶龙格四阶龙格- -库塔公式库塔公式(6.2),有,有由由(6.3)式可得式可得颈匿馒穗膨财垒腕煌悯摊惶起安横蔽渝电雅丘币帽鞋烙忌捎蓑反婴钝障硼第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页如果消去如果消去K1,K2,K3,K4,则上述格式可表示为,则上述格式可表示为这里这里胸泥宁濒榆冶诗墒撂宇残稀坯翼果杂踪径清本跨丧住砰夸短蠕暇涎地刀树第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页9.6.3 刚性方程组刚性方程组p378自学自学.汤罚稳彦林两溢霸咯佬纤院绚俩鞘钒冲睫醒包铀饱晃楼秘袱挟倒丽蜂卒叠第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法上页上页下页下页课程结束课程结束! !希望同学们好好复习希望同学们好好复习! !争取考个好成绩争取考个好成绩! !再见再见! !慎唯帕涂满烟赤新绳严蚀扬甚泛夜袋矽淑馆阂雏阉箱东晤彰邹仙厨淋锻虏第9章常微分方程初值问题数值解法第9章常微分方程初值问题数值解法