常微分方程出值问题得数值解法课件

,单击此处编辑母版标题样式,单击此处编辑母版文本样式,第二级,第三级,第四级,第五级,*,*,2024/9/15,1,常微分方程初值问题的数值解法,8.2,几种简单的单步法,8.3,Runge,Kutta,公式,8,.4,单步法的收敛性、,相容性与稳定性,8.5,线性多步法,第八章,8.,1 引言,数值算例,2024/9/15,2,本章着重讨论一阶常微分方程初值问题,在,区间,a,b,上的数值解法。,8.,1 引言,一、,问题,这些问题多数情况下求不出解析解，只能用近似的方法求解。常用的近似方法有两类。一类称为近似解析法，如级数解法，逐次逼近法等。另一类称为数值解法，它可以给出解在一些离散点上的近似值。,2024/9/15,3,若,f,（,x,y,）在区域,D,=,上连续,，,且,关于,y,满足,李普希兹（,Lipschitz,）,条件,即存在常,数,L,使,对,D,内任两个 均成立，其中,L,是与,x,，,y,无,关,的常数，,则上面的初值问题存在唯一解，且解是连续可微的。,二、解的存在唯一性,Remark,：在,f,(,x,y,),对,y,可微的情况下，若偏导数有界，则可取 ，此时有,此时,Lipschitz,连续条件显然成立。这是验证该条件的最简便的方法。,2024/9/15,4,解的,适定性,是指解的存在唯一性以及数值稳定性。此处主要是指解对于右端项以及初值扰动的适应性。关于适定性有如下的结论：,三、解的适定性,定理,：若,f,（,x,y,）在区域,D,=,上满足,Lipschitz,连续条件,，则初值问题是适定的。,四、等价的积分方程,若,y,(,x,),是初值问题的解，对方程两边同时积分，利用初始条件可得：,该方程为与初值问题同解的积分方程，我们可以从积分方程出发去构造初值问题的求解公式。,2024/9/15,5,常微分方程初值问题的数值解是求上述初值问题的解,y,(,x,),在区间,a,，,b,中的点,列,上的近似值 。 称为步长，,一般情况下我们取等步长，记为,h,。,五、数值解法,初值问题的,解析解（理论解）用 表示,数值解法的精确解用,表示，并记,f,n,=,f,(,x,n,y,n,),而 。,求初值问题的数值解一般是逐步进行的，即计算出,y,n,之后计算,y,n,+1,。,2024/9/15,6,数值解法一般分为：,（1）单步法,：在计算,时，只用到,和,，,即前一步的值。,（2）,多步法,：计算,时，除用到,和,以外，还要用 和,即前,k,步的值。,数值解法的分类,单步法和多步法都有,显,式,和隐式,方法之分。显,式,和隐式的单步法可以分别写成：,对多步法来说，显式和隐式方法具有相同的意义。,2024/9/15,7,六、离散化方法,建立常微分方程初值问题数值解法的过程，就是通过一定的离散化方法，将对连续性问题的求解转化为求解常微分方程在有限个离散节点上解的近似值的过程，这个过程通常称为,数值离散,常用的数值离散化方法有：,差商代替微商法,Taylor,展开法,数值积分法,2024/9/15,8,一、显式欧拉（,Euler）,公式,设节点为,，,欧拉方法的计算公式,8.2,几种简单的单步法,这是一种最简单的显式单步方法，该方法可以通过不同的途径获得。,1、差商法,用两点差商公式 代替导数,再用,表示,的近似值,则得到,2024/9/15,9,假设在,附近把,y,（,x,）,展成,Taylor,级数,取,h,的线性部分,并用,表示 的近似值,得,显式欧拉（,Euler）,方公式,（续）,2、,Taylor,展开法,3、数值积分法,对微分方程两端从,x,n,到,x,n,+1,积分，得等价的积分方程,对右端的积分部分采用左矩形公式近似即得,Euler,公式。,Remark,：,Taylor,展开法与数值积分法是构造微分方程数值解的两类主要的方法。,2024/9/15,10,欧拉（,Euler）,方法的几何意义,Euler,方法有明显的,几何意义,。如右图所示，一阶常微分方程初值问题的解曲线,y,(,x,),过点,P,0,(,x,0,y,0,),。从,P,0,出发以,f,(,x,0,y,0,),为斜率作一直线段，与,x,=,x,1,相交于点,P,1,(,x,1,y,1,),， 显然有,y,1,=,y,0,+,hf,(,x,0,y,0,),。,同理，再由,P,1,出发以,f,(,x,1,y,1,),为斜率作一直线段，与,x,=,x,2,相交于点,P,2,(,x,2,y,2,),，显然有,y,2,=,y,1,+,hf,(,x,1,y,1,),。这样一直做下去，得到一条折线,P,0,P,1,P,2,，作为,y,(,x,),的近似曲线。因此，显示,Euler,公式又称为,Euler,折线法,。,2024/9/15,11,二、隐式,E,u,ler,公式,中将积分用右矩形公式,代入，有,上式是一个隐式的单步方法，称为,隐式欧拉公式,或,后退的欧拉公式,。利用此公式，每一步都要把上式作为,y,n,+1,的一个方程来求解。从数值积分的误差分析，很难期望隐式欧拉法比显式欧拉法更精确。,若在等价积分方程,从而得到,2024/9/15,12,隐式,E,u,ler,公式（续）,通常情况下，,隐式欧拉公式,很难直接求出,y,n,+1,的值，故常用迭代法求解。,在实际计算时，该公式通常与显式,Euler,公式结合使用，并由显式,Euler,公式的结果作为迭代的初始值，从而有如下数值格式,对 循环计算，若,(,为给定的误差限,),，则取 作为,y,n,+1,的近似值,2024/9/15,13,隐式,E,u,ler,公式（续）,由于,f,(,x,y,),关于,y,满足,Lipschitz,条件，故有,故当 时该迭代法收敛到隐式,Euler,公式的解,y,n,+1,，其中,L,为,Lipschitz,常数,2024/9/15,14,三、梯形公式,为了得到更精确的方法，在等价的积分方程中用梯形公式,近似积分项，再分别用,y,n,，,y,n,+1,代替,y,(,x,n,),和,y,(,x,n,1,)，,即可得到,梯形公式,：,该方法也是一种隐式的单步方法。对该方法，,从,n,0,开始计算，每步都要求解,y,n,+1,的一个方程。一般来说，这是一非线性方程，可迭代计算如下：,2024/9/15,15,梯形公式（续）,使用上式时，先用第一式计算出,y,n,+1,的近似值 ，再用第二式反复进行迭代，得到数列 ，用 来控制是否继续进行迭代，其中,为允许误差。把满足要求的 作为,y,(,x,n,+1,),的近似值,y,n,+1,，类似地可得出,y,n,+2,，,y,n,+3,，。,2024/9/15,16,梯形公式（续）,当,f,(,x,y,),关于,y,满足,Lipschitz,条件时，且步长,h,满足 时，上述迭代过程是收敛的。这是因为：,2024/9/15,17,实用中，,h,取得较小时,为了简化计算，,梯形公式第二式只迭代一次就结束，得到,E,u,l,e,r,预,测,校正方法,(,改进的,E,u,l,e,r,方法,),：,其中第一式称为,预测,算式，第二式称为,校正,算式。,四、,Euler,预测校正公式,2024/9/15,18,E,u,ler,预测校正公式（续）,若将,E,u,l,e,r,预-校,方法中的第一式带入第二式，得,Remark,：,这是一种显示的单步方法。有时为了计算方便，常将上式改写成：,2024/9/15,19,五、单步法的局部截断误差和阶,设一般的单步法为：,显式公式：,隐式公式：,设,为数值方法的精确值，,y,（,）,为微分方程的精确解。,定义,1,：,为某一数值方法在 处的,整体截断误差,。,Remark,：整体截断误差不仅与 这步的计算有关，而且与前面所有点的计算的误差累计有关。为了简化误差分析，我们着重分析计算中的某一步。对一般的显式单步法，有如下定义：,2024/9/15,20,单步法的局部截断误差和阶,（续）,这就是上面定义中称,R,n+,1,为“局部”的含义，我们一般用该式作为定义。这里应该注意，,R,n+,1,和整体截断误差,e,n,+1,是不同的,。,若设 ，即第,n,步及以前各步都没有误差，,则由显示单步法计算一步所得之 与 之差为：,即在 的假设下， 。,定义,2,：对单步法，在,的假设下，称,为在 处的,局部截断误差,。,2024/9/15,21,单步法的局部截断误差和阶,（续）,Remark,：由前面的定义可知，若某个单步方法是一种,p,阶方法，则有,R,n+,1,=,O,(,h,p,+1,)，,即,p,阶方法的局部截断误差为,h,的,p,1,阶。我们往往比较关心,R,n,+1,按,h,展开式的第一项。,定义4,：若一个单步方法是一种,p,阶方法，其局部截断误差可以写成：,则,(,x,n,y,(,x,n,),h,p,+1,称为方法的,主局部截断误差,，或,局部截断误差的主项,。,定义,3,：若一个单步方法的局部截断误差为,O,(,h,p,+1,)，,即,则称该方法为,p,阶方法,（其中,p,为正整数）。,2024/9/15,22,例1,：求显示,Euler,公式的局部截断误差。,单步法的局部截断误差和阶,（续）,故,显示,Euler,公式,是一阶方法,局部截断误差为：,主局部截断误差为： 。,2024/9/15,23,例,2,：求,Euler,预测校正,公式,的局部截断误差。,单步法的局部截断误差和阶,（续）,2024/9/15,24,单步法的局部截断误差和阶,（续）,又由,故,故,Euler,预测,校正,方法为二阶方法。,2024/9/15,25,单步法的局部截断误差和阶,（续）,例,3,：求隐式,Euler,公式的局部截断误差。,故隐式,Euler,公式是一阶方法。,将上式中 、 均在 处做,Taylor,展开，整理得,2024/9/15,26,类似地可以证明，梯形公式的局部截断误差为：,单步法的局部截断误差和阶,（续）,即,梯形公式为二阶方法,。,2024/9/15,27,一、,Taylor,方法,假定初值问题的解,y,(,x,),以及函数,f,(,x,y,),是足够光滑的，有,8,.3,Runge,Kutta,公式,Taylor,展开法与数值积分法是推导高阶方法的常用手段。本节以,Taylor,展开法为基础，介绍如何推导高阶单步法。,当,h,充分小时，略去余项 ，将,y,(,k,),(,x,n,),用 来代替，则有,p,阶,Taylor,方法,：,其中,(,k,=1,2,p,),根据求导法则，其计算公式为：,2024/9/15,28,Taylor,方法（续）,Remark1,：显然，,p,阶,Taylor,方法的局部截断误差为,。当,p,=1,时，,Taylor,方法就是,Euler,方法。当,p,2,时，需要计算公式中的高阶导数。,Remark2,：显然，,Taylor,方法可以得到任意阶精度的方法。但在实际计算中，,Taylor,方法往往相当困难，因为,公式中的高阶导数会很复杂。故,Taylor,方法很少单独使用，但常用它来启发思路。,2024/9/15,29,二、,RungeKutta,方法,Euler,公式：,改进的,Euler,预测校正格式：,基本思想,：用不同点的函数值作线性组合，构造高阶单步的近似公式，把近似公式和解的,Taylor,展开式比较，使前面尽可能多的项完全相同。,这种方法间接应用,Taylor,展开的思想，避免了高阶导数计算的困难。,2024/9/15,30,一般的,Runge-Kutta,方法的形式为,其中， 为常数。选取这些常数的原则是使其截断误差阶尽可能高。,RungeKutta,方法（续）,2024/9/15,31,下面以二级,RungeKutta,公式为例进行具体推导。,对,要求适当选取系数 ，使当,时，上式的局部截断误差为 ，即成为二阶方法。,RungeKutta,方法（续）,2024/9/15,32,RungeKutta,方法（续）,按照局部截断误差的定义，有,2024/9/15,33,RungeKutta,方法（续）,2024/9/15,34,RungeKutta,方法（续）,要使上式等于 ，只需满足,在上式中有四个未知量，三个方程，故可以得到无穷多组解，也就是可以得到无穷多个二级二阶,Runge-Kutta,公式由于是二级方法，从而有 ，此时上式可以改写为,在上述方程中，选取不同的,c,2,，即可获得不同的二级二阶,Runge-Kutta,方法。常用的二阶,Runge-Kutta,方法有：,2024/9/15,35,Euler,预测,校正格式,中（间）点公式,RungeKutta,方法（续）,2024/9/15,36,二阶,Heun,（休恩）方法,此时可以论证：,此时 ，即二阶,Heun,方法是局部截断误差项数最少的方法且二级,Runge-Kutta,方法不可能达到三阶。,RungeKutta,方法（续）,2024/9/15,37,RungeKutta,方法（续）,用类似的方法可以研究其它级的,Runge-Kutta,方法。常用的其它级的,Runge-Kutta,方法有：,三,阶,Heun,（,休恩）公式,三,阶,Kutta,公式,2024/9/15,38,标准（经典）四阶,Runge-Kutta,公式,基尔（,Gill,）公式,RungeKutta,方法（续）,2024/9/15,39,定理（整体截断误差与局部截断误差的关系）,：若初值问题的一个单步方法之局部截断误差为,且单步法中函数 关于,y,满足,Lipschitz,条件，则,定义,：对于初值问题，如果一个单步显式方法产生近似解对于任一固定的,均有,则称该单步法是收敛的。,一、收敛性,8.4,单步法的收敛性、相容性与稳定性,Remark,：该定义可类似地应用于单步隐式方法以及后面的线性多步法。从定义可知，若格式收敛，整体截断误差,e,n,=,y,(,x,n,)-,y,n,必然趋于零。,2024/9/15,40,证明,：根据局部截断误差的定义，有,记,则有常数,c,0,，使得,整体截断误差与局部截断误差的关系（续）,由于,且 关于,满足,Lipschitz,条件，故有,2024/9/15,41,整体截断误差与局部截断误差的关系（续）,从而有,递推关系,：,2024/9/15,42,得：,故,当 固定时，,若 则 且,所以,整体截断误差与局部截断误差的关系（续）,证毕,2024/9/15,43,整体截断误差与局部截断误差的关系（续）,Remark2,：在该定理的条件下，欧拉方法是一阶方法，欧拉预测校正方法是二阶方法。当,f,(,x,y,),关于,y,满足利普希茨条件时，,r,级,Runge-Kutta,方法中的,关于,y,也满足利普希茨条件，所以定理中的条件得到满足，解的收敛性得到了保证。,Remark1,：该定理说明，数值方法的整体截断误差比局部截断误差低一阶。收敛的方法至少是一阶方法。要,判断一个单步法的收敛性，一是要知道其局部截断误差，二是要判断其增量函数,是否关于,y,满足,Lipschitz,条件。,2024/9/15,44,整体截断误差与局部截断误差的关系（续）,事实上，对于,Euler,预测校正方法，其,增量函数,故有,从而,Euler,预测校正方法的,增量函数,关于,y,满足,Lipschitz,条件 。,2024/9/15,45,二,.,相容性,故方法阶,p,1,的充要条件是 ，而 ，从而可以给出如下的定义：,由整体截断误差与局部截断误差的关系定理知，当,p,1,时单步方法收敛，且当,y,(,x,),是初值问题的准确解，显式单步方法有,p,阶精度时，若 关于,h,连续，则,其局部截断误差,2024/9/15,46,相容性（续）,定义,：如果单步方法的增量函数,(,x,y,h,),满足,相容的充要条件为,且,(,x,y,h,),关于,h,满足,Lipschitz,条件，则称单步方法,与初值问题,相容,。,Remark,：,以上讨论表明，若,(,x,y,h,),关于,h,满足,Lipschitz,条件，,p,阶方法当,p,1,时与初值问题相容，反之相容的方法至少是一阶方法,Euler,方法，,Euler,预测校正格式，,Runge-Kutta,方法等都与原微分方程相容。,2024/9/15,47,三、,稳定性,定义,1,*,用,一个数值方法求解微分方程初值问题时，对给定的步长,h,0,，,若在计算 时引入误差 （也称扰动），由此引起计算后面的 时的误差按绝对值均不增加，则称这个数值方法是,绝对稳定的,。,关于单步法收敛性的概念和收敛性定理都是在计算过程,中无任何舍入误差的前提条件下建立的。但实际计算时,通常会有舍入误差及其累积。数值求解微分方程的过程,是一个递推公式，必须考虑误差积累是否得到控制。,为保证微分方程本身的稳定性，这里假定 。,设,f,(,x,y,),关于,y,满足,Lipschitz,条件，我们总是针对模型,方程 进行讨论，其中 为复常数。,2024/9/15,48,定义,1,设,步长为,h,0,的单步法用于求解 时，,中由 引起的误差 满足 ，则称,单步法对于所用步长,h,和复数 是,绝对稳定,的。,若在计算 时有误差,但在计算后面的,稳定性（续）,Remark1,：在上面的定义中，可以取小于或等于关系符。取小于号是为了和线性多步法相一致。,Remark2,：从上面的定义可知，单步法是否稳定，与模型方程中的复数,以及所用步长,h,有关。,若对复平面上的某个区域,G,，当 时，,单步法绝对稳定，则称,G,为单步法的,绝对稳定区域,，,G,与实轴的交集为,绝对稳定区间,。,2024/9/15,49,Euler,显式公式,将,Euler,显式公式用于模型方程 ，则有,稳定性（续）,下面考察几个常用公式的稳定性。,设 有误差 ，,参与运算的量为 ，,误差为 ，则实际得到近似 的量为,由此引起的,要求误差不增加，即 ，,，即 。,与前面公式,相减，有 。,必须 。,是保证绝对稳定性对,步长,h,加的限制。当,为实数时，可以得到用,h,表示的绝对稳定的区间,(-2,0),。,2024/9/15,50,隐式,Euler,公式,后退的,Euler,公式用到 上，有 ，,故,即,令,得绝对稳定区域为 。,稳定性（续）,可见，若取,为实数，则对于任意的,h,都,是绝对稳定的。,2024/9/15,51,梯形公式,梯形公式用到 上，有,故,即,设,稳定性（续）,2024/9/15,52,当,x,=Re( )0,上式右端总小于,1,，故梯形法的绝对稳,定区域为,Re( )0,，即左半平面。,由于,稳定性（续）,2024/9/15,53,四阶经典,R-K,方法,四阶经典,R-K,方法用到 上，有,因此,稳定性（续）,2024/9/15,54,扰动满足,令,得四阶经典,R,K,方法的,绝对稳定区域为,稳定性（续）,2024/9/15,55,步的信息来预测 ，则可以期望获得较高的精度。,线性多步法的,基本思想,：如果充分利用前面多,8.5,线性多步法,一、线性多步法的一般公式,记 的近似值为,并记 ，,则,k,步线性多步方法一般形式为,其中 为常数， 不全为零。,2024/9/15,56,若 则为隐式方法，若 则为显式方法。,上式称为线性,k,步法,。,线性多步法的一般理论（续）,或写成,Remark,：,前面的,R,K,方法是增加一些非节点处的函数值的计算来提高单步法的精度的，这样使计算量增加了许多。而,线性多步法每步只需要计算一个函数值。,2024/9/15,57,线性多步法的一般理论（续）,对于隐式情形（ ）的公式，由于,f,(,x,y,),一般是非线性函数，故难以求解到,y,n,+1,的显示表达式，故常用迭代法求解：,其中 任意给出，,s,0,1,2,迭代到满足给定的精度要求为止。,可以证明，当,f,(,x,y,),满足,Lipschitz,条件（或 ）,时，只要 ，迭代关系式就是收敛的。,2024/9/15,58,定义 处的局部截断误差为,定义,：对于线性多步法,线性多步法的一般理论（续）,若,称 为线性多步法的,主局部截断误差,，并称该线性多步法为,p,阶方法。,2024/9/15,59,设 ，则有,线性多步法的一般理论（续）,将上式与线性多步法公式相减，有,其中,n,+1,介于,y,(,x,n,+1,),与,y,n,+1,之间，从而得到,2024/9/15,60,线性多步法的一般理论（续）,若线性多步法为显式方法（ ），则,这与显式单步方法的含义是一致的。如果线性多步法为隐式方法（ ），且该方法为,p,阶方法，则当,y,（,x,）充分可微时,即,R,n,+1,与,y,(,x,n,+1,)-,y,n,+1,的局部截断误差主项相同,。,2024/9/15,61,线性多步法的一般理论（续）,Remark,：可以证明，显式线性多步法的整体截断误差比局部截断误差低一阶。,2024/9/15,62,二、基于,数值积分的构造方法,将 方程两端从 积分得,（,*,）,记 ，利用,Lagrang,插值多项式来近似代替（,*,）式中的被积函数 ，则可以导出不同的线性多步法公式。下面导出常用的四阶,Adams,外插、内插公式。,2024/9/15,63,其插值余项为,基于,数值积分的构造方法（续）,在（,*,）式中取,k,=0,，并选择等距节点 作为插值节点，作函数,F,（,x,）的三次插值多项式,将 代入（,*,）式，得,2024/9/15,64,基于,数值积分的构造方法（续）,略去上式右端第三项，得,对上式积分部分作变量代换 ，并注意到,则,2024/9/15,65,线性四步,Adams,显式公式,基于,数值积分的构造方法（续）,其局部截断误差为,2024/9/15,66,基于,数值积分的构造方法（续）,利用第二积分中值定理，得,由于插值多项式 是在 上作出的，而积分区间是 ，故上述线性四步,Adams,显式公式称为,Adams,外插公式,。,2024/9/15,67,基于,数值积分的构造方法（续）,若在（,*,）式中取,k,=0,，但选择等距节点 作为插值节点，作函数,F,（,x,）的三次插值多项式，,仿照上面的做法，可以得到线性三步的,Adams,隐式公式及其局部截断误差,上式称为,Adams,内插公式,。,由于,Adams,内插公式是隐式公式，故用它计算时也需用迭代法。通常把,Adams,外插公式与内插公式联合起来交替使用，先由前者提供初始值，再由后者进行修正，即,2024/9/15,68,基于,数值积分的构造方法（续）,可以看出，无论四步,Adams,显式公式还是三步,Adams,隐式公式，均为四阶公式。类似于以上讨论，还可以采用数值积分方法得到其它,p,阶的线性多步法公式。,可以证明，当 时，迭代式收敛若上式中的第二式只迭代一次，则可以得到,Adams,预测校正公式,2024/9/15,69,三、,基于,Taylor,展开的构造方法,Taylor,展开法更具一般性。下面构造线性两步方法,利用线性多步法局部截断误差的定义，有,将上式在,x,n,处作,Taylor,展开，并按,h,的升幂整理排列，得到,2024/9/15,70,基于,Taylor,展开的构造方法（续）,2024/9/15,71,基于,Taylor,展开的构造方法（续）,令,求解上述方程组，得出 ，所得到的公式的局部截断误差至少为 。,2024/9/15,72,基于,Taylor,展开的构造方法（续）,可以只要求前面的几个方程成立，例如要求前面的四个方程成立时，所得公式的局部截断误差至少为 。由于此时方程个数少于未知量个数，故此种情形下方程组有无穷多组解。此时方程的解可以写成,若取 ，则 ，从而得到一个线性两步隐式公式,其局部截断误差为,上式称为,Simpson,公式,。该公式为四阶方法，它也可以由数值积分方法得到。,2024/9/15,73,基于,Taylor,展开的构造方法（续）,若取 ，则 ，从而得到一个线性两步显式公式,其局部截断误差为,需要指出的是，上式不能用数值积分方法得到,。,可以按照上述方法构造其它的计算格式，也可以用类似的方法构造其它的线性多步方法，如前述的,Adams,公式等。,2024/9/15,74,基于,Taylor,展开的构造方法（续）,Remark1,：用,Taylor,展开法比用数值积分法推导线性多步法更加灵活，推导出的方法也比积分法的更多。用数值积分法可以得到的方法，用,Taylor,展开法都可以得到，但用,Taylor,展开法得到的方法，用数值积分法却不一定能够得到。,Remark2,：应用线性多步法求解初值问题时，开头几点处的数值要用别的方法先计算出来。一般选用与多步法同阶或更高阶的单步法，如,Runge-Kutta,法、,Taylor,展开法等。,Remark3,：对线性隐式多步方法，除开头几点的函数值需要单独计算外，还需要迭代求解或采用预测,-,校正方法求解。,2024/9/15,75,数值算例,求解常微分方程初值问题,在,0,1,上的解，取步长,h,0.1,。,计算结果如下图所示：,该问题的理论解为：,2024/9/15,76,2024/9/15,77,2024/9/15,78,