数值分析第九章.pptx

《数值分析第九章.pptx》由会员分享，可在线阅读，更多相关《数值分析第九章.pptx（69页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、第九章 常微分方程的数值解一、Euler方法三、单步法的收敛性和稳定性二、Runge-Kutta方法四、线性多步法第1页/共69页很多科学技术和工程问题常用常微分方程的形式建立数学模型.但是对于绝大多数的微分方程问题，很难或者根本不可能得到它的解析解.本章重点考察一阶方程的初值问题的数值解法，就是寻求解y(x)在一系列离散点处的近似值 的方法.相邻两个节点间的距离 称为步长.第2页/共69页一、Euler方法1 欧拉公式由初值条件 表示积分曲线从出发，并在 处的切线斜率为因此可以设想积分曲线在 x=x0 附近可以用切线近似的代替曲线.切线方程为当x=x1时，代入有这样得到y(x1)的近似值y1

2、的方法.第3页/共69页重复上述方法，当 x=x2 时依次可以计算出x3,x4,处的近似值y3,y4,由此得到Euler公式：由于用折线近似代替方程的解析解，所以Euler方法也称为Euler折线法.例 用Euler法计算初值问题的解在x=0.3时的近似值，取步长h=0.1.第4页/共69页解：Euler公式的截断误差局部截断误差：一步Euler公式产生的误差;总体截断误差：Euler公式的累积总误差;第5页/共69页在假设yn=y(xn)，即第i步计算是精确的前提下，考虑的截断误差Rn=y(xn+1)yn+1 称为局部截断误差.定义欧拉法的局部截断误差：所以欧拉法具有 1 阶精度.若某算法的

3、局部截断误差为O(hp+1)，则称该算法有p 阶精度.定义第6页/共69页Lipschitiz条件：若存在正数L,使得对一切x,y1,y2有则称f(x,y)满足Lipschitiz条件.欧拉法的总体截断误差：那么设为局部截断误差，所以第7页/共69页第8页/共69页特别当n=m-1时，有总体误差与h是同阶的.上式还说明，当 时，有 即也就是说，ym收敛到方程的准确解第9页/共69页后退Euler公式（隐式欧拉法）（隐式欧拉公式）利用向后差商近似导数第10页/共69页 由于未知数yn+1同时出现在等式的两边，不能直接得到，故称为隐式欧拉公式，而前者称为显式欧拉公式.一般先用显式计算一个初值，再迭

4、代求解.隐式欧拉法的局部截断误差：即隐式欧拉公式具有1阶精度.第11页/共69页2 梯形公式和改进Euler方法梯形公式设y=y(x)是 的解,故由此得到用梯形公式近似第12页/共69页用yn来近似y(xn),用yn+1来近似y(xn+1),得梯形公式梯形公式是隐式的，可以用迭代法求解.第13页/共69页具有2阶精度.梯形公式的局部截断误差第14页/共69页中点欧拉公式中心差商近似导数假设 ,则可以导出即中点公式具有 2 阶精度.需要2个初值 y0和 y1来启动递推过程，这样的算法称为双步法，而前面的三种算法都是单步法.第15页/共69页方 法显式欧拉隐式欧拉梯形公式中点公式简单精度低稳定性最

5、好精度低,计算量大精度提高计算量大精度提高,显式多一个初值,可能影响精度有没有一种方法，既有这些方法的优点，而没有它们的缺点？第16页/共69页改进欧拉法(1)先用显式欧拉公式作预测，算出(2)再将 代入梯形公式的右边作校正,得到注：此法亦称为预测-校正法.可以证明该算法具有2阶精度，同时可以看到它是个单步递推格式，比隐式公式的迭代求解过程简单.第17页/共69页例 用梯形公式求解初值问题(步长h=0.2)解：梯形公式为于是整理得由 y(1)=y0=2 依次可得 y1,y2,y3,y4,y5.第18页/共69页例 用改进欧拉法求解初值问题要求步长h=0.2，并计算y(1.2)和y(1.4)解：

6、改进欧拉法公式为即第19页/共69页由 y(1)=y0=1 计算得第20页/共69页二、Runge-Kutta方法建立高精度的单步递推格式.单步递推法的基本思想是从(xn,yn)点出发，以某一斜率沿直线达到(xn+1,yn+1)点.欧拉法及其各种变形所能达到的最高精度为2阶.考察改进的欧拉法，可以将其改写为：斜率一定取K1 K2 的平均值吗？步长一定是一个h 吗？第21页/共69页首先希望能确定系数 1、2、p，使得到的算法格式有2阶精度，即在 的前提假设下，使得 将改进欧拉法推广为：(1)将K2在(xn,yn)点作 Taylor 展开1 二阶Runge-Kutta方法第22页/共69页(2)

7、将 K2 代入第1式，得到第23页/共69页(3)将yn+1与y(xn+1)在xn点的泰勒展开作比较要求 ,则必须有：这里有 个未知数，个方程。32所以存在无穷多个解！第24页/共69页所有满足上式的统称为2阶Runge-Kutta格式.若 则改进的欧拉方法若 则中点公式第25页/共69页2 四阶Runge-Kutta方法其中i(i=1,m)，i(i=2,m)和ij(i=2,m;j=1,i1)均为待定系数，确定这些系数的步骤与前面相似.第26页/共69页由于方程的个数少于未知量的个数，所以方程有无穷多个解，可以根据情况得到几种常用的解，即得到相应的四阶公式.最常用为四阶经典龙格-库塔法也称为标

8、准四阶龙格-库塔公式第27页/共69页Gill公式第28页/共69页753可达到的最高精度642每步须算Ki 的个数(2)龙格-库塔法的导出基于泰勒展开，故精度主要受解函数的光滑性影响.对于光滑性不太好的解,最好采用低阶算法而将步长h 取小.注：(1)龙格-库塔法的主要运算在于计算Ki 的值,即计算f的值.Butcher于1965年给出了计算量与可达到的最高精度阶数的关系：第29页/共69页例 用标准四阶Runge-Kutta法求初值问题在x=0.1处的近似值，取步长为h=0.1.解：所以第30页/共69页那么例 用标准四阶Runge-Kutta法求初值问题在x=0.4处的近似值，取步长为h=

9、0.2.第31页/共69页解：所以而所以第32页/共69页 若某算法对于任意固定的x=xi=x0+i h,当 h0(同时i )时有yi y(xi)，则称该算法是收敛的.定义1 单步法的收敛性三、单步法的收敛性和稳定性单步法是在计算yn+1时只用到前一步的信息yn.显式单步法的共同特征是它们都是将yn加上某种形式的增量，得出yn+1,计算公式如下：增量函数第33页/共69页Euler方法的增量函数改进Euler方法的增量函数 设y(x)是微分方程初值问题的准确解，定义则称为显式单步法在xn+1处的局部截断误差.第34页/共69页例：考察欧拉显式格式的收敛性:解：该问题的精确解为 欧拉公式为对任意

10、固定的 x=xi=ih，有第35页/共69页 设y(x)是微分方程初值问题的准确解，定义若存在最大整数p，使显式单步法的局部截断误差满足则称该方法具有p阶精度，或称为p方法.Tn+1按h展开的第一项，又称为主项.若局部截断误差的展开式写成则称 为局部截断误差的主项第36页/共69页单步法的收敛定理设单步法 具有p阶精度其增量函数 关于y满足Lipschitz条件即存在常数L，使对任何的 及任意的x有又设初值y0是准确的，即 则总体截断误差是p阶的，也就是特别的当 时,不论n为何值,总有即方法收敛.第37页/共69页在f(x,y)对y满足Lipschitz条件下,Euler法，改进Euler法和

11、Runge-Kutta法的增量函数 都对y满足Lipschitz条件，所以上述结论对这些方法都成立.例 设是求解微分方程的单步法，试求其局部截断误差的主项，并说出它具有几阶精度.解:第38页/共69页考虑 在xn处的Taylor展式所以该方法的局部截断误差的主项是具有一阶精度.第39页/共69页例 设试求出它具有几阶精度.解:考虑 在xn处的Taylor展式第40页/共69页所以该方法的局部截断误差的主项是具有二阶精度.2 单步法的稳定性收敛性是在假定每一步计算都准确的前提下,讨论步长 时，方法的总体截断误差是否趋于零的问题.稳定性是讨论舍入误差的积累能否对计算结果有严重的影响.第41页/共6

12、9页例：考察初值问题 在区间0.00.10.20.30.40.5精确解改进欧拉法 欧拉隐式欧拉显式 节点 xi 1.00002.0000 4.00008.0000 1.6000101 3.2000101 1.00002.5000101 6.25001021.56251023.90631039.76561041.00002.50006.25001.56261013.90631019.76561011.00004.97871022.47881031.23411046.14421063.05901070,0.5上的解，分别用欧拉显、隐式格式和改进的欧拉格式计算数值解.第42页/共69页 若一种数值解

13、法仅在节点值yn上有大小为的扰动，于以后各节点值ym(mn)上，仅由所引起的扰动都不超过时，称该方法是稳定的.定义一般分析时为简单起见，只考虑试验方程 为复数且Re()0设 在节点值yn处有扰动 令 那么于是第43页/共69页反复应用可得为使则可得如下定义 若 中 ，则称单步法是绝对稳定的.在复平面上，h满足 的区域称为方法的绝对稳定区域,它与实轴的交称为绝对稳定区间定义下面讨论已知的几种方法的绝对稳定区间和绝对稳定区域.第44页/共69页显式欧拉法:0-1-2ReImg在复平面上的绝对稳定区域是即以-1为中心，1为半径的圆域所以相应的绝对稳定区间是第45页/共69页隐式欧拉法(后退欧拉法)：

14、210ReImg在复平面上的绝对稳定区域是是以1为中心，1为半径的圆的外域所以相应的绝对稳定区间是即如果只考虑 0为实数)所以绝对稳定区域是所以因此是条件稳定的.第50页/共69页四、线性多步法在逐步推进的求解过程中，计算yn+1之前已经求出了一系列的近似值y0,y1,yn,如果充分利用前面信息来预测yn+1，则可期望会获得较高的精度，这就是线性多步法的基本思想.1 线性多步法的一般公式最常用的线性多步法公式为其中 为常数，yn-k为y(xn-k)的近似值fn-k=f(xn-k,yn-k)第51页/共69页特别的当 时,上式为显式,否则是隐式.设y(x)是微分方程初值问题的准确解，定义称为线性

15、多步法在xn+1上的局部截断误差.若 则称该方法具有p阶精度.若 则称局部截断误差的主项为 为误差常数.第52页/共69页例 设 yn+1=yn-1+2hf(xn,yn)为求解常微分初值问题的线性二步法，试求该二步公式的局部截断误差主项，和精度.解:由局部截断误差的定义可知考虑 在xn处的Taylor展式第53页/共69页代入可得所以局部截断误差的主项为具有二阶精度.例 试建立求解为微分方程初值问题具有如下形式的线性二步法，并使该方法具有二阶精度，同时求其局部截断误差的主项.解:局部截断误差为第54页/共69页考虑 在xn处的Taylor展式于是为使方法具有二阶精度则第55页/共69页解得因此

16、该方法为局部截断误差的主项为例 试建立求解为微分方程初值问题具有如下形式的线性二步法，并使该方法具有三阶精度，同时求其局部截断误差的主项.解:局部截断误差为第56页/共69页考虑 在xn处的Taylor展式第57页/共69页所以解得局部截断误差的主项为第58页/共69页2 Adams外推公式考虑用r+1个点(xn-k,f(xn-k,yn-k)构造一个r次多项式来近似的被积函数f(x,y(x),这里用yn-k作为y(xn-k)的近似值，令fn-k=f(xn-k,yn-k),用Newton向后插值公式来构造r次多项式,即第59页/共69页这里第60页/共69页而将yn+1作为y(xn+1)的近似值

17、，因此得到Adams外推公式显式格式当r=0时,为Euler公式,最常用的是r=3情况第61页/共69页i0123Adams外推系数第62页/共69页3 Adams内插公式用xn+1,xn,xn-r+1为插值节点,构造f(x,y(x)的r次多项式，得到内插公式隐式格式最常用的是r=3情况第63页/共69页i0123Adams外推系数第64页/共69页4 预报-校正公式不论单步法或多步法，隐式公式比显式公式稳定性好，但在实际使用隐式公式时，都会遇到两个问题：一个是隐式公式如何能方便地进行计算；另一个是实际计算步长取多大.如隐式梯形公式，每往前推进一步，不必进行多次迭代，而是采用一阶显式Euler

18、公式预测，二阶隐式梯形公式校正一次，构成显式改进Euler公式，能达到与梯形公式同阶的精度，即二阶精度.第65页/共69页预测-校正技术即保证了计算精度，又使隐式计算显式化，克服了隐式公式要反复迭代的困难.类似改进Euler方法，将Adams外推公式和内插公式结合起来，得到实用的预报-校正公式预报校正该方法的精度是四阶的第66页/共69页例 用Adams外推公式和预报-校正公式求初值问题在区间0,1上的数值解，取步长h=0.1解：为运用Adams方法进行计算，先用四阶Runge-Kutta法计算然后用Adams外推公式和内插公式分别计算结果如下表第67页/共69页kxkynR-KA-预报预报A-校正校正00.01.0000010.11.1103420.21.2428130.31.3997240.41.583641.5836550.51.797431.7974460.62.044232.0442470.72.327492.3275180.82.651072.6510890.93.010193.01921191.03.436553.43657第68页/共69页感谢您的观看。第69页/共69页