Matlab第六章课文.pptVIP

MATLAB数学实验 第六章 常微分方程 第六章 常微分方程 6.1 预备知识：常微分方程 6.2 解常微分方程的MATLAB指令 6.3 计算实验：Euler法和刚性方程组 6.4 建模实验: 导弹系统的改进 6.1 预备知识：常微分方程 数值解 准确解y(xn)未知, 数值解yn作为y(xn)的近似 6.2解常微分方程MATLAB指令 odefun= inline(y-2*t/y,t,y); [t,y]=ode45(odefun,[0,4],1); [t,y] plot(t,y,o-) [t,y]=ode45(odefun,0:1:4,1);[t,y] %M函数eg6_3fun.m function f=eg6_3fun(t,x) f(1)=-x(1)^3-x(2); f(2)=x(1)-x(2)^3; f=f(:); [t,x]=ode45(@eg6_3fun,[0 30],[1;0.5]) 作图 subplot(1,2,1);plot(t,x(:,1),t,x(:,2),:); subplot(1,2,2);plot(x(:,1),x(:,2)); eg6_5.m clear;close; sinit=bvpinit(0:4,[1;0]) %注意sinit的域名 odefun=inline([y(2);-abs(y(1))],t,y); bcfun=inline([ya(1);yb(1)+2],ya,yb); sol=bvp4c(odefun,bcfun,sinit) %注意sol的域名 t=linspace(0,4,101); y=deval(sol,t); plot(t,y(1,:),sol.x,sol.y(1,:),o,sinit.x,sinit.y(1,:),s) legend(解曲线,解点,粗略解) 6.3 计算：Euler法和刚性方程组 M函数euler.m function [t, y] = euler(odefun, tspan, y0, h) t=tspan(1):h:tspan(2); y(1)=y0; for i=1:length(t)-1, y(i+1)=y(i)+h*feval(odefun,t(i),y(i)); end t=t’;y=y’; 6.4 建模实验: 导弹系统的改进 例7 程序 close; fplot(‘exp(0.9*x)’, [0,10]); %模型1解析解 hold on; [t,x]=ode45(inline(‘0.9*x*(1-x/1000)’,‘t’,‘x’),[0 10],1); %模型2数值解 plot(t,x); axis([0 10 0 1500]); hold off; 模型比较 短期预报二个模型相近, 但作为长期预报,后者较前者合理。 第六章 习题 ex1(1)(2)(勘误: (2) -x’=2x+3y ) ex1(3)(6) ex2 ex7 ex8 例8 在导弹系统中设 a=90km/h, b=450km/h, T=0.1h. 求d, ?的有效范围? 解 有两个极端情形容易算出。 若? =0, 即敌舰正好背向行驶, 即x轴正向。 那么导弹直线飞行, 击中时间 t=d/(b-a)T 得d=T(b-a)=36km。 若? =?, 即迎面驶来, 类似地有d=T(a+b)=54km 一般地, 有36d54。 在线算法：对于测定的d 和?，可用 计算出t。如d=50, ? =?/2，写出M函数eg6_8fun，再用ode45求解 在线算法：eg6_8a 离线算法： 首先对于所有可能的d和?，计算击中所需时间，从而对不同? ，得d的临界值。具体应用时直接查表判断。 计算机模拟： 对于任意选定的参数a,b,d,T和?，利用MATLAB动画制作指令getframe和动画播放指令movie制作一个导弹追击敌舰演示工具 eg6_8b.m eg6_8c(90,450,30,0.3*pi,0.1) * * 1.微分方程的概念 常微分方程: f (t,y,y’,y’’,…,y(n))=0 微分方程组: 联系一些未知函数的一组微 分方程 线性常微分方程: y(n) + a1 (t) y(n-1) + … + an-1 (t) y’ + an (t) y = b(t) 若ai (t) (i =1, …,n) 与t无关, 称为常系数的 若b(t)=0,称为齐次的 2.初等积分法 分离变量法等 3.常系数线性微分方程 线性常微分方程的解为一个特解和相应 的齐次微分方程通解的叠加。 齐次微分方程的解可用特征根法求得 例1 求x’’+ 0.2 x’+3.92x = 0的通解 解 特征方程为?2 + 0.2? +3.