第12次常微分方程初值问题数值解法计算方法计算方法(Numerical Analysis)2021/7/11内容1.常微分方程初值问题解的存在性定理2.Euler公式3.梯形公式4.两步Euler公式5.欧拉法的局部截断误差6.改进型Euler公式7.龙格-库塔法8.算法实现2021/7/12常微分方程初值问题解的存在性定理2021/7/13第9章常微分方程初值问题数值解法包含自变量、未知函数及未知函数的导数的方程称为微分方程。微分方程中出现的未知函数最高阶导数的阶数称为微分方程的阶数。都是一次的,则称其为线性的,否则称为非线性的。如果未知函数y及其各阶导数9.1引言自变量个数只有一个的微分方程称为常微分方程。2021/7/14如下是一些典型方程求解析解的基本方法可分离变量法、常系数齐次线性方程的解法、常系数非齐次线性方程的解法等。的解就不能用初等函数及其积分来表达。但能求解的常微分方程仍然是很少的，大多数的常微分方程是不可能给出解析解。例如，一阶微分方程2021/7/15从实际问题当中归纳出来的微分方程，通常主要依靠数值解法来解决。(9.1)在区间axb上的数值解法。本章主要讨论一阶常微分方程初值问题2021/7/16定理1：如果函数f(x,y)在带形区域则方程(9.1)在 a,b 上存在唯一的连续可微分的解的解y=y(x)。内连续，且关于y满足李普希兹(Lipschitz)条件，即存在常数L(它与x,y无关)使2021/7/17推论：如果函数f(x,y)对y的偏导数在带形区域对R内的所有x，y都成立。即存在常数L(它与x,y无关)使则方程(9.1)在 a,b 上存在唯一的连续可微解y=y(x)。内有界。Home2021/7/18Euler公式本章假设微分方程初值问题(9.1)有解2021/7/19常微分方程初值问题(9.1)的数值解法的基本思想：算出精确解y(x)在区间 a,b 上的一系列离散节点的近似值处的函数值y=y(x)a=x0xn=bx1x2x3(未知)2021/7/110相邻两个节点的间距称为步长，步长可以相等，也可以不等。数值解法需要把连续性的问题加以离散化，从而求出离散节点的数值解。a=x0xn=bx1x2x3xn-1x4本章总是假定h为定数，称为定步长，这时节点可表示为2021/7/111常微分方程数值解法的基本出发点：离散化。采用“步进式”，即求解过程顺着节点排列的次序逐步向前推进。中的导数进行离散化处理。以便对初值问题计算的递推公式。算法：要求给出用已知信息2021/7/112欧拉（Euler）方法是解初值问题的最简单的数值方法。9.2简单的数值方法与基本概念的解y=y(x)代表通过点的一条称之为微分方程的积分曲线。积分曲线上每一点的切线的斜率等于函数在这点的函数值。9.2.1Euler公式初值问题2021/7/113Euler法的求解过程:从初始点P0(即点(x0,y0)出发，作积分曲线y=y(x)在P0点上切线，其斜率为y=y(x)x0xix1yx2P1(x1 , y1)P0Pnxi+1xnP2(x2 , y2)Pi(xi , yi) Pi+1(xi+1 , yi+1)y(x1)y(x2)y(xi)y(xi+1)y(xn)y(x0)2021/7/114这样就获得了P1点的坐标：(x1,y1)。将y1作为y(x1)的近似值(想象(x1,y1)在积分曲线y=y(x)上)当时，得过点P1(x1,y1)，作积分曲线y=y(x)的切线交直线x=x2于P2点。注意切线的斜率(近似)为直线方程为:当时,得由此获得了P2的坐标。直线的方程为：2021/7/115当时,得重复以上过程,对已求得点,以为（近似）斜率作直线这样,从x0逐个算出对应的数值解从图形上看,就获得了一条近似于曲线y=y(x)的折线。就获得了一系列的点:P1,P1,Pn。2021/7/116y=y(x)x0xix1yx2y1P0Pnxi+1xny2yiyi+1y(x1)y(x2)y(xi)y(xi+1)y(xn)y(x0)yn微分方程(9.1)的精确解y=y(x)的近似解为：y1，y2，,yn2021/7/117注：还可用数值积分法和泰勒展开法推导Euler公式（略）。Euler公式Euler法的计算公式可以表达为： ( 9.2 ) 其中，为常数，i=0,1,n2021/7/118解：取h=0.1，根据Euler公式，得例9.1：利用Euler公式求解微分方程的初值问题初值问题有解：由x0=0,y0=1，代入以上公式，得y1=1.1*y0-0.2*x0/y0=1.12021/7/119课堂练习：计算出x2,y2；x3,y3x0=0,y0=1x1=0.1,y1=1.12021/7/120xnyny(x n)0.11.10001.09540.21.19181.18320.31.27741.26490.41.35821.34160.51.43511.41420.61.50901.48320.71.58031.54920.81.64981.61250.91.71781.67331.01.78481.7321计算结果比较：初值问题有解：可以由此公式计算出准确解：y(xn)欧拉法准确值2021/7/121y=y(x)的近似解010.1 0.20.30.40.50.60.70.80.9Home11.522021/7/122梯形公式2021/7/1239.2.2梯形公式(9.4)改用梯形方法计算其积分项，即为了提高精度,对方程的两端在区间上积分得，2021/7/124(9.5)式的右端含有未知的yi+1,它是一个关于yi+1的函数方程,这类数值方法称为隐式方法。相反地,欧拉法是显式方法。代入(7.4)式,并用近似代替式中即可得到梯形公式(9.5)由于数值积分的梯形公式比矩形公式的精度高，所以梯形公式（9.5）比欧拉公式(9.2)的精度高。求解困难Home2021/7/125两步Euler公式2021/7/126对方程两端在区间上积分得(9.6)改用中矩形公式计算其积分项，即代入上式,并用yi近似代替式中y(xi)即可得到(9.7)9.2.3两步欧拉公式两步欧拉公式2个区间2021/7/127【注】欧拉方法和梯形方法，都是单步法,其特点是在计算yi+1时只用到前一步的信息yi;而两步欧拉公式(9.7)中除了yi外,还用到更前一步的信息yi-1,即调用了前两步的信息。Home2021/7/128欧拉法的局部截断误差2021/7/1299.2.4欧拉法的局部截断误差定义9.1在yi准确的前提下,即时,用数值方法计算yi+1的误差：衡量求解公式好坏的一个主要标准是求解公式的精度,因此引入局部截断误差和阶数的概念。称为该数值方法计算时yi+1的局部截断误差。2021/7/130(b)-(a），得在欧拉公式中，假定(近似地相等)，则有(a)而将真解y(x)在xi处按二阶泰勒展开，得(b)欧拉公式的截断误差推导：2021/7/131定 义 9.2 若 数 值 方 法 的 局 部 截 断 误 差 为 ，则称这种数值方法的精度阶数是P。步长(hN结束。9.2.6改进欧拉法算法实现2021/7/157（ 2 ）改改进进欧欧拉拉法法的的流流程程图图 (3)程序实现(改进欧拉法计算常微分方程初值问题)2021/7/158clearx=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程序2021/7/1599.3.5四阶龙格库塔法算法实现(1)计算步骤输入使用龙格库塔公式（9.20）计算出y1输出，并使转到直至nN结束。2021/7/160(2)四阶龙格库塔算法流程图(3)程序实现(4阶龙格-库塔法计算常微分方程初值问题)2021/7/161作业作业习题九习题九 P3161、2、5(1)Home2021/7/162 结束语结束语若有不当之处，请指正，谢谢！若有不当之处，请指正，谢谢！