东北农业大学《数值分析》课件-第九章（隋老师）.pptx

东北农业大学 数值分析 第九章 第九章 常微分方程数值解法 1 1、引言、引言 2yxy2y4xy4x2yf x y4y 13xdy2yf x y4dxxy 13:-(,)()-(,)()例例 方方程程令令：且且给给出出初初值值就就得得到到一一阶阶常常微微分分方方程程的的初初值值问问题题：(,)(,)(,)f x yyLipschitzLf x yf x yL yy只要函数适当光滑连续，且关于 满足条件，即存在常数，使得由常微分方程理论知，初值问题的解必存在且唯一。微分方程的数值解微分方程的数值解：设方程问题的解：设方程问题的解y(x)的存在区间是的存在区间是a,b，令，令a=x0 x1 xn=b,其中其中hk=xk+1-xk,如是等距节点如是等距节点h=(b-a)/n ,h h称为步长。称为步长。y(x)的解析表达式不容易得到或根本无法得的解析表达式不容易得到或根本无法得到，我们用数值方法求得到，我们用数值方法求得y(x)y(x)在每个节点在每个节点xk上上y(xk)的近似值，用的近似值，用yk表示，即表示，即yky(xk)，这这样样y0,y1,.,yn称为微分方程的数值解。称为微分方程的数值解。主要问题 如何将微分方程离散化，并建立求其数值解的递推公式；递推公式的局部截断误差，数值解与精确解的误差估计；递推公式的稳定性与收敛性。用差商代替微商用差商代替微商 数值积分数值积分 Taylor展开展开 微分方程离散化常用方法微分方程离散化常用方法 1 解常微分方程初值问题的解常微分方程初值问题的EulerEuler方法方法 EulerEuler方法方法 EulerEuler方法的误差分析方法的误差分析 向前向前EulerEuler公式公式（EulerEuler折线法或显格式折线法或显格式）向后向后EulerEuler公式公式（后退（后退EulerEuler公式）公式）梯形公式梯形公式（改进的（改进的EulerEuler公式）公式）EulerEuler预估校正格式预估校正格式 一、一、EulerEuler方法方法 1 1、向前、向前EulerEuler公式公式 用分段的折线用分段的折线逼近逼近函数逼近逼近函数 2、向后（后退的）、向后（后退的）Euler 方法方法 3、梯形公式 4、改进的尤拉公式、改进的尤拉公式 梯形公式虽然提高了精度，但使算法复杂。而在实际计算中只迭代一次，这样建立的预测校正系统称作改进的尤拉公式。1111 (,);(,)(,),2nnnnnnnnnnyyhf xyhyyf xyf xy预测校正11(,);(,);()/2.pnnncnnpnpcyyhf xyyyhf xyyyy二、二、Euler方法的误差分析方法的误差分析 2）总体方法误差总体方法误差 总体截断误差与局部截断误差的关系是：误差分析表误差分析表 EulerEuler方法方法 局部截局部截断误差断误差 总体截总体截断误差断误差 迭代收敛迭代收敛条件条件 向前向前EulerEuler方法方法 O(h2)O(h)向后向后EulerEuler方法方法 O(h2)O(h)0hL1 梯形公式梯形公式 O(h3)O(h2)0hL2(L为为Lip常数）常数）向后向后EulerEuler 方法收敛条件与截方法收敛条件与截断误差断误差 梯形公式的收敛性梯形公式的收敛性 2(01);(0)1.xyyxyy例：用尤拉公式和改进的尤拉公式解初值问题10.12().nnnnnhxEuleryyh yy解：取步长，公式为：112();2();1().2npnnnncnppnpcxyyh yyxyyh yyyyy改进的尤拉公式为：20.2(00.6);(0)1.hyyxyxy 例1：取步长，用欧拉法解初值问题122 ,0.2 0.80.2nnnnnnnnnnnEuleryyhfxyyyx yyx y解：格式为：0122 1 0.2 0.8,0.4 0.6144 0.6 0.461321yyyyyyy由计算得0.283(12);(1)2.5hyyxy例2：取步长，用梯形解初值问题小数点后至少保留 位。11111 ,20.2 83832nnnnnnnnnnhyyfxyfxyyyyy解：梯形公式为：1012345716 131312,1.22.307691.42.473371.62.562581.82.610622.02.63649nnyyyyyyyyyyyyyy故由计算得 0.23(12);(1)2.5hyxyxy例2：取步长，用梯形解初值问题小数点后至少保留 位。111111 ,2 0.1 33nnnnnnnnnnnnhyyf xyf xyyyx yxy解：梯形公式为：0111110011234511111110.30.1 332,2.6*,&,%,?nnnnkknnnnnnkyyx yyyx yxyyyyyyyyyy迭代格式：021234522222,*,&,%,yyyyyy20.2sin0 (1)1.1.21.45hyyyxyyy例3：取步长，用欧拉预校方法解初值问题计算及的近似值，小数点后至少保留 位。111121212111,2 0.2sin0.1sin +sinnnnnnnnnnnnnnnnnnnnnnnnyyhfxyhyyf xyf xyyyyyxyyyyxyyx解：欧拉预校格式为：0112211,0.63171 1.20.7154880.476961.40.52611yyyyyyyy故由计算得 50 (0)1.0.10.210yyyhy例4：写出用反复迭代的欧拉预校法解初值问题的计算公式，并取步长，计算，要求迭代误差不超过。01111101111,2 0,1,2,0,1,2,0.90.950.05nnnnkknnnnnnnnkknnnyyhfxyhyyfxyfxyknfx yyyyyyy 解：欧拉预校的迭代格式为：取 0012111341143751141101,0.9,0.905,0.90475,0.9047625,0.9047618756.25 1010,0.10.904761875yyyyyyyyyyyy故由计算得由于是取 0012232452254652252201,0.814286,0.818809,0.818583,0.818595,0.8185941010,0.20.818594yyyyyyyyyyyy故由计算得由于是取2.2.龙格龙格库塔方法库塔方法 基本思想基本思想 二阶二阶R-K方法方法 三阶三阶R-K方法方法 四阶四阶R-K方法方法 变步长变步长R-K方法方法 1112()()()()()()()2!()(,),()(,)(,)(,),nnnppnnnnxypTaylory xyy xhhy xhy xy xyxPy xf x yy xfx yfx y f x y若用 阶多项式近似函数有：其中。但由于公式中各阶偏导数计算复杂，不实用。一、基本思想一、基本思想(0)(0)(0)(1)(1)(1)(2)(2)(2)()(1);2,3,jjjjyffffyffxyffyffxyffyffjxy一般地有111112121(,)11()22 (,)(,)nnnnnnnnnnEulerEulerEuleryyhKKf xyEuleryyhKKKf xyKf xh yhK如果将公式与改进公式写成下列形式：公式改进公式11(,)()(,)(,)nnf x yy xyf x yf x y以上两组公式都使用函数在某些点上的值的线性组合来计算的近似值。Euler公式：每步计算一次的值，为一阶方法。改进Euler公式：需计算两次的值，二阶方法。(,)(,)()-nnnf x yxyTaylory xxTaylorR K于是可考虑用函数在若干点上的函数值的线性组合来构造近似公式，构造是要求近似公式在处的展开式与解在 处的展开式的前面几项重合，从而使近似公式达到所需要的阶数。即避免求偏导，又提高了方法的精度，此为方法的基本思想。11111-(,)(,)(2,3,),pnniiinniininijjjiijiR Kyyhc KKf xyKf xa h yhb Kipab c一般地，方法设近似公式为其中，都是参数.,(,)()iijinnnab cxyTaylory xxTaylor确定参数，的原则是:使近似公式在处的展开式与在 处的展开式的前面项尽可能多地重合。二、二阶龙格库塔方法二、二阶龙格库塔方法 111221222112()(,)(,)nnnnnnpyyh c Kc KKf xyKf xa h yhb K当时，近似公式为 112221123221(,)(,)(,(,)(,)(,)(,)(,)(,)()nnnnnnnnnnnnnnnxnnynnnnnxyTayloryyh c f xyc f xa h yhb f xyyh c f xycf xya hfxyhb fxyf xyO hy上式在处的展开式为12222321()(,)(,)+(,)(,)()nnxnnynnnnccf xy hc a fxyb fxyf xyhO h123123()()()()()()2(,)(,)2 (,)(,)()nnnnnnnnnxnnynnnny xxTaylorhy xy xhy xy xO hhyf xy hfxyfxyf xyO h在 处的展开式为122232 211 1/2 1/2(),ccc ac bO h有无穷多组解，每一组解得一近似公式，局部截断误差均为这些方法统称二阶方法。43()()O hO h可以证明，无论这四个参数如何选择，都不能使局部截断误差达到，也即在计算两次函数值的情况下，局部截断误差的阶最高为。122211121211,1,2()/2(,)(,)nnnnnnccabEuleryyh KKKf xyKf xh yhK取此为改进公式。近似公式为 122211212110,1,2(,)(2,2)nnnnnnccabyyhKKf xyKf xhyhK取此为常用的二阶公式，称为中点公式 三、三阶龙格库塔方法三、三阶龙格库塔方法 1123121312(4)6(,)(,)22(,2)nnnnnnnnRKhyyKKKKf xyhhKf xyKKf xh yhKhK常用的三阶公式为：四、四阶龙格库塔方法四、四阶龙格库塔方法 112341213243(22)6(,)(,)22(,)22(,)nnnnnnnnnnRKhyyKKKKKf xyhhKf xyKhhKf xyKKf xh yhK常用的四阶公式为：0.2,-83;(0)2.0.44 hR Kyyyy例：设取步长写出用经典（标准的）四阶方法求解初值问题 的计算公式，计算的近似值，小数点后至少保留 位。12,8-3,=0.2,0.2,0.4fx yyhyyyy解：112341122343(22);6,83;,5.62.1;22,6.322.37;22,4.2081.578.nnnnnnnnnnnnnnhyyKKKKKf xyyKhKfxyyKhKfxyyKf xh yKy由经典的四阶龙格库塔公式得 10121.20160.549402,0.22.30040.42.4654nnyyyyyyyy由于11,2,3,4,454652pRKppRKpppRKRKTaylor两点说明：）当时，公式的最高阶数恰好是 当时，公式的最高阶数不是，如时仍为，时 公式的最高阶数为。）方法的导出基于展开，故要求所求问题的解具有较高的光滑度。RKEulerRKRKEuler当解充分光滑时，四阶 方法确实优于改进法。对一般实际问题，四阶方法一般可达到精度要求。如果解的光滑性差，则用四阶 方法解的效果不如改进法。五、变步长的龙格五、变步长的龙格库塔库塔方法方法()1()51115(2)1,(),2,2nhnhnnnnhnxhyy xychhxxhyc以经典四阶龙格库塔公式为例。从节点 出发，以 为步长求一近似值将步长折半，即取为步长从 跨两步到，求一近似值每跨一步的截断误差是5(2)11(2)11()11(2)(2)()1111()