劉 林湯靖師侯錫云
(1南京大學(xué)天文與空間科學(xué)學(xué)院南京210046)
(2南京大學(xué)空間環(huán)境與航天動(dòng)力學(xué)研究所南京210093)
關(guān)于定點(diǎn)在地-月系共線平動(dòng)點(diǎn)附近的探測(cè)器軌道外推(或稱預(yù)報(bào))問(wèn)題,并不是一個(gè)新問(wèn)題,該類型的飛行器,實(shí)際上就是一顆帶有一定軌道特征的遠(yuǎn)地衛(wèi)星.但是,盡管理論上繞地運(yùn)行周期與月球一致,但它受到的月球引力影響卻不是一個(gè)小擾動(dòng),而是一個(gè)幾乎與地球引力作用同等重要的外力源.因此,該問(wèn)題也無(wú)法處理成一個(gè)簡(jiǎn)單的受攝二體問(wèn)題.如果仍然像構(gòu)造平動(dòng)點(diǎn)特殊軌道(如暈軌道)問(wèn)題那樣,在地-月系中來(lái)處理該問(wèn)題,也無(wú)法像衛(wèi)星運(yùn)動(dòng)那樣處理成一個(gè)簡(jiǎn)單的受攝二體問(wèn)題,其難點(diǎn)有如下兩個(gè)方面:
(1)即使基本模型采用簡(jiǎn)單的圓型限制性三體問(wèn)題[1?3],其基本解也較復(fù)雜.如果考慮月球軌道偏心率,處理成橢圓型限制性三體問(wèn)題,不僅基本解更復(fù)雜,而且無(wú)濟(jì)于事.因?yàn)樘?yáng)引力攝動(dòng)影響與月球軌道偏心率的效應(yīng)相當(dāng),無(wú)論采用上述哪類基本模型,要在此基礎(chǔ)上構(gòu)造攝動(dòng)解,都不可能像一般受攝二體問(wèn)題那樣簡(jiǎn)單,無(wú)法滿足實(shí)際應(yīng)用的需求.
(2)即使形式上處理成一般的受攝二體問(wèn)題,但由于探測(cè)器到月球的距離與探測(cè)器到地球的距離相當(dāng),這是處理第三體攝動(dòng)問(wèn)題的一個(gè)難點(diǎn).
鑒于上述分析,對(duì)于地面測(cè)控而言,宜在J2000地心天球坐標(biāo)系中采用數(shù)值方法實(shí)現(xiàn)地-月系平動(dòng)點(diǎn)探測(cè)器的軌道外推,其數(shù)學(xué)模型就是一個(gè)表面形式上的受攝二體問(wèn)題,狀態(tài)運(yùn)動(dòng)方程如下:
對(duì)于航天器的軌道外推問(wèn)題,實(shí)為一個(gè)定量問(wèn)題,因涉及的各物理量有不同的量綱和大小,不便于問(wèn)題的分析和表達(dá),擬采用無(wú)量綱單位(可簡(jiǎn)稱為歸一化單位)來(lái)處理問(wèn)題.以地球衛(wèi)星的軌道運(yùn)動(dòng)為背景,包括月球之類的遠(yuǎn)地衛(wèi)星(如地-月系的L1或L2點(diǎn)探測(cè)器)軌道,通常采用的歸一化單位系統(tǒng),即長(zhǎng)度單位[L]、質(zhì)量單位[M]和時(shí)間單位[T]如下:
其中,GE是地心引力常數(shù),ae是地球參考橢球體的赤道半徑,時(shí)間單位[T]是導(dǎo)出單位,目的是使該計(jì)算單位系統(tǒng)中,引力常數(shù)G=1和中心天體引力常數(shù)μ=GE=1.如果地球引力模型采用當(dāng)今地固坐標(biāo)系所對(duì)應(yīng)的WGS-84系統(tǒng),則相應(yīng)的地球動(dòng)力學(xué)扁率為J2=1.082636022×10?3.
在采用上述歸一化單位后,狀態(tài)運(yùn)動(dòng)方程(1)式變?yōu)橄铝行问?
該式中的是月球、太陽(yáng)等第三體的無(wú)量綱質(zhì)量,即
其中,GM和GS分別為月心和日心引力常數(shù).
其中,是航天器到第j個(gè)天體的位置矢量,是航天器到第j個(gè)天體的距離,是地心天球坐標(biāo)系中月球、太陽(yáng)等第三體坐標(biāo)矢量.
地球的引力常數(shù)即GE=398600.4418 km3/s2,而月球、太陽(yáng)、水星、金星、火星、木星和土星的引力常數(shù)相對(duì)地球的大小依次分別為0.0123000383、332946.050895、0.055273598、0.814998108、0.107446732、317.8942053和95.1574041.
攝動(dòng)量級(jí)的近似估計(jì)式為
該式中的即各大天體相對(duì)地球的質(zhì)量比,其值即上面給出的引力常數(shù)比,r和分別為探測(cè)器和攝動(dòng)天體到地球的距離,估計(jì)中可取的平均值,但因水星軌道的偏心率較大(e>0.2),對(duì)應(yīng)的r′值將分別按近日點(diǎn)和遠(yuǎn)日點(diǎn)兩種狀態(tài)設(shè)定平均距離.
考慮定點(diǎn)在地-月系共線平動(dòng)點(diǎn)L1和L2的探測(cè)器到地球的距離r分別為
其中是地-月平均距離.由(6)式給出上述各大天體對(duì)L1和L2點(diǎn)探測(cè)器軌道的攝動(dòng)量級(jí)分別如下:
主項(xiàng)(扁率J2)的攝動(dòng)量級(jí)估計(jì)
定點(diǎn)在地-月系共線平動(dòng)點(diǎn)L1和L2的一個(gè)通常尺度(包括質(zhì)量和承受光壓的等效截面)的探測(cè)器,太陽(yáng)光壓攝動(dòng)的量級(jí)估計(jì)如下:
其中,κ=1.44,面質(zhì)比S/m=109,ρ⊙為1 au處的太陽(yáng)輻射壓強(qiáng),ρ⊙=0.3169×10?17.
根據(jù)上述對(duì)外力因素的量級(jí)估計(jì),對(duì)于定點(diǎn)在地-月系共線平動(dòng)點(diǎn)L1和L2或其附近的探測(cè)器軌道,若考慮10?6以上的攝動(dòng)因素,相應(yīng)的力模型中只需要考慮如下攝動(dòng)源:月球、太陽(yáng)和金星的質(zhì)點(diǎn)引力,地球非球形引力位的扁率J2和太陽(yáng)光壓,其中最主要的是月球和太陽(yáng)的質(zhì)點(diǎn)引力.
這里給出地-月系圓型限制性三體問(wèn)題模型下兩個(gè)簡(jiǎn)單的算例,初始?xì)v元為2016-09-30UTC0:00:0.0,對(duì)應(yīng)的TDT(地球力學(xué)時(shí))是57661.0007891667(MJD),探測(cè)器定點(diǎn)在地-月系L1和L2點(diǎn)的各一條軌道上.經(jīng)簡(jiǎn)單坐標(biāo)轉(zhuǎn)換,即可獲得J2000.0地心天球坐標(biāo)系中相應(yīng)的兩條軌道的初始位置(x,y,z)、速度(˙x,˙y,˙z)和相應(yīng)的軌道根數(shù),分別列于表1和表2.
表1L 1點(diǎn)和L 2點(diǎn)的位置和速度Table 1 The positions and velocities of the points L1 and L2
表2L 1點(diǎn)和L 2點(diǎn)的軌道根數(shù)Table 2 The orbital elements of the points L1 and L2
兩條軌道的圖像見(jiàn)圖1–2,這表明在初始時(shí)刻瞬間,實(shí)際上都是一條偏心率較大的環(huán)繞地球的橢圓軌道,探測(cè)器均處于該軌道的遠(yuǎn)地點(diǎn)和近地點(diǎn)(對(duì)讀者而言,這一特點(diǎn)是容易理解的,無(wú)需做過(guò)多解釋),兩圖中的坐標(biāo)單位ae是地球參考橢球體的赤道半徑.就地-月+探測(cè)器系統(tǒng)而言,這都是初始瞬時(shí)軌道,而在月球的引力作用下,探測(cè)器與月球軌道“同步”做相同的圓軌道運(yùn)動(dòng).
第1節(jié)引言中已指出,應(yīng)在J2000地心天球坐標(biāo)系中處理其軌道運(yùn)動(dòng)問(wèn)題,并采用數(shù)值方法實(shí)現(xiàn)相應(yīng)的軌道外推.為了定量顯示這類軌道外推中誤差傳播狀態(tài)的主要特征,顯然應(yīng)選擇地-月-日+探測(cè)器的質(zhì)點(diǎn)引力系統(tǒng),相應(yīng)的狀態(tài)運(yùn)動(dòng)方程即
其中是地球引力加速度:而和分別為月球和太陽(yáng)的無(wú)量綱質(zhì)量,見(jiàn)前面的(4)式.相應(yīng)的月球和太陽(yáng)引力攝動(dòng)加速度、的具體形式分別為
圖1L1點(diǎn)初始軌道在J2000地心天球坐標(biāo)系中(赤道面內(nèi))的圖像Fig.1 The figure(in the equatorial plane)of point L1initial orbit in the J2000 reference system
圖2L 2點(diǎn)初始軌道在J2000地心天球坐標(biāo)系中(赤道面內(nèi))的圖像Fig.2 The figure(in the equatorial plane)of pointL 2initial orbit in the J2000 reference system
這里考查的平動(dòng)點(diǎn)軌道,包括如下3種類型:
(1)初始時(shí)刻定點(diǎn)在地-月系的L1點(diǎn)或L2點(diǎn)處的平動(dòng)點(diǎn)軌道,以下簡(jiǎn)稱該類型軌道為L(zhǎng)1點(diǎn)軌道或L2點(diǎn)軌道;
(2)初始時(shí)刻定點(diǎn)在地-月系的L1點(diǎn)或L2點(diǎn)附近的halo軌道上,以下簡(jiǎn)稱該類型軌道為L(zhǎng)1點(diǎn)暈軌道或L2點(diǎn)暈軌道;
(3)初始時(shí)刻定點(diǎn)在地-月系的L1點(diǎn)或L2點(diǎn)附近的Lissajous軌道上,以下簡(jiǎn)稱該類型軌道為L(zhǎng)1點(diǎn)Lissajous軌道或L2點(diǎn)Lissajous軌道.
經(jīng)初步設(shè)計(jì)(對(duì)應(yīng)所采用的質(zhì)點(diǎn)引力系統(tǒng))分別給出6條軌道,在J2000地心天球坐標(biāo)系中各對(duì)應(yīng)的軌道初值如下:所有初始時(shí)刻對(duì)應(yīng)歷元為2016-09-30UTC0:00:0.0(相應(yīng)TDT的MJD為57661.0007891667),位置、速度和相應(yīng)的軌道根數(shù)分別列于表3和表4.表中的軌道類型1、2、···、6依次為L(zhǎng)1點(diǎn)軌道、L1點(diǎn)暈軌道、L1點(diǎn)Lissajous軌道、L2點(diǎn)軌道、L2點(diǎn)暈軌道、L2點(diǎn)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
下面首先對(duì)上述6類軌道作7d和27d的軌道外推,給出一個(gè)誤差傳播的定量輪廓,在此基礎(chǔ)上再作定性分析.
為簡(jiǎn)單起見(jiàn)(也不失一般性),在考查誤差傳播中,將初始誤差全部集中在最重要的軌道半長(zhǎng)徑上(根據(jù)目前定軌的實(shí)際狀況,軌道半長(zhǎng)徑的精度為10 m量級(jí)),7d的軌道外推結(jié)果列于表5–6.
表5 平動(dòng)點(diǎn)軌道外推7d的軌道狀態(tài)Table 5 The states of the libration point orbits propagated for 7days
表6 平動(dòng)點(diǎn)軌道外推7d的空間位置狀態(tài)Table 6 The positions of the libration point orbits prop a gated for 7days
盡管探測(cè)器的定點(diǎn)只是近似的,實(shí)際運(yùn)行過(guò)程中必須通過(guò)不斷軌控才能保持這類軌道,外推弧段增長(zhǎng)至27d,只是為了進(jìn)一步了解這類特殊軌道的動(dòng)力學(xué)特征及其相應(yīng)的誤差傳播規(guī)律.初始誤差仍全部集中在軌道半長(zhǎng)徑上,27d的軌道外推結(jié)果列于表7–8,其中,第2和3兩類軌道(即L1點(diǎn)暈軌道和L1點(diǎn)Lissajous軌道)只外推了22 d,其原因?qū)⒃?.5小節(jié)中具體說(shuō)明.
表7 平動(dòng)點(diǎn)軌道外推27d的軌道狀態(tài)Table 7 The states of the libration point orbits propagated for 27days
對(duì)于探測(cè)器的軌道運(yùn)動(dòng)而言,通常所說(shuō)的長(zhǎng)期位置預(yù)報(bào)和短期位置預(yù)報(bào)中的長(zhǎng)期或短期并不是簡(jiǎn)單的時(shí)間間隔,而是運(yùn)行弧段的長(zhǎng)短.因此,為了比較上述各條軌道之間外推誤差傳播的定量大小,需要了解它們的軌道運(yùn)行周期,這6條特殊軌道的初始運(yùn)行周期TS值依次為
L1類:TS=13.708768 d,12.221658 d,13.113512 d;
L2類:TS=126.58784 d,58.847196 d,93.089267d.
由此便于了解軌道外推7d和27d對(duì)上述6條軌道所對(duì)應(yīng)的弧段長(zhǎng)短.
在已知6條軌道自身的動(dòng)力學(xué)特性和軌道外推弧段長(zhǎng)短的前提下,不難看出表5所列出的外推位置誤差所反映的一些動(dòng)力學(xué)規(guī)律,基本上可歸納如下:
(1)軌道外推7d均為短弧,位置誤差都在1 km以內(nèi);
(2)對(duì)于L1點(diǎn)軌道和L2點(diǎn)軌道,外推7d或27d,位置誤差的累積仍不嚴(yán)重,其誤差傳播的特征就是一般Kepler運(yùn)動(dòng)特征的反映.由于初值誤差(?a0=10 m)較小,實(shí)為小擾動(dòng),既不會(huì)激發(fā)其初值不穩(wěn)定的固有特征,又不會(huì)明顯改變短弧段誤差累積的效果,反而周期性的效果比長(zhǎng)期累積效應(yīng)更明顯,見(jiàn)表5中軌道1(即L1點(diǎn)軌道)的誤差定量狀態(tài).
(3)對(duì)于L1點(diǎn)和L2點(diǎn)的halo軌道和Lissajous軌道而言,由于其位置已經(jīng)“遠(yuǎn)離”不穩(wěn)定平動(dòng)點(diǎn)L1和L2,而嚴(yán)格的halo軌道和Lissajous軌道設(shè)計(jì)又無(wú)法實(shí)現(xiàn),探測(cè)器的定點(diǎn)只是近似的,在同樣是?a0=10 m的初值誤差情況下,已不能再簡(jiǎn)單地只看作對(duì)halo軌道和Lissajous軌道的小擾動(dòng),而更重要的是對(duì)平動(dòng)點(diǎn)的大擾動(dòng)起作用,在不太長(zhǎng)的外推弧段內(nèi),其平動(dòng)點(diǎn)本身的不穩(wěn)定特征即顯現(xiàn)無(wú)遺.見(jiàn)表5中的軌道2、軌道3、軌道5和軌道6,特別是軌道2和3,相對(duì)而言,27d的弧段顯得更長(zhǎng),外推弧段超過(guò)22 d后,其軌道偏心率很快就會(huì)達(dá)到e≈1.0的狀態(tài).
表8 平動(dòng)點(diǎn)軌道外推27d的空間位置狀態(tài)Table 8 The positions of the libration point orbits propagated for 27days
從上述6條軌道的外推計(jì)算結(jié)果已能看出拉格朗日點(diǎn)軌道外推中位置誤差傳播的基本特征,首先將由?a0=10 m導(dǎo)致的位置誤差集中列于表9.
表9 拉格朗日點(diǎn)軌道外推位置誤差的定量狀態(tài)(單位:km)Table 9 The quantitative state of position error s of the Lagrange point orbit propagation(unit:km)
綜上幾小節(jié),就平動(dòng)點(diǎn)探測(cè)器軌道運(yùn)行誤差傳播狀態(tài)的簡(jiǎn)單計(jì)算和分析,可以表明:在地-月系L1和L2平動(dòng)點(diǎn)軌道設(shè)計(jì)中,確實(shí)很難實(shí)現(xiàn)較長(zhǎng)弧段的無(wú)動(dòng)力控制運(yùn)行軌道,而不是設(shè)計(jì)者本身的問(wèn)題.那么,在這樣較短的弧段內(nèi),就軌道預(yù)報(bào)而言,要達(dá)到較高的精度是容易實(shí)現(xiàn)的,本節(jié)最后一小節(jié)將會(huì)給出具體算例,并有實(shí)測(cè)結(jié)果的檢驗(yàn).
為了進(jìn)一步揭示平動(dòng)點(diǎn)飛行器軌道外推中誤差傳播的動(dòng)力學(xué)機(jī)制,這里再給出一條與上述L1點(diǎn)暈軌道相近的環(huán)繞地球運(yùn)行的大橢圓逆行軌道,同樣對(duì)應(yīng)一遠(yuǎn)地衛(wèi)星.初始時(shí)刻仍為2016-09-30UTC0:00:0.0,初始位置和速度及相應(yīng)的軌道根數(shù)分別列于表10和表11.該軌道的初始運(yùn)行周期為12.208831 d.對(duì)此軌道同作27d的外推,計(jì)算結(jié)果列于表12–13.
表10L1點(diǎn)暈軌道的位置和速度Table 10 The position and velocity of the halo orbit at pointL 1
表11L 1點(diǎn)暈軌道的軌道根數(shù)Table 11 The orbital elements of the halo or bit at pointL 1
表12L 1點(diǎn)暈軌道外推27d的軌道狀態(tài)Table 12 The states of the halo orbit prop a gated for 27days at pointL 1
表13L 1點(diǎn)暈軌道外推27d的空間位置狀態(tài)Table 13 The positions of the halo orbit prop a gated for 27days at pointL 1
這樣一條有別于L1點(diǎn)暈軌道繞地運(yùn)行的大橢圓逆行軌道,本質(zhì)上就是一條普通的Kepler軌道,其固有的初值不穩(wěn)定性在繞地運(yùn)行27d(僅2圈)的短弧內(nèi)不會(huì)有明顯體現(xiàn),這是一個(gè)常規(guī)問(wèn)題,本文引進(jìn)這一算例的目的,是從另一角度體現(xiàn)地-月系平動(dòng)點(diǎn)軌道誤差快速傳播的固有不穩(wěn)定性特征.至于逆行軌道自身的運(yùn)動(dòng)規(guī)律及其動(dòng)力學(xué)特征,已超出本文論述的范疇,不再介紹,如有需要,可見(jiàn)文獻(xiàn)[4]及其有關(guān)作者的研究工作.
既然嚴(yán)格的拉格朗日點(diǎn)軌道設(shè)計(jì)無(wú)法實(shí)現(xiàn),探測(cè)器的定點(diǎn)只是近似的,運(yùn)行過(guò)程中必須通過(guò)不斷的軌控才能保持,那么,對(duì)于地面測(cè)控和星上控制,只能從短弧角度來(lái)考慮問(wèn)題.
對(duì)于平動(dòng)點(diǎn)軌道的短弧定軌和軌道預(yù)報(bào)而言,相對(duì)地球低軌或高軌衛(wèi)星的同類問(wèn)題,實(shí)無(wú)任何特殊困難和特別需要處理的難題.采用南京大學(xué)空間環(huán)境與航天動(dòng)力學(xué)研究所自主編寫(xiě)的定軌軟件和利用國(guó)內(nèi)的USB(Unified S-band)測(cè)量數(shù)據(jù),在沒(méi)有任何其他輔助信息的前提下,對(duì)嫦娥3號(hào)的相關(guān)任務(wù)探測(cè)器進(jìn)行了定軌,并與北京航天飛行控制中心的事后定軌結(jié)果作了對(duì)比.在此定軌的基礎(chǔ)上,采用非常簡(jiǎn)單的數(shù)值外推方法(只考慮地、月、日三體的質(zhì)點(diǎn)引力和簡(jiǎn)單的光壓模型,外推中的6個(gè)軌道初值采用相關(guān)任務(wù)的定軌結(jié)果)進(jìn)行了軌道預(yù)報(bào),毫無(wú)困難地達(dá)到了較高精度.略去不必要的細(xì)節(jié)說(shuō)明,將有關(guān)結(jié)果一并列于表14–15.
表14L 2點(diǎn)暈軌道短弧外推3 d與事后精密定軌結(jié)果的比較Table 14 The comparison of two precise orbit determination method s for the pointL 2orbit after propagating 3 days
表15L 2點(diǎn)暈軌道短弧外推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分別對(duì)應(yīng)北京航天飛行控制中心和南京大學(xué)空間環(huán)境與航天動(dòng)力學(xué)研究所的結(jié)果,表中的結(jié)果基本上已能說(shuō)明問(wèn)題,但為了讓讀者對(duì)這類探測(cè)器的定軌和外推精度有更清晰的了解,下面進(jìn)一步作些必要的說(shuō)明.
(1)關(guān)于光壓模型,在不了解探測(cè)器的具體細(xì)節(jié)情況下,作者們根據(jù)獨(dú)立定軌中獲得的有關(guān)估計(jì)值,獲得了包括衛(wèi)星表面熱性能在內(nèi)的等效面質(zhì)比,從而給出了相應(yīng)的經(jīng)驗(yàn)?zāi)P?一個(gè)等效的平面模型.
(2)盡管沒(méi)有具體給出兩個(gè)單位的定軌(包括測(cè)量數(shù)據(jù))細(xì)節(jié),但表14–15所給出的計(jì)算結(jié)果,已能說(shuō)明本文要體現(xiàn)的這類特殊軌道的定軌和外推精度了.因?yàn)楸碇械慕Y(jié)果是外推3 d和7d與事后精密定軌結(jié)果的比較,且兩個(gè)單位的定軌結(jié)果之差基本上在500 m之內(nèi),這樣的比較更能體現(xiàn)兩個(gè)單位定軌結(jié)果的真實(shí)性以及本文所采用的力模型的合理性.
上述計(jì)算結(jié)果和兩點(diǎn)補(bǔ)充說(shuō)明充分表明:盡管這類探測(cè)器的軌道特殊,初值誤差的傳播程度遠(yuǎn)比一般的環(huán)繞型探測(cè)器的軌道顯著,但相應(yīng)的短弧定軌和高精度軌道預(yù)報(bào)并無(wú)特殊困難.
就定軌和預(yù)報(bào)的需求,顯然是在J2000地心天球坐標(biāo)系中進(jìn)行相關(guān)問(wèn)題的處理,而對(duì)這類具有特殊性質(zhì)的軌道,軌控又必須通過(guò)相應(yīng)的地-月系旋轉(zhuǎn)坐標(biāo)系來(lái)處理.這就涉及到兩種坐標(biāo)系之間的轉(zhuǎn)換問(wèn)題,其本身是容易實(shí)現(xiàn)的,而在具體的航天任務(wù)中,各有關(guān)部門(mén)根據(jù)實(shí)際需求,對(duì)相應(yīng)的地-月系旋轉(zhuǎn)坐標(biāo)系實(shí)有不同取法,故這里不再做相應(yīng)討論.