《数值分析常微分方程数值解法.pptx》由会员分享,可在线阅读,更多相关《数值分析常微分方程数值解法.pptx(105页珍藏版)》请在taowenge.com淘文阁网|工程机械CAD图纸|机械工程制图|CAD装配图下载|SolidWorks_CaTia_CAD_UG_PROE_设计图分享下载上搜索。
1、1近似求解第1页/共105页2一阶常微分方程初值问题为保证上述问题解的存在唯一性,对函数 f 做如下假设:f 关于y满足Lipschitz条件,即存在常数 L0,使得当(L为正常数)时,f 关于y满足Lipschitz条件.常微分方程解的存在唯一性第2页/共105页3数值解法:就是寻求解 y(x)在一系列离散点上的近似值相邻两个节点的间距 称为步长.本章均假设此时节点 常微分方程数值解法第3页/共105页4 初值问题的各种差分方法都采用“步进式”,即求解过程顺着节点排列的次序一步一步地向前推进.只要给出从已知信息 计算 的递推公式,此类计算格式统称为差分格式.常微分方程差分格式第4页/共105
2、页5数值求解一阶常微分方程初值问题难点:如何离散 y?常见离散方法 差商近似导数 数值积分方法 Taylor展开方法第5页/共105页6结果有设用 y(xn)的近似值 yn代入上式右端,记所求结果为yn+1,这样导出的计算公式将方程 y=f(x,y)中的 y在 xn 点 用向前差商近似 差商近似导数(Euler公式)Euler公式,单步法,显式公式第6页/共105页7 Euler公式计算总结设的精确解为将区间a,b N等分,分点 xn=x0+nh(n=1,2,N)步长用Euler公式求得 yn(n=1,2,N.)第7页/共105页8结果有由此得差分格式计算公式将方程 y=f(x,y)中的 y在
3、 xn+1 点用向后差商近似向后Euler公式,单步法,隐式公式 差商近似导数(向后Euler公式)第8页/共105页9设将方程 y=f(x,y)的两端从 xn 到xn+1 求积分,得用不同的数值积分方法近似上式右端积分,可以得到计算 y(xn+1)的不同的差分格式.用左矩形公式计算积分项再离散化,得Euler公式 数值积分方法(Euler公式)第9页/共105页10用右矩形公式计算积分项再离散化,得向后Euler公式 数值积分方法(向后Euler公式)第10页/共105页11用梯形法计算积分项再离散化,即可得如下计算公式梯形公式,单步法,隐式公式 数值积分方法(梯形公式)第11页/共105页
4、12 Taylor展开方法(Euler公式)由此亦可得Euler公式将 y(xn+1)在 xn 点作Taylor展开,得第12页/共105页131 Euler方法用差分格式来近似微分方程初值问题第13页/共105页14例 用Euler公式求解初值问题解取步长h=0.1,计算结果见下表计算公式为第14页/共105页15 y(xn)1.09541.18321.26491.34161.41421.48321.54921.61251.67331.7321xn0.10.20.30.40.50.60.70.80.91.0 yn1.10001.19181.27741.35821.43511.50901.58
5、031.64981.71781.7848 两者比较可以看出Euler法的精度很差.第15页/共105页16h h h Straight line approximationEuler法的几何意义x0 x1x2x3y y0第16页/共105页17 向后Euler法这是一个隐式公式.已知 yn,上式是关于yn+1的非线性方程,不能从中直接求出yn+1,需要用非线性方程求根的方法来求解,通常用迭代法.第17页/共105页18 已知 yn 求yn+1的迭代法则当hL1时,迭代函数 向后Euler法计算量大,但数值稳定性好!第18页/共105页19 Euler法的局部截断误差局部截断误差:用Euler法
6、计算一步所产生的误差,即假设第n步的计算是精确的 yn=y(xn),y(xn+1)与 yn+1 的差,即这种误差称为局部截断误差,简称截断误差.第19页/共105页20当 yn=y(xn)时,局部截断误差主项 Euler法的局部截断误差第20页/共105页21 向后Euler法的局部截断误差定义其局部截断误差为向后Euler法的计算公式第21页/共105页22 向后Euler法的局部截断误差局部截断误差主项第22页/共105页23在不考虑舍入误差的情况下,称 y(xn+1)与 yn+1之差为整体截断误差,记为令则 Euler法的整体截断误差第23页/共105页24 Euler法的整体截断误差第
7、24页/共105页25注意e0=0,故第25页/共105页26所以 这表明,Euler方法的整体截断误差与h同阶.当h0时,eN 0.如果某种数值方法的局部截断误差为 O(hp+1),则称它的精度是 p阶的,或称之为 p 阶方法.由此可知Euler方法仅为一阶方法.第26页/共105页27 梯形公式2 改进的Euler法从直观上看,用梯形公式计算数值积分要比矩形公式精度高.相应地求解常微分方程的梯形公式的精度也要比Euler公式精度高.第27页/共105页28 梯形公式的局部截断误差局部截断误差主项 故梯形公式是二阶方法第28页/共105页29则当hL2时,梯形公式是一个隐式公式.已知 yn,
8、用迭代法计算yn+1 用梯形法计算,从 yn 到yn+1要进行迭代计算,计算量较大.实际计算常把Euler法和梯形法结合起来,这就是改进的Euler法.第29页/共105页30 先用Euler法求得一个初步的近似值,记为 ,称之为预测值,然后用它替代梯形法右端的yn+1 再直接计算 fn+1,得到校正值 yn+1,这样建立的预测校正系统称为改进的Euler法:预测校正 改进的Euler法第30页/共105页31 实践表明,改进Euler法的精度明显高于Euler法.利用Taylor展开可以证明,改进Euler的局部截断误差为O(h3),故它是二阶方法.它有下列平均化形式:第31页/共105页3
9、2例 用改进Euler法求解解取步长h=0.1,计算结果见下表改进Euler法公式为第32页/共105页33 y(xn)1.09541.18321.26491.34161.41421.48321.54921.61251.67331.7321xn0.10.20.30.40.50.60.70.80.91.0 yn1.09591.18411.26621.34341.41641.48601.55251.61651.67821.7379 同Euler法的计算结果比较,改进Euler法明显改善了精度.第33页/共105页34 改进Euler法的几何意义(K1+K2)/2K2K1xnxn+1第34页/共10
10、5页35考察差商 根据微分中值定理,存在点,利用所给方程 y=f(x,y)得称 K*=f(,y()为区间xn,xn+1上的平均斜率,这样只要对平均斜率 K*提供一种算法,相应地便导出一种计算格式.3 Runge-Kutta法 Runge-Kutta法的基本思想第35页/共105页36梯形法:改进Euler法:其中Euler法:向后Euler法:第36页/共105页37RungeKutta法的基本思想就是设法在xn,xn+1内多预报几个点的斜率值,然后把它们加权平均作为平均斜率,以期望构造出更高精度的计算格式.Runge-Kutta法的基本思想第37页/共105页38考察区间xn,xn+1内一点
11、xn+p=xn+ph,0 n)上产生的偏差值均不超过,则称该方法是绝对稳定的.用某种数值方法求解上述模型方程,对固定的 h,使得该方法绝对稳定的h和 的全体称为该方法的绝对稳定区域.第92页/共105页93 Euler方法的绝对稳定区域模型问题计算格式设实际计算yn时产生误差n,它使节点值 yn+1产生误差n+1,则从而为使误差不扩大,需第93页/共105页94 向后Euler方法的绝对稳定区域模型问题计算格式设实际计算yn时产生误差n,它使节点值yn+1产生误差n+1,则向后Euler方法的绝对稳定区域第94页/共105页95 梯形公式的绝对稳定区域模型问题计算格式第95页/共105页96类
12、似于Euler方法的讨论,得梯形法的绝对稳定区域为即 梯形公式的绝对稳定区域第96页/共105页97 RK方法的绝对稳定区域 二阶RK方法的绝对稳定区域对于模型问题改进Euler法第97页/共105页98当 是实数时,2.51 h 0当 是实数时,2.78 h 0 二阶RK方法的绝对稳定区域当 是实数时,2 h 0 三阶和四阶RK方法的绝对稳定区域分别为第98页/共105页99例 y=20y (0 x 1),y(0)=1,分别取 h=0.1及 h=0.2 用经典的四阶RK方法计算.h=0.24.9825.0125.0625.03125.0 xn0.20.40.60.81.0h=0.10.931
13、010.121010.141020.151030.17104解用四解RK方法计算,其误差见下表 本例中=20,h分别为2及4,前者在绝对稳定区间内,后者则不在.从以上结果看到,如果步长h不满足绝对稳定条件,误差增长很快.第99页/共105页100前面研究了单个方程 y=f 的数值解法,只要把 y 和 f 理解为向量,所提供的各种算法即可推广应用到一阶方程组的情形.6 一阶微分方程组与高阶微分方程的数值解法 一阶方程组的数值解法第100页/共105页101以 yn,zn表示节点 xn 上 y(x),z(x)的近似解,则其改进的Euler法具有形式例 对于方程组令原方程组变为第101页/共105页102预报即第102页/共105页103校正即第103页/共105页104关于高阶微分方程(或方程组)的初值问题,原则上总可以归结为一阶微分方程组来求解.引进新变量 z=y,即可化为一阶方程组的初值问题从而可以运用前面介绍的算法来求解.高阶微分方程的数值解法例 对于下列二阶方程的初值问题第104页/共105页105感谢您的观看。第105页/共105页