资源描述：
《matlab在科学计算中的应用ppt课件.ppt》由会员上传分享，免费在线阅读，更多相关内容在教育资源-天天文库。
1、5.2.4微分方程转换5.2.4.1单个高阶常微分方程处理方法12例：函数描述为：functiony=vdp_eq(t,x,flag,mu)y=[x(2);-mu*(x(1)^2-1)*x(2)-x(1)];>>x0=[-0.2;-0.7];t_final=20;>>mu=1;[t1,y1]=ode45('vdp_eq',[0,t_final],x0,[],mu);>>mu=2;[t2,y2]=ode45('vdp_eq',[0,t_final],x0,[],mu);>>plot(t1,y1,t2,y2,':')>>figure;plot(y1(:
2、,1),y1(:,2),y2(:,1),y2(:,2),':')3>>x0=[2;0];t_final=3000;>>mu=1000;[t,y]=ode45('vdp_eq',[0,t_final],x0,[],mu);由于变步长所采用的步长过小，所需时间较长，导致输出的y矩阵过大，超出计算机存储空间容量。所以不适合采用ode45()来求解，可用刚性方程求解算法ode15s()。46.2.4.2高阶常微分方程组的变换方法5例：67描述函数：functiondx=apolloeq(t,x)mu=1/82.45;mu1=1-mu;r1=sqrt((x(
3、1)+mu)^2+x(3)^2);r2=sqrt((x(1)-mu1)^2+x(3)^2);dx=[x(2);2*x(4)+x(1)-mu1*(x(1)+mu)/r1^3-mu*(x(1)-mu1)/r2^3;x(4);-2*x(2)+x(3)-mu1*x(3)/r1^3-mu*x(3)/r2^3];8求解：>>x0=[1.2;0;0;-1.04935751];%输入初值>>tic,[t,y]=ode45('apolloeq',[0,20],x0);tocelapsed_time=0.5000>>length(t),>>plot(y(:,1),y(
4、:,3))ans=689得出的轨道不正确，默认精度RelTol设置得太大，从而导致的误差传递，可减小该值。9改变精度：>>options=odeset;options.RelTol=1e-6;>>tic,[t1,y1]=ode45('apolloeq',[0,20],x0,options);tocelapsed_time=0.2970>>length(t1),>>plot(y1(:,1),y1(:,3)),ans=187310>>min(diff(t1))ans=1.8927e-004>>plot(t1(1:end-1),…diff(t1))11例
5、：12>>x0=[1.2;0;0;-1.04935751];>>tic,[t1,y1]=rk_4('apolloeq',[0,20,0.01],x0);tocelapsed_time=0.8590>>plot(y1(:,1),y1(:,3))%绘制出轨迹曲线显而易见，这样求解是错误的，应该采用更小的步长。13>>tic,[t2,y2]=rk_4('apolloeq',[0,20,0.001],x0);tocelapsed_time=12.4380％计算时间过长>>plot(y2(:,1),y2(:,3))%绘制出轨迹曲线严格说来某些点仍不满足10－
6、6的误差限，所以求解常微分方程组时建议采用变步长算法，而不是定步长算法。14例：试将二元方程组：若两个高阶微分方程同时含有最高价导数项，需要先进行相应处理，再应用上述变换方法先消去其中一个高阶导数求解：1516用MATLAB符号工具箱求解，令％>>symsx1x2x3x4>>[dx,dy]=solve(‘dx+2*x4*x1=2*dy’,‘dx*x4+…3*x2*dy+x1*x4-x3=5’,‘dx,dy’)%dx,dy为指定变量dx=-2*(3*x4*x1*x2-5+x4*x1-x3)/(3*x2+2*x4)dy=(2*x4^2*x1+5-x4*
7、x1+x3)/(3*x2+2*x4)对于更复杂的问题来说，手工变换的难度将很大，所以如有可能，可采用计算机去求解有关方程，获得解析解。如不能得到解析解，也需要在描写一阶常微分方程组时列写出式子，得出问题的数值解。175.3特殊微分方程的数值解 5.3.1刚性微分方程的求解刚性微分方程－刚性过程：微分方程描述的变化过程中，若包含着多个相互作用但变化速度相差十分悬殊的子过程，这样一类过程就认为具有“刚性”。－刚性微分方程：描述“刚性”过程的微分方程称为刚性微分方程，相应的初值问题称为“刚性问题”。MATLAB采用求解函数ode15s()，该函数的调用格
8、式和ode45()完全一致。[t,x]=ode15s(Fun,[t0,tf],x0,options,p1,p2,…)18例
