地球?火星引力輔助下的地?火無動力循環軌道

臧詩慧 武 迪

(清華大學航天航空學院，北京100084)

2020 年是火星探測大年，阿聯酋的“希望號”火星探測器、中國的“天問一號”火星探測器以及美國的“毅力號”火星車均成功發射，掀起了火星探測的熱潮?；鹦潜砻孑d人科考站是人類走向火星的必由之路，此類科考站需定期更換人員、補充物資，并要求能夠往返于地球與火星之間的航天器執行相應的客貨運任務。這些航天器需要配備大質量的生命維持系統以及足夠大的船員活動空間，因此每次變軌都會消耗大量的燃料。地?火無動力循環軌道是一種往返于地球與火星之間的無動力軌道，運行在這種軌道上，航天器能夠消耗較少燃料執行地球與火星之間的客貨運任務。由此飛船結構強度的要求能夠有所降低，飛船的使用壽命也可相應延長，可以用于執行地球與火星之間定期的貨運或者客運任務。

目前，火星短期居住站的設計已經涉及到地?火循環軌道，利用運行在這種軌道上的航天器以及地球與火星各自的接駁航天器，可以實現定期的人員更換與物資補充[1-2]。地火循環軌道器的設計與結構也在近年來被提出[3]。Hollister[4]在 1969 年對地球?金星往返軌道的研究中，首次提出循環軌道的概念，并計算得到地球?金星循環軌道,同時首次提出了地?火循環軌道。

1985 年 Aldrin 提出 “Aldrin 循環軌道”，使用兩艘飛船進行地?火與火?地間的運輸，但是航天器與地球、火星交會時的速度比較大，實用性受到限制[5]。McConaghy 等[6]在理想的太陽系模型下，忽略火星引力的作用，利用自由返回軌道的組合，得到了若干地?火循環軌道。Ocampo 基于上述結果，總結出地火循環軌道的系統性搜索方法，使用p,h,s,i四個參數確定一條地火循環軌道，其中，p表示一個循環周期包含的地球火星會合周期數，h表示一個循環周期包含的整圈轉移軌道或半圈轉移軌道數，s表示一個循環周期包含的一般轉移軌道數，i表示一個循環周期內飛船繞日的圈數。優化甩擺角的大小以獲得一條唯一確定的循環軌道。遍歷四個參數，并以有解性、甩擺角比為標準，得到了較為實用的地?火循環軌道集合[7]。一些研究在地火循環軌道中引入推力，Chen 等[8]在地火循環軌道中引入小推力，通過減少航天器與火星交會時的剩余速度以優化循環軌道，Vergaaij 等[9]提出了基于太陽帆動力的地火循環軌道。

2015 年，Rogers 等[10]基于速度杠桿變軌以及小推力過程提出了航天器進入循環軌道的優化方案。2016 年，Anderson[11]研究了在與地球交會時，接駁航天器與循環軌道器交會的軌跡規劃與燃料消耗問題，進一步完善了地火循環軌道的應用研究。

然而，這些模型均忽略了火星引力的影響，默認航天器在火星引力影響球外飛掠。同時部分軌道的效率較低，航天器每運行十幾年才能夠完成一次地球?火星之間的往返。本文將在此基礎上引入地球與火星的引力輔助。多圈轉移軌道會顯著地增加地?火轉移的時間，拉低軌道器的效率，故本文放棄使用整圈轉移軌道，僅使用半圈轉移軌道，尋找更加高效的地?火無動力循環軌道。

1 地?火無動力循環軌道的動力學模型

1.1 動力學模型

本文采用理想化的太陽系動力學模型進行計算，僅考慮太陽的中心引力作用，地球和火星位于共面圓軌道上。在實際軌道設計中，可以將本文模型的計算結果作為初值，采用攝動法求解高精度模型下的無動力循環軌道。因此，本文采用的坐標系和動力學模型定義如下：

(1)太陽位于慣性坐標系原點，僅考慮太陽的中心引力作用，地球與火星位于共面圓軌道上；

(2) 地球軌道半徑為1 AU，周期為1 a，火星軌道周期為1.875 a;

(3)地球的引力常數μe=3.985 00×105km3/s2，火星的引力常數μm=4.281 73×104km3/s2。

根據文獻[6]，火星的軌道周期設為1.875 a。在此情況下，地球與火星的相對幾何關系具有周期性，周期為 15 a。因此，在本文中，地火轉移的窗口搜索設置為轉移初始時地球與末端時火星間的相角差，每隔15 a 均會有一次相同的發射窗口。實際應用中，應將此相角差轉換為真實的發射時刻。

1.2 地火間的 Lambert 轉移

進行地火轉移時，探測器需要在特定的時間Te?m內從地球軌道上的特定位置re轉移到火星軌道上的特定位置rm，這個問題可以用Lambert 轉移的解法求解。本文使用Lambert 問題已有的求解程序[12]，僅考慮飛行器轉移軌道公轉小于一周期的情況。初末位置矢量與轉移時間確定，則Lambert 轉移軌道唯一確定。

以圖1 為例，航天器從re轉移到rm，轉移時間為Te?m，可以使用求解程序得到航天器的初始速度矢量、末端速度矢量，從而得到相應的Lambert 轉移軌道(黑色曲線)。由于地球和火星的軌道為共面圓軌道，因此初始的地火Lambert 轉移軌道由轉移時間Te?m和相角差完全確定。此外，火地Lambert轉移可以進行類似求解，轉移時間記為Tm?e。

圖 1 地火 Lambert 轉移

1.3 半圈轉移軌道

本文以地火半圈轉移軌道作為輔助軌道，從而實現地?火無動力循環軌道的設計與搜索。傳統的自由返回軌道需要在地火Lambert 轉移后，經過火星引力輔助進入Lambert 轉移返回，但此種軌道存在剩余速度過大、發射窗口較少的問題。因此，本文引入了半圈轉移軌道作為輔助軌道，即在進行引力輔助時可以選擇進入和中心天體周期相同、具有半周期自由返回特性的圓軌道，以選擇合適的返回時機降低剩余速度，實現地?火無動力循環。

設在引力輔助過程中，航天器僅受行星引力的影響。航天器在飛掠行星前的速度可以分解為行星公轉速度Vp與航天器相對行星的剩余速度V∞的矢量和。由于航天器飛掠時間較短，飛掠前后Vp近似不變，V∞大小不變，方向改變。調整航天器近星點的赤經赤緯以及高度，可以控制飛掠后V∞的指向。如圖2，由于飛掠后V∞大小不變，方向可以改變，所以飛掠后航天器的速度分布在一個以Vp終點為球心，V∞大小為半徑的球面A上。

圖2 引力輔助中的速度疊加

球面A與以Vp起點為球心、以Vp大小為半徑的球面B相交，形成交線圓，分布在交線圓上的速度矢量大小與Vp相等。如果引力輔助后的航天器的速度矢量分布在這個交線圓上，同時航天器的日心距與行星的日心距相同，根據以下軌道能量的公式，航天器半長軸與行星半徑相同，周期與行星的周期相同。

式中a為軌道的半長軸，μp為中心天體的引力常數，v為航天器的速度大小，r為航天器與中心天體的距離。

由于航天器位置矢量垂直于Vp，過Vp并垂直于位置矢量的平面一定與交線圓有兩個交點，即交線圓上一定存在兩點C、D與位置矢量垂直。則引力輔助后的航天器將進入與行星軌道半徑相同的圓軌道，這個軌道與行星的圓軌道有一定的夾角，其夾角γ可以根據剩余速度求得。

運行在這個軌道上的航天器半年后會再次回到行星公轉的平面，與原行星交會。假如微調軌道，使得航天器在再次交會時與原天體保持一定的距離，不受其引力影響，則航天器可以每半年與行星交會一次。這個軌道被稱為半圈轉移軌道，如圖3。

圖3 半圈轉移軌道

2 地?火無動力循環軌道的設計

文獻[2-7]將循環軌道定義為一個周期后天體和航天器相對幾何位置復原的軌道，其中相對幾何關系包括航天器相對地球的位置以及火星相對地球的位置。如果在這個定義下同時考慮地球與火星的引力輔助作用，則軌道既要滿足時間約束，又要滿足地?火、火?地轉移Lambert 轉移的約束，導致問題約束較多，無法尋找到具有使用價值的循環軌道。

事實上，相對幾何關系復原并不是地?火循環軌道功能的本質要求，只要航天器在這個軌道上能不斷無動力往返于地球與火星之間，那么這條軌道就具有同樣的功能。本文考慮放寬地火無動力循環軌道的定義，不再要求這個軌道具有周期性，能夠不斷無動力往返于地球與火星之間即可。

2.1 約束條件

2.1.1 甩擺角(引力輔助高度約束)

圖 4 中δ為航天器進行引力輔助時的甩擺角。計算公式為

式中，μp為中心天體的引力常數，rp為雙曲線軌道近星點半徑(也稱甩擺半徑)[13]。

圖4 引力輔助中的甩擺角

當剩余速度固定時，甩擺角越大，則近星點高度rp越小。由于現實中的地球與火星并不是一個質點，航天器在飛掠的時候近星點必須在大氣層外，因此甩擺角不能太大。本文要求航天器近星點高度的最小值為200 km。相應的甩擺角為最大甩擺角δmax。

定義甩擺角比Tr=δmax/δ，其中δ為實際甩擺角。如果Tr>1，則引力輔助可行，如果Tr<1，則引力輔助不可行。

2.1.2 剩余速度約束

剩余速度決定了航天器在天體引力場內的軌道能量，從天體表面或者近地軌道出發的航天器要與循環軌道器進行對接，需要加速使軌道能量近似相等，航天器之間才有可能進行對接。循環軌道器的剩余速度越高，則接駁航天器要有更強的變軌能力，如此將會消耗更多的能量。本文考慮將循環軌道器在地球與火星引力場內的剩余速度大小v∞e與v∞m限制在 10 km/s 以內，如果進行一次地?火或火?地轉移后的剩余速度大于10 km/s，則判定此次轉移不符合限制條件，航天器進入半圈轉移軌道等待下一次轉移，直到轉移之后的剩余速度小于10 km/s。

2.2 計算步驟

地?火往返過程可以分為以下5 個步驟：

(1) 從地球轉移到火星，初始轉移角度為θ0，轉移時間為t0；

(2)在火星軌道上進行半圈循環，等待合適的時機進行火?地轉移；

(3) 從火星轉移到地球；

(4)在地球軌道上進行半圈循環，等待合適的時機進行地?火轉移；

(5) 進行下一步地?火轉移。

2.2.1 從地球轉移到火星

設航天器第一步地?火轉移在地球與火星公轉平面內進行，轉移具有兩個自由度，因此定義初始轉移角度θ0與初始轉移時間t0為確定第一步地?火轉移的兩個自由變量。同時定義兩個參數θe(t)與θm(t)，作為描述地球與火星位置的參數，不妨定義θe(0)=0，其中t為航天器在軌道上的運行時間，單位為a。地球與火星在時刻t時的坐標為

式中，rm為火星公轉半徑，單位為AU。地球與火星的公轉速度矢量為

式中，ve和vm分別為地球與火星公轉的線速度大小。見圖5，可得航天器第一次轉移的終態位置R?為

圖 5 第一次地?火轉移 (θ0 =1.3π，t0 =1.62 a)

當t= 0 時地球的坐標R+= (1,0,0) AU，則第一步地?火轉移的Lambert 方程確定，有唯一解，可以由式(9) 解出轉移前后的航天器速度。

式中，L為Lambert 轉移的求解函數。

同時可得第一次地?火轉移開始時，火星與地球的相位差θm(0)

計算得到航天器在地球引力場內剩余速度大小v∞e和航天器在火星引力場內剩余速度大小v∞m

2.2.2 在火星軌道上進行半圈循環

首先計算火?地轉移，以航天器轉移前后位置矢量之間的夾角θa為未知數，θa∈(0,2π)，則第二次轉移的初始位置矢量R+為

終態位置矢量R?為

轉移時間ta為

相應的Lambert 轉移有唯一解

以θ0=1.3π，t0=1.62 a 時為例，畫出關于θa的函數圖像，見圖6。由于航天器在火星引力場內進行的是無動力的引力輔助，因此引力輔助前后航天器的剩余速度不變，要求由圖可知，這個方程存在兩個解，為盡可能減少航天器轉移的時間，本文均取轉移時間較短的解。使用擬牛頓迭代法可以在10 步左右得到θa較為精確的數值解。計算得到火 ? 地轉移的甩擺角δ=?Vpm+Vpm??以及地球引力場內的剩余速度v∞e=|Vpe+?Ve(t0+ta)|。如果δ > δmax即Tr< Trmin或者v∞>10 km/s，則判定此次轉移不可行，進入火星半圈轉移軌道等待下一次窗口(圖7)。如果轉移可行，則直接進入第三步。

圖 7 進入半圈轉移軌道 (θ0 =1.3π，t0 =0.8 a)

探測器進入半圈轉移軌道的速度為

Tr=?Vm+Vh?，如果Tr<1，則此次循環軌道終止，無法進行下一步轉移，退出計算。若Tr>1，則航天器進入半圈轉移軌道。

2.2.3 進行火?地轉移

航天器運行半個周期后，重復火?地轉移的計算，如果可行，則進行火?地轉移(圖8)。如果不可行，則繼續重復第二步的半圈循環，定義Nm為航天器在火星半圈轉移軌道上等待的次數。

圖 8 火地轉移 (θ0 =1.3π，t0 =1.62 a)

為保證航天器等待的時間不過長，限制航天器每次在半圈轉移軌道上的等待次數N不超過10 次。

第四步的計算過程與第二步類似，航天器在地球半圈轉移軌道上循環，直到存在合適的窗口，定義Ne為航天器在地球半圈轉移軌道上的等待次數。之后進行第五步的計算，與第三步計算過程類似，進行下一步地?火轉移。

完成五步計算之后，從第二步開始進行相同步驟的循環，計算下一個地?火循環，直到參數滿足算法中的退出條件，停止計算。

圖9 為計算的流程圖。

3 地?火無動力循環軌道計算結果

3.1 初始參數對范圍

并非所有初始參數對(θ0,t0)都能滿足初始的限制要求。令初始地?火轉移中地球與火星引力場內航天器的剩余速度小于10 km/s，可以得到圖10。圖中實線區域內的初始參數對滿足v∞m<10 km/s，虛線區域內的初始參數滿足v∞e<10 km/s，取交集，得到圖 11。使用圖 11 區域中的初始參數對 (θ0,t0)進行遍歷，尋找域內的地?火無動力循環軌道。

圖9 計算流程圖

圖10 使得剩余速度滿足條件的區域

圖11 初始參數對范圍

3.2 n 次地?火無動力循環軌道

根據流程圖所示 (圖 9) 進行計算 (使用 CPU AMD Ryzen 5 4600H，主頻 3.0 GHz，內存 8.00 GB配置的筆記本電腦，Matlab 程序運行耗時約 7 h)，得到一系列地?火無動力循環軌道。挑選其中循環次數較多的幾個軌道，參數與相關數據列于表1。

以θ0= 1.357π，t0= 0.832 a 為例說明，此軌道為計算結果中除Aldrin 循環軌道外，循環次數最多的軌道(圖12)。計算第一步，可得航天器初始時刻從地球離開時速度矢量為(?6.20,31.77,0)km/s，到達火星時速度矢量為(19.10,?8.50,0)km/s，地球引力場內剩余速度大小為6.48 km/s，火星引力場內的剩余速度大小為3.39 km/s，轉移前θe=0，θm=1.48，轉移后θe=5.23，θm=4.26。

計算第二步，首先判斷直接火?地轉移的可能性，計算發現不存在能夠在一個軌道周期內進行火?地Lambert 轉移的軌道，于是轉入火星半圈轉移軌道，等待合適的窗口。根據式(2)，半循環軌道與火星運行平面的夾角γ= 0.1。循環結束后再次嘗試計算火?地轉移。最終在循環六次 (三個火星年) 后發現合適的窗口。此時θe= 2.87,θm= 4.26，探測器在火星引力場內剩余速度大小未變，仍然是3.39 km/s。

進入第三步，計算火?地轉移的參數，直接代入相應函數，得到航天器在離開火星時的速度矢量為 (18.49,?9.90, 0) km/s，到達地球時速度矢量為 (?20.55,?25.07, 0) km/s，火星引力場內剩余速度大小為3.39 km/s，地球引力場內剩余速度大小為6.29 km/s，轉移結束時θe= 0.25，θm= 6.22，轉移時間為212 d。

表1 地?火無動力循環軌道計算

圖12 地?火無動力循環軌道算例

進入第四步，首先判斷直接地?火轉移的可能性，計算發現可以直接轉移，跳過第四步進入第五步。代入相應的Lambert 轉移函數，得到航天器在離開地球時速度矢量為(?3.13,33.57,0)km/s，到達火星時速度矢量為(12.36,?19.3,0)km/s，地球引力場內剩余速度為6.29 km/s，火星引力場內剩余速度為 7.19 km/s，轉移結束時θe=0.95，θm=3.66，轉移時間為406 d。

第五步結束后，重新開始第二步計算，進行下一次循環，直到計算無法進行時，退出程序。

最終該軌道能夠使航天器在地球與火星之間無動力往返四次，轉移最長時間406 d，最短時間107 d。四次往返中，有兩次火?地轉移需要航天器在火星軌道上待機六個半圈轉移軌道，即三個火星年，期間航天器可以在軌道上閑置，待最后一次火星半圈轉移軌道結束后，再將人員物資轉移到航天器上，進行下一次火?地轉移。整個軌道耗時17.76 a，平均每次循環耗時4.44 a，其中第二、第三次循環各自僅耗時 1.95 a，2.16 a。

根據文獻[7]，循環軌道周期一般以地火會合周期(約2.14 a) 為單位，通常有兩個及以上的會合周期，一個軌道周期中，航天器只完成一次地?火軌道相同，優于三個及以上會合周期的循環軌道，并且本文提出的地?火循環軌道中的部分循環能夠在短時間內進行，進行地球與火星間的快速轉移。

在高精度模型中，本文所提到的地?火轉移軌道具有穩定性。舉例而言，上文中初始參數θ0=1.357π，t0= 0.832 a 的軌道在高精度軌道中漂移至θ0=1.352π，t0=0.822 a，誤差較小。

航天器完成該軌道時，可以通過變軌進入新的循環軌道或停泊進入地球環繞軌道，執行下一個任務。

4 總結

本文在前人基礎上，考慮地球與火星引力輔助效應，引入半圈轉移軌道，在剩余速度與引力輔助近星點半徑的限制下，給出一種新的地?火無動力循環軌道算法。利用初始轉移角度與初始轉移時間這兩個自由參數進行遍歷計算。算例表明，存在能夠在18 年間于地球與火星之間無動力往返4 次的循環軌道，效率較高。這些軌道能夠作為未來地球與火星之間客運、貨運的備選軌道，從而減少航天器燃料消耗，降低運輸成本。