cy52常微分方程初值问题(二).docx

cy52常微分方程初值问题(二) cy52常微分方程初值问题(二) 第五章常微分方程初值问题数值解法(二)第五节线性多步法 第六节 方程组和高阶方程 上页下页返回 cy52常微分方程初值问题(二) §5线性多步法 当b-10时，为隐式公式;b-1=0则为显式公式.用若干节点处的y及y值的线性组合来近似fj=f(xj,yj) y(xi+1).其通式可写为： yi+1=a0yi+a1yi-1+.+akyi-k+h(b-1fi+1+b0fi+b1fi-1+.+bkfi-k) 基于数值积分的构造法将y=f(x,y)在上积分，得到 y(xi+1)-y(xi)= xi+1xi f(x,y(x)dxxi+1i 只要近似地算出右边的积分Ikf(x,y(x)dx,则可通x过yi+1=yi+Ik近似y(xi+1).而选用不同近似式Ik,可得到不同的计算公式. 上页下页返回 cy52常微分方程初值问题(二) 亚当姆斯显式公式利用k+1个节点上的被积函数值fi,fi-1,.,fi-k构造k阶牛顿后插多项式Nk(xi+th),t,有 xi+1xi f(x,y(x)dx=Nk(xi+th)hdt+Rk(xi+th)hdt00 1 1 yi+1=yi+hNk(xi+th)dt0 1 /*显式计算公式*/Newton1 局部截断误差为：Ri=y(xi+1)-yi+1=hRk(xi+th)dt0 插值余项 例：k=1时有N1(xi+th)=fi+tfi=fi+t(fi-fi-1)hyi+1=yi+hdt=yi+(3fi-fi-1)0221df(,y()153xxRi=hth(t+1)hdt=hy(i)201 上页下页返回 dx 2! 12 51238251720 常用的是k=3的4阶亚当姆斯显式公式yi+1h=yi+(55fi-59fi-1+37fi-2-9fi-3)24上页下页返回 cy52常微分方程初值问题(二) 亚当姆斯隐式公式利用k+1个节点上的被积函数值fi+1,fi,fi-k+1构造k阶牛顿后插多项式.与显式多项式完全类似地可得到一系列隐k+2(k+2)式公式，并有Ri=Bkhy(i),其中Bk与fi+1,fi,fi-k+1的系数亦可查表得到.k01fi+1112512924128121924 fi fi-1 fi-2 Bk12 小于Bk112524124 - 112 23 12419720 常用的是k=3的4阶亚当姆斯隐式公式yi+1=yi+h(9fi+1+19fi-5fi-1+fi-2)24 较同阶显式稳定上页下页返回 cy52常微分方程初值问题(二) 亚当姆斯预估-校正系统Predicted留意：三步所用公valuepi+1Step1:用Runge-Kutta法计算前k个初值；式的精度必需相同。通常用经典RungeStep2:用Adams显式计算预估值；Kutta法协作4阶Adams公式。Step3:用同阶Adams隐式计算校正值。2515(5)hy(i)4阶Adams显式公式的截断误差为7 20XX5(5)y(x)y=hy(i)4阶Adams隐式公式的截断误差为i+1i+1720y(x)-y251当h充分小时，可近似认为ii,则：y(xi+1)-yi+1-19i+1i+1y(xi+1)-yi+1= 251(yi+1-yi+1)27019y(xi+1)yi+1(yi+1-yi+1)270y(xi+1)yi+1+ 外推技术ModifiedfinalCorrectedModified/*extrapolation*/valuevalueyci+1valuemi+1i+1 上页下页返回 cy52常微分方程初值问题(二) 基于泰勒绽开的构造法yi+1=a0yi+a1yi-1+.+akyi-k+h(b-1fi+1+b0fi+b1fi-1+.+bkfi-k) 将通式中的右端各项yi-1,yi-k;fi+1,fi-1,fi-k分别在xi点作泰勒绽开，与精确解y(xi+1)在xi点的泰勒绽开作比较.通过令同类项系数相等，得到足以确定待定系数a0,ak;b-1,b0,bk的等式，则可构造出线性多步法的公式.上页下页返回 cy52常微分方程初值问题(二) 例：设yi+1=a0yi+a1yi-1+a2yi-2+h(b0y+b1y+b2y+b3yii-1i-2i-3)确定式中待定系数a0,a1,a2,b0,b1,b2,b3,使得公式具有4阶精度.211314(4)+O(h5)i+2hyi-6hyi+24hyi解：yi-1=yi-hy24324(4)+O(h5)yi-2=yi-2hyi+2hyi-3hyi+3hyi 1213(4)+O(h4)yi-1=yi-hyi+2hyi-6hyi243(4)+O(h4)yi-2=yi-2hyi+2hyi-3hyi 9293(4)+O(h4)yi-3=yi-3hyi+2hyi-2hyi1213y(xi+1)=yi+hyi+2hyi+6hyi+124 h4yi(4)+O(h5)/*y(xi)=yi*/ a0+a1+a2=1h(-a1-2a2+b0+b1+b2)=h7个未知数21h2(1a+2ab2b3b)=h12123225个方程94113h3(-1aa+b+2b+b)=h上页1212363226491141下页h4(24a1+2abbb)=h2123363224返回 cy52常微分方程初值问题(二) 令a1=a2=0Adams显式公式以yi+1取代yi-1,并取a1=a2=0Adams隐式公式以yi-3取代yi-3,则可导出另一组4阶显式算法，其中包含了闻名的米尔尼/*Milne*/公式4h取a1=1,a2=0yi+1=yi-3+(2yyii-1+2yi-2)3得到辛甫生公式与Milne公式匹配使用145(5)=hy(i),i(xi,xi+1)其局部截断误差为在区间上积分，并用-1ii+145Simpson数值积分公式来近似积y=f(x,y)在区间分项，亦可得此Simpson公式。注：上式也可通过数值积分导出，即将xi+1上积分，得到y(xi+1)=y(xi-3)+f(x,y(x)dx,xi-3f,f,f再过ii-1i-2做f的插值多项式即可. 辛甫生公式 hyi+1=yi-1+(yi+1+4yi+yi-1)3 上页下页返回 cy52常微分方程初值问题(二) Milne-Simpson系统的缺点是稳定性差，为改善稳定性，考虑 另一种隐式校正公式：yi+1=a0yi+a1yi-1+a2yi-2+h(b-1yi+1+b0yi+b1yi-1) 要求公式具有4阶精度。通过泰勒绽开，可得到5个等式，从中解出6个未知数，则有1个自由度.哈明用a1的不同数值进行试验，发觉当a1=0时，公式的取a1=1得稳定性较好，即：Simpson公式13hyi+1=(9yi-yi-2)+(yi+1+2yi-yi-1)88 15(5)其局部截断误差为Ri=-hy(i),i(xi,xi+1)40 注：哈明公式不能用数值积分方法推导出来. 上页下页返回 cy52常微分方程初值问题(二) §6微分方程组与高阶方程一阶微分方程组(x)=f1(x,y1(x),.,ym(x)y1.=fm(x,y1(x),.,ym(x)ym(x)000初值y1(x0)=y1,y2(x0)=y2,.,ym(x0)=ym IVP的一般形式为： 0y1y1f1.,f=.,y0=.将问题记作向量形式，令：y=.yfy0mmm y(x)=f(x,y)y(x0)=y0 前述全部公式皆适用于向量形式. 上页下页返回 cy52常微分方程初值问题(二) 高阶微分方程y(n)=f(x,y,y,.,y(n1)y(x0)=a0,y(x0)=a1,.,y(n1)(x0)=an-1 化作一阶微分方程组求解.引入新变量y1=y,y2=y,.,yn=y(n-1)y1yn-1ny=y2.=yn=f(x,y1,.,yn) 初值条件为：y1(x0)=a0y2(x0)=a1.yn(x0)=an-1上页下页返回 cy52常微分方程初值问题(二) §6边值问题的数值解2阶常微分方程边值问题x(a,b)y=f(x,y,y)y(a)=a,y(b)=b 打靶法先猜测一个初始斜率y(a)=s,通过解初值问题y=f(x,y,y)y(a)=ay(a)=sy 每计算一个(s)都必需解一个ODE.(s0) 斜率=s0y(x) b(s1) y(b)=(s)斜率=s1 找出s*使得(s*)=b,即把问题转化为求方程(s)-b=0的根。 上页0ab 下页x返回 cy52常微分方程初值问题(二) 有限差分法将求解区间等分为N份，取节点xi=a+ih(i=0,N),在每一个节点处将y和y离散化。 泰勒绽开y(x+h)-y(x)y(x)-y(x-h)2hhhy(x)=y(4)()h12y(x+h)-2y(x)+y(x-h)2=+O(h)2hy(x)=y(x+h)-y(x-h)+O(h2)2h yi+1-yi-1yi+1-2yi+yi-1=f(xi,yi,)2h2hy0=a,yN=b i=1,.,N-1 上页下页返回 - 9 -