好文档就是一把金锄头!
欢迎来到金锄头文库![会员中心]
电子文档交易市场
安卓APP | ios版本
电子文档交易市场
安卓APP | ios版本

MATLAB作业5.docx

6页
  • 卖家[上传人]:拖***
  • 文档编号:291070598
  • 上传时间:2022-05-11
  • 文档格式:DOCX
  • 文档大小:17.67KB
  • / 6 举报 版权申诉 马上下载
  • 文本预览
  • 下载提示
  • 常见问题
    • 本文格式为Word版,下载可任意编辑MATLAB作业5 MATLAB 1、 试求出下面线性微分方程的通解 dy(t)dt55作业5 ?13dy(t)dt44?64dy(t)dt33?152dy(t)dt22?176dy(t)dt?80y(t)?e?2t[sin(2t??3)?cos(3t)] y=dsolve(['D5y+13*D4y+64*D3y+152*D2y+176*Dy+80*y=','exp(-2*t)*(sin(2*t+pi/3)+cos(3*t))'],'y(0)=1','y(1)=3','y(pi)=2','Dy(0)=1','Dy(1)=2'); vpa(y,20) ans = .20576131687242798354e-2*exp(-2.*t)*cos(3.*t)+.15538705805619602373e-1*exp(-2.*t)*sin(2.*t)+.76830587084294035587e-2*exp(-2.*t)*cos(2.*t)+98.159206062620455336*exp(-2.*t)*t+59.405044899367325899*exp(-2.*t)*t^3-106.24422608844727795*exp(-2.*t)*t^2-30.741892776456442810*exp(-2.*t)+.20576131687242798354e-2*exp(-2.*t)*sin(3.*t)+31.732152104579289128*exp(-5.*t) 2、 试求解下面微分方程的通解以及得志x(0)?1,x(?)?2,y(0)?0条件下的解析解。

      ????(t)?4x(t)?3y(t)?ex(t)?5xsin(4t) ? ?6t??2y(t)?y(t)?4x(t)?6x(t)?ecos(4t)??6t [x,y]=dsolve('D2x+5*Dx+4*x+3*y=exp(-6*t)*sin(4*t)','2*Dy+y+4*Dx+6*x=exp(-6*t)*cos(4*t)','x(0)=1','x(pi)=2','y(0)=0'); >> vpa(x,10) ans = 0.08586166859*exp(t) - 0.057658489325149275828152894973755/(exp(7.549834435*t)^(1/4)*exp(t)^(13/4)) + (0.9469805542*exp(7.549834435*t)^(1/4))/exp(t)^(13/4) + (0.02481626654*cos(4.0*t))/exp(6.0*t) - (0.016682998530139342028763560499272*sin(4.0*t))/exp(6.0*t) >> vpa(y,10) ans = - 0.28620556196983670815825462341309*exp(t) + 0.09045056185/(exp(7.549834435*t)^(1/4)*exp(t)^(13/4)) + (0.3018304533*exp(7.549834435*t)^(1/4))/exp(t)^(13/4) - (0.10607545320921207832043364760466*cos(4.0*t))/exp(6.0*t) + (0.0683488486*sin(4.0*t))/exp(6.0*t) 3、 试求出微分方程??y(x)?(2?1x?(x)?(1?)y1x)y(x)?xe2?5x的解析解通解,并求出得志 边界条件y(1)??,y(?)?1的解析解。

      >> syms x y; y=dsolve('D2y-(2-1/x)*Dy+(1-1/x)*y=x^2*exp(-5*x)') y = C3*exp(t) + C2*exp((t*(x - 1))/x) + x^3/(exp(5*x)*(x - 1)) >> y=dsolve('D2y-(2-1/x)*Dy+(1-1/x)*y=x^2*exp(-5*x)','y(1)=pi','y(pi)=1') y = (exp(t)*(exp((x - 1)/x) - x*exp((x - 1)/x) - pi*exp((pi*(x - 1))/x) - (x^3*exp((pi*(x - 1))/x))/exp(5*x) + (x^3*exp((x - 1)/x))/exp(5*x) + pi*x*exp((pi*(x - 1))/x)))/(exp(pi)*exp((x - 1)/x) - exp(1)*exp((pi*(x - 1))/x) - x*exp(pi)*exp((x - 1)/x) + x*exp(1)*exp((pi*(x - 1))/x)) - (exp((t*(x - 1))/x)*(exp(1) - x*exp(1) - pi*exp(pi) - (x^3*exp(pi))/exp(5*x) + (x^3*exp(1))/exp(5*x) + pi*x*exp(pi)))/(exp(pi)*exp((x - 1)/x) - exp(1)*exp((pi*(x - 1))/x) - x*exp(pi)*exp((x - 1)/x) + x*exp(1)*exp((pi*(x - 1))/x)) + x^3/(exp(5*x)*(x - 1)) 4、 Lotka-Volterra扑食模型方程为??(t)?4x(t)?2x(t)y(t)?x?(t)?x(t)y(t)?3y(t)?y,且初值为x(0)?2,y(0)?3, 试求解该微分方程,并绘制相应的曲线。

      >> syms x y t; >> f=inline('[4*x(1)-2*x(1)*x(2); x(1)*x(2)-3*x(2)]','t','x'); >> [t,x]=ode45(f,[0,10],[2;3]);plot(t,x) 5.554.543.532.521.51012345678910 5、 是给出求解下面微分方程的MATLAB命令, y(3)?y?ty??y?ty22?e?ty,?(0)???y(0)?2,yy(0)?0 并绘制出y(t)曲线试问该方程存在解析解吗?选择四阶定步长Runge-Kutta算法求解该方程时,步长选择多少可以得出较好的精度,MATLAB语言给出的现成函数在速度、精度上 举行对比 该方程的解析解不存在 >> f=inline('[x(2); x(3); -t^2*x(1)*x(2)-t^2*x(2)*x(1)^2+exp(-t*x(1))]','t','x'); [t,x]=ode45(f,[0,10],[2;0;0]); >> plot(t,x) 2.521.510.50-0.5-1-1.5-2022345678910 6、 试用解析解和数值解的方法求解下面的微分方程组 ????(t)?e,x(t)??2x(t)?3x??(t)?4y?(t)?sint,y(t)?2x(t)?3y(t)?4x????5t?(0)?2x(0)?1,x?(0)?4y(0)?3,y 解析解: >> syms t x y >> [x,y]=dsolve('D2x=-2*x-3*Dx+exp(-5*t)','D2y=2*x-3*y-4*Dx-4*Dy-sin(t)','x(0)=1','Dx(0)=2','y(0)=3','Dy(0)=4') x = 17/(4*exp(t)) - 10/(3*exp(2*t)) + 1/(12*exp(5*t)) y = 100/(3*exp(2*t)) - 265/(16*exp(t)) - 71/(5*exp(3*t)) + 11/(48*exp(5*t)) + cos(t)/5 - sin(t)/10 + (51*t)/(4*exp(t)) 数值解: function dx=apolloeq(t,x) dx=[x(2);-2*x(1)-3*x(2)+exp(-5*t);x(4);2*x(1)-3*x(3)-4*x(2)-4*x(4)-sin(t)]; >> x0=[1;2;3;4]; >> [t,y]=ode45('apolloeq',[0,20],x0); >> plot(y(:,1),y(:,3)) 3.532.521.510.50-0.500.20.40.60.811.21.4 7、 下面的方程在传统微分方程教程中经常被认为是刚性微分方程。

      使用常规微分方程解法 和刚性微分方程解法分别求解这两个微分方程的数值解,并求出解析解,用状态变量曲线对比数值求解的精度 1??1?9y1?24y2?5cost?sint,y??3?1?y?2??24y1?51y2?9cost?sint,?3?y1(0)?y2(0)?1323(1) 解:function ydot = lorenzeq(t,y) ydot=[9*y(1)+24*y(2)+5*cos(t)-1/3*sin(t);-24*y(1)-51*y(2)-9*cos(t)+1/3*sin(t)]; >> t_final=100; y0=[1/3;2/3]; >> [t,y]=ode45('lorenzeq',[0,t_final],y0); >> plot(t,y) 10.80.60.40.20-0.2-0.4-0.60102030405060708090100 >> opt=odeset; opt.RelTol=1e-6; >> [t,y]=ode15s('lorenzeq',[0,t_final],y0,opt); >> plot(t,y) %刚性解法 10.80.60.40.20-0.2-0.4-0.60102030405060708090100 ?1??0.1y1?49.9y2,y1(0)?1?y??2??50y2,y2(0)?2 ?y?yy3(0)?1??3?70y2?120y3,(2)解:function ydot = lorenzeq(t,y) ydot=[-0.1*y(1)-49.9*y(2);-50*y(2); 70*y(2)-120*y(3)]; >> t_final=100; y0=[1;2;1]; [t,y]=ode45('lorenzeq',[0,t_final],y0); >> plot(t,y) — 6 —。

      点击阅读更多内容
      关于金锄头网 - 版权申诉 - 免责声明 - 诚邀英才 - 联系我们
      手机版 | 川公网安备 51140202000112号 | 经营许可证(蜀ICP备13022795号)
      ©2008-2016 by Sichuan Goldhoe Inc. All Rights Reserved.