MATLAB作业5

MATLAB作业51、试求出下面线性微分方程的通解。543225432()()()()()136415217680()[sin(2)cos(3)]3tdytdytdytdytdytytettdtdtdtdtdt假设上述微分方程满足已知条件，(0)1,(1)3,()2,(0)1,(1)2yyyyy&&试求出满足该条件的微分方程的解析解。解：symsty;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、试求解下面微分方程的通解以及满足(0)1,()2,(0)0xxy条件下的解析解。66()5()4()3()sin(4)2()()4()6()cos(4)ttxtxtxtytetytytxtxtet[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、试求出微分方程2511()(2)()(1)()xyxyxyxxexx的解析解通解，并求出满足边界条件(1),()1yy的解析解。symsxy;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扑食模型方程为()4()2()()()()()3()xtxtxtytytxtytyt，且初值为(0)2,(0)3xy，试求解该微分方程，并绘制相应的曲线。symsxyt;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)01234567891011.522.533.544.555.55、是给出求解下面微分方程的MATLAB命令，(3)22,(0)2,(0)(0)0tyytyytyyeyyy并绘制出()yt曲线。试问该方程存在解析解吗？选择四阶定步长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)012345678910-2-1.5-1-0.500.511.522.56、试用解析解和数值解的方法求解下面的微分方程组5()2()3(),(0)1,(0)2()2()3()4()4()sin,(0)3,(0)4txtxtxtexxytxtytxtyttyy解析解：symstxy[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))数值解：functiondx=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))00.20.40.60.811.21.4-0.500.511.522.533.57、下面的方程在传统微分方程教程中经常被认为是刚性微分方程。使用常规微分方程解法和刚性微分方程解法分别求解这两个微分方程的数值解，并求出解析解，用状态变量曲线比较数值求解的精度。11212122119245cossin,(0)33(1)1224519cossin,(0)33yyyttyyyytty解：functionydot=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)0102030405060708090100-0.6-0.4-0.200.20.40.60.81opt=odeset;opt.RelTol=1e-6;[t,y]=ode15s('lorenzeq',[0,t_final],y0,opt);plot(t,y)%刚性解法0102030405060708090100-0.6-0.4-0.200.20.40.60.81112122232330.149.9,(0)1(2)50,(0)270120,(0)1yyyyyyyyyyy解：functionydot=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)0102030405060708090100-1-0.500.511.52opt=odeset;opt.RelTol=1e-6;[t,y]=ode45('lorenzeq',[0,t_final],y0,opt);plot(t,y)%刚性解法opt=odeset;opt.RelTol=1e-6;[t,y]=ode45('lorenzeq',[0,t_final],y0,opt);plot(t,y)0102030405060708090100-1-0.500.511.528、试求出习题3中给出的微分方程边值问题数值解，绘制出()yt曲线，并和该习题得出的解析解比较精度。9、试用数值方法求解偏微分方程22220,00,001,00,0xyyxuuxyuuxy，并绘制出u函数曲面。pdetool