计算方法常微分方程初值问题数值解法.ppt

《计算方法常微分方程初值问题数值解法.ppt》由会员分享，可在线阅读，更多相关《计算方法常微分方程初值问题数值解法.ppt（62页珍藏版）》请在淘文阁 - 分享文档赚钱的网站上搜索。

1、第12次常微分方程初值问题数值解法计算方法(Numerical Analysis)内容1.常微分方程初值问题解的存在性定理2.Euler公式3.梯形公式4.两步Euler公式5.欧拉法的局部截断误差6.改进型Euler公式7.龙格-库塔法8.算法实现常微分方程初值问题解的存在性定理第9章常微分方程初值问题数值解法 包含自变量、未知函数及未知函数的导数的方程称为微分方程。微分方程中出现的未知函数最高阶导数的阶数称为微分方程的阶数。都是一次的,则称其为线性的,否则称为非线性的。如果未知函数y及其各阶导数9.1引言 自变量个数只有一个的微分方程称为常微分方程。如下是一些典型方程求解析解的基本方法 可

2、分离变量法、常系数齐次线性方程的解法、常系数非齐次线性方程的解法等。的解就不能用初等函数及其积分来表达。但能求解的常微分方程仍然是很少的，大多数的常微分方程是不可能给出解析解。例如，一阶微分方程 从实际问题当中归纳出来的微分方程，通常主要依靠数值解法来解决。(9.1)在区间axb上的数值解法。本章主要讨论一阶常微分方程初值问题定理1：如果函数f(x,y)在带形区域则方程(9.1)在a,b 上存在唯一的连续可微分的解的解y=y(x)。内连续，且关于y满足李普希兹(Lipschitz)条件，即存在常数L(它与x,y无关)使推论：如果函数f(x,y)对y的偏导数在带形区域对R内的所有x，y都成立。即

3、存在常数L(它与x,y无关)使则方程(9.1)在a,b上存在唯一的连续可微解y=y(x)。内有界。HomeEuler公式本章假设微分方程初值问题(9.1)有解常微分方程初值问题(9.1)的数值解法的基本思想：算出精确解y(x)在区间a,b 上的一系列离散节点的近似值处的函数值y=y(x)a=x0 xn=b x1x2x3(未知)相邻两个节点的间距称为步长，步长可以相等，也可以不等。数值解法需要把连续性的问题加以离散化，从而求出离散节点的数值解。a=x0 xn=b x1x2x3xn-1x4 本章总是假定h为定数，称为定步长，这时节点可表示为 常微分方程数值解法的基本出发点：离散化。采用“步进式”，

4、即求解过程顺着节点排列的次序逐步向前推进。中的导数进行离散化处理。以便对初值问题计算的递推公式。算法：要求给出用已知信息欧拉（Euler）方法是解初值问题的最简单的数值方法。9.2简单的数值方法与基本概念的解y=y(x)代表通过点的一条称之为微分方程的积分曲线。积分曲线上每一点的切线的斜率等于函数在这点的函数值。9.2.1Euler公式初值问题Euler法的求解过程:从初始点P0(即点(x0,y0)出发，作积分曲线y=y(x)在P0点上切线，其斜率为y=y(x)x0 xix1yx2P1(x1,y1)P0Pnxi+1xnP2(x2,y2)Pi(xi,yi)Pi+1(xi+1,yi+1)y(x1)

5、y(x2)y(xi)y(xi+1)y(xn)y(x0)这样就获得了P1点的坐标：(x1,y1)。将y1作为y(x1)的近似值(想象(x1,y1)在积分曲线y=y(x)上)当时，得过点P1(x1,y1)，作积分曲线y=y(x)的切线交直线x=x2于P2点。注意切线的斜率(近似)为直线方程为:当时,得由此获得了P2的坐标。直线的方程为：当时,得重复以上过程,对已求得点,以为（近似）斜率作直线这样,从x0逐个算出对应的数值解从图形上看,就获得了一条近似于曲线y=y(x)的折线。就获得了一系列的点:P1,P1,Pn。y=y(x)x0 xix1yx2y1P0Pnxi+1xny2yiyi+1y(x1)y(

6、x2)y(xi)y(xi+1)y(xn)y(x0)yn微分方程(9.1)的精确解y=y(x)的近似解为：y1，y2，,yn注：还可用数值积分法和泰勒展开法推导Euler公式（略）。Euler公式Euler法的计算公式可以表达为：(9.2)其中，为常数，i=0,1,n解：取h=0.1，根据Euler公式，得例9.1：利用Euler公式求解微分方程的初值问题初值问题有解：由x0=0,y0=1，代入以上公式，得y1=1.1*y0-0.2*x0/y0=1.1课堂练习：计算出x2,y2；x3,y3x0=0,y0=1x1=0.1,y1=1.1xnyny(x n)0.1 1.1000 1.09540.2 1

7、.1918 1.18320.3 1.2774 1.26490.4 1.3582 1.34160.5 1.4351 1.41420.6 1.5090 1.48320.7 1.5803 1.54920.8 1.6498 1.61250.9 1.7178 1.67331.0 1.7848 1.7321计算结果比较：初值问题有解：可以由此公式计算出准确解：y(xn)欧拉法准确值y=y(x)的近似解010.1 0.20.30.4 0.5 0.60.70.8 0.9Home11.52梯形公式9.2.2梯形公式(9.4)改用梯形方法计算其积分项，即为了提高精度,对方程的两端在区间上积分得，(9.5)式的右端

8、含有未知的yi+1,它是一个关于yi+1的函数方程,这类数值方法称为隐式方法。相反地,欧拉法是显式方法。代入(7.4)式,并用近似代替式中即可得到梯形公式(9.5)由于数值积分的梯形公式比矩形公式的精度高，所以梯形公式（9.5）比欧拉公式(9.2)的精度高。求解困难Home两步Euler公式对方程两端在区间上积分得(9.6)改用中矩形公式计算其积分项，即代入上式,并用yi近似代替式中y(xi)即可得到(9.7)9.2.3两步欧拉公式两步欧拉公式2个区间【注】欧拉方法和梯形方法，都是单步法,其特点是在计算yi+1时只用到前一步的信息yi;而两步欧拉公式(9.7)中除了yi外,还用到更前一步的信息

9、yi-1,即调用了前两步的信息。Home欧拉法的局部截断误差9.2.4欧拉法的局部截断误差定义9.1在yi准确的前提下,即时,用数值方法计算yi+1的误差：衡量求解公式好坏的一个主要标准是求解公式的精度,因此引入局部截断误差和阶数的概念。称为该数值方法计算时yi+1的局部截断误差。(b)-(a），得 在欧拉公式中，假定(近似地相等)，则有(a)而将真解y(x)在xi处按二阶泰勒展开，得(b)欧拉公式的截断误差推导：定 义 9.2 若 数 值 方 法 的 局 部 截 断 误 差 为，则称这种数值方法的精度阶数是P。步长(h1)越小，P越高，则局部截断误差越小。计算精度越高。评论：欧拉公式的精度讨

10、论两步欧拉公式的局部截断误差为：从而两步欧拉公式的阶数是2.推导过程省略。欧拉公式的局部截断误差为：y(xi+1)yi+1=O(h2)欧拉方法仅为一阶方法。Home改进的Euler公式9.2.5改进的欧拉公式欲综合欧拉公式和梯形公式，得到改进的欧拉公式。计算工作量小，但精度低。显式欧拉公式(9.2)梯形公式精度高些,但为隐式公式,需用迭代法求解,计算量大。(9.5)(9.10)预测：校正：先用欧拉公式(9.2)求出一个初步的近似值,称为预测值,它的精度不高；改进的思路：这种预测-校正方法称为改进的欧拉公式：再用梯形公式(9.5)对校正一次，即迭代一次，求得yi+1,称为校正值。可以证明,公式(

11、9.10)的精度为二阶。(9.11)通常表示成下列平均化形式(9.12)(9.10)可以改写为如下的一步显式格式：例9.2用改进欧拉法解初值问题区间为0,1，取步长h=0.1解:改进欧拉法公式课堂练习：x0=0,y0=1计算出xnynyny(x n)0.1 1.1000 1.0959 1.09540.2 1.1918 1.1841 1.18320.3 1.2774 1.2662 1.26490.4 1.3582 1.3434 1.34160.5 1.4351 1.4164 1.41420.6 1.5090 1.4860 1.48320.7 1.5803 1.5525 1.54920.8 1.6

12、498 1.6165 1.61250.9 1.7178 1.6782 1.67331.0 1.7848 1.7379 1.7321计算结果比较(Matlab)：初值问题有解：按照此公式计算出：y(xn)a)并且和由欧拉方法计算结果进行比较。b)与改进的欧拉方法计算结果进行比较对比结果：改进欧拉方法精确度更高改进欧拉法 欧拉法准确值Home龙格-库塔法9.3.1龙格-库塔(Runge-Kutta)法的基本思想9.3龙格-库塔（Runge-Kutta）法Euler公式可改写成(*)Euler公式的误差公式为：Euler公式是1阶方法。改进的Euler公式又可改写成 欧拉公式与改进欧拉公式在形式上有

13、一个共同点:都是用f(x,y)在某些点上值的线性组合得出y(xi+1)的近似值yi+1,而且增加计算f(x,y)的次数,可提高截断误差的阶。欧拉公式:每步计算一次f(x,y)的值,局部截断误差为 于是可考虑用函数f(x,y)在若干点上的函数值的线性组合来构造近似公式。如此，可以构造出更高精度的计算格式，这就是龙格-库塔法的基本思想。改进欧拉公式：需计算两次f(x,y)的值，它是2阶方法，局部截断误差为。事实上，将方程的两端在区间上积分得，针对右端的积分可以使用不同的积分公式进行近似求解。(*)和号中f(x)的取值节点越多，就越精确。表达式(*)可近似地表达为（和号代表积分公式）：Simpson

14、公式：积分区间a,b划分为2等分，3个节点牛顿柯特斯：将积分区间a,b划分为n等分,n+1个节点。在上取两点，经过精心构造，得到如下的2阶龙格-库塔格式(有2阶精度)9.3.2二阶龙格-库塔法（9.14）其中（9.17）满足条件(9.17)的格式(9.14)的局部截断误差为此条件保证了(9.14)近似效果最好。式(9.17)中有三个未知量,但只有两个方程,因而有无穷多解。若取,则p=1，将以上所解的值代入式(9.14)并改写可得 不难发现，上面的格式就是改进的欧拉格式。9.3.4四阶龙格库塔法如果需要再提高精度，用类似上述的处理方法，只需在区间上用四个点处的斜率加权平均作为平均斜率k*的近似值

15、，构成一系列四阶龙格-库塔公式。具有四阶精度，即局部截断误差是。（9.20）经过精心的推导与构造（过程从略）得到最常用的一种4阶经典龙格-库塔公式。局部截断误差是。例9.3取步长h=0.2，用经典龙格-库塔公式求解初值问题。解:由四阶龙格-库塔公式可得：课堂：取h=0.2,计算y1取h=0.2,计算y1;x0=0,y0=1,取n=0取h=0.2,计算y2;x1=0.2,y1=1.1832,取n=1xnynynyny(x n)0.1 1.1000 1.0959 1.09540.2 1.1918 1.1841 1.1832 1.18320.3 1.2774 1.2662 1.26490.4 1.3

16、582 1.3434 1.34162 1.34160.5 1.4351 1.4164 1.41420.6 1.5090 1.4860 1.4833 1.48320.7 1.5803 1.5525 1.54920.8 1.6498 1.6165 1.6125 1.61250.9 1.7178 1.6782 1.67331.0 1.7848 1.7379 1.7321 1.7321改进欧拉法欧拉法准确值四阶龙格-库塔法 龙格-库塔方法的推导基于Taylor展开方法，因而它要求所求的解具有较好的光滑性。如果解的光滑性差，那么，使用四阶龙格库塔方法求得的数值解，其精度可能反而不如改进的欧拉方法。在实际

17、计算时，应当针对问题的具体特点选择合适的算法。龙格-库塔方法总结：算法实现这部分不会出现在笔试卷子中；只作为Matlab编程的参考使用。（1）计算步骤输入使用以下改进欧拉法公式进行计算输出，并使转直至nN结束。9.2.6改进欧拉法算法实现（2）改进欧拉法的流程图(3)程序实现(改进欧拉法计算常微分方程初值问题)clearx=0,yn=1%初始化forn=1:10yp=yn+0.1*(yn-2*x/yn);%预测x=x+0.1;yc=yn+0.1*(yp-2*x/yp);%校正yn=(yp+yc)/2%平均（再矫正）end本题的精确解为改进欧拉法的Matlab程序9.3.5四阶龙格库塔法算法实现(1)计算步骤输入使用龙格库塔公式（9.20）计算出y1输出，并使转到直至nN 结束。(2)四阶龙格 库塔算法流程图(3)程序实现(4阶龙格-库塔法计算常微分方程初值问题)作业习题九 P3161、2、5(1)Home