《常微分方程的数值解法实用教案》由会员分享，可在线阅读，更多相关《常微分方程的数值解法实用教案（36页珍藏版）》请在金锄头文库上搜索。
1、1第一节引言(ynyn)在常微分方程课程里面讨论的是一些典型方程求解解析解的基本方法。然而在生产实践和科学研究中遇到的微分方程往往比较复杂，在很多情况下，不能给出解的解析表达式；有时候即时能用解析表达式来表示，又因为计算量太大而不实用，有时候一些是已经有了求解的基本方法的典型方程，但实际使用时也是有困难的。以上情况说明用求解解析解的基本方法来求微分方程的解往往是不适宜的，甚至很难办到。实际问题中，对于求解微分方程，一般只要求得到解的若干个点上的近似值或者解的便于(biny)计算的近似表达式。本章研究微分方程的数值解法，而且着重讨论微分方程中最简单的一类问题一阶方程的初值问题。第1页/共35页第
2、一页，共36页。2第一节引言(ynyn)1、一阶方程的初值问题假定上式在区间a,b上存在唯一且足够光滑的解y(x)。所谓数值解法就是寻求解y(x)在一系列离散(lsn)点，也称为节点处的值：要计算出解函数y(x)在一系列节点a=x0x1xn=b处的近似值第2页/共35页第二页，共36页。3第一节引言(ynyn)节点间距，即步长为： 通常采用等距节点，即hi = h (常数)等间距节点在这些节点上采用离散化方法（通常用数值积分、微分、泰勒展开等）将上述初值问题化成关于(guny)离散变量的相应问题。把这个相应问题的解yn作为y(xn)的近似值。这样求得的yn就是上述初值问题在节点xn上的数值解。
3、一般说来，不同的离散化导致不同的方法。第3页/共35页第三页，共36页。4第二节欧拉方法(fngf)一、欧拉法Euler1、向前(xinqin)差商近似导数记为记为第4页/共35页第四页，共36页。5第二节欧拉方法(fngf)2、举例例1用欧拉法求初值问题当h=0.02时在区间0,0.10上的数值解。解：根据欧拉公式可以得到：此外(cwi)，可以得到方程的真解：第5页/共35页第五页，共36页。6第二节欧拉方法(fngf)求解过程(guchng)如下：nxnyny(xn) n = y(xn) - yn001.00001.0000010.020.98200.98250.000520.040.96
4、500.96600.000530.060.94890.95030.001440.080.93360.93540.001850.100.91920.9230.0021第6页/共35页第六页，共36页。7第二节欧拉方法(fngf)3、欧拉方法(fngf)的几何意义根据已知条件：曲线y(x)上的点（x0,y0）及该点处曲线的导数(do sh)f(x0,y0)，则可以得到过该点的直线：该直线与xx1的交点P1，则P1的纵坐标y1为：就用y1作为y(x1)的近似值逐次进行后可以得到一条折线P0P1Pn，该折线看作是初值问题的积分曲线的近似，因此欧拉方法也称为欧拉折线法第7页/共35页第七页，共36页。8
5、第二节欧拉方法(fngf)从上述几何意义上得知，由Euler法所得的折线明显偏离了积分曲线(qxin)，可见此方法非常粗糙即误差太大。4、欧拉法的局部截断误差（1）截断误差定义在假设yi=y(xi)，即第i步计算是精确的前提下，考虑的截断误差Ri+1=y(xi+1)yi+1，称为局部截断误差如图所示：APi+1即为欧拉方法在xi+1点的截断误差（2）如果某种方法的局部截断误差是则称该方法具有p阶精度A第8页/共35页第八页，共36页。9第二节欧拉方法(fngf)（3）则截断误差的大小？写出y(xn+1)的泰勒(til)展开式：由欧拉方法可以得到：则上面两个公式相减得到：具有(jyu)1阶精度第
6、9页/共35页第九页，共36页。10第二节欧拉方法(fngf)二、改进的欧拉法一阶方程(fngchng)的初值问题与如下积分方程(fngchng)是等价的：当x=x1时可以借助于数值积分，求y(x1)的值1、用矩形公式第10页/共35页第十页，共36页。11第二节欧拉方法(fngf)可以推导出：用矩形法计算( jsun)右端的积分与用欧拉法计出的结果完全相同2、用梯形公式则可以推导出：第11页/共35页第十一页，共36页。12第二节欧拉方法(fngf)梯形(txng)公式的截断误差：梯形(txng)公式具有二阶精度，比欧拉方法有了进步 第12页/共35页第十二页，共36页。13第二节欧拉方法(
7、fngf)和欧拉公式相比较(bjio)，梯形公式在计算yi+1时候也只用到前一步的值yi，但是若yi已知，将yi带入公式求解时候，一般不能直接得到yi+1，而需要通过其他方法（比如迭代法）求解，所以梯形公式被称为隐式公式。3、改进的欧拉方法梯形公式是隐式的，一般用迭代法求解，计算量较大。实际中常将欧拉公式和梯形公式联合使用，先用欧拉公式得出一个y(xi+1)的近似值称为预估值，然后对预估值使用梯形公式对它进行精确化，得到较为精确的近似值yi+1，称之为校正值，计算公式为：这样的预估校正系统称为改进的欧拉方法。第13页/共35页第十三页，共36页。14第二节欧拉方法(fngf)为了(wile)便
8、于编写程序，常将上面的公式改写为如下式：第14页/共35页第十四页，共36页。15第二节欧拉方法(fngf)4、举例P90，例题5-1在区间0,1.5上，取h=0.1，求解。解：（1）用欧拉法计算公式如下(rxi)：（2）用改进欧拉法计算公式如下(rxi)：第15页/共35页第十五页，共36页。16本题的精确解为，可用来检验近似( jns)解的精确程度。计算结果如表：xn欧拉法欧拉法yn迭代一次迭代一次改进欧拉法改进欧拉法yn准确解准确解01110.11.11.0959091.0954450.21.1918181.1840961.1832160.31.2774381.2602011.26491
9、10.41.3582131.3433601.3416410.51.4351331.4161021.4142140.61.5089661.4829561.4832400.71.5803381.5525151.5491930.81.6497831.6164761.6124520.91.7177791.6781681.6733201.01.7847701.7378691.7320511.11.851181.7958221.7888541.21.9174641.8522421.8439091.31.9840461.9073231.8973671.42.0514041.9612531.9493591.5
10、2.1200522.0142072.000000第16页/共35页第十六页，共36页。17第二节欧拉方法(fngf)P108习题(xt)5-1第17页/共35页第十七页，共36页。18第二节欧拉方法(fngf)5、欧拉两步公式(gngsh)中心差商第18页/共35页第十八页，共36页。19第三节龙格库塔方法(fngf)一、龙格库塔法的基本思想1、平均斜率考察差商：根据微分中值定理：在闭区间a,b上连续，开区间(a,b)上可导，则至少存在(a,b)上的一点(ydin)，使得下式成立：根据上面公式可以得到：称为(chn wi)区间xi,xi+1上的平均斜率第19页/共35页第十九页，共36页。20
11、第三节龙格库塔方法(fngf)因此只要对K*提供一种算法，就可以求得数值解，根据该观点对欧拉法及改进的欧拉进行分析。2、基于平均斜率对欧拉法和改进的欧拉法进行分析（1）在欧拉公式(gngsh)中，是取了一个点xi上的斜率值f(xi,yi)作为平均斜率K*的近似值的，已经知道其精度较低。（2）对于改进的欧拉公式(gngsh)：可以看出(kn ch)它用两个点xi 和 xi+1上的斜率K1和K2 的算术平均值作为平均斜率 K* 的近似值的，而xi+1处的斜率是由已知信息预测得到的。第20页/共35页第二十页，共36页。21第三节龙格库塔方法(fngf)根据上面(shngmin)的分析得到：如果在区
12、间xi,xi+1上多测几个点的斜率值，然后取其加权平均作为平均斜率的近似值，有可能构造出具有更高精度的计算公式，此即为龙格库塔算法的基本思想。二、二阶龙格库塔方法1、推广改进的欧拉方法，考察区间xi,xi+1上的一点：用xi和xp两个点上的斜率值K1和K2的加权平均作为平均斜率K*的近似值：即取：与改进的欧拉法类似，有：如何得到xp的斜率？第21页/共35页第二十一页，共36页。22第三节龙格库塔方法(fngf)如何得到xp的斜率？根据改进的欧拉法，可以利用欧拉法预测(yc)的值：则可以得到点xi+p斜率K2：则可以得到算法的具体表达式：第22页/共35页第二十二页，共36页。23第三节龙格库
13、塔方法(fngf)2、选择参数使得算法具有( jyu)2阶精度计算上面公式的局部截断误差：根据泰勒公式有：另外：第23页/共35页第二十三页，共36页。24第24页/共35页第二十四页，共36页。25第三节龙格库塔方法(fngf)考虑(kol)到：则有：二元泰勒(ti l)展开两式相减可得到：第25页/共35页第二十五页，共36页。26第三节龙格库塔方法(fngf)两式相减可得到：则要使得上式满足二阶精度，即Ri+1=O(h3)，只需要：共有3个参数，但只需要满足2个条件，因此满足该式的参数不止一组，而是一簇，所有(suyu)满足条件的公式通称为二阶龙格库塔公式。第26页/共35页第二十六页，
14、共36页。27第三节龙格库塔方法(fngf)该公式的特别情况：有：此时(csh)，二阶龙格库塔公式就是改进的欧拉公式。第27页/共35页第二十七页，共36页。28第三节龙格库塔方法(fngf)如果p=1/2，则有：此时(csh)二阶龙格库塔公式变成：为变形的欧拉公式。第28页/共35页第二十八页，共36页。29第三节龙格库塔方法(fngf)3、高阶龙格库塔公式为了进一步提高精度，可以考虑在区间xi,xi+1上，除了xi，xp外，再增加几个( j)点。如再增加一个点：xi+m利用该三点处的斜率的加权平均值作为K*近似值：此时计算公式为：式中K1和K2的取法与前相同，如何获得K3？第29页/共35
15、页第二十九页，共36页。30第三节龙格库塔方法(fngf)为了得到xi+m处的斜率K3，需要(xyo)确定该点的函数值y(xi+m)，可以用二阶龙格库塔公式估计y(xi+m)的近似值yi+m：然后通过计算函数值，就可以得到K3：由此得到的3阶龙格库塔公式为：第30页/共35页第三十页，共36页。31第三节龙格库塔方法(fngf)通过选择未知参数使得公式具有3阶精度，采用前面相同的方法，可以得到(ddo)参数需满足如下条件：第31页/共35页第三十一页，共36页。32第三节龙格库塔方法(fngf)常用的三阶(snji)龙格库塔公式为：其中(qzhng)第32页/共35页第三十二页，共36页。33
16、第三节龙格库塔方法(fngf)经典( jngdin)的四阶龙格库塔公式为：第33页/共35页第三十三页，共36页。34第三节龙格库塔方法(fngf)第34页/共35页第三十四页，共36页。35谢谢大家(dji)观赏！第35页/共35页第三十五页，共36页。内容(nirng)总结1。第1页/共35页。假定上式在区间a, b上存在唯一且足够 光滑的解y(x)。第2页/共35页。在这些节点上采用离散化方法（通常用数值积分、微分、泰勒展开等）将上述初值问题化成关于离散变量的相应问题。从上述几何意义上得知(d zh)，由Euler法所得的折线明显偏离了积分曲线，可见此方法非常粗糙即误差太大。用矩形法计算右端的积分与用欧拉法计出的结果完全相同。称为区间xi,xi+1上的平均斜率。第34页/共35页第三十六页，共36页。
