精通MATLAB的请帮忙

jonwang · 发表于 2004-12-30 21:42:41

<

>行星和太阳的质量，行星到太阳的平均距离，行星的平均速度
<

>采用牛顿万有引力定律和加速度定律模拟太阳系的运动轨迹。
<

>运行过程中观察4个行星的运动。修改该程序：
 % total number of big objects in solar system
Nsphere = 10;
 % color of planets
Col = [1:Nsphere];
 % Gravitational constant in kg^-2*m^3*s^-2
Grav = 6.67259*(10^(-11));

% mass of planets in kg
Mass = zeros(Nsphere,1);
Mass(1) = 0.33*(10^24); % Mercury
Mass(2) = 4.87*(10^24); % Venus
Mass(3) = 5.97*(10^24); % Earth
Mass(4) = 0.642*(10^24); % Mars
Mass(5) = 1899*(10^24); % Jupiter
Mass(6) = 568*(10^24); % Saturn
Mass(7) = 86.8*(10^24); % Uranus
Mass(8) = 102*(10^24); % Neptune
Mass(9) = 0.0125*(10^24); % Pluto
Mass(10) = 1047.3486*Mass(5); % The sun 
 % initial position in m
X = zeros(Nsphere,1);
Y = zeros(Nsphere,1);
X(1) = 57.9*(10^9); % Mercury
X(2) = 108.2*(10^9); % Venus
X(3) = 149.6*(10^9); % Earth
X(4) = 227.9*(10^9); % Mars
X(5) = 778.6*(10^9); % Jupiter
X(6) = 1433.5*(10^9); % Saturn
X(7) = 2872.5*(10^9); % Uranus
X(8) = 4495.1*(10^9); % Neptune
X(9) = 5870.0*(10^9); % Pluto
X(10) = 0; % The sun
 scatter(X,Y,'or');
 % initial velocity in m*s^-1 (using mean orbital velocity)
VX = zeros(Nsphere,1);
VY = zeros(Nsphere,1);
 VY(1) = 47.89*(10^3); % Mercury
VY(2) = 35.03*(10^3); % Venus
VY(3) = 29.79*(10^3); % Earth
VY(4) = 24.13*(10^3); % Mars
VY(5) = 13.06*(10^3); % Jupiter
VY(6) = 9.64*(10^3); % Saturn
VY(7) = 6.80*(10^3); % Uranus
VY(8) = 5.40*(10^3); % Neptune
VY(9) = 4.70*(10^3); % Pluto
VY(10) = 0; % The sun 
 % acceleration of spheres
AX = zeros(Nsphere,1);
AY = zeros(Nsphere,1);
 % Use Newton's law of gravitation to move the planets in space
dt = 300;
day = 24*60*60/dt;

for tstep=1:100000000
 for sphere=1:Nsphere
 rest = setdiff( (1:Nsphere),sphere);
 dists = (X(rest)-X(sphere)).^2+(Y(rest)-Y(sphere)).^2;
 dists = dists.^(-3/2);
 AX(sphere) = Grav*Mass(rest)'*((X(rest)-X(sphere)).*dists);
 AY(sphere) = Grav*Mass(rest)'*((Y(rest)-Y(sphere)).*dists);
end

X = X + dt*VX;
Y = Y + dt*VY;

VX = VX + dt*AX;
VY = VY + dt*AY;
 if(~mod(tstep,day))
 title(sprintf('day=%d', tstep/day));
 % only plot sun and first four planets
 RANGE = max(max(abs(X(1:4))),max(abs(Y(1:4))));
 hold on; scatter(X,Y,'o',Col);
 axis([-RANGE RANGE -RANGE RANGE]);
 hold off;

 drawnow;
end
 end
A 添加一个彗星（需要增加Nsphere 数组等）.
假设这个和彗星以一个初速度从Jupiter附近开始运动。
B 在地球附近添加月球
C 能够演示出彗星撞地球或其他行星；找到月球对地球应力最大的时刻
请问怎么修改?

		自动登录	找回密码
密码			注-册-帐-号