常微分方程数值解法 2

第六章,常微分方程数值解法,本章研究常微分方程初值问题(6.1.1)的数值解法.并且假定f(x,y)满足解的存在唯一性定理及相当光滑等条件.初值问题(6.1.1)的精确解记为y(x).,建立数值解法的基本思想,本章的数值解法，它不是求(6.1.1)的解y(x)的解析表达式或近似表达式，而是通过某种离散化方法，将连续变量的初值问题转化为关于离散量的差分方程的初值问题来求一系列离散点上的解值y(xi)的近似值yi.利用计算机解微分方程主要使用数值方法.,取一系列点x0,x1,xn,y(x0)=y0,y(x1)y1,y(xn)yn,y0,y1,yn,称为数值解用离散化方法建立求y(xn)的近似值yn的递推格式(差分方程)，求得ynh=xnxn-1称为步长本章都取定步长.,解初值问题(6.1.1)的数值解法，其特点都是采取步进式的方法，即求解过程顺着节点排列的次序一步一步向前推进.这种数值解法分为两大类：(1)单步法：若求yn+1，只需利用它前一步的信息yn，则称这种方法为单步法。它由y0出发，可求得y1，y2，y3(2)多步法：若求yn+1，需利用它前面至少两个点的息，则称这种方法为多步法.,数值解法研究的主要问题：,(1)方法的推导：采用的离散化手段，精度准则.(2)收敛性：差分方程的解是否充分逼近初值问题的解.(3)稳定性：初始数据、计算过程中每步产生的误差对以后各步解的影响，这种误差传播是否在衰减.,具体的数值方法还应考虑(1)误差估计(2)解的起动方法(3)步长如何选取(4)隐式方法的如何计算,建立数值解法的基本途径,常用的离散化方法(1)Taylor展开(2)化导数为差商(3)数值积分,6.1单步法及基本概念,Euler折线法,利用Taylor展开法将y(xn+1)在xn处Taylor展开得差分方程yn+1=yn+hf(xn,yn)（n=0,1,2,）称此方法为Euler折线法或矩形法,利用化导数为差商的方法得差分方程yn+1=yn+hf(xn,yn)（n=0,1,2,）(6.1.2),利用数值积分的方法在xn,xn+1上对y(x)=f(x,y(x)积分得用左矩形求积公式计算定积分有y(xn+1)y(xn)+hf(xn,y(xn)以此得差分方程yn+1=yn+hf(xn,yn)（n=0,1,2,）,梯形法,用不同的近似公式计算定积分的值,就得到解初值问题的不同数值解法.用梯形求积公式计算积分得(6.1.3)这个方法称为梯形法.它是隐格式.运用它常采用下面的迭代格式,改进的Euler方法,若梯形法只迭代一次,便得改进的Euler方法(6.1.4)也可写成形式(6.1.4),数值方法精度的衡量准则,定义设yn=y(xn)，则称Tn=y(xn+1)-yn+1为方法的从xn到xn+1这一步的局部截断误差.定义若差分方程对所有y(x)Mr都精确成立，而对于某个r+1次多项式不能精确成立，则称这个数值方法是r阶的.等价定义若数值方法的局部截断误差为O(hr+1)，则称这种方法为r阶方法，这里r为非负整数.方法的阶数越高，逼近效果越好.,由Taylor展开式知(6.1.5)对于Euler法Tn=y(xn+1)-y(xn)-hf(xn,y(xn)=y(xn+1)-y(xn)-hy(xn)Euler法是1阶方法,对于梯形法由Taylor展开式知(6.1.6)(6.1.5)与(6.1.6)相减并利用y(xn+1)=y(xn)+O(h)梯形法法是2阶方法.,对于改进Euler法，设yn=y(xn)，利用Taylor展开式K1=f(xn,y(xn)=y(xn),将其代入中Tn有改进Euler法是2阶方法.,Runge-Kutta方法,Taylor展开式理论上讲，只要解y(x)充分光滑，通过保留Taylor展开式的若干项就可得到任意阶的近似公式，但计算y(x)的各阶导数很麻烦。可间接利用这种思想,Euler法也可写成形式其局部截断误差为O(h2)，是一阶方法.每步计算f的值一次,Euler预报-校正公式也可写成形式局部截断误差为O(h3)，是二阶方法每步计算f的值二次可以通过增加计算f的值的次数，提高公式的阶数（精度）,以f在不同点上的函数值的线性组合来代替yn+1yn，其中有一些可待定选取的待定参数，通过Taylor展开确定这些待定参数使建立的数值方法按要求达到一定的阶数，这种思想就是Runge-Kutta方法的思想,Runge-Kutta法的一般形式其中Ri,ai,bij都是常数，a0=0，b1j=0，(j=1,2,i-1),二阶Runge-Kutta法的一般形式其中R1,R2,a,b为待定常数其局部截断误差为O(h3)，是二阶方法每步计算f的值二次,设yn=y(xn)K1=f(xn,yn)=f(xn,y(xn)把K2中f在(xn,y(xn)处Taylor展开,再将K1,K2代入yn+1中，将其与y(xn+1)泰勒展开式比较，要使y(xn+1)-yn+1=O(h3)，含h0,h1,h2的项相同即有,个未知数，个方程满足条件的解不止一组取R1=R2=1/2,a=b=1,就是改进Euler法；取R1=0,R2=1,a=1/2,b=1/2的中点方法.,四阶Runge-Kutta法的一般形式其中有13个待定常数局部截断误差为O(h)，是四阶方法每步计算f的值四次,设yn=y(xn)把K2,K3,K4中f在(xn,y(xn)处泰勒展开后，将K1,K2,K3,K4代入yn+1中，再将yn+1按h的幂重新整理后与y(xn+1)泰勒展开式比较，要使y(xn+1)-yn+1=O(h5)，含h0,h1,h2,h3,h4的项相同从而确定各参数,13个未知数，11个方程满足条件的解不止一组最常用的是标准四阶Runge-Kutta公式(6.1.7),单步法的一般形式为yn+1=yn+h(xn,yn,h)(6.1.8)其中多元函数依赖于f.,单步法的收敛性,en=y(xn)yn称为整体截断误差.收敛性就是讨论当x=xn固定且h=(xn-x0)/n0时en0的问题.定义如果对f(x,y)满足解的存在唯一性条件的初值问题(6.1.1)，差分方程(6.1.8)的解对每个确定xa,b满足则称单步法(6.1.8)是收敛的.,定义若增量函数(x,y,h)使(x,y,0)=f(x,y)成立，则称单步法(6.1.8)与(6.1.1)相容.满足相容条件的方法至少是1阶的.定理设增量函数(x,y,h)在区域axb,-y,0hh0中连续，并且对变量y满足Lipschitz条件,即(x,y1,h)-(x,y2,h)Ly1-y2在这个前提下单步法(6.1.8)收敛的充分必要条件是相容性条件成立.,单步法的绝对稳定性,当步长取定后，计算中的误差随着步数的增加会不会积累到超出我们许可的范围，这就是稳定性问题.单步法(6.1.8)应用于模型方程y=y(0)设得到的解为yn+1=R(h)yn当0时,实验方程的精确解y(x)=e(x-a)按模递减的,这要求(6.1.8)的解yn是递减的，误差也是递减的，即要求满足R(h)1,定义单步法(6.1.8)应用于模型方程y=y(0)，若得到的解yn+1=R(h)yn满足R(h)1,则称单步法(6.1.8)是绝对稳定的.使R(h)1,(6.3.3)的解当n时无界；误差也是无界的.这种情况下，如果误差相对于真解是小的，就说方法是相对稳定的.当h充分小,初始误差较小时,有k0d0,kik0这时要求ri(h)r0(h)，i=1,p.,当0时,实验方程的精确解y(x)=e(x-a)按模递减的;要求(6.3.3)的解当n时递减到0，误差也是递减的。这种情况下，就需要ri(h)1,i=0,1,p来保证.这时就称方法是绝对稳定的.方法(6.2.1)的稳定性取决于特征多项式(r;h)零点的性质，故又称(r;h)为(6.2.1)的稳定多项式.,稳定性的定义,记=h定义设方法(6.2.1)是收敛的.ri()(i=0,1,p)是稳定多项式(r;)=0的根,r0()是满足(6.3.5)的根.(1)若对任意的,R有ri()r0()，i=1,p.且当ri()=r0()时,ri()是单根,则称方法(6.2.1)在,上是相对稳定的，,称为此方法的相对稳定区间.(2)若对任意的(,)R有ri()0时,相对稳定性重要；当0时,绝对稳定重要，才有意义.定义若一个方法的绝对稳定区间是(-,0),则称该方法是A稳定的.,(1)收敛的方法的相对稳定区间不空。因为方法稳定的前提是收敛，故满足根条件和相容条件,0,.(2)只有0的情况讨论绝对收敛性才有意义.(3)从误差分析的角度看，绝对稳定的方法是理想的.(3)相对稳定区间和绝对稳定区间越大越好，可使尽可能大的一类微分方程是稳定的.,例讨论梯形方法的相对和绝对稳定区间.,解梯形方法为对于试验方程y=y，梯形方法成为即其稳定多形式为它只有一个根，记为且r0(0)=r0=1,由于由于其只有一个根，所以梯形方法的相对稳定区间是(-，+)当0时,r0()1当0时,r0()1这样，梯形方法的绝对稳定区间是(-，0).即梯形方法是A稳定的.,例讨论Simpson方法的相对和绝对稳定区间,Simpson方法为利用该格式求解试验方程y=y得其稳定多形式为它有两个根,考察当0时，上式小于等于1；当0时,当0时,因此，Simpson方法不存在绝对稳定区间.,例线性多步法（1）确定法方法中的局部截断误差主项，并指出方法的阶数（2）讨论该方法的收敛性和绝对稳定性,解p=1,a0=1/2,a1=1/2,b-1=0,b0=7/4,b1=-1/4,步长的选取:稳定性与步长有关，步长的选取一定要保证方法是稳定的，即h属于稳定区间.收敛性是反映递推公式本身的整体截断误差对计算结果的影响；稳定性反映某一计算步骤中出现的误差对计算结果的影。,6.5预测-校正方法,线性多步法(6.2.1)当b-10时，是隐格式估计出yn+1的一个初值y(0)n+1，利用迭代y(j)n+1是yn+1的第j次近似值.当hb-11时，此迭代收敛.,用显式方法来作预测值y(0)n+1，用隐式方法迭代校正一次得y(1)n+1值，这种显式与隐式联合使用构成的方法称为预报-校正法.作为预报的显式公式称为预测式，用于校正的隐式公式称为校正式.例如，改进的Euler方法.,预测公式与校正公式选取同阶方法，可使截断误差用预测值和校正值表示，再用截断误差来修正，得到提高精度的方