DRO計算及其在地月系中的攝動力研究

吳小婧，曾凌川，鞏應奎

（中國科學院空天信息創(chuàng)新研究院，北京100094）

隨著深空探測任務的開展，有許多不同的關鍵技術需要研究，首先需要解決的關鍵問題就是航天器的軌道設計問題。相較于近地空間，深空探測器所處的引力環(huán)境具有多樣性，不再局限于經典的二體開普勒軌道，其基本動力學模型可簡化為一個受攝的圓形限制性三體問題（Circular Restricted Three-Body Problem，CRTBP）［1］。在CRTBP中，以地月（Earth-Moon，EM）系為例，由于地月系相較日地系，第二主天體的軌道偏心率更大、太陽作為第三引力體對地月系的引力影響更強烈，因此地月系三體軌道設計將面臨更大挑戰(zhàn)［2］。本文嘗試以地月系為背景對遠距逆行軌道（Distant Retrograde Orbit，DRO）設 計 做 初 步研究。

DRO屬于圓形限制性三體問題中一類特殊的平面對稱軌道［3］，其圍繞較小主天體運動，在會合坐標系中，運行方向與主天體的公轉方向相反，因此，軌道是逆行的。理想的周期DRO只存在于圓形限制性三體問題中，在實際力模型中，由于攝動力的影響，周期DRO變?yōu)閿M周期DRO［4-5］。當這些擬周期軌道受到比其他三體軌道族更大的擾動時，仍能保持穩(wěn)定［6-7］，因此非常適于深空導航組網、中繼通信以及科學數據采集等任務。

NASA 曾計劃在其木星冰月軌道任務（NASA’s Jupiter Icy Moon Orbiter，JIMO）中使用DRO來平衡木星和其衛(wèi)星的引力攝動，該任務原計劃于2015年左右發(fā)射［8］。NASA的小行星重定向任務（ARM）曾計劃捕獲一顆近地小行星［9-10］，并使用低推力航天器將其拖拽到DRO上［11］，由于DRO的長期穩(wěn)定性，當受到如非對稱引力場、太陽輻射壓以及來自太陽系其他天體的引力等攝動時，在未來的幾十年甚至更長時間，小行星都不需要位置保持機動［10，12］。學者們在基于DRO的深空探測方面也做了大量的研究工作。Ocampo和Rosborough提出利用日地系DRO構建太陽風暴預警系統(tǒng)的應用設想［13］，介紹了日地系中DRO的幾種類型，并針對日地系DRO設計了脈沖轉移軌道。針對此應用設想，Demeyer和Gurfil研究了從地球到日地DRO的轉移軌道的設計方法［14］；Scott和Spencer利用微分校正法計算了轉移到日地DRO的轉移軌道族［15］。Stramacchia等提出利用日地DRO部署天基望遠鏡網絡探測近地危險小行星的應用設想［16］。徐明與徐世杰提出了將中繼衛(wèi)星布置在地月DRO的設想［17］，研究了空間雙圓模型下DRO的軌道穩(wěn)定性，并研究了轉移軌道的設計方法。Murakami和Yamanaka基于未來在月球附近建造一個新的空間站或燃料補給港的設想，對在地月DRO上的交會對接技術進行了研究［18］。Conte等對直接從近地軌道（Low Earth Orbit，LEO）轉移到火星和采用月球DRO作為中轉站轉移到火星2種地火轉移方式開展了比較研究［19］。這些研究均是在CRTBP模型下，并沒有考慮月球的引力攝動。

此外，為了設計DRO，Lara用攝動分析方法得到了DRO的低階解析解和高階解析解［20］，并利 用 高 階 解 析 解 計 算 了 DRO 周 期 軌 道［21］。Bezrouk和Parker研究了在長達數萬年的時間里，地月系中幾個不同大小的DRO的演化問題［22］。

以上研究主要針對需要完全放手、長時間隔離軌道的任務，例如火星樣本返回任務或小行星重定向任務等對軌道精度要求不高的任務，而且高階模型存在形式復雜和計算繁瑣等不足。本文針對需要高軌道精度的深空導航和通信服務，著重高精度動力學模型的建立，通過數值方法研究了DRO的設計問題。

考慮到設計符合動力學環(huán)境需求的標稱任務軌道不僅能有效地節(jié)省軌道保持的燃料［23］，而且可以為在軌導航制導提供理想的參考軌道。因此本文針對實際工程應用中標稱任務軌道的設計問題，以地月系為背景，研究了DRO的設計方法，分析了在實際力環(huán)境下地月DRO擬周期軌道的軌道特性及主要攝動因素。首先，針對地月DRO的動力學特征，利用構造流函數法求解其周期軌道，并通過周期軌道在地月會合坐標系和地月慣性坐標系中的對應關系，從共振角度對DRO軌道特性進行了分析和總結；然后，在實際力模型下，研究地月DRO的主要攝動因素，在考慮太陽引力攝動和月球軌道偏心率攝動的星歷模型中，分析地月DRO的非開普勒特征與穩(wěn)定性。仿真研究表明，構造流函數法可以有效應用于DRO的數值求解，從而擺脫了經典微分校正法迭代初值對近似解析解的依賴，可以在近似解析解失效的更廣大區(qū)域求得DRO周期軌道的數值解。不同于以往的研究，本文通過不同攝動力模型下DRO軌道特性的分析對比，得到影響地月DRO的主要攝動因素及其影響程度，可以為設計符合動力學環(huán)境需求的標稱任務軌道提供理論依據。

1 動力學模型

1.1 坐標系定義

1）地月會合坐標系O-xyz。在圓形限制性三體條件下的地月會合坐標系O-xyz，其原點在地月系的質心O，x軸由質心指向月心方向，z軸指向系統(tǒng)角速度方向，y軸與x軸、z軸構成右手坐標系［24］。

2）地月慣性坐標系OE-XYZ。其原點在地球質心OE，月球在地月參考平面上繞地球作圓周運動，t＝0時刻，地月慣性坐標系和地月會合坐標系的軸線對齊，月球位于＋X軸上［24］。

1.2 圓形限制性三體問題

根據開普勒定律，假設地球和月球這2個主天體都可看成1個質點，它們圍繞彼此作圓周運動。第3個小天體航天器，受主天體的引力影響，但不影響主天體的軌道。對各單位進行無量綱化處理，用兩主天體的質量和作為單位質量，較小主天體的質量為m′。對于地月系，m′＝0.012 150 5。兩主天體間的平均距離無量綱化為1。本文中采用的地月距離為384 400 km。

CRTBP地月會合坐標系會隨著兩主天體旋轉，由于它們之間的距離是恒定的，所以這2個主天體在會合坐標系中看起來是靜止的。無量綱化處理后，原點位于2個主天體的質心，2個主天體位于x軸上，較大主天體位于-μ處，較小主天體位于1-μ處，其中，μ＝m′＝0.012 150 5。z軸沿兩主天體軌道的法向，y軸與x軸和z軸構成右手坐標系。對其進行無量綱化處理使得兩主天體的軌道周期為2π。對于地月系，CRTBP的時間間隔2π對應于27.28 d。CRTBP的運動方程 為［24］

式中：

其中：d1、d2分別為航天器到兩主天體的距離；r＝（x，y，z）為航天器的位置向量；v＝（vx，vy，vz）為速度向量。

CRTBP方程（1）有一個運動積分，被稱為雅可比積分：

式中：

雅可比積分實際上是質點總能量的量度，雅可比積分越大，則質點能量越小，反之亦然。雅可比積分雖然不能確定質點的全部運動規(guī)律，但通過研究雅可比積分，仍能得到有關質點運動的豐富信息。

2 利用構造流函數法計算DRO

2.1 構造流函數法

由于限制性三體問題的非線性，目前只有依靠微分方程的數值解法才能得到滿足實際任務精度要求的DRO周期解。經典的微分校正法是利用軌道的狀態(tài)轉移矩陣構造牛頓迭代過程的，迭代初值可以由CRTBP的線性化模型給出，也可利用Lindstedt-Poincaré攝動分析方法得到更高階的初值。本文針對地月DRO的動力學特征，利用構造流函數法求解其周期軌道。該方法的優(yōu)點是不需要以近似解作為迭代初值，因此不受非線性的影響，可以在近似解析解失效的更廣大區(qū)域求取周期軌道的數值解。

下面給出流函數的概念［25-26］。設系統(tǒng)式（1）以X（t0）＝X0為初始狀態(tài)的解軌跡為φ（t，t0），則φ（t，t0）：X（t0）→X（t）定義了動力學系統(tǒng)式（1）的流映射，它將t0時刻的初始狀態(tài)映射到t時刻狀態(tài)。這里的φ被稱為流函數。本文使用一種特殊的流函數，它能將初始狀態(tài)X0映射到未來某時刻位于x軸上的另一狀態(tài)，將其表示為φ（k，X0），其中k為流函數與x軸相交的次數，相對于φ以給定時間作為積分終止條件，則φ以與x軸的第k次相交作為終止條件。由于DRO屬于平面圓形限制性三體問題［3］，因此系統(tǒng)狀態(tài)向量可表示為X＝［x y vxvy］T，為了計算φ（k，X0）的值，引入如下判別函數：

由于DRO軌道具有沿X軸的對稱性，其與X軸相交時必須滿足垂直相交條件，即y方向速度為零。根據這個特征，判別函數具有簡單形式，只是取狀態(tài)X的y分量，且初始積分狀態(tài)可以取為X0＝［x00 0 vy0］T。若限定雅可比積分為C，則根據雅可比積分公式可得

對初始狀態(tài)X0逐步積分，并在每一步積分中計算判別函數的值，當其符號第k次發(fā)生改變時，表明流函數第k次跨越了x軸，這時以Crit（X）符號改變前后的t值作為初值，恰好可以應用Brent方法來解方程：

設解為t*，則

顯然，φ（k，X0）的解算包含了一個數值積分和一個數值求根過程，雖然無法得到其解析表達式，卻仍然可以用數值方法研究φ（k，X0）的性質。若以φvx表示φ的vx分量，則計算雅可比積分C下的DRO軌道等價于求解方程：

式中：φ*（x）為以x為自變量的軌跡函數在vx方向的分量?？紤]到DRO屬于圓形限制性三體問題中的平面對稱軌道［3］，X表達式為

式（9）同樣可以應用數值求根方法解出，迭代初值當然也可以由線性化解析解給出，但本文采用更方便的方法，即觀察φ*（x）的曲線（見圖1（a）），尋找其與x軸的交點，由圖中得到交點P的近似坐標作為初值。圖1（b）中L1和L2分別表示地月系的拉格朗日點L1、L2。

構造流函數法利用了CRTBP的對稱性，它不但可以計算DRO軌道，而且可以系統(tǒng)地計算一大類具有x軸對稱性的周期軌道，如圖2所示。即利用構造流函數法計算了1 ∶3地球共振軌道在地月會合坐標系和地月慣性坐標系中的運動軌跡，在地月會合坐標系下的一個軌道周期里，航天器在地月慣性坐標系下環(huán)繞地球3圈，Re為地球半徑。從圖2可以明顯看出，共振軌道在地月會合坐標系中是閉合的，其拱線在地月慣性坐標系中旋轉。這類軌道可以看做地月系中繞地球旋轉的DRO軌道。

圖1 地月系的φ*（x）曲線和以點P為初值迭代得到的DRO（C＝2.93）Fig.1 φ*（x）curve of the Earth-Moon system and DRO obtained by iteration with point P as initial value（C＝2.93）

2.2 DRO計算及其非開普勒特征分析

本文通過在圓形限制性三體模型下延拓雅可比積分C計算地月DRO周期軌道族。DRO存在較大的幅值范圍，當幅值較小時，DRO完全可以看作低軌的環(huán)月軌道，此時，DRO具有較高的雅可比常數，接近L1（CL1＝3.18834）和L2（CL2＝3.17216）的雅可比常數。隨著雅可比常數的減小，DRO幅值將逐漸增大，越靠近地球時對應的雅可比常數越小。圖3在很大的范圍內繪制了地月CRTBP中的DRO周期軌道族。每個軌道的顏色表示雅可比常數，由右邊的顏色欄指定。很明顯，一些軌道已經延伸到離月球很遠的地方，因此被稱為DRO。越靠近月球，DRO的周期越短，隨著軌道周期通過與月球1∶4、1∶3和1∶2的共振，DRO的大小增加，并逐漸接近與月球1 ∶1的共振。

圖2 地球DRO（1∶3）共振軌道Fig.2 The Earth resonant orbit of DRO（1∶3）

圖3 DRO周期軌道族與雅可比常數CFig.3 DRO periodic family and Jacobi constant C

以上是在地月會合坐標系下對DRO進行了計算。然而，在航天器任務分析中，了解地月慣性坐標系解的性質是至關重要的。為此，本文定義了一個地月慣性參考框架，其具有以下性質：參考框架的坐標原點是地球，月球在地月參考平面上繞地球作圓周運動。在t＝0時刻，地月慣性坐標系和地月會合坐標系的軸線對齊。因此，在t＝0時刻，月球位于地月慣性坐標系的＋X軸上。由CRTBP的對稱性可知，判別式（7）的解t*即為1／2軌道周期，聯合式（9）和式（10）通過數值方法得到雅可比積分C與DRO軌道周期的對應關系，從而可找到1∶2、1 ∶3及1 ∶4的月球共振軌道，圖4（a）、圖4（b）、圖4（c）分別展示了它們在地月慣性坐標系下的運動軌跡。從圖4可以看出，這3種共振軌道均為非開普勒軌道，它們在地月慣性坐標系中是周期閉合的。圖4（d）是共振比非整數的DRO在地月慣性坐標系中的運行情況，其軌跡不再是封閉的曲線，而看起來像是在月球軌道附近打水漂。

圖4 DRO在地月慣性坐標系中的運動軌跡Fig.4 Motion trajectory of DRO in the Earth-Moon inertial coordinate system

3 攝動力模型下的DRO仿真分析

第2節(jié)DRO研究是基于圓形限制性三體模型的。本節(jié)研究地月攝動力模型下，DRO的動力學特性。地月系中，攝動因素主要包括太陽引力、太陽光壓、月球非對稱引力、其他行星的引力等。表1顯示了以月球為中心的DRO所受加速度的種類及其近似量級。在實際力模型下，CRTBP中理想的周期DRO不再呈現周期性，運行軌跡為擬周期DRO。

為了衡量每種攝動對DRO動力學模型精度的影響，選用STK中的地月擬周期DRO軌道初值［x，y，z，vx，vy，vz］＝［155 932 km，0 km，0 km，0 km／s，-0.869 709 km／s，0 km／s］外 推1年，每 次 從高精度模型（考慮了月球非對稱引力、木星引力、金星引力、太陽光壓、太陽引力）中去掉一種攝動項，然后和高精度模型下的DRO擬周期軌道形態(tài)進行比較。具體模型如下所示。

表1 攝動平均量級Tab le 1 Average m agnitude of perturbation

模型1：高精度模型。

模型2：高精度模型中不考慮月球非對稱引力。

模型3：高精度模型中不考慮木星引力和金星引力。

模型4：高精度模型中不考慮太陽光壓。

模型5：高精度模型中不考慮太陽引力。

模型1～模型5的外推結果如圖5所示。從圖5可以看出，模型1到模型4，擬周期DRO形態(tài)無顯著變化，到模型5，由于模型中去掉了太陽引力，擬周期軌道形態(tài)發(fā)生了顯著變化。因此圖5的仿真結果表明，太陽引力是影響擬周期DRO形態(tài)最主要的攝動因素，在實際工程建模中必須考慮太陽引力攝動的影響。

接下來通過模型1的高精度模型和模型2～模型5的各種攝動力模型的比較來分析各種攝動模型的軌道外推精度。具體的，對地月會坐標系中的x坐標、y坐標以及到月球距離在各種攝動模型與高精度模型之間做差，而后統(tǒng)計軌道差的均方差作為各種攝動模型軌道外推精度的估計值。圖6為模型1和模型2的差值，即不考慮月球非對稱引力攝動的動力學模型軌道外推精度；圖7為模型1和模型3的差值，即不考慮木星、金星引力攝動的動力學模型軌道外推精度；圖8為模型1和模型4的差值，即不考慮太陽光壓的動力學模型軌道外推精度；圖9為模型1和模型5的差值，即不考慮太陽引力的動力學模型軌道外推精度。

圖5 攝動力模型下地月擬周期DROFig.5 The Earth-Moon quasi-periodic DRO under perturbative forcemodel

圖6 模型1和模型2在x、y坐標及到月球距離的偏差Fig.6 Deviation of x，y coordinate and distance to the Moon between Model 1 and Model 2

圖7 模型1和模型3在x、y坐標及到月球距離的偏差Fig.7 Deviation of x，y coordinate and distance to the Moon between Model 1 and Model 3

圖8 模型1和模型4在x、y坐標及到月球距離的偏差Fig.8 Deviation of x，y coordinate and distance to the Moon between Model 1 and Model 4

圖9 模型1和模型5在x、y坐標及到月球距離的偏差Fig.9 Deviation of x，y coordinate and distance to the Moon between Model 1 and Model 5

從圖6～圖9可以看出，曲線隨著積分時間的推進逐漸發(fā)散，下面分別對以上積分的前1、2、3個月的均方差進行統(tǒng)計，統(tǒng)計結果如表2所示。

從表2可以看出，不考慮太陽引力的模型誤差最大，量級為幾千km，其次是不考慮太陽光壓的模型，誤差量級為幾十到幾百km，再次是不考慮月球非對稱引力的模型，誤差為1 km左右，最后為不考慮金星和木星引力的模型，誤差量級不到1 km。在實際工程應用時，在一定精度要求下，一些次要的攝動因素可以忽略不計，本文提出了一種結構簡單的高精度地月DRO動力學模型，下面就該模型的建立和精度進行分析。

對于地月DRO飛行器而言，地月會合坐標系相對地月慣性坐標系之間的旋轉直接反映了太陽引力和月球軌道偏心率對其運行軌跡的影響。因此可以使用標準星歷數據來表示太陽和月球的運動狀態(tài)，獲得星歷模型下的擬周期DRO，從而實現在動力學模型中考慮太陽引力和月球軌道偏心率等攝動因素的目的。

表2 不同攝動模型軌道外推精度的均方差統(tǒng)計結果Tab le 2 Statistical resu lts of m ean square error of trajectory extrapolation p recision for different pertu rbation m odels km

為了檢驗星歷模型的精確性。本文在全力模型（考慮了月球軌道偏心率、太陽引力、金星引力、木星引力以及太陽光壓5種攝動）和只考慮太陽引力和月球軌道偏心率的星歷模型下對STK中的地月擬周期DRO初值進行了積分，積分器為龍格庫塔7～8階變步長積分器，行星星歷表使用JPL的DE430。圖10展示了全力模型和星歷模型下，航天器坐標y、坐標x以及到月球距離在3個月積分時間內2種模型差值的均方差隨時間的變化情況。

圖10的仿真結果表明，在前10天的積分時間內，星歷模型的誤差在km量級，在前1個月的積分時間內，模型誤差在幾十km量級，隨著積分時間的增加，模型誤差逐漸增大。雖然在軌道外推時，模型誤差很快達到km量級，但在地月系這樣大尺度的空間范圍內，仍然可以利用星歷模型來分析DRO在實際力環(huán)境中的運動特性，為任務軌道設計提供理論依據。在工程實用中，導航測量弧段和軌道保持策略的確定要和模型精度統(tǒng)一進行閉環(huán)考慮。

圖10 星歷模型和全力模型x、y坐標及到月球距離偏差的均方差Fig.10 Mean square error of deviation between ephemeris model and full forcemodel of x，y coordinate and distance to the Moon

4 結 論

通過流函數法在圓形限制性三體模型下求解得到DRO周期軌道族，并分析了軌道的非開普勒特性以及軌道形狀在地月會合坐標系與地月慣性坐標系中隨共振比的變化情況。其后分析了在攝動力模型下DRO航天器的運動特性以及影響DRO的主要攝動因素及其影響程度。仿真分析表明：

1）利用流函數法可以在近似解析解失效的更廣大區(qū)域求取周期DRO的數值解。

2）整數共振比DRO（例如1∶2、1∶3、1∶4等）在地月慣性坐標系中是閉合的非開普勒周期軌道，而非整數共振比DRO在地月慣性坐標系中則是不閉合的，看起來像是在月球軌道附近打水漂。這是由于在CRTBP模型下，只有整數比的共振軌道航天器在地月慣性坐標系下與月球的相對位置變化呈周期性，即如果共振比為1 ∶n（n為正整數），DRO航天器在旋轉系中逆行繞月n圈后，在地月慣性坐標系中與月球同時完成繞地球一圈，回到起點位置，從而在地月慣性坐標系中形成一條閉合的軌道。

3）攝動力模型下，DRO將演變?yōu)閿M周期軌道，擬周期DRO仍可保持很好的穩(wěn)定性，并且在地月系中，影響DRO軌道穩(wěn)定性的主要攝動力是太陽引力，其次是太陽光壓，接下來是金星、木星等行星引力攝動，分別比前幾項攝動小3～5個數量級，在實際力模型中可以不予考慮。

4）在地月系中，只考慮太陽引力和月球軌道偏心率的星歷模型能夠近似地反映DRO在實際力模型中的運動狀態(tài)，因此在實際工程軌道設計中可以使用星歷模型來完成一些任務特性分析，為軌道設計奠定理論基礎。

本文不同于以往的研究，針對深空導航和通信應用，得到了不同攝動項對DRO動力學模型精度的影響結果，提出了一種結構簡單的高精度DRO動力學模型，可以用于高精度軌道方案設計，或者作為自主導航和軌道保持的在軌動力學使用。但是，本文最后建立的動力學模型沒有計入太陽光壓攝動，在后續(xù)的研究中有待進一步研究探討，以提高模型的精度。