水星的運(yùn)動(dòng)規(guī)律

1、水星的運(yùn)動(dòng)規(guī)律摘要本文主要在已知水星的遠(yuǎn)日點(diǎn)和繞日運(yùn)行的線速度的條件下，通過建立微分方程模型，使用解析法和數(shù)值方法求解水星的軌道方程與位置。解析法的求解的過程中，結(jié)合了開普勒三大定律， 準(zhǔn)確的給出了微分方程的精確解， 求得水星到太陽的最近距離 rm定 4.6016x101(m),水星繞太陽運(yùn)行的周期約為 88 天。數(shù)值計(jì)算求解水星自遠(yuǎn)日點(diǎn)運(yùn)行50 天后的位置時(shí)，本文分別采用了 Simpson 求積法，基于壓縮映射的求根方法以及經(jīng)典的四階龍格一庫塔法，使用 matlab 數(shù)學(xué)軟件編程，得到了較為合理的行星運(yùn)行模型的近似解，三種方法所得結(jié)果對(duì)應(yīng)分用=3.791,1球 4.767 父 1010,也=

2、3.791,r2歸 4.767 父 1010及司=3.802,r3劃 4379M1010。關(guān)鍵詞行星軌道微分方程 Simpson 法四階龍格一庫塔法 matlab一.問題重述水星到太陽的最遠(yuǎn)距離為0.6982父1011m 此時(shí)水星繞太陽運(yùn)行的線速度為3.886父104m/s。試求問題一水星到太陽的最近距離問題二水星繞太陽運(yùn)行的周期問題三從遠(yuǎn)日點(diǎn)開始的第 50天(地球大)結(jié)束時(shí)水星的位置并畫出軌道曲線二.問題分析求水星到太陽的最近距離以及水星繞太陽運(yùn)行的周期等，需要先將水星軌道方程2r求出，因此可以根據(jù) Newton 第二定律及萬有引力定律-嚶生七二m二|，建立微rdt分方程模型，將原問題轉(zhuǎn)化為

3、求解帶有初值條件的微分方程問題，進(jìn)而采用解析法或數(shù)值方法求解遠(yuǎn)日點(diǎn)和周期。三.模型假設(shè)1.水星運(yùn)行的軌道是以太陽為一個(gè)焦點(diǎn)的橢圓2 .從太陽指向水星的線段在單位時(shí)間內(nèi)掃過的面積相等3 .水星運(yùn)行周期的平方與其運(yùn)行軌道橢圓長(zhǎng)軸的立方之比為常量四.符號(hào)系統(tǒng)1.Vo水星在遠(yuǎn)日點(diǎn)的線速度2.M太陽的質(zhì)量3.m 水星的質(zhì)量4.ro水星在遠(yuǎn)日點(diǎn)的距離5.T周期五.建立模型與求解模型一水星的軌跡方程設(shè)太陽中心所在的位置為復(fù)平面的原點(diǎn) O,在時(shí)刻 t,水星位于 Z(t)=rei0所表示的點(diǎn) P。這里r=r(t),6=t)均為 t 的函數(shù)，分別表示Z(t)的模和輻角。于是水星的速度為烏=曳 89+舟9電一雙曳+

4、ir吧),加速度為dtdt出dt出咚=9(曰-r(西)2)+i(r曰+2支妙J(1.1),而太陽對(duì)行星的引力依萬有引dt:dtdtdtdtdt_力定律，大小為mMG,方向由行星位置 P 指向太陽的中心 O,故為-mMGe，其rr中 M=1.989M1030(kg)為太陽的質(zhì)量， m 為水星的質(zhì)量， G=6.672101(Nm2/kg2)為萬有引力常數(shù)。行星的速度為 v0。那么就有初值條件:drdt由dt因此問題轉(zhuǎn)化為求解帶初值問題的微分方程組ro_d日drd8rr+2dt2dtdtir,er。8 匕=0_V0t-0一0又將rd|+2如此=0兩邊同乘以 r,即得色“2吧）=0,從而r23=c,d

5、t2dtdtdtdtdttt1cdic.t其中 C1=rv,這樣有向線段 OP 在時(shí)間事內(nèi)掃過的面積等于 1r2tdt=E,t2dt2依 Newton 定律，我們得到一嗎。rd27=md4(1.2),將 1.1)代入(1.2),dt然后比較實(shí)部與虛部，就有dSr-dtd2r羽cdrd?,2 二 0dtdt,du、2MG-r()=dtr2這是兩個(gè)未知函數(shù)的二階微分方程組。條件。假設(shè)當(dāng) t=0 時(shí)，行星正處于遠(yuǎn)日點(diǎn),在確定某一行星軌道時(shí)，需要加上定解而遠(yuǎn)日點(diǎn)位于正實(shí)軸上，距原點(diǎn)。為 r0,t=0tz0二r0二0tZ0=0_V0t=0=一d2rdt2-注)dtdrdtt=0=0史dt這個(gè)正是 Kep

6、ler 的第二定律，從太陽指向水星的線段在單位時(shí)間內(nèi)掃過的面積相等。,222將（1.3）代入d4r（變）2=-學(xué)得町-c3=-嚶，于是我們可以得到水dt2出r2dt2r3r2星運(yùn)行的較為簡(jiǎn)單形式的數(shù)學(xué)模型：22drG2 一dtrddtdrd,du、d2ud122d2u二C|()=-C1TV=-C1u77dtdtd?ddtd12/2,與，化簡(jiǎn)后為 d4+u=-(1.4),其中 p=u,引進(jìn) u=u-1，立r2d12pMGp即可以求出u-，=u=Acos（8-日0）,這里 A 和日。是待定的常數(shù)P以寫為p1-ecos(1-2)這個(gè)就是水星的軌道方程，是一條平面二次曲線。由于水星繞太陽運(yùn)行，故必有0

7、e1。由于 r 在 t=0 時(shí)取道最大值 r0（遠(yuǎn)日點(diǎn)），這個(gè)就意味著此時(shí)函數(shù)cos。-備）取道最大值 1.于是就有為=0,e=1-2，從而軌跡方程為r0r=p0對(duì)于水星而言，r=0.698M1011（m）,v。=3.886M104（m/s）,又水星的1-ecos近日點(diǎn)到太陽的距離rm=p=p。依據(jù)已知數(shù)據(jù)，可知1-ecos二1eMGdrdtt=o=0為了求得行星的軌跡方程，要消去變量t,令 r=1,那么變=與可以改寫為udtrd2一二c1udt,22drG2dtr將上式代入記e=Ap,上式可152G210G=r0v0ft2.713x10(m2/s),p5.547=10(m),MG模型二水星的

8、運(yùn)行周期(1.5)將有關(guān)數(shù)據(jù)代入，易得 T：7.6025106(s):87.9919(d)模型三水星的位置由于水星的運(yùn)行滿足 Kepler 第二定律，則該式可改寫為 Er2d8=CQt,從而可得p2dl-t)C1(1-ecos?)如果我們要求 t=TI時(shí)相應(yīng)的日和 r,則意味著首先要解方程，F(xiàn)R)=CT1,p其中在求出了 t=T1 時(shí)的日后，立即可以由 r=-J 得到相應(yīng)的 r。1 -ecos-卜面用數(shù)值方法求解水星的位置1 .Simpson 法e=1-fc0.2055,從而計(jì)算水星到太陽的最近距離為rm4.60161010(m)設(shè)水星的周期為Kepler第二定律，我們有T1r2d%tJCIT

9、02dt2(1.4)上式左端為水星軌跡橢圓所圍的面積，記為 S,由于橢圓的半長(zhǎng)軸a，半短1一e2軸 b=,p,從而有 S=nab=3231e(1-e2)2將上式代入式(1.4),解得 T=27Tp23C1(1-e2)2C1T1CIT,2abp3(1-e2)23.5703由被積函數(shù)12 的恒正性可知F(e)單調(diào)， 從而方程 F(e)=C3 的根必(1-ecosi)2p存在且唯一。取 Ae=h,ek=kh(k=1,2,.),記 Fk=F(%)。若Fn3 理，七書之理，PP那么 3 位于斗與斗卅之間，在 h 適當(dāng)小時(shí)，可取 6&en0nn-in計(jì)算F(6)可采用不同的數(shù)值積分法，本文采用 Simps

10、on 法，取步長(zhǎng) h=0.001,具體求解過程見附錄一，最后結(jié)果為三=3.791,r，4.76710102 .基于壓縮映像的求根方法我們引入水星軌道橢圓的參數(shù)方程，由于橢圓的半長(zhǎng)軸a=P2,半短軸1-eb=,從而中心到焦點(diǎn)的距離為 la2-b2=ae。因左焦點(diǎn)為原點(diǎn)，故橢圓中1 -e2心位于(ae,0),于是得到參數(shù)方程,Lx;a(ecos)y-bsin它們與r,e的關(guān)系為x2y2=r2,-=tanx一.、_W44 一“4此式可改與成 C/t=(xy-yx)d=abesin(：小)-esinabG當(dāng) t=TI時(shí)解方程esin-CarCT記九=C-,g*=esin、那么上式即邛=gF),就是說要

11、去求函 ab數(shù)gT)的不動(dòng)點(diǎn)，求解方程不動(dòng)點(diǎn)可以采用簡(jiǎn)單迭代法，對(duì)于水星，我們已計(jì)算出e*0.2055,由于 e 很小，因此迭代收斂理論上可以很快，當(dāng)時(shí)間從遠(yuǎn)日點(diǎn)開始的第 50 天結(jié)束時(shí)，意味著 T1=0.432 父 107(s),從而不妨取中0=0,于是1=-esin0=3.57032-esin1=3.65575=-esin4=3.67476=-esin5=3.6747故一二 3.6747由式x=a(e+cosy=bsin中，x2+y2=r2,y=tan6,可以計(jì)算出相應(yīng)的 0,x即由一 bsin:tarn=0.75849a(ecos)得 8=0.64891,而 6=8+n=3.791此時(shí)的

12、距離 r 為r=Ja(e+cos)2十bsin/2=4.7668父1010(m3.經(jīng)典四階 Runge-Kutte 法由我們將由最初的微分方程組求解水星的位置，方程組見下亞上 MGdt2-r3一一 r2dc,dt-r2drt田=0drI 八出t=0=0令qd,那么我們可以得到一階微分方程組:dt_2_dqCiMGdrqdtd-Cidtr2t_0drdtt=0日,t衛(wèi)二0若記這個(gè)微分方程組中方程的右端依次為 Q(t,q,r,0),R(t,q,r,8)和 S(t,q,r,6),則相應(yīng)的四階 Runge-Kutte 迭代格式法為h,qk1=qk(KI2K22K3K4)6和久2h八1一人一(N12叫2

13、N3N4)6h,這里對(duì)于qk=qk+(K1+2K2+2K3+K4),有6初值為 q。=0,r。=0,=0,則對(duì)于給定的步長(zhǎng)值 h,類似可以逐步計(jì)算一系列的 qkJk,由于行星繞著太陽運(yùn)行，只需取W2 兀,5 書之 2 兀，而取得行星軌道系列點(diǎn)的近似坐標(biāo)(rk,ek),再通過極坐標(biāo)與直角坐標(biāo)的轉(zhuǎn)換，繼而可以繪出軌道曲線。通過 matlab 編程求解得日=3.802,r比4.779父101,軌道曲線如下KI=Q(tk,qk,r-k)K2hhKi=Q(tk-,qk-/k22.hLi-T,A.嗎2K4h二Q(tk-,qkhK22，rk2)-Q(tkh,qkhK3,rk也,，h)程序見附錄二。六.模型推

14、廣本文建立的微分方程模型對(duì)于求解行星繞日運(yùn)行軌道具有廣泛的應(yīng)用空間，只需給出行星的遠(yuǎn)日點(diǎn)和在遠(yuǎn)日點(diǎn)的運(yùn)行線速度即可計(jì)算出軌道方程，用數(shù)學(xué)軟件繪出近似的軌道曲線，對(duì)于研究天體運(yùn)行有所幫助。止匕外，本文采用的求解微分方程的數(shù)值方法，具有較為快速且準(zhǔn)確的收斂效果，可以用來求解其他類似的微分方程模型。七.參考文獻(xiàn)11 樂經(jīng)良，數(shù)學(xué)實(shí)驗(yàn)，北京，高等教育出版社，1999 年 10 月2周品，matlab 數(shù)值分析，北京，機(jī)械工業(yè)出版社，2009 年 1 月附錄附錄一functionq1=y2(x)q1=(1-0.2055*cos(x)A-2;h=0.001;k=1;x=h*k;f=quad(y2,0,x)

15、whilef(k)3.8091k=k+1;x=k*h;f(k)=quad(y2,0,x);endx附錄二formatlongc1=2.7132e15;M=1.989e30;G=6.672e-11;Q=inline(2.7132e15A2/(rA3)-1.989e30*6.672e-11/(rA2);R=inline(q);S=inline(2.7132e15/(rA2);q=0;r=0.6982e11;theta=0;t=0;k=1;h=0.001e7;whiletheta=2*piK1=Q(r);L1=R(q);N1=S(r);K2=Q(r+h/2*L1);L2=R(q+h/2*K1);N2

16、=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);%水星任意位置的橫坐標(biāo)yy(k)=rr(k)*sin(ee(k);%水星任意位置的縱坐標(biāo)k=k+1;end;plot(xx,yy)%畫

17、出水星的軌道曲線text(0,0,太陽)text(0.6982e11,0,遠(yuǎn)日點(diǎn))text(-4.6078e+010,0,近日點(diǎn))holdon;plot(0,0,r.,MarkerSize,20);holdoffholdon;plot(0.6982e11,0,r，MarkerSize,20);holdoffholdon;plot(-4.6078e+010,0,r.,MarkerSize,20);holdofftitle(水星繞太陽運(yùn)行的軌道曲線)clcq=0;r=0.6982e11;theta=0;t=0;k=1;h=0.001e7;whilet=50*24*3600%求水星自遠(yuǎn)日點(diǎn)開始第 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