第五章 常微分方程数值解.pptVIP

§5.3 初始问题的线性多步法 线性多步法的好处就是：少量的函数值计算可得到高阶收敛的方法！！！ 一、一些常用的多步法 1、Adams方法 Adams显式方法 绝对稳定区间 （-1，0） （-6/11，0） （-3/10，0） * 第五章 常微分方程数值解 §5.1　问题的提出 本章只介绍初始问题的解法，原因是初始问题系统科学中研究的最重要的问题。 §5.2 初始问题的单步法 以上说明显式欧拉法是一阶收敛的算法 2、隐式欧拉法 二、龙格-库塔法（RK) 我们的目的是构造更高阶的单步法！ 称此结构为m+1级龙格-库塔法。它最多是m+1阶的。 一些常用的显式RK方法 1、二级二阶RK法 改进的欧拉法（p277),中点公式（p278),休恩方法（p278)等。 2、三级三阶RK法 休恩三阶方法（P279),库塔三阶方法（P279)等 3、四级四阶RK法 经典四级四阶RK法（P280),基尔四级四阶RK法(P280) 三、单步法的数值稳定性 这样做的根据是：1）对实验模型数值不稳定的方法，不可用； 2）一般的初始问题在其解的存在区域内，可局部线性化转化为试验方程。 几个常用算法的数值稳定性 实例分析 *分别用显式Euler法、改进的Euler法、经典4级4阶RK法计算* 考虑数值稳定性对步长h的限制： 四、变步长和误差控制算法 从前面的分析可知，求解常微分方程初始问题步长的选择是关键。然而，事先确定一个既保证收敛又数值稳定的步长很困难，另外方程的解可能在求解区间的某些部分变化平缓，可用大步长，而在其他部分变化剧烈，要用小步长。所以讨论变步长的算法是必要的。 对于一个选定的算法，我们可以事先分析确定保证数值稳定的步长h的一个可选范围[h_min,h_max]。接下来，我们 讨论变步长的误差控制算法。 function [t,x]=weifenfchqiujie(f_name,t0,x0,h_min,h_max,epsl,t1) %解dx/dt=f(x,t),x(t0)=x0,终点t1.用变步长的算法。 t(1)=t0; x(1)=x0;%初始值 h=(h_max+h_min)/2;%h选在绝对稳定区间内 n=1; while 1 if t(n)-t1==0,break,end; if t1-t(n)h h=t1-t(n); end k1=feval(f_name,x(n),t(n)); k2=feval(f_name,x(n)+h/4*k1,t(n)+h/4); k3=feval(f_name,x(n)+3/32*h*k1+9/32*h*k2,t(n)+3/8*h); k4=feval(f_name,x(n)+1932/2197*h*k1- 7200/2197*h*k2+7296/2197*h*k3,t(n)+12/13*h); k5=feval(f_name,x(n)+439/216*h*k1-8*h*k2+3680/513*h*k3-845/4104*h*k4,t(n)+h); k6=feval(f_name,x(n)-8/27*h*k1+2*h*k2-3544/2565*h*k3+1859/4104*h*k4-11/40*h*k5,t(n)+1/2*h); x(n+1)=x(n)+h*(16/135*k1+6656/12825*k3+28561/56430*k4-9/50*k5+2/55*k6); %5阶RK-Fehlberg法 x1(n+1)=x(n)+h*(25/216*k1+1408/2565*k3+2197/4104*k4-1/5*k5);%4阶RK-Fehlberg法 %d=h*(1/360*k1-128/4275*k3-2197/75240*k4+1/50*k5+2/55*k6);%一步误差 d=x1(n+1)-x(n+1); if (abs(d)/h)epsl seita=0.84*sqrt(sqrt(epsl*h/abs(d)));%调整步长参数 h=seita*h; else t(n+1)=t(n)+h; n=n+1; end end 实例1 function f=ffff(x,t) f=-50*x; [t,x]=weifenfchqiujie(ffff,0,1,0.1,0.5,0.00001,1); plot(t,x) 实例2 function