《计算方法-常微分方程课件.ppt》由会员分享，可在线阅读，更多相关《计算方法-常微分方程课件.ppt（18页珍藏版）》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、第八章第八章 常微分方程数值解常微分方程数值解/*Numerical Methods for Ordinary Differential Equations*/y 为人口数量，为人口数量，r为人口增长率，为人口增长率，y0为时刻为时刻x0的人口数量。的人口数量。求常微分方程求常微分方程解析解解析解的方法有多种多样，但是利用这些方法，的方法有多种多样，但是利用这些方法，我们只能对极少数特殊类型的常微分方程求解；科学研究和我们只能对极少数特殊类型的常微分方程求解；科学研究和工程技术上的大量常微分方程的求解需借助于工程技术上的大量常微分方程的求解需借助于数值计算方法数值计算方法。比如，描述人口增长的
2、著名比如，描述人口增长的著名人口模型：人口模型：常微分方程常微分方程是常用的数学模型。是常用的数学模型。该方程有该方程有解析解解析解y(x)=y0er(x-xo)。考虑考虑一阶一阶常常微分方程的微分方程的初值问题初值问题/*Initial-Value Problem*/:只要只要 f(x,y)在在a,b R1 上连续，且关于上连续，且关于 y 满足满足 Lipschitz 条条件件，即存在与，即存在与 x,y 无关的常数无关的常数 L 使使对任意定义在对任意定义在 a,b 上的上的 y1(x)和和 y2(x)都成立，则上述都成立，则上述IVP存在唯一解存在唯一解。要计算出解函数要计算出解函数
3、y(x)在一系列节点在一系列节点 a=x0 x1 xn=b 处的近似值处的近似值节点间距节点间距hi=xi+1-xi称为称为步长步长，当，当 hi=h 为常数时称为为常数时称为等步长等步长。利用利用y0求节点求节点x1处的近似值处的近似值y1,再从再从y1来求出来求出y2,直直至求出所有的至求出所有的yn.称之为称之为步进法步进法.1 欧拉方法欧拉方法/*Eulers Method*/欧拉公式：欧拉公式：xixi+1向前差商近似导数向前差商近似导数对于对于xia,b,有有 则则 有有 将上式中的函数值将上式中的函数值y(xi)都用近似值都用近似值yi来表示来表示,则有数值计算则有数值计算格式格
4、式7.1 Eulers Method例例7.1 求解初值问题求解初值问题 取步长分别为取步长分别为h=0.1和和h=0.05,进行计算进行计算,结果有结果有解解:该方程的解析解是该方程的解析解是y=(1+2x)1/2.欧拉格式是欧拉格式是k xi y(xi)yi (h=0.1)yi (h=0.05)012345.0.00.10.20.30.40.5.1.01.10001.19181.27741.35821.4351.1.01.09771.18761.27131.35021.4174.1.01.09541.18321.26491.34161.4142.在在假假设设 yi=y(xi)，即即第第 i
5、 步步计计算算是是精精确确的的前前提提下下，考考虑虑的的截截断断误误差差 Ri=y(xi+1)yi+1 称称为为局局部部截截断断误误差差/*local truncation error*/。若若某某算算法法的的局局部部截截断断误误差差为为O(hp+1)，则则称称该该算算法法有有p 阶阶精度。精度。欧拉法的局部截断误差：欧拉法的局部截断误差：欧拉法具有欧拉法具有 1 阶精度。阶精度。Ri 的的主项主项/*leading term*/局部截断误差：局部截断误差：7.1 Eulers Method 欧拉公式的改进：欧拉公式的改进：隐式欧拉法隐式欧拉法/*implicit Euler method*/
6、向后差商近似导数向后差商近似导数xixi+1y(xi+1)yi+hf(xi+1,y(xi+1)1,.,0(),(111=+=+niyxfhyyiiii由于未知数由于未知数 yi+1 同时出现在等式的两边，不能直接得到，故同时出现在等式的两边，不能直接得到，故称为称为隐式隐式/*implicit*/欧拉公式，而前者称为欧拉公式，而前者称为显式显式/*explicit*/欧拉公式。欧拉公式。隐式隐式欧拉法的局部截断误差：欧拉法的局部截断误差：即即隐式欧拉公式具有隐式欧拉公式具有 1 阶精度。阶精度。7.1 Eulers Method 梯形公式梯形公式/*trapezoid formula*/显、隐
7、式两种算法的显、隐式两种算法的平均平均注：注：的确有局部截断误差的确有局部截断误差 ，即梯形公式具有即梯形公式具有2 阶精度，比欧拉方法有了进步。阶精度，比欧拉方法有了进步。但注意到该公式是但注意到该公式是隐式隐式公式，不便于实际计算。公式，不便于实际计算。)1,.,0(),(111=+=+niyxfhyyiiii7.1 Eulers Method 改进欧拉法改进欧拉法/*modified Eulers method*/Step 1:先用先用显式显式欧拉公式作欧拉公式作预测预测，算出，算出),(1iiiiyxfhyy+=+Step 2:再将再将 代入代入隐式隐式梯形公式的右边作梯形公式的右边作
8、校正校正，得到，得到1+iy),(),(2111+=iiiiiiyxfyxfhyy7.1 Eulers Method或者或者 例例7.2 求解初值问题求解初值问题 取步长分别为取步长分别为h=0.1和和h=0.05,进行计算进行计算,结果有结果有解解:该方程的解析解是该方程的解析解是y=(1+2x)1/2.改进欧拉格式改进欧拉格式是是k xi y(xi)yi(h=0.1)yi(h=0.05)012345.0.00.10.20.30.40.5.1.01.09591.18411.26621.34321.4164.1.01.09561.18351.26531.34211.4148.1.01.0954
9、1.18321.26491.34161.4142.7.1 Eulers Method2 龙格龙格-库塔法库塔法/*Runge-Kutta Method*/建立高精度的单步递推格式。建立高精度的单步递推格式。xnxn+1xn+hK*考察欧拉公式考察欧拉公式 K*K1=f(xi,yi)改进欧拉公式改进欧拉公式 K*(K1+K2)/2斜率斜率一定取一定取K1 K2 的的平均值平均值吗吗？步长一定是步长一定是一个一个h 吗吗？这里称为函这里称为函数数y(x)在区间在区间xn,xn+1上的上的平均斜平均斜率率.几何意义是几何意义是y(x)在区间在区间xn,xn+1上点上点xn+h处的处的斜率斜率.根据微
10、分中值定根据微分中值定理理,存在存在01 2 Runge-Kutta Method首先希望能确定系数首先希望能确定系数 1、2、p，使得到的算法格式有使得到的算法格式有2阶阶精度，即在精度，即在 的前提假设下，使得的前提假设下，使得 Step 1:将将 K2 在在(xi,yi)点作点作 Taylor 展开展开将将改进欧拉法推广为：改进欧拉法推广为：),(),(12122111phKyphxfKyxfKKKhyyiiiiii+=+=+Step 2:将将 K2 代入公式，得到代入公式，得到2 Runge-Kutta MethodStep 3:将将 yi+1 与与 y(xi+1)在在 xi 点点的的
11、泰勒泰勒展开作比较展开作比较要求要求 ，则必须有：，则必须有：这里有这里有 个未知个未知数，数，个方程。个方程。32存在存在无穷多个解无穷多个解。所有满足上式的格式统称为。所有满足上式的格式统称为2阶龙格阶龙格-库库塔格式塔格式。注意到，注意到，就是改进的欧拉法。就是改进的欧拉法。Q:为获得更高的精度，应该如何进一步推广？为获得更高的精度，应该如何进一步推广？其中其中 i (i=1,m)，i (i=2,m)和和 ij(i=2,m;j=1,i 1)均为待均为待定系数，确定这些系数定系数，确定这些系数的步骤与前面相似。的步骤与前面相似。2 Runge-Kutta Method).,(.),(),(
12、),(.1122112321313312122122111 +=+=+=+=mm mmmmimiiiiiimmiihKhKhKyhxfKhKhKyhxfKhKyhxfKyxfKKKKhyy 常用的三阶常用的三阶龙格龙格-库塔法库塔法/*Classical Runge-Kutta Method*/：2 Runge-Kutta Method 常用的四阶常用的四阶龙格龙格-库塔法库塔法：例例7.3 用四阶龙格库塔公式求解例题用四阶龙格库塔公式求解例题7.1 取取h=0.1计算计算,结果有结果有解解:对于例题对于例题7.1有以下四阶龙格库塔计算格式有以下四阶龙格库塔计算格式k xi y(xi)yi K
13、1 K2 K3 K4 012345.0.00.10.20.30.40.5.1.00.912870.845120.790570.74536.0.954760.878180.817470.767810.72622.1.01.09541.18321.26491.34161.4142.0.952290.876050.815620.766200.72481.0.912620.844940.790390.745200.70697.1.01.09541.18321.26491.34161.4142.2 Runge-Kutta Method3 微分方程组与高阶方程微分方程组与高阶方程 /*Systems of
14、 Differential Equations and Higher-Order Equations*/一阶微分方程组一阶微分方程组IVP的一般形式为：的一般形式为：=)(,.),(,()(.)(,.),(,()(1111xyxyxfxyxyxyxfxymmmm初值初值0002020101)(,.,)(,)(mmyxyyxyyxy=将将问题记作问题记作向量形式向量形式，令：，令：前述所有公式皆适前述所有公式皆适用于向量形式。用于向量形式。改进欧拉公式的形式改进欧拉公式的形式 3 Systems of DEs and Higher-Order Equations 四阶龙格库塔公式四阶龙格库塔公式 预报预报：校正校正：例例 高阶微分方程高阶微分方程3 Systems of DEs and Higher-Order Equations=10)1(1000)1()()(,.,)(,)(),.,(nnnnaxyaxyaxyyyyxfy化作化作一阶微分方程组一阶微分方程组求解。求解。引入新变量引入新变量初值条件为：初值条件为：