数值分析线性多步法

1第五节线性多步法2©2009,HenanPolytechnicUniversity2§5线性多步法第五章常微分方程数值解法龙格-库塔方法是一类重要算法，但这类算法在每一步都需要先预报几个点上的斜率值，计算量比较大。考虑到计算yi+1之前已得出一系列节点上的斜率值，能否利用这些已知值来减少计算量呢？这就是线性多步法的设计思想。构造线性多步法有多种途径，这里介绍两种：基于数值积分的构造方法；基于Taylor展开的构造方法。)...(...110111101rnrnnnrnrnnnffffhyyyy--+---+++++++++=bbbbaaa线性多步法的通式可写为：),(jjjyxff=当b-10时，为隐式公式;b-1=0则为显式公式。3©2009,HenanPolytechnicUniversity3§5线性多步法第五章常微分方程数值解法5.5.1基于数值积分的构造法dxxyxfxyxykkxxkk+=-+1))(,()()(1将在上积分，得到),(yxfy=],[1+kkxx只要近似地算出右边的积分，则可通过近似y(xk+1)。而选用不同近似式I，可得到不同的计算公式。例如利用左矩形积分公式得到欧拉公式；梯形积分公式得到梯形公式。dxxyxfIkkxx+1))(,(Iyykk+=+14©2009,HenanPolytechnicUniversity4§5线性多步法第五章常微分方程数值解法一般地，利用插值原理所建立的一系列数值积分方法也可以导出解微分方程的一系列计算公式。运用插值方法的关键在于选取合适的插值节点。假设已构造出f(x,y(x))的插值多项式Pr(x),则dxxPdxxyxfkkkkxxrxx++11)())(,(dxxPyykkxxrkk++=+1)(15©2009,HenanPolytechnicUniversity5§5线性多步法第五章常微分方程数值解法5.5.2亚当姆斯显式公式利用r+1个节点上的被积函数值构造r阶牛顿后插多项式,有11,,,,,,nnnnnrnrxfxfxf----11100(,())()()nnnxxthxrnrnxfxyxdxNxthhdtRxthhdt+=+=+++Newton插值余项0()11,rjjrnnjjjntNxthfjxxtjh-=-+=--=，其中0表示阶向前差分，6©2009,HenanPolytechnicUniversity6§5线性多步法第五章常微分方程数值解法11100(,())()()nnnxxthxrnrnxfxyxdxNxthhdtRxthhdt+=+=+++1110001000()11rjjnnrnnnjjrrjjjnnjnjnjjjtyyhNxthdtyhfdtjtyhdtfyhfja+-=--==-=++=+--=+-=+Adams显式公式dtjtrnjjj--=10)1(aa的系数，和为不依赖于jaj0111/225/1233/87©2009,HenanPolytechnicUniversity7§5线性多步法第五章常微分方程数值解法实际计算时，将差分展开000000111jijnjniijrrrrriijjnjjnijnirinijjiijiijffijjffffiiaaab--=----=======-=-=-=10rnnriniiyyhfb+-==+jrijiriijab=-=)1(其中：8©2009,HenanPolytechnicUniversity8§5线性多步法第五章常微分方程数值解法局部截断误差为：1110(2)12(2)102(2)()()()1!()rnnrnrnrrrrnrrrnRyxyhRxthdtyhtdtBhyrBhyx++++++++=-=+==+例：1110,,1,(3)2nnnnnnnryyhfhryyff++-==+==+-35()12nRhyx常用的是r=3的4阶亚当姆斯显式公式)9375955(243211---+-+-+=iiiiiiffffhyy5(5)251()720nRhyx21()2nRhyx9©2009,HenanPolytechnicUniversity9§5线性多步法第五章常微分方程数值解法Adams显式公式用作为插值节点，在求积区间[xn,xn+1]上用插值函数Nr(x)近似f(x,y(x))，而xn+1不在插值节点内，因此是一个外推的过程。虽然公式是显式的，便于计算，但效果并不理想，比如稳定性较差等。因此改用通过r+1个节点的插值多项式Nr(x)近似f(x,y(x))，由于xn+1是其中一个插值点，因此是内插多项式，但导出的公式是隐式的。11,,,,,,nnnnnrnrxfxfxf----1111,,,,,,nnnnnrnrxfxfxf++-+-+10©2009,HenanPolytechnicUniversity10§5线性多步法第五章常微分方程数值解法5.5.3亚当姆斯隐式公式利用r+1个节点上的被积函数值fn+1,fn,…,fn-r+1构造r阶牛顿前插多项式。与显式多项式完全类似地可得到一系列隐式公式。0**1110***110,1,1rjjnnjnjjjrrinnrinirijijityyhfdtjjyyhfiaabba++--=++-==-=+=-=+=-差分展开，得j011-1/22-1/123-1/24*ja局部截断误差为：2(2)112(2)()()()rrrnnrnrrrnRyxyBhyBhyx++++++=-=11©2009,HenanPolytechnicUniversity11§5线性多步法第五章常微分方程数值解法常用的是k=3的4阶亚当姆斯隐式公式)5199(242111--+++-++=iiiiiiffffhyy较同阶显式稳定例：11110,,1,()2nnnnnnnryyhfhryyff++++==+==++31()12nRhyx-21()2nRhyx-5(5)19()720nRhyx-12©2009,HenanPolytechnicUniversity12§5线性多步法第五章常微分方程数值解法5.5.4亚当姆斯预测-校正系统Step1:用Runge-Kutta法计算前i个初值；Step2:用Adams显式计算预测值；Step3:用同阶Adams隐式计算校正值。注意：三步所用公式的精度必须相同。通常用经典Runge-Kutta法配合4阶Adams公式。+-++=-+-+=--+++---+]519),(9[24]9375955[24211113211iiiiiiiiiiiiifffyxfhyyffffhyy13©2009,HenanPolytechnicUniversity13§5线性多步法第五章常微分方程数值解法4阶Adams隐式公式的截断误差为4阶Adams显式公式的截断误差为)(720251)()5(5111fhyxyeiii=-=+++)(72019)()5(5111fhyxyeiii-=-=+++1111()19()251nnnnyxyyxy++++---14©2009,HenanPolytechnicUniversity14§5线性多步法第五章常微分方程数值解法1111251()()270nnnnyxyyy++++---111119()()270nnnnyxyyy++++--预测值与校正值的事后误差估计1111251()()270nnnnyxyyy++++--111119()()270nnnnyxyyy+++++-15©2009,HenanPolytechnicUniversity15§5线性多步法第五章常微分方程数值解法Adams显隐预测-校正-改进系统)9375955(243211---+-+-+=iiiiiiffffhyp预测：)(27025111iiiipcpm-+=++改进：)519),(9(2421111--++++-++=iiiiiiifffmxfhyc校正：)(270191111++++--=iiiipccy改进：16©2009,HenanPolytechnicUniversity16§5线性多步法第五章常微分方程数值解法5.5.5基于泰勒展开的构造法)...(...110111101rnrnnnrnrnnnffffhyyyy--+---+++++++++=bbbbaaa将通式中的右端各项yn-1,…,yn-k;fn+1,fn-1,…,fn-k(y’n+1,…,y’n-k)分别在xn点作泰勒展开，与精确解y(xn+1)在xn点的泰勒展开作比较。通过令同类项系数相等，得到足以确定待定系数a0,…,ar;b-1,b0,…,br的等式，则可构造出线性多步法的公式。),(jjjyxff=当b-10时，为隐式公式;b-1=0则为显式公式。17©2009,HenanPolytechnicUniversity17§5线性多步法第五章常微分方程数值解法++-+-+=-+=++==-)1(1)(1)(1)!1()(!)(!)(mimjimjjijijjikiymkhyjkhyyjkhyy+-+--=--=-=+=-=-+=-)1()(11)(11)1(0!)()!1()()!1()(!)(mimjimjjjijjjijjkiymkhyjkhyjkhyjkhy=-=--++=pkpkkikkikiybhyay01118©2009,HenanPolytechnicUniversity18§5线性多步法第五章常微分方程数值解法=-=--++=pkpkkikkikiybhyay011-=+=-==++++-+--+++-+-+=pkmimjmjijkpkmjmimjijikiymkhyjkhbhymkhyjkhyay1)1(1)(101)1(1)(1...!)()!1()(...)!1()(!)(...])()([)!1(])()([!)1(10111)(1100+-+-++-+-+=+-==++=-=-==mipkmkpkmkmmjjipkjkpkjkjpkikykbkamhykbjkajhya19©2009,HenanPolytechnicUniversity19§5线性多步法第五章常微分方程数值解法)()(!)(!2)()()(1)(21+++++++=mimmiiiihOxymhxyhxyhxyxy==-+-==-=-=pkpkjkjkpkkmjkbjkaa0110),,2,1(,1)()(1=-=++++++++-+---+=-=pkpkmimmkmkmiiihOxykbmkamhyxye012)1(11111)()(])()1()(1[)!1(~)(...])()([)!1(])()([!)1(10111)(11001+-+-++-+-+=+-==++=-=-==+mipkmkpkmkmmjjipkjkpkjkjpkikiykbkamhykbjkajhyay20©2009,HenanPolytechnicUniversity20§5线性多步法第五章常微分方程数值解法)(33221102211020301-----==--+++++++=+=iiiiiiikkkikkikiybybybybhyayayaybhyay=---+=+++--=---+=++++--=++11083241612712381642412132121321213212132102121