15-线性多步法解析课件.ppt

1、机动 目录 上页 下页 返回 结束 1.5 1.5 线性多步法线性多步法 0010010110Newton(),(),()()()=(0,1,)nnniNxf xf x xxxf x xxxxxxxxxxihin插值多项式：当节点等距分布时 00000(1)10Newton=+,()(+)()()()(1)(),(,)(1)!(1)(2).(1),10!nknnknnnnxxtx xthNxNxthf xkfRxt ttn hx xnttt tttkkk 前插公式（一般当 靠近 时用）设则一、预备知识预备知识 1.5 1.5 线性多步法线性多步法 111011Newton,(),(),()()

2、()nnnnnnnnnnxxNxf xf xxxxf xxxxxxxxx后插公式（一般当 靠近 时用）将节点顺序倒置 则 0(1)10,()(+)()()()(1)(),(,)(1)!nkknnnnnknnnntxxthNxNxthf xkfRxt ttn hx xn设则（-1）二、线性多步法欧拉法和梯形法的回顾求未知函数在 的近似值 ，基本思想：对于积分表达式1ntt1nu11()()(,()nntnntu tu tf t u t dt将被积函数 用水平直线或连接 ，两点的直线 代替。1(,),nnf t u tt t(,)(,()nnf t uf t u t(,(,()nnntf t u

3、t111(,(,()nnntf tu t111(,)(,()(,()nnnnnnttttf t uf t u tf tu thh010nnnttnhuuu-=+单步法：计算节点的近似值 仅用到前一节点的值，故从初值 可算出以后各节点的值。111222(,(),(,(,(),(,(,()nnnnnnnnt u ttf tu ttf tu t2(,()L t u t的曲线为122221122212()()()()(,()(,()(,()2()()(,()2nnnnnnnnnnnnttttttttL t u tf tu tf tu thhttttf t u th112()()(,()nntnntu

4、tu tL t u t dt221151623()(,()(,()(,()121212nnnnnnnu thf tu tf tu tf t u t像欧拉格式一样，我们可得近似求解格式如下：1221151623(,)(,)(,)121212nnnnnnnnuuhf tuf tuf t u曲线，如我们利用 Lagrange Lagrange 插值可得经过1,nnt t(,()f t u t然而为了近似上的曲线也可用多点插值1,(,()nnt tf t u t+用它近似中的，则得到积分近似值2211(,),(,),(,)nnnnnntutut u为提高精度，构造线性多步法，一般形式为其中00(,),

5、0,0njnjnjjjkff tu和是常数，且和不同时为(1.41)1100()()()kkjnjkn kkn kknjnjjjun un un uhf=11122(,)(,),(,),nnnnnnnnuuhf t uf tuf tu+-=+上面的格式与欧拉格式不同之处在于增加了包括的两项即由()111.5(,)1nnnnut uu+计算，与欧拉格式仅由前面一点计算的这种单步法不同，格式称为多步法。用(1.41)计算 需用到前 k 个节点的值 故称为多步法或 k 步法。又因(1.41)关于 是线性的，故称为线性多步法。11,nnn ku uu，n kunjf注意：需附加初值121,;ku uu

6、00kk为显式；为隐式00(,)|t tduf t udtuutu0121ktttt0121kuuuu现在，设已给出常微分方程初值问题011011(),.,.,kku tt ttu uu-的解在处的近似值，或者说给出表头(一)数值积分法适当取k+1个节点，用的k次Lagrange插值多项式近似11()()(,()nntnntu tu tf t u t dt将方程 写成积分形式，比如在 上积分得,nnt th(,)uf t u(1.42)(,()f t u t(,()f t u t,()n kLt 线性多步法不同的插值节点导出不同的多步法(1)Adams外插法 （显式多步法）,(),n kfLa

7、grangeLt的插值多项式则,(,()()()n kn kf t u tLtrt1,nnn kt tt-取为节点，构造其中 是插值余项，代入(1.42)得knr,111,()()()()nnnnttnnn kn kttu tu tLt dtrt dt舍去余项1,()nntn kn ktRrt dt(),jjuu t用代替即得11,()nntnnn ktuuLt dt,n kR 为局部截断误差(1.43)(1.44)下面给出(1.44)的具体形式。由插值节点等距，且被插值点 靠近最后一个节点 ，故用牛顿向后插值公式 1,nntt tnt2,12(1)()()1!2!n kn knnnnLtLt

8、hfff(1).(1)!kn kkfk 引进记号10,!)1).(2)(1(sjjssssjs则,0()(1)()kjjn knjjLtfj因此10kjnnjnjjuuhaf(1.45)式中10(1)0,1,2,.jjadtjj(,).jnjnjnjjff tu式中表示 阶向前差分，式(1.45)就是著名的 外插公式,即为Euler法Adams0k 0|1,|1jjjjata t显然，故当时级数收敛，且0()jjjG ta tja系数 的计算：定义生成函数11000010()(1)()(1)(1)ln(1)jjjjjG tt dtdjjttdtt ln(1)1()1tG ttt从而nt比较 的

9、系数，则有22201211(1)()1+23ttaa ta ttt 两端展开成幂级数有12011110,1,231nnnaaaann2.1Adamsja由此可依次算出系数，见表。最终得到外插公式例如0100:(,)nnnnnnnnkuuhafuhfuhf t uj101111111:)231(,)(,)22nnnnnnnnnnnnnkuuh a fafuh fffuhf t uf tu（2101122112111211222:()15()()21215()()()21223165(,)(,)(,)121212nnnnnnnnnnnnnnnnnnnnnnnnnnkuuh a fafafuh ff

10、fffuh fffffffuhf t uf tuf tu 111223355593:(,)(,)2424379(,)(,)2424nnnnnnnnnnkuuhf t uf tuf tuf tu0(1)jjlnjn lljffl 则Adams外插公式可改写为10()knnkln lluuhb f显式其中(1)klkljj ljbal 的值见表2.2。kibAdams1,()nntn kn ktRrt dt11(2),()(1)(),1kkkn kn knn knrrthhuttk 外插公式的进一步分析：将(1.45)右端的差分表为 的线性组合，为此由差分公式njf 外插格式的局部截断误差为则12

11、1(2)2(2),10(1)()d()1kkkkkn kkRhuahuk2()kn knO htt因此，外插格式的局部截断阶为 。外插法是个k+1步k+1阶显示格式。Adams)(2khOAdams(2)Adams内插法 （隐式多步法）(1),(1),()(,()(),n kn ku tf t u tkLagrangeLtrt构造或的+1次插值多项式插值余项为(),则(1)(1),(,()()()n kn kf t u tLtrt从而11(1)(1)1,()()()()nnnnttnnn kn kttu tu tLt dtrt dt舍去余项1(1)(1),()nntn kn ktRrt dt(

12、),Adamsjjuu t并用代替即得内插法1(1)1,()nntnnn ktuuLt dt(1),0Eulern kkR 时即为改进的法，为内插法的局部截断误差(1.46)111,()nnnn knttttt+-+取为节点，比外插多一点内插法(1.46)的具体化1*110kjnnjnjjuuhaf 0*1(1)0,1,.,1jjadjkj 其中，由牛顿向后插值公式 (1)(1)(1)(1)(1)!jjjjjj 1,0,其中二项系数1(1)(),+110()()(1)ktjjn kn knnjjLtLthfj(1.47)(1.47)代人(1.46)得Adams内插公式用生成函数法可导出系数 的

13、递推公式，见表2.3*ja(1.48)0(1)jjlnjn lljffl 则Adams内插公式可改写为+1*11,10()knnkln lluuhbf 隐式其中1*1,(1)klkljj ljbal 的值见表2.4。*1,klbAdams11(1)(1)(1),102(3)11()()(1)(),2nnnnttn kn kn knttkkkn knRrt dtrth dthuhdttk 内插公式的进一步分析：将(1.48)右端的差分表为 的线性组合，为此由差分公式njf 内插格式的局部截断误差为则(1)*3(3),213(),()kkn kkn knkRahuttO h举例：*0*1100 1

14、10:nnnnkuuh afaf 111()2nnnnuhfff12nnnhuff 梯形格式*1*21011211:nnnnnkuuh a fafaf111111()()()()212nnnnnnnnuhfffffff 11581121212nnnnuhfff11111112:()(2)212nnnnnnnnkuuhffffff1121(33)24nnnnffff1129195124242424nnnnnuhffffAdams内插法与外插法的比较达到相同的误差阶，内插法比外插法可少用一个初始已知量*1,1.|,|jjklklaabb内插法的舍入误差影响比外插法小322.(),()kkO hO

15、h内插法的局部截断误差阶为外插法的局部截断误差阶为3.外插法显式格式；内插法隐式格式Adams内插法是k+1步k+2阶的内插法的迭代求解(0)11(1)*()*11,0111,111(,)(,)nkmmnnknnkin ln lluuuhbf tuhbf tu 初始值0,1,m 收敛条件*111,01,0(,)|1nnkkf tuhbh bLu取h充分小，可满足收敛条件实际上，由(1)*()111,01111(,)(,)mmnnknnnnuuhbf tuf tu(1)*()1111,011(,)|mmnnnknnf tuuh buuu故当时迭代法收敛。*111,01,0(,)|1nnkkf t

16、uhbh bLu当 充分小时,迭代程序收敛，越小，收敛速度越快。可以用相应的外插公式作为预报格式，然后用内插格式进行迭代：hh(0)101(1)*()*11,0111,111(,)(,)knnkln llkmmnnknnkin ln lluuhb fuuhbf tuhbf tu 预报格式迭代格式初始近似：(用外插法给出)(0)10knnkln lluuhb f即内插法的求解格式如下由于预报格式已是很好的近似，故为了达到较高精度，仅需迭代很少次数，一般二、三次迭代就已足够，特别如只进行一次迭代，就得到如下预报修正格式：例如，上式中令 ，则有3k(0)1123(0)11112(5559379)24

17、9(,)19524nnnnnnnnnnnnnhuuffffhuuf uufff预报格式修正格式(0)101*(0)*11,0111,111(,)(,)knnkln llknnknnkin ln lluuhb fuuhbf tuhbf tu 预报格式修正格式注：1.Adams内外插法思路：先将初值问题改写为积分形式，再用适当的数值积分离散；2.也可将初值问题写成其他积分形式，例如22()()(,()nntnntu tu tf t u t dt再用适当的数值积分代替右端积分，例如，用Simpson公式，可得212(4)3nnnnnhuufff为线性二步法。3.数值积分法只能构造一类特殊的多步法，其

18、系数1,1,0(,),kk mllkm k 当(二)待定系数法0();()().kjjjL u t hu tjhhu tjh令()()()Tayloru tu tjhu tjht设为初值问题的解，将和在点 用公式展开，按h的同次幂合并同类项，得(1)()01();()()()qqqL u t hc u tc hutc h ut其中00111201121112,2()1(2)!1(2),2,3,(1)!kkkqqqkqqkcckckqkqq01+1()+2,0,0,jjppjju tpcccc若有次连续微商，则可选取k(足够大)和使而，即选满足011201111212=0,2()=011(2)(

19、2)=0,!(1)!kkkppppkkkkkpp1(1)21();()()ppppL u t hchutO h此时()=(,(),u tf t u t由于则,01(1)2,1()(,()()()kjnjnnn kjpppn kpu tjhhf tjh u tjhRRchutO h,1()(,)(),()n knjnnjnjnjppRuu tjhff tuO hO hpk略去余项，并用代替，用记，即得如下线性多步法，局部截断误差为可以证明其整体截断误差为，故称为 阶 步法。11000()()()kkjnjkn kkn knjnjjjun un un uhf=kkk不妨设=1，=0为显式格式；0为

20、隐式格式jj待定系数法构造多步法的思路：选取，使局部截断误差的阶尽可能高。例：讨论一般的二步法.20101201232,.0kcccc此时=1记，其余四个系数，由确定，即有011101221123112+1=0,2()=01(4)(2)=0211(8)(4)=062cccc2121(1)(5)8(1)(15)12nnnnnnhuuufff=41125112111(16)(8)(1)24624111(32)(16)(1713)12024360cc 而且一般的二步法为10121(1),(15),1221(1),(5)312 22143nnnnnhuufff+为4阶二步法，是具有最高阶的二步法，称为

21、Milne法。若取 为二步Adams内插法；若取 ，则为显式格式05 44510,=1=0,0ccc 当时为三阶2步法；当时但，格式变为(三)预估校正算法将隐式k步法 改写为,0,1,1,njn kujku若已求出则上式为关于的非线性方程，采用如下迭代法求解(,)njnjnjff tu其中00kkjnjjnjjjuhf1100(,)kkn kjnjjn kn kjnjjjuuhf tuhf11100(,),0,1,n kn kkkmmjnjjn kjnjjjuuhf tuhfm001/(|),Lipschitzn kn kkuhLLfuu其中为给定的迭代初值。若为 关于的常数，初值选择合适，迭

22、代法收敛。校正算式(C算式)0n ku可用显式多步法计算，:(,)njnjmmnjEff tu预估算式和校正算式统称为预估校正算法，简称预校算法或PC算法。计算过程如下（一）110*00,n kjjkknjnjjjuuhfP预估算式(算式)110,1,1(,)n kn kn kMMMn kmMuff tu其中。得到的结果和为下一步预估计算所用110*100:,n kjnjjnjkkMMjjP uuhf111100:,n knjn knjkkmMmMjkjjjC uuhfhf()MP EC:(,)njnjmmnjEff tu（二）110*100:,n kjnjjnjkkMMjjP uuhf111

23、100:,01n knjn knjkkmMmMjkjjjC uuhfhfmM:(,)njnjMMnjEff tu()MP ECE任一显式多步法和隐式多步法都可搭配成预校算法(四)多步法计算的问题1.需附加初值为保持精度，计算附加初值时要节点加密或采用同样精度的Runge-Kutta法011,ku uu2.如何选择阶p(或步数k):要考虑光滑性和稳定性及总的工作量；高阶多步法的绝对稳定域小，3.如何选择步长h:可根据精度要求由粗到细逐渐调整，计算过程中可以改变步长，但计算过程复杂 212(;)()ppu t hu tAhA hA h(五)Richardson Richardson 外推法外推法利

24、用 为计算步长，则得2h212(;)()()()2222pphhhhu tu tAAA()(;)(RichardsonO htu t h对于多步法而言，表头的值式非常重要的，如果表头的值精度低，则用精度再高的多步格式也无济于事，要得到高精度的表头值，除了用高阶单步法外，也可从较低阶单步法如欧拉法，配合缩短步长以及外推技巧外推法 等计算得高精度的开始值，例如设我们利用欧拉法计算得开始点 处的值为，由于欧拉法整体截断误差为，故可写成则23232112(;)(;)()(1)222hu tu t hu tA hA h因此，取 作为 的近似值，误差将提高一阶为 ，如果用 作为的近似值，整体截断误差为 ，

25、则有2(;)(;)2hu tu t h()u t)(2hO()u t)(phO112(;)()ppu t hu tAhA h利用 为计算步长，有2h1121(;)()222pppphhhu tu tAA消去 项，则得ph112(;)(;)()()212ppphu tu t hu tO h()1()()nO hRichardsou tn+以它作为的近似值的截断误差，精度提高一阶，这就是利用外推技巧外推法提高格式精度的方法。谢谢 谢谢作业：作业：作业要求：写出程序，列表或用图形显示结果，并作业要求：写出程序，列表或用图形显示结果，并给出图或表所说明的结果给出图或表所说明的结果P20-22P20-22，例，例2.1,2.22.1,2.2和习题和习题 1 1，3 3