卫星轨道和位置.docx
《卫星轨道和位置.docx》由会员分享，可在线阅读，更多相关《卫星轨道和位置.docx（11页珍藏版）》请在冰豆网上搜索。
卫星轨道和位置
卫星轨道和位置
水星的轨道和位置
摘要
本文主要在已知水星的远日点和绕日运行的线速度的条件下，通过建立微分方程模型，使用解析法和数值方法求解水星的轨道方程与位置。
解析法的求解的过程中，
，求得水星到太阳的最近结合了开普勒三大定律，准确的给出了微分方程的精确解
10距离，水星绕太阳运行的周期约为88天。
数值计算求解水星自r,4.6016，10（m）m
远日点运行50天后的位置时，本文分别采用了Simpson求积法，基于压缩映射的求根方法以及经典的四阶龙格—库塔法，使用matlab数学软件编程，得到了较为合
10,,3.791理的行星运行模型的近似解，三种方法所得结果对应分，，r,，4.7671011
1010,,3.802，及，。
r,，4.76710r,，4.77910,,3.7912323
关键词行星轨道微分方程Simpson法四阶龙格—库塔法matlab
一（问题重述
1140.698210，3.88610，水星到太阳的最远距离为m，此时水星绕太阳运行的线速度为m,s。
试求
问题一水星到太阳的最近距离
问题二水星绕太阳运行的周期
问题三从远日点开始的第50天（地球天）结束时水星的位置并画出轨道曲线
二（问题分析
求水星到太阳的最近距离以及水星绕太阳运行的周期等，需要先将水星轨道方程
2mMGdZ,i,,em求出，因此可以根据Newton第二定律及万有引力定律，建立微22rdt
PAGE1
分方程模型，将原问题转化为求解带有初值条件的微分方程问题，进而采用解析法或数值方法求解远日点和周期。
三（模型假设
1（水星运行的轨道是以太阳为一个焦点的椭圆
2（从太阳指向水星的线段在单位时间内扫过的面积相等
3（水星运行周期的平方与其运行轨道椭圆长轴的立方之比为常量
四（符号系统
1（水星在远日点的线速度v0
M2.太阳的质量
3.水星的质量m
4.水星在远日点的距离ro
T5.周期
五（建立模型与求解
模型一水星的轨迹方程
i,设太阳中心所在的位置为复平面的原点O，在时刻t，水星位于Ztre（）,
rrtt,,（）,（）,,Zt（）所表示的点P。
这里均为t的函数，分别表示的模和辐角。
于是
dZdrddrd,,iii,,,水星的速度为，加速度为,，,，（）eireeirdtdtdtdtdt
222，，dZdrdddrd,,,i,2（1.1），而太阳对行星的引力依万有引（（））
（2）erir,,，，,,222dtdtdtdtdtdt，，
mMGmMGi,力定律，大小为，方向由行星位置P指向太阳的中心O,故为，其,e22rr
30,1122中为太阳的质量，m为水星的质量，Mkg,，1.98910（）GNmkg,，,6.67210（/）为万有引力常数。
PAGE2
2mMGdZ,i依Newton定律，我们得到,,（1.2），将（1.1）代入（1.2），em22rdt
然后比较实部与虚部，就有
2,ddrd,,r20，,,2,dtdtdt,2drdMG,2,r（）,,,22,dtdtr,
这是两个未知函数的二阶微分方程组。
在确定某一行星轨道时，需要加上定解条件。
假设当t=0时，行星正处于远日点，而远日点位于正实轴上，距原点O为，r0行星的速度为。
那么就有初值条件:
v0
,rrt,00,,,0,t,0
,dr,0,t,0dt,
vd0,,,t,0dtr,0,
因此问题转化为求解带初值问题的微分方程组
2,ddrd,,r，,20,2dtdtdt,2drdMG,,2,,,r（）22,dtdtr,
rr,,00t,,,,0t,0,,
dr,0,t,0dt,
vd0,,,0t,dtr,0,
2ddrd,,dd,d,22r，,20又将两边同乘以r，即得，从而（1.3）,,（）0r,rc12dtdtdtdtdtdt
tt，,1d,ct,21,t其中，这样有向线段在时间内扫过的面积等于，这crv,OPrdt,100,22dtt个正是Kepler的第二定律，从太阳指向水星的线段在单位时间内扫过的面积相等。
PAGE3
222drdMG,drMGc21将（1.3）代入,,,r（）得，于是我们可以得到水,,,22232dtdtrdtrr
星运行的较为简单形式的数学模型:
22,drMGc1,,,,232dtrr,cd,,1,2,dtr,rr,,00t,
dr,,00t,dt,
,00t,,,
c1d,1为了求得行星的轨迹方程，要消去变量t，令，那么可以改写为,r,2udtr
22drddududdu,d,222,,,,,,cccu（）从而将上式代入,cu111122dtdtdddtddt,,,
2222du11drMGcc11，,uuu,,,,,p,，化简后为（1.4），其中，引进，立2232,dppdtrrMG
1eAp,,,,,,,uuAcos（）即可以求出，这里A和是待定的常数。
记，上式可,00p
以写为
pr,1cos（）e,,,,0
这个就是水星的轨道方程，是一条平面二次曲线。
由于水星绕太阳运行，故必
有01,,e。
由于r在t=0时取道最大值（远日点），这个就意味着此时函数r0
p,,,,0,1e取道最大值1.于是就有，从而轨迹方程为cos（）,,,00r0
p114rrmvms,，,，0.69810（）,3.88610（/）,。
对于水星而言，，又水星的001cose,,
pp近日点到太阳的距离。
依据已知数据，可知,,rm1cos1,,，ee
PAGE4
2cp101152pm,,，5.54710（），，，从而计e,,,10.2055crvms,,，2.71310（/）100MGr0
10算水星到太阳的最近距离为r,4.6016，10（m）m
模型二水星的运行周期
设水星的周期为T，那么利用Kepler第二定律，我们有
T11d,2（1.4）,rdtCT1,022dt
p上式左端为水星轨迹椭圆所围的面积，记为S，由于椭圆的半长轴，半短a,21,e
22pp,2p,b,S,ab,轴，从而有将上式代入式（1.4）,解得T,,3321,e2222C（1,e）（1,e）1
6（1.5）将有关数据代入，易得T,7.6025，10（s）,87.9919（d）
模型三水星的位置
，,,2由于水星的运行满足Kepler第二定律，则该式可改写为，从rd,,C,t1,,而可得
2,pdt,,2,0C（1ecos）,,1
CT11F（,）,如果我们要求时相应的和，则意味着首先要解方程，，,t,Tr12p其中
1（）,F,d,2,0（1,cos）e,
pr在求出了时的后，立即可以由得到相应的r。
t,T,,,11cose,,
下面用数值方法求解水星的位置
1.Simpson法
PAGE5
CT111F（）,由被积函数的恒正性可知单调，从而方程F（,）,的根必22,,（1cos）ep
CTCT1111,,,存在且唯一。
取，记。
若FF，,,,,,,hkhk,（1,2,...）FF,（）,nn，1kkk22pp那么位于与之间，在h适当小时，可取。
,,,,,nnn，1
F（）,计算可采用不同的数值积分法，本文采用Simpson法，取步长h=0.001，
10r,，4.76710具体求解过程见附录一，最后结果为，,,3.791
2.基于压缩映像的求根方法
p我们引入水星轨道椭圆的参数方程，由于椭圆的半长轴，半短轴a,21,e
p22b,，从而中心到焦点的距离为。
因左焦点为原点，故椭圆中a,b,ae21,e
心位于（ae，0），于是得到参数方程
xae,，（cos）,,,yb,sin,,
r,,它们与的关系为
y222x，y,r,,tan,x
,，,Ctxyyxdabeeab,,,,，,,，,（''）[sin（）sin],,,,,此式可改写成1,,
当时解方程t,T1
CT11esin,，,,abCT
11g（,）,,,esin,,,,g（）记，，那么上式即，就是说要去求函,,ab
g（）,数的不动点，求解方程不动点可以采用简单迭代法，对于水星，我们已计算出e,0.2055，由于e很小，因此迭代收敛理论上可以很快，当时间从远日点开始的第
750天结束时，意味着,从而T,0.432，10（s）1
3CTCT21112,,,（1,e）,3.57032abp
PAGE6
不妨取，于是,,00
,,,esin,,3.570310
,,,esin,,3.655721
......
,,,,,esin3.674754
,,,,,esin3.674765
故,,3.6747
y222x,a（e，cos,）,y,bsin,由式，，可以计算出相应的，x，y,r,,tan,,x
即由
bsin,tan,,0.75849
ae,（，cos）
得0.64891，而,,,,,,，,3.791
2210r,[a（e，cos,）]，[bsin,],4.7668，10此时的距离为（m）r
3.经典四阶Runge-Kutte法
由我们将由最初的微分方程组求解水星的位置，方程组见下
22,cdrMG1,,,,232dtrr,cd,,1,2,dtr,rr,,00t,
dr,,00t,dt,
,00t,,,
drq,令，那么我们可以得到一阶微分方程组:
dt
PAGE7
2,CdqMG1,,,32dtrr,
dr,,q,dt,,Cd,1,,2dtr,r,r,t,00,dr,,0t,0,dt
,,0t,0,
若记这个微分方程组中方程的右端依次为，Q（t,q,r,,）,R（t,q,r,,）和S（t,q,r,,）
则相应的四阶Runge-Kutte迭代格式法为
hq,q，（K，2K，2K，K）k，1k12346
hr,r，（L，2L，2L，L）k，1k12346
h,,,，（N，2N，2N，N）k，1k12346
h这里对于，有,，，，，（22）qqKKKKkk，112346
K,Q（t,q,r,,）1kkkk
hKhLhNh111,（，,，,，,,，）KQtqrkkkk22222
hKhLhNh222,（，,，,，,,，）KQtqrkkkk32222
K,Q（t，h,q，hK,r，hL,,，hN）4kk3k3k3初值为，则对于给定的步长值h，类似可以逐步计算一系列的q,0,r,r,,,00000
，由于行星绕着太阳运行，只需取，而取得行星轨道上一qr,,,,,,,,,2,2kkknn，1系列点的近似坐标（r,,），再通过极坐标与直角坐标的转换，继而可以绘出轨道kk
10r,，4.77910曲线。
通过matlab编程求解得，，轨道曲线如下,,3.802
PAGE8
程序见附录二。
六（模型推广
本文建立的微分方程模型对于求解行星绕日运行轨道具有广泛的应用空间，只需给出行星的远日点和在远日点的运行线速度即可计算出轨道方程，用数学软件绘出近似的轨道曲线，对于研究天体运行有所帮助。
此外，本文采用的求解微分方程的数值方法，具有较为快速且准确的收敛效果，可以用来求解其他类似的微分方程模型。
七（参考文献
【1】乐经良，数学实验，北京，高等教育出版社，1999年10月【2】周品，matlab数值分析，北京，机械工业出版社，2009年1月
PAGE9
八（附录
附录一
functionq1=y2（x）
q1=（1-0.2055*cos（x））.^-2;
h=0.001;
k=1;
x=h*k;
f=quad（'y2',0,x）
whilef（k）<3.8091
k=k+1;
x=k*h;
f（k）=quad（'y2',0,x）;
end
x
附录二
formatlong
c1=2.7132e15;M=1.989e30;G=6.672e-11;Q=inline（'2.7132e15^2/（r^3）-1.989e30*6.672e-11/（r^2）'）;
R=inline（'q'）;
S=inline（'2.7132e15/（r^2）'）;q=0;r=0.6982e11;theta=0;t=0;k=1;h=0.001e7;
whiletheta<=2*pi
K1=Q（r）;L1=R（q）;N1=S（r）;
K2=Q（r+h/2*L1）;L2=R（q+h/2*K1）;N2=S（r+h/2*L1）;
K3=Q（r+h/2*L2）;L3=R（q+h/2*K2）;N3=S（r+h/2*L2）;
K4=Q（r+h*L3）;L4=R（q+h*K3）;N4=S（r+h*L3）;
t=t+h;
PAGE10
q=q+h/6*（K1+2*K2+2*K3+K4）;
r=r+h/6*（L1+2*L2+2*L3+L4）;
theta=theta+h/6*（N1+2*N2+2*N3+N4）;
rr（k）=r;
ee（k）=theta;
xx（k）=rr（k）*cos（ee（k））;%水星任意位置的横坐标
yy（k）=rr（k）*sin（ee（k））;%水星任意位置的纵坐标
k=k+1;
end;
plot（xx,yy）%画出水星的轨道曲线
text（0,0,'太阳'）
text（0.6982e11,0,'远日点'）
text（-4.6078e+010,0,'近日点'）
holdon;plot（0,0,'r.','MarkerSize',20）;holdoffholdon;plot（0.6982e11,0,'r.','MarkerSize',20）;holdoff
holdon;plot（-4.6078e+010,0,'r.','MarkerSize',20）;holdoff
title（'水星绕太阳运行的轨道曲线'）
clc
q=0;r=0.6982e11;theta=0;t=0;k=1;
h=0.001e7;
whilet<=50*24*3600%求水星自远日点开始第50天的位置
K1=Q（r）;L1=R（q）;N1=S（r）;
K2=Q（r+h/2*L1）;L2=R（q+h/2*K1）;N2=S（r+h/2*L1）;
K3=Q（r+h/2*L2）;L3=R（q+h/2*K2）;N3=S（r+h/2*L2）;
K4=Q（r+h*L3）;L4=R（q+h*K3）;N4=S（r+h*L3）;
t=t+h;
q=q+h/6*（K1+2*K2+2*K3+K4）;
r=r+h/6*（L1+2*L2+2*L3+L4）;
theta=theta+h/6*（N1+2*N2+2*N3+N4）;
rr（k）=r;
ee（k）=theta;
xx（k）=rr（k）*cos（ee（k））;
PAGE11
yy（k）=rr（k）*sin（ee（k））;
k=k+1;
end
r=rr（k）
theta=theta%水星自远日点开始第50天的位置（半径及角度）text（xx（k）,yy（k）,'自远日点运行50天后的位置'）
holdon
plot（xx（k）,yy（k）,'g.','MarkerSize',20）;holdoff
PAGE12
升级会员 