11.2Runge-Kutta法 11.3线性多步法

龙格--库塔法欧拉法实际问题导入第十一章常微分方程数值解/*OrdinaryDifferentialEquationsNumericalSolution*/线性多步法11.2Runge-Kutta方法一、Runge-Kutta法的基本思想1、基本思想Runge和Kutta不是用求微商的方法，而是计算不同点上的函数值，然后对这些函数值作线性组合，建立近似公式，把近似公式和Taylor展开式作比较，使得前面若干项相吻合，从而使近似公式达到一定的阶数。对于Euler公式，写成平均化的形式为每步计算f一次总体误差O(h)001)(),(yxyyxhfyynnnn2、具体做法),(111nnnnyxhfKKyy对于改进的Euler公式，写成平均化的形式为),(),(2121121211KyhxhfKyxhfKKKyynnnnnn),(),(2),(1111nnnnnnnnnnyxfyxfhyyyxhfyy每步计算f两次总体误差O(h2)注：两组公式的相同点是f(x,y)在某点处的函数值的线性组合，得到近似值，而且增加f(x,y)的计算次数，可提高总体截断误差的阶数，按此想法，若设法在[xn,xn+1]内多预估几个点，这些点处函数值做加权平均值，建立高精度公式。一般的Runge－Kutta法构造，设计计算格式如下适当的选择参数，使其具有m阶精度。),(......),(),(),(11,22112321313312122122111mmmmmnmnmnnnnnnmmnnKbKbKbyhaxhfKKbKbyhaxhfKKbyhaxhfKyxhfKKcKcKcyy11miic二、二阶Runge-Kutta公式适当的选择参数c1,c2,a,b，使其具有二阶精度。),(),(12122111bKyahxhfKyxhfKKcKcyynnnnnn按照二元函数泰勒级数展开：将K1,K2,代入yn+1得：),(12bKyahxhfKnn......),(1yxnnfbKahfyxfh......)()(22211fbfafchfcchyyyxnn另一方面：)()(1hxyxynn...)(!3)(!2)()(32nnnnxyhxyhxyhxy...)2(6)(!232fffffffffhfffhhfyyxyyyxyxxyxn22111KcKcyynn......2yxbffafhhf欲使设计的格式具有二阶精度，只需满足上式条件的参数，无论如何取法，都要计算f两次，局部截断误差O(h3)，都是二阶方法，统称为二阶Runge-Kutta法。改进的Euler格式212112221bcaccc1,2121bacc)32,32(),()3(41121211KyhxhfKyxhfKKKyynnnnnn32,43,4121baccHeun公式三、三阶Runge-Kutta公式问题：若不增加函数值的计算次数，能否适当选取4个参数，使近似公式的局部截断误差的阶再提高，达到O(h4)？),(12bKyahxhfKnn试算一下：)(22)(422232hOffbfabffahfbfafhhfyyxyxxyx)(22)()(42222322211hOffbfabffachfbfafhcfcchyyyyxyxxyxnn)()(!3)(!2)()()(4321hOxyhxyhxyhxyxynnnnn)()2(6)(!2)(4232hOfffffffffhfffhhfxyyxyyyxyxxyxn注：如果不增加函数值的计算，无法提高公式的精度！三阶Runge-Kutta法8个参数满足常见的三阶Runge-Kutta公式),(),(),(232131331212213322111KbKbyhaxhfKKbyhaxhfKyxhfKKcKcKcyynnnnnnnn61312113322323222332232313212321cbacacacacabbabaccc)2,()2,2(),()4(612131213211KKyhxhfKKyhxhfKyxhfKKKKyynnnnnnnn四、四阶Runge-Kutta法13个参数满足11个方程可得到四阶Runge-Kutta公式。),(),(),(),(3432421414423213133121221443322111KbKbKbyhaxhfKKbKbyhaxhfKKbyhaxhfKyxhfKKcKcKcKcyynnnnnnnnnn常见的四阶Runge-Kutta法),()2,2()2,2(),()22(61342312143211KyhxhfKKyhxhfKKyhxhfKyxhfKKKKKyynnnnnnnnnn优点：1、一步法，给定初值后，可逐步计算；2、精度高；3、计算过程便于改变步长。缺点：计算量大，一步计算四个函数值1)0(yyxy例：112341213243(22)62222nnnnnnnnnnhyyKKKKKxyhhKxyKhhKxyKKxhyhK四阶Runge-Kutta法格式为取步长h=0.5,计算结果为1,0,),(00yxyxyxf11.51.6252.3131.7962.303.123.3264.4593.4334.435.796.1307.9986.4567.9610.1210.7513.8311.7613.7617.4518.3723.4420.83K1K2K3K4yn}yxi()11.3线性多步法线性多步法的建立常用方法•利用数值积分的方法•基于泰勒展开式的待定系数法单步法----计算yn+1时只使用yn的值。多步法的概念多步法----计算yn+1时使用前面的k个yi值，即由yn-k+1,yn-k+2,…,yn-1,yn计算yn+1。(k=1,2,……)线性多步法----计算yn+1的公式由yn-k+1,yn-k+2,…,yn-1,yn的线性组合表达。基本思想求初值问题等价于求积分方程1))(,()()(1nnxxnndxxyxfxyxy00)(],[),,(yxybaxyxfy取不同的数值方法计算积分，可导出不同的计算公式：1、用矩形法计算积分项：))(,())(,()()(11nnxxnnxyxhfdxxyxfxyxynn得到：(欧拉公式)),(1nnnnyxhfyy基本思想2、用梯形法计算积分项：))(,())(,(2))(,()()(1111nnnnxxnnxyxfxyxfhdxxyxfxyxynn得到：(改进Euler公式)),(),(2111nnnnnnyxfyxfhyy3、推广：由插值法可建立一系列数值积分，运用这些方法可导出求初值问题的系列计算公式：具体做法：利用插值多项式近似取代，)(xk)())(,(xFxyxf用近似取代得到公式：1)(nnxxkdxx1))(,(nnxxdxxyxf1)(1nnxxknndxxyy四阶阿达姆斯(Adams)外插（显式）公式设xn-3,xn-2,xn-1,xn处的函数值分别为yn-3,yn-2,yn-1,yn如何利用这些值计算yn+1？令：)())(,(xFxyxf已知：))(,()),(,()),(,()),(,(112233nnnnnnnnxFxxFxxFxxFx利用插值多项式原理，建立三次Lagrange插值多项式：)())()(())()(()())()(())()(()())()(())()(()())()(())()(()(3231332123212231131211323213213nnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnnxFxxxxxxxxxxxxxFxxxxxxxxxxxxxFxxxxxxxxxxxxxFxxxxxxxxxxxxx余项：))()()((!4)()(321)4(3nnnnxxxxxxxxFxR代入可得：1))(,()()(1nnxxnndxxyxfxyxy1)()()()(331nnxxnndxxRxxyxy易算：)(9)(37)(59)(5524)(32131nnnnxxxFxFxFxFhdxxnn略去余项，得到四阶阿达姆斯外插公式：)9375955(243211nnnnnnffffhyy),(),(720251)(3)5(531nnxxxxyhdxxRnn余项：这是显式公式。公式的局部截断误差为O(h5)。四阶阿达姆斯(Adams)内插（隐式）公式设xn-2,xn-1,xn,xn+1处的函数值分别为yn-2,yn-1,yn,yn+1利用这些值计算yn+1：令：)())(,(xFxyxf已知：))(,()),(,()),(,()),(,(111122nnnnnnnnxFxxFxxFxxFx利用同样的方法可得：)(72019)()(5)(19)(924)()(*)5(52111yhxFxFxFxFhxyxynnnnnn略去余项，得到四阶阿达姆斯外插公式：)5199(242111nnnnnnffffhyy),(),(7201912**)5(5nnxxyh余项：这是隐式公式。公式的局部截断误差为O(h5)。阿达姆斯预测-校正公式阿达姆斯内插公式是隐式公式，与改进的欧拉公式用法类似需要用迭代法计算yn+1。为了减小迭代次数，可与阿达姆斯外插公式联合使用，构造只迭代1次的预测-校正格式，即阿达姆斯预测-校正公式22111113322111,,5,19,924,9,37,59,5524nnnnnnnnnnnnnnnnnnnnyxfyxfyxfyxfhyyyxfyxfyxfyxfhyy公式仍有四阶精度。例用阿达姆斯预测-校正公式求初值问题1.0,10,0)0(2hxyyxy的数值解。解,),(2yxyxf0,000yx使用四阶龙格－库塔公式求出发值，计算得04488.0,3.0,01998.0,2.0,00500.0,1.0332211yxyxyx再用阿达姆斯预测－校正公式，07949.0,,5,19,92407951.0,9,37,59,552411223344340011223334yxfyxfyxfyxfhyyyxfyxfyxfyxfhyy类似可计算.,,......,,,,10106655yyyyyy计算结果如表：0.40.50.60.70.80.91.00.07950.12350.17620.23690.30460.37790.4555xiyi}yi