劉 林湯靖師侯錫云

(1南京大學(xué)天文與空間科學(xué)學(xué)院南京210046)

(2南京大學(xué)空間環(huán)境與航天動力學(xué)研究所南京210093)

1 引言

關(guān)于定點在地-月系共線平動點附近的探測器軌道外推(或稱預(yù)報)問題,并不是一個新問題,該類型的飛行器,實際上就是一顆帶有一定軌道特征的遠地衛(wèi)星.但是,盡管理論上繞地運行周期與月球一致,但它受到的月球引力影響卻不是一個小擾動,而是一個幾乎與地球引力作用同等重要的外力源.因此,該問題也無法處理成一個簡單的受攝二體問題.如果仍然像構(gòu)造平動點特殊軌道(如暈軌道)問題那樣,在地-月系中來處理該問題,也無法像衛(wèi)星運動那樣處理成一個簡單的受攝二體問題,其難點有如下兩個方面:

(1)即使基本模型采用簡單的圓型限制性三體問題[1?3],其基本解也較復(fù)雜.如果考慮月球軌道偏心率,處理成橢圓型限制性三體問題,不僅基本解更復(fù)雜,而且無濟于事.因為太陽引力攝動影響與月球軌道偏心率的效應(yīng)相當(dāng),無論采用上述哪類基本模型,要在此基礎(chǔ)上構(gòu)造攝動解,都不可能像一般受攝二體問題那樣簡單,無法滿足實際應(yīng)用的需求.

(2)即使形式上處理成一般的受攝二體問題,但由于探測器到月球的距離與探測器到地球的距離相當(dāng),這是處理第三體攝動問題的一個難點.

鑒于上述分析,對于地面測控而言,宜在J2000地心天球坐標(biāo)系中采用數(shù)值方法實現(xiàn)地-月系平動點探測器的軌道外推,其數(shù)學(xué)模型就是一個表面形式上的受攝二體問題,狀態(tài)運動方程如下:

2 外力分析

對于航天器的軌道外推問題,實為一個定量問題,因涉及的各物理量有不同的量綱和大小,不便于問題的分析和表達,擬采用無量綱單位(可簡稱為歸一化單位)來處理問題.以地球衛(wèi)星的軌道運動為背景,包括月球之類的遠地衛(wèi)星(如地-月系的L1或L2點探測器)軌道,通常采用的歸一化單位系統(tǒng),即長度單位[L]、質(zhì)量單位[M]和時間單位[T]如下:

其中,GE是地心引力常數(shù),ae是地球參考橢球體的赤道半徑,時間單位[T]是導(dǎo)出單位,目的是使該計算單位系統(tǒng)中,引力常數(shù)G=1和中心天體引力常數(shù)μ=GE=1.如果地球引力模型采用當(dāng)今地固坐標(biāo)系所對應(yīng)的WGS-84系統(tǒng),則相應(yīng)的地球動力學(xué)扁率為J2=1.082636022×10?3.

在采用上述歸一化單位后,狀態(tài)運動方程(1)式變?yōu)橄铝行问?

該式中的是月球、太陽等第三體的無量綱質(zhì)量,即

其中,GM和GS分別為月心和日心引力常數(shù).

其中,是航天器到第j個天體的位置矢量,是航天器到第j個天體的距離,是地心天球坐標(biāo)系中月球、太陽等第三體坐標(biāo)矢量.

2.1 各大天體的引力攝動量級估計

地球的引力常數(shù)即GE=398600.4418 km3/s2,而月球、太陽、水星、金星、火星、木星和土星的引力常數(shù)相對地球的大小依次分別為0.0123000383、332946.050895、0.055273598、0.814998108、0.107446732、317.8942053和95.1574041.

攝動量級的近似估計式為

該式中的即各大天體相對地球的質(zhì)量比,其值即上面給出的引力常數(shù)比,r和分別為探測器和攝動天體到地球的距離,估計中可取的平均值,但因水星軌道的偏心率較大(e>0.2),對應(yīng)的r′值將分別按近日點和遠日點兩種狀態(tài)設(shè)定平均距離.

考慮定點在地-月系共線平動點L1和L2的探測器到地球的距離r分別為

其中是地-月平均距離.由(6)式給出上述各大天體對L1和L2點探測器軌道的攝動量級分別如下:

2.2 地球非球形引力攝動量級

主項(扁率J2)的攝動量級估計

2.3 太陽光壓的攝動量級

定點在地-月系共線平動點L1和L2的一個通常尺度(包括質(zhì)量和承受光壓的等效截面)的探測器,太陽光壓攝動的量級估計如下:

其中,κ=1.44,面質(zhì)比S/m=109,ρ⊙為1 au處的太陽輻射壓強,ρ⊙=0.3169×10?17.

根據(jù)上述對外力因素的量級估計,對于定點在地-月系共線平動點L1和L2或其附近的探測器軌道,若考慮10?6以上的攝動因素,相應(yīng)的力模型中只需要考慮如下攝動源:月球、太陽和金星的質(zhì)點引力,地球非球形引力位的扁率J2和太陽光壓,其中最主要的是月球和太陽的質(zhì)點引力.

3 定點在地-月系L1和L2點的探測器軌道概貌

這里給出地-月系圓型限制性三體問題模型下兩個簡單的算例,初始歷元為2016-09-30UTC0:00:0.0,對應(yīng)的TDT(地球力學(xué)時)是57661.0007891667(MJD),探測器定點在地-月系L1和L2點的各一條軌道上.經(jīng)簡單坐標(biāo)轉(zhuǎn)換,即可獲得J2000.0地心天球坐標(biāo)系中相應(yīng)的兩條軌道的初始位置(x,y,z)、速度(˙x,˙y,˙z)和相應(yīng)的軌道根數(shù),分別列于表1和表2.

表1L 1點和L 2點的位置和速度Table 1 The positions and velocities of the points L1 and L2

表2L 1點和L 2點的軌道根數(shù)Table 2 The orbital elements of the points L1 and L2

兩條軌道的圖像見圖1–2,這表明在初始時刻瞬間,實際上都是一條偏心率較大的環(huán)繞地球的橢圓軌道,探測器均處于該軌道的遠地點和近地點(對讀者而言,這一特點是容易理解的,無需做過多解釋),兩圖中的坐標(biāo)單位ae是地球參考橢球體的赤道半徑.就地-月+探測器系統(tǒng)而言,這都是初始瞬時軌道,而在月球的引力作用下,探測器與月球軌道“同步”做相同的圓軌道運動.

4 平動點軌道外推問題的相關(guān)分析

4.1 平動點軌道分析采用的數(shù)學(xué)模型

第1節(jié)引言中已指出,應(yīng)在J2000地心天球坐標(biāo)系中處理其軌道運動問題,并采用數(shù)值方法實現(xiàn)相應(yīng)的軌道外推.為了定量顯示這類軌道外推中誤差傳播狀態(tài)的主要特征,顯然應(yīng)選擇地-月-日+探測器的質(zhì)點引力系統(tǒng),相應(yīng)的狀態(tài)運動方程即

其中是地球引力加速度:而和分別為月球和太陽的無量綱質(zhì)量,見前面的(4)式.相應(yīng)的月球和太陽引力攝動加速度、的具體形式分別為

圖1L1點初始軌道在J2000地心天球坐標(biāo)系中(赤道面內(nèi))的圖像Fig.1 The figure(in the equatorial plane)of point L1initial orbit in the J2000 reference system

圖2L 2點初始軌道在J2000地心天球坐標(biāo)系中(赤道面內(nèi))的圖像Fig.2 The figure(in the equatorial plane)of pointL 2initial orbit in the J2000 reference system

4.2 平動點軌道類型的選擇

這里考查的平動點軌道,包括如下3種類型:

(1)初始時刻定點在地-月系的L1點或L2點處的平動點軌道,以下簡稱該類型軌道為L1點軌道或L2點軌道;

(2)初始時刻定點在地-月系的L1點或L2點附近的halo軌道上,以下簡稱該類型軌道為L1點暈軌道或L2點暈軌道;

(3)初始時刻定點在地-月系的L1點或L2點附近的Lissajous軌道上,以下簡稱該類型軌道為L1點Lissajous軌道或L2點Lissajous軌道.

經(jīng)初步設(shè)計(對應(yīng)所采用的質(zhì)點引力系統(tǒng))分別給出6條軌道,在J2000地心天球坐標(biāo)系中各對應(yīng)的軌道初值如下:所有初始時刻對應(yīng)歷元為2016-09-30UTC0:00:0.0(相應(yīng)TDT的MJD為57661.0007891667),位置、速度和相應(yīng)的軌道根數(shù)分別列于表3和表4.表中的軌道類型1、2、···、6依次為L1點軌道、L1點暈軌道、L1點Lissajous軌道、L2點軌道、L2點暈軌道、L2點Lissajous軌道,表4–8類同.

表3 6條軌道的位置和速度Table 3 The positions and velocities of 6 oribits

表4 6條軌道的軌道根數(shù)Table 4 The orbital elements of 6 oribits

下面首先對上述6類軌道作7d和27d的軌道外推,給出一個誤差傳播的定量輪廓,在此基礎(chǔ)上再作定性分析.

4.3 6條平動點軌道外推7d的狀態(tài)

為簡單起見(也不失一般性),在考查誤差傳播中,將初始誤差全部集中在最重要的軌道半長徑上(根據(jù)目前定軌的實際狀況,軌道半長徑的精度為10 m量級),7d的軌道外推結(jié)果列于表5–6.

表5 平動點軌道外推7d的軌道狀態(tài)Table 5 The states of the libration point orbits propagated for 7days

表6 平動點軌道外推7d的空間位置狀態(tài)Table 6 The positions of the libration point orbits prop a gated for 7days

4.4 6條軌道外推27d的狀態(tài)

盡管探測器的定點只是近似的,實際運行過程中必須通過不斷軌控才能保持這類軌道,外推弧段增長至27d,只是為了進一步了解這類特殊軌道的動力學(xué)特征及其相應(yīng)的誤差傳播規(guī)律.初始誤差仍全部集中在軌道半長徑上,27d的軌道外推結(jié)果列于表7–8,其中,第2和3兩類軌道(即L1點暈軌道和L1點Lissajous軌道)只外推了22 d,其原因?qū)⒃?.5小節(jié)中具體說明.

表7 平動點軌道外推27d的軌道狀態(tài)Table 7 The states of the libration point orbits propagated for 27days

對于探測器的軌道運動而言,通常所說的長期位置預(yù)報和短期位置預(yù)報中的長期或短期并不是簡單的時間間隔,而是運行弧段的長短.因此,為了比較上述各條軌道之間外推誤差傳播的定量大小,需要了解它們的軌道運行周期,這6條特殊軌道的初始運行周期TS值依次為

L1類:TS=13.708768 d,12.221658 d,13.113512 d;

L2類:TS=126.58784 d,58.847196 d,93.089267d.

由此便于了解軌道外推7d和27d對上述6條軌道所對應(yīng)的弧段長短.

在已知6條軌道自身的動力學(xué)特性和軌道外推弧段長短的前提下,不難看出表5所列出的外推位置誤差所反映的一些動力學(xué)規(guī)律,基本上可歸納如下:

(1)軌道外推7d均為短弧,位置誤差都在1 km以內(nèi);

(2)對于L1點軌道和L2點軌道,外推7d或27d,位置誤差的累積仍不嚴(yán)重,其誤差傳播的特征就是一般Kepler運動特征的反映.由于初值誤差(?a0=10 m)較小,實為小擾動,既不會激發(fā)其初值不穩(wěn)定的固有特征,又不會明顯改變短弧段誤差累積的效果,反而周期性的效果比長期累積效應(yīng)更明顯,見表5中軌道1(即L1點軌道)的誤差定量狀態(tài).

(3)對于L1點和L2點的halo軌道和Lissajous軌道而言,由于其位置已經(jīng)“遠離”不穩(wěn)定平動點L1和L2,而嚴(yán)格的halo軌道和Lissajous軌道設(shè)計又無法實現(xiàn),探測器的定點只是近似的,在同樣是?a0=10 m的初值誤差情況下,已不能再簡單地只看作對halo軌道和Lissajous軌道的小擾動,而更重要的是對平動點的大擾動起作用,在不太長的外推弧段內(nèi),其平動點本身的不穩(wěn)定特征即顯現(xiàn)無遺.見表5中的軌道2、軌道3、軌道5和軌道6,特別是軌道2和3,相對而言,27d的弧段顯得更長,外推弧段超過22 d后,其軌道偏心率很快就會達到e≈1.0的狀態(tài).

表8 平動點軌道外推27d的空間位置狀態(tài)Table 8 The positions of the libration point orbits propagated for 27days

4.5 平動點軌道運行中誤差傳播狀態(tài)的基本特征

從上述6條軌道的外推計算結(jié)果已能看出拉格朗日點軌道外推中位置誤差傳播的基本特征,首先將由?a0=10 m導(dǎo)致的位置誤差集中列于表9.

表9 拉格朗日點軌道外推位置誤差的定量狀態(tài)(單位:km)Table 9 The quantitative state of position error s of the Lagrange point orbit propagation(unit:km)

綜上幾小節(jié),就平動點探測器軌道運行誤差傳播狀態(tài)的簡單計算和分析,可以表明:在地-月系L1和L2平動點軌道設(shè)計中,確實很難實現(xiàn)較長弧段的無動力控制運行軌道,而不是設(shè)計者本身的問題.那么,在這樣較短的弧段內(nèi),就軌道預(yù)報而言,要達到較高的精度是容易實現(xiàn)的,本節(jié)最后一小節(jié)將會給出具體算例,并有實測結(jié)果的檢驗.

4.6L1點暈軌道的運行狀態(tài)

為了進一步揭示平動點飛行器軌道外推中誤差傳播的動力學(xué)機制,這里再給出一條與上述L1點暈軌道相近的環(huán)繞地球運行的大橢圓逆行軌道,同樣對應(yīng)一遠地衛(wèi)星.初始時刻仍為2016-09-30UTC0:00:0.0,初始位置和速度及相應(yīng)的軌道根數(shù)分別列于表10和表11.該軌道的初始運行周期為12.208831 d.對此軌道同作27d的外推,計算結(jié)果列于表12–13.

表10L1點暈軌道的位置和速度Table 10 The position and velocity of the halo orbit at pointL 1

表11L 1點暈軌道的軌道根數(shù)Table 11 The orbital elements of the halo or bit at pointL 1

表12L 1點暈軌道外推27d的軌道狀態(tài)Table 12 The states of the halo orbit prop a gated for 27days at pointL 1

表13L 1點暈軌道外推27d的空間位置狀態(tài)Table 13 The positions of the halo orbit prop a gated for 27days at pointL 1

這樣一條有別于L1點暈軌道繞地運行的大橢圓逆行軌道,本質(zhì)上就是一條普通的Kepler軌道,其固有的初值不穩(wěn)定性在繞地運行27d(僅2圈)的短弧內(nèi)不會有明顯體現(xiàn),這是一個常規(guī)問題,本文引進這一算例的目的,是從另一角度體現(xiàn)地-月系平動點軌道誤差快速傳播的固有不穩(wěn)定性特征.至于逆行軌道自身的運動規(guī)律及其動力學(xué)特征,已超出本文論述的范疇,不再介紹,如有需要,可見文獻[4]及其有關(guān)作者的研究工作.

4.7 平動點軌道的定軌和短弧預(yù)報精度的檢驗

既然嚴(yán)格的拉格朗日點軌道設(shè)計無法實現(xiàn),探測器的定點只是近似的,運行過程中必須通過不斷的軌控才能保持,那么,對于地面測控和星上控制,只能從短弧角度來考慮問題.

對于平動點軌道的短弧定軌和軌道預(yù)報而言,相對地球低軌或高軌衛(wèi)星的同類問題,實無任何特殊困難和特別需要處理的難題.采用南京大學(xué)空間環(huán)境與航天動力學(xué)研究所自主編寫的定軌軟件和利用國內(nèi)的USB(Unified S-band)測量數(shù)據(jù),在沒有任何其他輔助信息的前提下,對嫦娥3號的相關(guān)任務(wù)探測器進行了定軌,并與北京航天飛行控制中心的事后定軌結(jié)果作了對比.在此定軌的基礎(chǔ)上,采用非常簡單的數(shù)值外推方法(只考慮地、月、日三體的質(zhì)點引力和簡單的光壓模型,外推中的6個軌道初值采用相關(guān)任務(wù)的定軌結(jié)果)進行了軌道預(yù)報,毫無困難地達到了較高精度.略去不必要的細節(jié)說明,將有關(guān)結(jié)果一并列于表14–15.

表14L 2點暈軌道短弧外推3 d與事后精密定軌結(jié)果的比較Table 14 The comparison of two precise orbit determination method s for the pointL 2orbit after propagating 3 days

表15L 2點暈軌道短弧外推7d與事后精密定軌結(jié)果的比較Table 15 The comparison of two precise orbit determination methods for the pointL 2orbit after prop a gating 7days

表14–15中的精密定軌A和B分別對應(yīng)北京航天飛行控制中心和南京大學(xué)空間環(huán)境與航天動力學(xué)研究所的結(jié)果,表中的結(jié)果基本上已能說明問題,但為了讓讀者對這類探測器的定軌和外推精度有更清晰的了解,下面進一步作些必要的說明.

(1)關(guān)于光壓模型,在不了解探測器的具體細節(jié)情況下,作者們根據(jù)獨立定軌中獲得的有關(guān)估計值,獲得了包括衛(wèi)星表面熱性能在內(nèi)的等效面質(zhì)比,從而給出了相應(yīng)的經(jīng)驗?zāi)Ｐ?一個等效的平面模型.

(2)盡管沒有具體給出兩個單位的定軌(包括測量數(shù)據(jù))細節(jié),但表14–15所給出的計算結(jié)果,已能說明本文要體現(xiàn)的這類特殊軌道的定軌和外推精度了.因為表中的結(jié)果是外推3 d和7d與事后精密定軌結(jié)果的比較,且兩個單位的定軌結(jié)果之差基本上在500 m之內(nèi),這樣的比較更能體現(xiàn)兩個單位定軌結(jié)果的真實性以及本文所采用的力模型的合理性.

上述計算結(jié)果和兩點補充說明充分表明:盡管這類探測器的軌道特殊,初值誤差的傳播程度遠比一般的環(huán)繞型探測器的軌道顯著,但相應(yīng)的短弧定軌和高精度軌道預(yù)報并無特殊困難.

4.8 關(guān)于兩種坐標(biāo)系中平動點探測器軌道形式的轉(zhuǎn)換問題

就定軌和預(yù)報的需求,顯然是在J2000地心天球坐標(biāo)系中進行相關(guān)問題的處理,而對這類具有特殊性質(zhì)的軌道,軌控又必須通過相應(yīng)的地-月系旋轉(zhuǎn)坐標(biāo)系來處理.這就涉及到兩種坐標(biāo)系之間的轉(zhuǎn)換問題,其本身是容易實現(xiàn)的,而在具體的航天任務(wù)中,各有關(guān)部門根據(jù)實際需求,對相應(yīng)的地-月系旋轉(zhuǎn)坐標(biāo)系實有不同取法,故這里不再做相應(yīng)討論.