欧拉公式的改进.ppt

《欧拉公式的改进.ppt》由会员分享，可在线阅读，更多相关《欧拉公式的改进.ppt（22页珍藏版）》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。

1、1EulersMethod 欧拉公式的改进：欧拉公式的改进：隐式欧拉法隐式欧拉法/*implicit Euler method*/向后差商近似导数向后差商近似导数x0 x1)(,()(1101xyxfhyxy+)1,.,0(),(111=+=+niyxfhyyiiii由于未知数由于未知数yi+1同时出现在等式的两边，不能直接得到，故同时出现在等式的两边，不能直接得到，故称为称为隐式隐式/*implicit*/欧拉公式，而前者称为欧拉公式，而前者称为显式显式/*explicit*/欧拉公式。欧拉公式。一般先用显式计算一个初值，再一般先用显式计算一个初值，再迭代迭代求解。求解。隐式隐式欧拉法的局部

2、截断误差：欧拉法的局部截断误差：即隐式欧拉公式具有即隐式欧拉公式具有1阶精度。阶精度。Hey!Isnttheleading termofthelocal truncation errorofEulersmethod?Seemsthatwecanmakeagooduseofit1EulersMethod梯形公式梯形公式/*trapezoid formula*/显、隐式两种算法的显、隐式两种算法的平均平均注：注：的确有局部截断误差的确有局部截断误差，即梯形公式具有即梯形公式具有2阶精度，比欧拉方法有了进步。阶精度，比欧拉方法有了进步。但注意到该公式是但注意到该公式是隐式隐式公式，计算时不得不用到公

3、式，计算时不得不用到迭代法，其迭代收敛性与欧拉公式相似。迭代法，其迭代收敛性与欧拉公式相似。中点欧拉公式中点欧拉公式/*midpoint formula*/中心差商近似导数中心差商近似导数x0 x2x1假设假设，则可以导出，则可以导出即中点公式具有即中点公式具有2阶精度。阶精度。需要需要2个初值个初值y0和和y1来启动递推来启动递推过程，这样的算法称为过程，这样的算法称为双步法双步法/*double-step method*/，而前面的三种算法都是，而前面的三种算法都是单步法单步法/*single-step method*/。方方法法 1EulersMethod显式欧拉显式欧拉隐式欧拉隐式欧拉

4、梯形公式梯形公式中点公式中点公式简单简单精度低精度低稳定性最好稳定性最好精度低精度低,计算量大计算量大精度提高精度提高计算量大计算量大精度提高精度提高,显式显式多一个初值多一个初值,可能影响精度可能影响精度Cantyougivemeaformulawithalltheadvantagesyetwithoutanyofthedisadvantages?Doyouthinkitpossible?Well,callmegreedyOK,letsmakeitpossible.改进欧拉法改进欧拉法/*modified Eulers method*/Step 1:先用先用显式显式欧拉公式作欧拉公式作预测预

5、测，算出，算出),(1iiiiyxfhyy+=+Step 2:再将再将代入代入隐式隐式梯形公式的右边作梯形公式的右边作校正校正，得到，得到1+iy),(),(2111+=iiiiiiyxfyxfhyy注：注：此法亦称为此法亦称为预测预测-校正法校正法/*predictor-corrector method*/。可以证明该算法具有可以证明该算法具有2阶精度，同时可以看到它是个阶精度，同时可以看到它是个单单步步递推格式，比隐式公式的迭代求解过程递推格式，比隐式公式的迭代求解过程简单简单。后面将。后面将看到，它的看到，它的稳定性高稳定性高于显式欧拉法。于显式欧拉法。1EulersMethod2龙格龙

6、格-库塔法库塔法/*Runge-Kutta Method*/建立高精度的单步递推格式。建立高精度的单步递推格式。单步递推法的单步递推法的基本思想基本思想是从是从(xi,yi)点出发，以点出发，以某一斜某一斜率率沿直线达到沿直线达到(xi+1,yi+1)点。欧拉法及其各种变形所点。欧拉法及其各种变形所能达到的最高精度为能达到的最高精度为2阶阶。考察改进的欧拉法，可以将其改写为：考察改进的欧拉法，可以将其改写为：斜率斜率一定取一定取K1 K2的的平均值平均值吗？吗？步长一定是一个步长一定是一个h吗吗？2Runge-KuttaMethod首先希望能确定系数首先希望能确定系数 1、2、p，使得到的算法

7、格式有，使得到的算法格式有2阶阶精度，即在精度，即在的前提假设下，使得的前提假设下，使得 Step 1:将将K2在在(xi,yi)点作点作Taylor展开展开将改进欧拉法推广为：将改进欧拉法推广为：),(),(12122111phKyphxfKyxfKKKhyyiiiiii+=+=+Step 2:将将K2代入第代入第1式，得到式，得到2Runge-KuttaMethodStep 3:将将yi+1与与y(xi+1)在在xi 点的点的泰勒泰勒展开作比较展开作比较要求要求，则必须有：，则必须有：这里有这里有个未知个未知数，数，个方程。个方程。32存在存在无穷多个解无穷多个解。所有满足上式的格式统称为

8、。所有满足上式的格式统称为2阶龙格阶龙格-库库塔格式塔格式。注意到，注意到，就是改进的欧拉法。就是改进的欧拉法。Q:为获得更高的精度，应该如何进一步推广？为获得更高的精度，应该如何进一步推广？其中其中 i (i=1,m)，i (i=2,m)和和 ij(i=2,m;j=1,i 1)均为待均为待定系数，确定这些系数定系数，确定这些系数的步骤与前面相似。的步骤与前面相似。2Runge-KuttaMethod).,(.),(),(),(.1122112321313312122122111 +=+=+=+=mm mmmmimiiiiiimmiihKhKhKyhxfKhKhKyhxfKhKyhxfKyxf

9、KKKKhyy 最常用为四级最常用为四级4阶阶经典龙格经典龙格-库塔法库塔法/*Classical Runge-Kutta Method*/：2Runge-KuttaMethod注：注：龙格龙格-库塔法库塔法的主要运算在于计算的主要运算在于计算Ki的值，即计算的值，即计算f的的值。值。Butcher于于1965年给出了计算量与可达到的最高精年给出了计算量与可达到的最高精度阶数的关系：度阶数的关系：753可达到的最高精度可达到的最高精度642每步须算每步须算Ki 的个数的个数由于龙格由于龙格-库塔法的导出基于泰勒展开，故精度主要受库塔法的导出基于泰勒展开，故精度主要受解函数的光滑性影响。对于光滑

10、性不太好的解，最好解函数的光滑性影响。对于光滑性不太好的解，最好采用采用低阶算法低阶算法而将步长而将步长h取小取小。HW:p.202#1,22Runge-KuttaMethodLab15.Runge-Kutta(OrderFour)UseRunge-Kuttamethodoforderfourtoapproximatethesolutionoftheinitial-valueproblem,and(1)YouaresupposedtowriteafunctionvoidRK4(double(*f)(),doublea,doubleb,doubley0,intn,FILE*outfile)toa

11、pproximatethesolutionofProblem(1)withy=f,xina,b,andtheinitialvalueofybeingy0.Outputtheapproximatingvaluesofyonthen+1equallyspacedgridpointsfromatobtooutfile.InputThereisnoinputfile.Instead,youmusthandinyourfunctionina*.hfile.Theruleofnamingthe*.hfileisthesameasthatofnamingthe*.cor*.cppfiles.Output(r

12、epresentsaspace)Foreachtestcase,youaresupposedtoprintn+1lines,andeachlineisinthefollowingformat:fprintf(outfile,“%8.4f%12.8en”,x,y);2Runge-KuttaMethodSample Judge ProgramSample Judge Program#include#include#include98115001_15.hdoublef0(doublex,doubley)return(y x*x+1.0);voidmain()FILE*outfile=fopen(o

13、ut.txt,w);intn;doublea,b,y0;a=0.0;b=2.0;y0=0.50;n=10;RK4(f0,a,b,y0,n,outfile);fprintf(outfile,n);fclose(outfile);Sample Output Sample Output(represents a space)represents a space)0.00005.00000000e 0010.20008.29293333e 0010.40001.21407621e+0000.60001.64892202e+0000.80002.12720268e+0001.00002.64082269

14、e+0001.20003.17989417e+0001.40003.73234007e+0001.60004.28340950e+0001.80004.81508569e+0002.00005.30536300e+0003收敛性与稳定性收敛性与稳定性 /*Convergency and Stability*/收敛性收敛性/*Convergency*/若若某某算算法法对对于于任任意意固固定定的的x=xi=x0+i h，当当h0(同时同时i )时有时有yiy(xi)，则称该算法是，则称该算法是收敛收敛的。的。例：例：就初值问题就初值问题考察欧拉显式格式的收敛性。考察欧拉显式格式的收敛性。解：解：该

15、问题的精确解为该问题的精确解为欧拉公式为欧拉公式为对任意固定的对任意固定的x=xi=i h，有，有 3ConvergencyandStability 稳定性稳定性/*Stability*/例：例：考察初值问题考察初值问题在区间在区间0,0.5上的解。上的解。分别用欧拉显、隐式格式和改进的欧拉格式计算数值解。分别用欧拉显、隐式格式和改进的欧拉格式计算数值解。0.00.10.20.30.40.5精确解精确解改进欧拉法改进欧拉法 欧拉隐式欧拉隐式欧拉显式欧拉显式 节点节点xi1.0000 2.00004.0000 8.00001.6000 101 3.2000 1011.00002.5000 10

16、16.2500 10 21.5625 10 23.9063 10 39.7656 10 41.00002.50006.25001.5626 1013.9063 1019.7656 1011.00004.9787 10 22.4788 10 31.2341 10 46.1442 10 63.0590 10 7Whatiswrong?!AnEngineercomplains:Maththeoremsaresounstablethatasmallperturbationontheconditionswillcauseacrashontheconclusions!3ConvergencyandStab

17、ility若若某某算算法法在在计计算算过过程程中中任任一一步步产产生生的的误误差差在在以以后后的的计计算算中中都都逐逐步步衰衰减减，则则称称该该算算法法是是绝绝对对稳稳定定的的/*absolutely stable*/。一般分析时为简单起见，只考虑一般分析时为简单起见，只考虑试验方程试验方程/*test equation*/常数，可以常数，可以是复数是复数当步长取为当步长取为h时，将某算法应用于上式，并假设只在初值时，将某算法应用于上式，并假设只在初值产生误差产生误差，则若此误差以后逐步衰减，就称该算法，则若此误差以后逐步衰减，就称该算法相对于相对于绝对稳定绝对稳定，的全体构成的全体构成绝对稳

18、定区域绝对稳定区域。我们。我们称称算法算法A比算法比算法B稳定稳定，就是指，就是指A的绝对稳定区域比的绝对稳定区域比B的的大大。h h=h3ConvergencyandStability例：例：考察显式欧拉法考察显式欧拉法由此可见，要保证初始误差由此可见，要保证初始误差 0以后逐步衰减，以后逐步衰减，必须满足：必须满足：0-1-2ReImg例：例：考察隐式欧拉法考察隐式欧拉法可见绝对稳定区域为：可见绝对稳定区域为：210ReImg注：注：一般来说，隐式欧拉法的绝对稳定性比同阶的显式法一般来说，隐式欧拉法的绝对稳定性比同阶的显式法的好。的好。3ConvergencyandStability例：例

19、：隐式龙格隐式龙格-库塔法库塔法而而显式显式14阶方法的绝对稳定阶方法的绝对稳定区域为区域为其中其中2阶方法阶方法的绝对稳定区域为的绝对稳定区域为0ReImgk=1k=2k=3k=4-1-2-3-123ReImg无条件稳定无条件稳定HW:p.202#64线性多步法线性多步法 /*Multistep Method*/用用若干若干节点处的节点处的y及及y 值的值的线性组合线性组合来近似来近似y(xi+1)。).(.110111101kikiiikikiiiffffhyyyy +=其通式可写为：其通式可写为：当当 1 0时，为时，为隐式公式隐式公式;1=0则为则为显式公式显式公式。基于数值积分的构造

20、法基于数值积分的构造法将将在在上积分，得到上积分，得到只要只要近似地算出右边的积分近似地算出右边的积分，则可通过，则可通过近似近似y(xi+1)。而。而选用不同近似式选用不同近似式Ik，可得到不同的计算公式，可得到不同的计算公式。4MultistepMethod亚当姆斯显式公式亚当姆斯显式公式/*Adams explicit formulae*/利用利用k+1个节点上的被积函数值个节点上的被积函数值构造构造k阶牛顿阶牛顿后后插插多项式多项式,有有Newton插值余项插值余项/*显式计算公式显式计算公式*/局部截断误差为：局部截断误差为：例：例：k=1时有时有4MultistepMethod注：

21、注：一般有一般有，其中，其中Bk 与与yi+1计算公式中计算公式中fi,fi k各项的各项的系数系数均可查表得到均可查表得到。10123kfifi 1fi 2fi 3BkMisprintonp.204常用的是常用的是k=3的的4阶亚当姆斯显式公式阶亚当姆斯显式公式4MultistepMethod亚当姆斯隐式公式亚当姆斯隐式公式/*Adams implicit formulae*/利用利用k+1个节点上的被积函数值个节点上的被积函数值fi+1,fi,fi k+1 构造构造k阶阶牛顿牛顿前插前插多项式。与显式多项式完全类似地可得到一系列多项式。与显式多项式完全类似地可得到一系列隐式公式隐式公式，并

22、有，并有，其中，其中与与fi+1,fi,fi k+1的系数亦可查表得到。的系数亦可查表得到。10123kfi+1fifi 1fi 2Bk常用的是常用的是k=3的的4阶亚当姆斯隐式公式阶亚当姆斯隐式公式小于小于Bk较同阶显较同阶显式式稳定稳定4MultistepMethod亚当姆斯预测亚当姆斯预测-校正系统校正系统 /*Adams predictor-corrector system*/Step 1:用用Runge-Kutta法法计算前计算前k 个初值；个初值；Step 2:用用Adams 显式显式计算计算预测预测值；值；Step 3:用同阶用同阶Adams 隐式隐式计算计算校正校正值。值。注意

23、：注意：三步所用公三步所用公式的精度必须相同。式的精度必须相同。通常用通常用经典经典Runge-Kutta 法法配合配合4阶阶Adams 公式公式。Hey!LookatthelocaltruncationerroroftheexplicitandimplicitAdams methods:andDontyouthinktheressomethingyoucando?4阶阶Adams隐式公式的截断误差为隐式公式的截断误差为4阶阶Adams显式公式的截断误差为显式公式的截断误差为当当h充分小充分小时，可近似认为时，可近似认为 i i，则：，则：Predictedvaluepi+1Modifiedv

24、aluemi+1Correctedvalueci+1Modifiedfinalvalueyi+1外推技术外推技术/*extrapolation*/4MultistepMethodAdams 4th-Order predictor-corrector AlgorithmToapproximatethethesolutionoftheinitial-valueproblemAt(N+1)equallyspacednumbersintheintervala,b.Input:endpointsa,b;integerN;initialvaluey0.Output:approximationyatthe(

25、N+1)valuesofx.Step 1Seth=(b a)/N;x0=a;y0=y0;Output(x0,y0);Step 2Fori=1,2,3ComputeyiusingclassicalRunge-Kuttamethod;Output(xi,yi);Step 3Fori=4,Ndosteps4-10Step 5 ;/*predict*/Step 6 ;/*modify*/Step 7;/*correct*/Step 8 ;/*modifythefinalvalue*/Step 9Output(xi+1,yi+1);Step 10Forj=0,1,2,3Setxi=xi+1;yi=yi+1;/*Preparefornextiteration*/Step 11STOP.应为应为(ci+1 pi+1),但因但因ci+1尚未尚未算出，只好用算出，只好用(ci pi)取代之。取代之。