行星引力加速

行星引力加速

姓名: 刘志强

99级3班

实验说明：

旅行者Ⅱ

引力加速原理

现在说明引力加速原理. 如下页

图所示,当飞船沿椭圆轨道飞行到C

点时与某行星(设为木星)相遇, 在进

入木星引力为主的范围后, 飞船被木

星短时俘获, 这时, 飞船一面被木星携

带着以巨大的木星轨道速度运动, 一

面被木星较强的引力改变着运动方向.

在被俘获的较短时间内, 太阳对飞船

的引力作用相对木星引力作用很小,

而由于飞船质量相对于木星质量很小,

其几乎不改变木星的运动,在短时间

内木星视为匀速直线运动, 木星坐标

系(即原点建立在木星中心, 三坐标轴相对于太阳系不能转动的坐标系)可视为惯性

系. 飞船相对该惯性系的运动与电粒子在核场中的散射类似, 遵守机械能守恒定律,

因此, 飞船再次飞离木星引力场时的相对速度的大小Uf,与它进入木星引力场的速

度大小是相等的,木星引力场的作用仅使相对速度的方向发生偏转.

设木星的轨道速度Vt, 飞船与木星相互作

用前后相对日心系的速度分别为Vi和Vf, 则

:

Vi = Vj + Ui

Vf = Vj +Uf

飞船获得的动能增量为:

ΔT = mVj·(Uf-Ui)

通过改变瞄准距离, 速度大小的增量可接

近Vi, 飞船的轨道改为双曲线, 从而大大缩短飞

行时间.

程序原理

诸行星的椭圆轨道的离心率皆近似为0, 故近似

看成圆. 以太阳为原点, 可根据各个行星的半径画出其轨道. 以木星为例: 设木星轨道半径为Rjup, 周期为Tjup, 则:

木星角速度为:Wjup = 2π/Tjup.

令 t =0 时, 木星角坐标为 Θjup0, 则 t 时刻木星的角坐标为:

Θjup(t)=Θjup0 + tWjup =Θjup0 + 2πt/Tjup

令木星的公转轨道半径为rjup，则木星t时刻的位置可表示为：

x(t)=rjupcosΘjup(t)

y(t)=rjupsinΘjup(t)

于是木星的位置可表示出。同理可表示出其他诸星的位置。

可设飞船的位置为(x(t),y(t))，则飞船受力作用可表示如下：

∑∑

+−=+−=

23222322

)()()()(

iiiiyiiiix

yxyymGMFyxxxmGMF

于是，飞船的加速度在x方向及y方向的分量为：

∑∑

+−=+−=

23222322

)()()()(

iiiiyiiiix

yxyyGMayxxxGMa

至此，飞船的动力学方程可列出。

程序说明

开始 调用函数actual, 进行旅行者II号运行的仿真模

拟。

旅行者轨迹 调用orbit函数，展示旅行者II号运行

的

仿真模拟的轨迹。

加速模型 为了清晰地展示出飞船经加速后轨道由椭

圆变为双曲线，本按键调用的model_1函数只计入太阳及木

星的引力作用。 程序说明 即此论文。

离开 关闭程序。

actual函数

G=6.67e-11 万有引力常数

Msun=1.99e30 太阳质量

Mear=5.976e24 地球质量

Mjup=1.9e27 木星质量

Msat=5.688e26 土星质量

Mura=8.684e25 天王星质量

Mnep=1.024e26 海王星质量

unit=1.4959787e11 天文单位

rear=unit 地球与太阳距离

rjup=5.203*unit 木星与太阳距离

rsat=9.522*unit 土星与太阳距离

rura=19.201*unit 天王星与太阳距离

rnep=30.074*unit 海王星与太阳距离

%the original value of therax:tx0

tear0=0.1; 地球初始角坐标

tjup0=2.0601 木星初始角坐标

tsat0=3.3207; 土星初始角坐标

tura0=4.4899; 天王星初始角坐标

tnep0=4.88; 海王星初始角坐标

%draw orbits

ti=0:pi/60:2*pi;

axis off

plot(rear*cos(ti),rear*sin(ti),'c',rjup*cos(ti),rjup*sin(ti),'c',rsat*co

s(ti),rsat*sin(ti),'c', ...

rura*cos(ti),rura*sin(ti),'c',rnep*cos(ti),rnep*sin(ti),'c');

hold on

xo=(rear+10.37e6)*cos(tear0); 地球初始x坐标

xd=3.7e9*(-sin(tear0)); 地球初始x方向速度

yo=(rear+10.37e6)*sin(tear0); 地球初始y坐标

yd=3.7e9*(cos(tear0)); 地球初始y方向速度

[t,x]=ode23('v2',[0:1:4516],[xo,xd,yo,yd],...

[],tear0,tjup0,tsat0,tura0,tnep0,Msun,Mear,Mjup,Msat,Mura,Mnep,rear,rjup

,rsat,rura,rnep,G);

save data t , x; 保存数据：时间及飞船的位置及速度

以下为模拟动化

Sun=line(0,0,'color','r','marker','.','markersize',25,'erasemode','xor');

Earth=line(rear*cos(tear0),rear*sin(tear0),'color','b','marker','.','markersize',15,'erasemode','xor');

Jupiter=line(rjup*cos(tjup0),rjup*sin(tjup0),'color','g','marker','.','m

arkersize',25,'erasemode','xor');

Saturn=line(rsat*cos(tsat0),rsat*sin(tsat0),'color','r','marker','.','ma

rkersize',20,'erasemode','xor');

Uranus=line(rura*cos(tura0),rura*sin(tura0),'color','r','marker','.','ma

rkersize',15,'erasemode','xor');

Neptune=line(rnep*cos(tnep0),rnep*sin(tnep0),'color','k','marker','.','m

arkersize',20,'erasemode','xor');

voyager=line(rear*cos(tear0),rear*sin(tear0),'color','k','marker','*','m

arkersize',5,'erasemode','xor');

for t=1:1:4517

thetaear=tear0+t/365*2*pi;

thetajup=tjup0+t/(11.86*365)*2*pi;

thetasat=tsat0+t/(29.46*365)*2*pi;

thetaura=tura0+t/(84.01*365)*2*pi;

thetanep=tnep0+t/(164.79*365)*2*pi;

set(Earth,'xdata',rear*cos(thetaear),'ydata',rear*sin(thetaear));

set(Jupiter,'xdata',rjup*cos(thetajup),'ydata',rjup*sin(thetajup));

set(Saturn,'xdata',rsat*cos(thetasat),'ydata',rsat*sin(thetasat));

set(Uranus,'xdata',rura*cos(thetaura),'ydata',rura*sin(thetaura));

set(Neptune,'xdata',rnep*cos(thetanep),'ydata',rnep*sin(thetanep));

set(voyager,'xdata',x(t,1),'ydata',x(t,3));

drawnow;

end

hold off

图象如下

v2函数

function xdot=voyager(t,x,flag,tear0,tjup0,tsat0,tura0,tnep0,Msun,...

Mear,Mjup,Msat,Mura,Mnep,rear,rjup,rsat,rura,rnep,G)

xdot=[x(2);

3600*24*3600*24*(G*Msun*(-

x(1))/(x(1)^2+x(3)^2)^1.5+G*Mear*(rear*cos(tear0+t/365*2*pi)-

x(1))/((rear*cos(tear0+t/365*2*pi)-x(1))^2+(rear*sin(tear0+t/365*2*pi)-

x(3))^2)^1.5+G*Mjup*(rjup*cos(tjup0+t/(11.86*365)*2*pi)-

x(1))/((rjup*cos(tjup0+t/(11.86*365)*2*pi)-

x(1))^2+(rjup*sin(tjup0+t/(11.86*365)*2*pi)-

x(3))^2)^1.5+G*Msat*(rsat*cos(tsat0+t/(29.46*365)*2*pi)-

x(1))/((rsat*cos(tsat0+t/(29.46*365)*2*pi)-

x(1))^2+(rsat*sin(tsat0+t/(29.46*365)*2*pi)-

x(3))^2)^1.5+G*Mura*(rura*cos(tura0+t/(84.01*365)*2*pi)-

x(1))/((rura*cos(tura0+t/(84.01*365)*2*pi)-

x(1))^2+(rura*sin(tura0+t/(84.01*365)*2*pi)-

x(3))^2)^1.5+G*Mnep*(rnep*cos(tnep0+t/(164.79*365)*2*pi)-

x(1))/((rnep*cos(tnep0+t/(164.79*365)*2*pi)-

x(1))^2+(rnep*sin(tnep0+t/(164.79*365)*2*pi)-x(3))^2)^1.5);

x(4);

3600*24*3600*24*(G*Msun*(-

x(3))/(x(1)^2+x(3)^2)^1.5+G*Mear*(rear*sin(tear0+t/365*2*pi)-

x(3))/((rear*cos(tear0+t/365*2*pi)-x(1))^2+(rear*sin(tear0+t/365*2*pi)-

x(3))^2)^1.5+G*Mjup*(rjup*sin(tjup0+t/(11.86*365)*2*pi)-

x(3))/((rjup*cos(tjup0+t/(11.86*365)*2*pi)-

x(1))^2+(rjup*sin(tjup0+t/(11.86*365)*2*pi)-