《计算方法 常微分方程初值问题数值解法-euler公式-龙格-库塔法》由会员分享，可在线阅读，更多相关《计算方法 常微分方程初值问题数值解法-euler公式-龙格-库塔法（62页珍藏版）》请在金锄头文库上搜索。
1、第12次 常微分方程初值问题数值解法,计算方法 (Numerical Analysis),内容,常微分方程初值问题解的存在性定理 Euler公式 梯形公式 两步Euler公式 欧拉法的局部截断误差 改进型Euler公式 龙格-库塔法 算法实现,常微分方程初值问题 解的存在性定理,第9章 常微分方程初值问题数值解法,包含自变量、未知函数及未知函数的导数的方程称为微分方程。,微分方程中出现的未知函数最高阶导数的阶数称为微分方程的阶数。,9.1 引言,自变量个数只有一个的微分方程称为常微分方程。,如下是一些典型方程求解析解的基本方法 可分离变量法、 常系数齐次线性方程的解法、 常系数非齐次线性方程的
2、解法等。,从实际问题当中归纳出来的微分方程，通常主要依靠数值解法来解决。,定理1：如果函数f(x, y)在带形区域,则方程( 9.1 ) 在a, b上存在唯一的连续可微分的解的解 y=y(x) 。,内连续，且关于y满足李普希兹(Lipschitz)条件，即存在常数L(它与x, y无关)使,推论：如果函数f(x,y)对y的偏导数 在带形区域,对R内的所有x，y 都成立。,即存在常数L(它与x,y无关)使,则方程( 9.1 ) 在a , b上存在唯一的连续可微解y=y(x) 。,内有界。,Home,Euler公式,本章假设微分方程初值问题(9.1)有解,常微分方程初值问题(9.1)的数值解法的基本
3、思想： 算出精确解y(x)在区间a,b上的一系列离散节点,的近似值,处的函数值,y=y(x),a=x0,xn=b,x1,x2,x3,(未知),相邻两个节点的间距 称为步长，步长可以相等，也可以不等。,数值解法需要把连续性的问题加以离散化，从而求出离散节点的数值解。,常微分方程数值解法的基本出发点：离散化。采用“步进式”，即求解过程顺着节点排列的次序逐步向前推进。,中的导数 进行离散化处理。,以便对初值问题,欧拉（Euler）方法是解初值问题的最简单的数值方法。,9.2 简单的数值方法与基本概念,的解y=y(x)代表通过点 的一条称之为微分方程的积分曲线。,积分曲线上每一点 的切线的斜率 等于函
4、数 在这点的函数值。,9.2.1 Euler公式,初值问题,Euler法的求解过程: 从初始点P0( 即点(x0,y0)出发，作积分曲线y=y(x)在P0点上切线 ，其斜率为,y=y(x),x0,xi,x1,y,x2,P1(x1 , y1),P0,Pn,xi+1,xn,P2(x2 , y2),Pi(xi , yi),Pi+1(xi+1 , yi+1),y(x1),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)，作积分
5、曲线y=y(x)的切线交直线x=x2于P2点。注意切线 的斜率(近似)为,直线 方程为:,当 时,得,由此获得了P2的坐标。,直线 的方程为：,当 时,得,重复以上过程, 对已求得点 ,以 为（近似）斜率作直线,y=y(x),x0,xi,x1,y,x2,y1,P0,Pn,xi+1,xn,y2,yi,yi+1,y(x1),y(x2),y(xi),y(xi+1),y(xn),y(x0),yn,微分方程(9.1)的精确解y=y(x)的近似解为： y1， y2， , yn,注：还可用数值积分法和泰勒展开法推导 Euler公式（略）。,Euler公式,解：取h=0.1，根据Euler公式，得,例9.1：
6、利用Euler公式求解微分方程的初值问题,由x0 = 0, y0 =1，代入以上公式，得 y1 = 1.1*y0 - 0.2* x0/y0 = 1.1,课堂练习：计算出 x2 , y2；x3 , y3,x0 = 0, y0 =1 x1 = 0.1 , y1= 1.1,计算结果比较：,初值问题有解：,可以由此公式计算出准确解：y(xn),欧拉法,准确值,y=y(x)的近似解,0,1,0.1,0.2,0.3,0.4,0.5,0.6,0.7,0.8,0.9,Home,1,1.5,2,梯形公式,9.2.2 梯形公式,(9.4 ),改用梯形方法计算其积分项，即,为了提高精度,对方程 的两端在区间 上积分
7、得，,(9.5)式的右端含有未知的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时只用到前
8、一步的信息yi;,而两步欧拉公式 (9.7)中除了yi外,还用到更前一步的信息yi-1,即调用了前两步的信息。,Home,欧拉法的局部截断误差,9.2.4 欧拉法的局部截断误差,定义9.1 在yi准确的前提下, 即 时, 用数值方法计算yi+1的误差：,衡量求解公式好坏的一个主要标准是求解公式的精度, 因此引入局部截断误差和阶数的概念。,称为该数值方法计算时yi+1的局部截断误差。,(b)-(a），得,欧拉公式的截断误差推导：,定义9.2 若数值方法的局部截断误差为 ，则称这种数值方法的精度阶数是P。,步长(h1) 越小，P越高，则局部截断误差越小。计算精度越高。,评论：,欧拉公式的精度讨论,
9、两步欧拉公式的局部截断误差为：,从而两步欧拉公式的阶数是2.推导过程省略。,欧拉公式的局部截断误差为 ： y(xi+1 ) yi+1 = O(h2) 欧拉方法仅为一阶方法。,Home,改进的Euler公式,9.2.5 改进的欧拉公式,欲综合欧拉公式和梯形公式，得到改进的欧拉公式。,计算工作量小，但精度低。,显式欧拉公式,(9.2),(9.10),预测：,校正：,先用欧拉公式(9.2)求出一个初步的近似值, 称为预测值, 它的精度不高；,改进的思路：,这种预测-校正方法称为改进的欧拉公式：,再用梯形公式(9.5)对 校正一次，即迭代一次，求得yi+1,称为校正值。,可以证明,公式(9.10)的精
10、度为二阶。,(9.11),通常表示成下列平均化形式,(9.12),(9.10)可以改写为如下的一步显式格式：,例9.2 用改进欧拉法解初值问题,区间为0,1，取步长h=0.1,解: 改进欧拉法公式,课堂练习： x0 = 0, y0 = 1 计算出,计算结果比较(Matlab)：,初值问题有解：,按照此公式计算出：y(xn) 并且和由欧拉方法计算结果进行比较。 与改进的欧拉方法计算结果进行比较 对比结果： 改进欧拉方法精确度更高,改进欧拉法,欧拉法,准确值,Home,龙格-库塔法,9.3.1 龙格-库塔(Runge-Kutta)法的基本思想,9.3 龙格-库塔（Runge-Kutta）法,Eul
11、er公式可改写成,(*),改进的Euler公式又可改写成,欧拉公式与改进欧拉公式在形式上有一个共同点: 都是用f(x,y)在某些点上值的线性组合得出y(xi+1)的近似值yi+1,而且增加计算f(x,y)的次数,可提高截断误差的阶。,欧拉公式: 每步计算一次f(x,y)的值,局部截断误差为,于是可考虑用函数f(x,y)在若干点上的函数值的线性组合来构造近似公式。,如此，可以构造出更高精度的计算格式，这就是龙格-库塔法的基本思想。,改进欧拉公式：需计算两次f(x,y)的值，它是2阶方法，局部截断误差为 。,事实上，将方程 的两端在区间 上积分得，,针对右端的积分可以使用不同的积分公式进行近似求解
12、。,(*),和号中f(x)的取值节点越多，就越精确。,表达式(*)可近似地表达为（和号代表积分公式）：,Simpson公式：积分区间 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)并改写可得,不难发现，
13、上面的格式就是改进的欧拉格式。,9.3.4 四阶龙格库塔法,如果需要再提高精度，用类似上述的处理方法，只需在区间 上用四个点处的斜率加权平均作为平均斜率k*的近似值，构成一系列四阶龙格-库塔公式。具有四阶精度，即局部截断误差是 。,（9.20）,经过精心的推导与构造（过程从略）得到最常用的一种4阶经典龙格-库塔公式。,局部截断误差是 。,例9.3 取步长h=0.2，用经典龙格-库塔公式求解 初值问题。,解:,由四阶龙格-库塔公式可得：,课堂：取h=0.2,计算y1,取h=0.2,计算y1 ; x0 =0, y 0 = 1,取n=0,取h=0.2,计算y2 ; x1 =0.2, y1 = 1.1
14、832, 取n=1,改进 欧拉法,欧拉法,准确值,四阶龙格 -库塔法,龙格-库塔方法的推导基于Taylor展开方法，因而它要求所求的解具有较好的光滑性。 如果解的光滑性差，那么，使用四阶龙格库塔方法求得的数值解，其精度可能反而不如改进的欧拉方法。 在实际计算时，应当针对问题的具体特点选择合适的算法。,龙格-库塔方法总结：,算法实现,这部分不会出现在笔试卷子中；只作为Matlab编程的参考使用。,（1）计算步骤 输入 使用以下改进欧拉法公式进行计算, 输出 ，并使 转 直至n N 结束。,9.2.6 改进欧拉法算法实现,（2）改进欧拉法的流程图,(3) 程序实现(改进欧拉法计算常微分方程初值问题
15、 ),clear x=0, yn=1 %初始化 for n = 1:10 yp = 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 输出 ，并使 转到 直至n N 结束。,(2) 四阶龙格库塔算法流程图,程序实现(4阶龙格-库塔法计算常微分方程初值问题),作业 习题九 P316 1、2、5(1),Home,
