dengyue.zip_2014数学建模

clc clear um = 4.902*10^12; % 月球引力常数 R = 1738*10^3; % 月球半径 m= 2400 ; % 飞行器质量--------需要改 %------------ x0 = 1753*10^3; y0 = 0; z0 = 0; vx0 =0; vy0 = 1700; vz0=0; %---------- theta = 19.51; %　着陆的月心角 ----- 需要改 x1= R*cos(theta/57.3); y1= R*sin(theta/57.3); z1= 0; vx1= 0; vy1=0; vz1=0; %----------- r0 = sqrt(x0^2+y0^2+z0^2); a0 = -um*x0/r0^3; b0 = -um*y0/r0^3; c0 = -um*z0/r0^3; t1= 750; % 　飞行的总时间------------需要改 a1 =-6*(a0*t1^2+(vx1+3*vx0)*t1-4*(x1-x0))/t1^3; b1 = -6*(b0*t1^2+(vy1+3*vy0)*t1-4*(y1-y0))/t1^3; c1 = -6*(c0*t1^2+(vz1+3*vz0)*t1-4*(z1-z0))/t1^3; a2= 6*(a0*t1^2+2*(vx1+2*vx0)*t1-6*(x1-x0))/t1^4; b2= 6*(b0*t1^2+2*(vy1+2*vy0)*t1-6*(y1-y0))/t1^4; c2=6*(c0*t1^2+2*(vz1+2*vz0)*t1-6*(z1-z0))/t1^4; t= 1:t1; % 需要修改 ax = a0+a1*t+a2*t.^2; ay = b0+b1*t+b2*t.^2; az = c0+c1*t+c2*t.^2; vx = vx0+a0*t+a1*t.^2/2+a2*t.^3/3; vy = vy0+b0*t+b1*t.^2/2+b2*t.^3/3; vz = vz0+c0*t+c1*t.^2/2+c2*t.^3/3; x = x0+vx0*t+a0*t.^2/2+a1*t.^3/6+a2*t.^4/12; y = y0+vy0*t+b0*t.^2/2+b1*t.^3/6+b2*t.^4/12; z = z0+vz0*t+c0*t.^2/2+c1*t.^3/6+c2*t.^4/12; r = sqrt(x.^2+y.^2+z.^2); gx = -um*x/r.^3; gy = -um*y/r.^3; gz = -um*z/r.^3; apx = ax-gx; apy = ay-gy; apz =az-gz; ap = sqrt(apx.^2+apy.^2+apz.^2); %推力加速度随时间的函数 p = ap*m; % p是推力 plot(t,p) % 画推力随时间变化的图 title('推力大小随时间过程的变化'); %------------------------------------ % 下面来画积分图 syms t ; f =