一種月球表面飛躍轉移軌跡設計方法

王浩帆，張洪華，王澤國，關軼峰

北京控制工程研究所，北京 100190

20世紀70年代以來，巡視器一直是行星表面無人探測的主要工具[1]。然而，這種探測方式存在固有的缺點：著陸作業時，需要一個著陸器將巡視器部署到行星表面。這既增加了任務飛行器的總質量，也增加了到達地面后巡視器與著陸器分離的復雜性。從功能的角度來看，傳統的巡視器可以探測的地形類型是有限的，適合在相對平坦的地面上活動，但卻不能爬上陡峭的懸崖，不能越過裂縫，也不能從地面的危險中逃生。如果目標是探索更危險、更多樣的行星、衛星和小行星，目前的巡視器技術可能會限制實現任務的能力。

飛躍器是一種在星體表面上跳躍飛行的交通工具，與現有的巡視器相比具有顯著的優勢[2]。飛躍器可以快速移動至待探測區域，并能克服地形障礙，如巨石和火山口。因此，飛躍器作為一種可以為行星探索提供重要用途的技術，在未來具有重大的發展潛力。2007年，南加州大學的信息科學研究所和航天與空間技術部門開發了月球進入與接近的地面研究平臺(lunar entry and approach platform for research on ground, LEAPFROG)，將其作為月球著陸技術的開發和測試相關挑戰的解決方案[3]。2009年，為了發展一種可以在月球表面著陸，穿越指定距離，并將數據發回地球的飛行器，谷歌舉辦了X-Prize競賽。麻省理工學院的下一個巨大飛躍團隊(next giant leap, NGL)將飛躍器作為主要研究項目，并開發了陸地人造月球與簡化重力模擬器(terrestrial artificial lunar and reduced gravity simulator, TALARIS)試驗臺[4]，用來研究飛躍器的制導、導航與控制算法。文獻[5]中設計了一種新型高原地形飛躍器(highland terrain hopper, HOPTER)，該飛躍器結構簡單，機動性強，成本低廉，可以編隊組網飛行對復雜地形進行探測。文獻[6]中設計了一種基于氣囊的飛躍器，具有魯棒性好、成本低等優勢。

標準軌跡設計是實施制導和控制的重要指標，要求能夠預定著陸場景，同時滿足多種約束條件，從而提高任務的安全性和可靠性。從產生標準軌跡的計算復雜性和實時應用考慮，將標準軌跡生成可以劃分為離線軌跡優化和在線軌跡優化兩種。采用離線生成軌跡作為標準軌跡，是在任務之前進行計算，然后將數據存儲到星載計算機上，對計算成本要求較低，更多是考慮如何在多約束的復雜條件下尋找到最優解。離線軌跡優化的本質是解決復雜約束下的最優控制問題。間接法和直接法是最常用的兩類方法[7]。間接法是利用極大值原理將最優控制問題轉化成兩點或者多點邊值問題，需要引入協態變量，但是協態變量沒有明確的物理意義，算法對協態變量的初值猜測很敏感。直接法是將變量參數化，將連續動力學方程轉化為一系列代數方程，原最優控制問題就變成了非線性規劃問題，然后通過已有的數值算法求得最優軌跡。直接法相對間接法應用更為廣泛，而且有很多不同類型的直接法。偽譜法是一種求解非線性規劃問題很有效的直接法[8]，通過使用非均勻網格，以及計算可行解時需要的較少的節點，從而得到更平滑的解。文獻[9]針對臨近空間滑翔式高超聲速的特點，采用高斯偽譜法對固定翼飛行器和折疊翼飛行器進行了軌跡優化，提出了一種綜合目標與約束的軌跡優化思想。但是偽譜法由于需要對控制變量參數化， 無法保證全局最優性，而且計算成本高，應用受到了一定的限制[10]。近年來，全局智能優化算法的發展為解決軌跡優化問題提供了一種新的直接方法。文獻[11]針對升力式火星探測器進入火星大氣段，采用遺傳算法對進入軌道進行優化設計，提出不同的推力發動機制動方案。文獻[12]提出一種基于隨機森林的模型對月球著陸過程軌跡重規劃技術，利用離線訓練好的模型根據航天器狀態在線計算其控制量，降低了航天器位置速度誤差。但是，由于智能方法目前仍缺乏完善的理論證明，可靠性較低，限制了其實際應用。

另一種近年來發展迅速的直接法是凸優化方法，對于一個優化問題而言，得到全局最優解的條件取決于系統的凸性而不是線性，極大擴展了使用優化方法求解的問題[13]。在此條件下，問題可以實時解決，而且如果問題是可行的，計算出的解就是全局最優解。文獻[14]中將火星動力下降問題轉化為求解二階錐問題(second order cone problem, SOCP)，通過內點法使問題在多項式時間內得到求解，從而使在線軌跡優化成為可能。隨后采用凸優化方法求解飛行器最優軌跡的研究逐漸增多。文獻[15]中提出了一種連續凸化技術，解決了動力下降段中存在氣動阻力和非凸約束的問題。文獻[16]中考慮了導航及障礙檢測敏感器視場約束及制動發動機推力大小，研究了凸優化方法處理含有復雜約束的著陸軌跡優化問題。文獻[17]針對傳統內點法求解效率不高的問題，提出了一種機載使用的實時算法。文獻[18]中，美國國家航空航天局(national aeronautics and space administration, NASA)基于凸優化理論開發了大轉向燃料最優制導(guidance for fuel-optimal large diverts, G-FOLD)軟件，并在Xombie火箭上驗證了其有效性。文獻[19]中考慮了小天體的不規則地形和弱引力場，提出了一種基于凸優化方法的小天體飛躍的制導算法，將飛躍器在初始段和著陸段受到的約束均考慮為二階錐約束，具有較高的魯棒性和精度。這些方法可以適用于飛躍器單個階段的最優燃料軌跡規劃，但不能解決全任務階段的軌跡規劃問題。這是由于這些研究都是在整個運動過程中約束相同的條件下進行的，而在實際工程中，飛躍器在起飛時受發動機指向限制且不能與周圍環境相撞，因此需要垂直上升。上升到一定高度后發動機轉向開始飛躍。到達目標點時，為了保證飛躍器安全穩定地著陸，加入了垂直下降約束。這三個階段受到不同的高度和推力約束，因此需要對凸優化方法進行擴展。

本文針對飛躍器轉移過程中存在垂直起飛、飛躍、垂直下降階段約束不同的問題，從凸優化求解算法出發，通過假設上升和著陸時間相同，簡化了分析模型，然后利用黃金分割法搜索最優燃耗對應的上升、著陸時間以及飛躍全程的時間，在整個凸優化框架內對不同的階段采用不同的推力和位置約束，采用內點法得到了最優燃耗解，并通過兩個仿真算例校驗了算法的合理性和有效性。

1 飛躍過程建模

飛躍過程是指飛躍器從起始點以靜止狀態垂直起飛，利用變推力發動機到達目標點上方，然后垂直著陸至目標點。

1.1 坐標系建立

忽略月球表面曲率影響，建立如圖1所示的月球表面飛躍坐標系O-XYZ，原點O為目標著陸點，OX軸指向從著陸點到月球質心的反方向，OY軸指向初始飛躍點，OZ軸由右手法則確定。

圖1 飛躍坐標系Fig.1 Hop trajectory frame

1.2 動力學方程建模

行星飛躍主要有兩種方式——巡航式和彈道式，飛躍方式如圖2所示。由于受到姿態要求和推力指向要求，彈道式飛躍在飛躍器上升到一定高度后，增加推力進行彈道式飛躍，達到一定的高度和速度后減小發動機推力或關閉發動機。當飛躍器即將接近目標時進行減速制動，飛躍器到達目標上空一定高度時下降，實現定點軟著陸。通常還需要滿足半矩形約束，即除了上升段和著陸段外，飛行軌跡不能與虛線相交。這種方式只需要發動機在發射段和著陸段提供推力，可以節省燃料。

巡航式飛躍是指飛躍器先上升到一定高度，然后水平向目標移動，運動過程中，垂直方向上必須始終由發動機施加與月球重力加速度相等的控制加速度，水平方向上采用先加速，然后滑行，最后減速的方式接近目標點。由于發動機在整個飛行過程中都是開機狀態，因此預計飛行燃耗會更高[20]。

圖2 兩種類型的飛躍軌跡Fig.2 Two types of hop trajectories

不考慮月球自轉的影響，并且將月球加速度視為常量處理，飛躍器的動力學方程為：

(1)

(2)

式中：Isp為發動機比沖；ge為地球表面加速度大小。

1.3 優化問題建模

對于深空探測任務，燃耗是重要的性能指標，節約燃料越多，越有利于有效載荷的增加和成本的節約。同時，燃耗少也能節省整個飛行器的質量，更加有效地進行探測任務，因此把燃料消耗最少即燃料最優作為著陸軌跡優化目標的研究可以為著陸之后的探測任務提供更好的基礎。

以燃料最優為目標的模型可以描述如下：

(1)性能指標

(3)

式中：tf表示到達目標點的時間；J為燃料最優的性能指標。

(2)約束條件

對于彈道式飛躍器來說，首先由月面飛行到一定高度，然后進行彈道式飛躍，飛躍至目標點上方一定距離后，緩慢到達預定著陸點，狀態約束條件如下：

(4)

式中：r0，r1，r2，rf，v0和vf∈R3為常值。

此外，因為初始上升段和末端著陸段要求是垂直上升和下降的，因此這兩個階段的推力方向約束為:

(5)

式中：e1=[1,0,0]T表示OX方向的單位列向量。

考慮一定的工程實際，由于推進系統的能力有限，飛躍器的控制推力幅值約束為:

0≤Tmin≤‖T‖≤Tmax

(6)

式中：Tmin和Tmax分別為推進系統的最小推力幅值和最大推力幅值。

飛躍過程中，為了保證飛躍器不與障礙碰撞，需要加入高度約束，即:

rx(t)≥hd,t∈[t1,t2]

(7)

式中：rx為r在OX方向的分量；hd為根據任務軌跡附近月球表面的起伏程度所確定的高度。

2 飛躍轉移過程的凸優化方法

2.1 控制約束的凸化處理

由以上分析，燃料最優飛躍模型是求在滿足式(4)到式(7)的約束條件和式(1)的運動方程下使得式(3)的性能指標最小的最優控制變量和最優狀態變量。由于式(6)描述的約束條件中含有幅值下界從而成為一個非凸約束問題，直接求解可能會陷入局部最優，所以需要將非凸約束問題進行一定的轉化后再求解。由文獻[14]可知，可以通過引入松弛變量Γ而將式(6)轉化為如下的凸約束：

推力方向約束為

(8)

約束(8)導致無法利用傳統的凸優化算法進行計算，因為在飛躍的過程中，不同階段約束不同，對于這種情況，可以采取兩種處理方式——分段凸優化和全程凸優化。分段凸優化指將飛躍過程分為上升段、飛躍段、下降段，每個階段分別進行計算，但是這種方式需要給定上升段和飛躍段末端速度約束，否則無法得到全局最優解。而如果采用全程凸優化，一個關鍵的難點是無法確定推力約束的切換時間點。因此，本文考慮采用線搜索的方法得到推力切換點。

2.2 推力切換點的計算

對于一般工程情況，上升段和下降段的時間遠小于飛躍時間，不妨先假設上升段和下降段的時間相同，通過下述方法來確定飛行時間的上界和下界。

考慮如果用最短的時間上升，即采用最大推力上升，所需時間為:

(9)

式中：av為飛躍器的垂向加速度大??；gl的定義見式(1)；hd的定義見式(7)。式(9)在計算過程中未考慮質量變化，由于采用最大推力時，飛躍器質量變化較大，因此實際上升時間更少。

對于時間上界的計算，如果按照最小推力計算，即：

(10)

(11)

2.3 變量等效與離散化

為了方便采用數值算法來求解，考慮如下的變量等效：

(12)

式中：σ為松弛變量的等效變量；u為等效推力矢量；z為等效質量。

將式(12)代入式(1)，則動力學方程式(1)的第二、三個等式變為：

(13)

(14)

從式(14)可以推出：

(15)

式中：m0是飛躍器的初始質量。

因為α>0，則燃耗最小性能指標函數等價于：

(16)

用等效變量表示的控制約束如下：

‖u(t)‖≤σ(t),?t∈[0,tf]

(17)

(18)

將式(12)代入式(14)和(18)得：

(19)

Tmine-z(t)≤σ(t)≤Tmaxe-z(t),?t∈[0,tf]

(20)

不等式(20)的第一部分定義了可行解的凸區域，但是第二部分不是。因此考慮采用二階錐和線性近似來處理這些不等式約束。將第一部分用二階錐來近似(通過用e-z(t)泰勒展開的前三項來近似)：

Tmine-z0(t){1-[z(t)-z0(t)]+

式中：z0(t)=ln(m0-αTmaxt)。

式(20)的第二部分用線性近似(通過用e-z(t)泰勒展開的前兩項來近似)：

σ(t)≤Tmaxe-z0(t){1-[z(t)-z0(t)]}

(21)

從而得到式(20)的二階錐和線性近似:

將整個著陸過程的時間[0,tf]均勻分為N-1段，每段間隔為Δt，則對應的時間節點為

tk=kΔt,k=0,…N

式中：k表示離散時間序列的第k個時間節點，終端時間為第N個時間節點，即tf=NΔt。

然后將各離散時間節點的控制量u和σ通過M個基函數來表示，即:

式中，φj表示第j個離散點的控制量基函數，pj為待優化的常系數列向量。

本文選取分段常值基函數，且取M=N，則最優控制問題轉化為參數優化問題:

則動力學方程(1)可以表示為：

(22)

式中：

式中：

式中：Ac和Bc為連續系統狀態轉移矩陣和控制輸入矩陣;A和B為對應的離散系統狀態轉移矩陣和控制輸入矩陣。

定義如下常值矩陣:

(23)

定義垂直上升段與飛躍段的切換點為k1，飛躍段與垂直下降段的切換點為k2，滿足:

1≤k1

則飛躍軌跡燃耗最優問題最終轉化為如下問題。

(1)性能指標

(24)

(2)約束條件

(25)

(26)

(27)

(28)

(29)

對于任意給定的N，式(24)～(29)定義了一個有限維二階錐問題，通過使用現有的SOCP算法，可以非常有效地求解，并保證收斂到全局最優解。需要注意到N也是一個解參數，對應運動時間上界tf，且滿足0

(30)

式中：mwet為飛躍器總質量;mfuel為燃料質量;α定義見式(2)。

式中：δ為需要調節的參數。

則參數N滿足

Nmin≤N≤Nmax

式中：

式中：int為向下取整函數。

2.4 飛躍全程的凸優化算法

為了便于全程凸優化，還需對垂直起飛和著陸段進行約束。利用式(26)對式(25)進行近似處理，首先將2.2節計算出的tas換算為離散時間點k1，在此時間節點前，滿足：

(31)

在k1點以后的飛躍段，施加如下約束：

(32)

垂直下降階段約束與垂直上升段約束相同：

(33)

最終即可把飛躍問題轉化為求解SOCP問題，具體的算法步驟見式(34)至式(39)，算法流程圖見圖3。本文算法通過固定上升、下降時間，將分段凸優化問題轉化為全程凸優化問題進行求解，因此得到的結果優于分段優化的結果。

(34)

令i=0，轉步驟2。

(35)

(36)

(37)

(38)

(39)

圖3 算法流程Fig.3 Algorithm flowchart

3 仿真校驗

對于已建立的SOCP問題，采用MATLAB求解二階錐規劃問題的工具包YALMIP[22]并借助SDPT3[23]解算器來求解，仿真條件見表1?？紤]平地飛躍與向坑內飛躍兩種情況進行仿真。

表1 仿真條件

由于本文坐標系OY軸指向起始點，且初始和終端速度都為0，因此任何3維飛躍問題都可以轉化到2維平面。若不采用本文的坐標系定義方法，對于3維空間，只需對OZ軸方向進行與OY軸同樣的約束，即將式(23)中eu向量改為[0,1,1,0]T，并適當增加飛行時間上界即可。

3.1 平地飛躍

考慮如圖4所示工況，其中虛線為半矩形約束，飛行軌跡不能落在該區域內。為進行燃耗的比較，考慮了分段凸優化和全程凸優化兩種情況。

圖4 月球飛躍軌跡Fig.4 Hop trajectory on the Moon

(1)分段凸優化

將運動軌跡分為三段，分別進行凸優化求解燃料最優軌跡。起始點位置為r0=[0,15000,0]T，速度為v0=0，上升段為從月面垂直起飛到30 m高處懸停。飛躍段的初始位置為r1=[30,15 000,0]T，初始速度為v1=0，飛躍過程保持高度大于30 m，橫向移動15 km，即末端位置為r2=[30,0,0]T，末端速度為v2=0。著陸段為從高度30 m垂直著陸到月面，即rf=0，vf=0。仿真結果如圖5所示，其中虛線為30 m高度半矩形約束。

從圖5可以看出，飛行軌跡滿足半矩形約束；飛躍器先加速到一定速度后，發動機關機進行滑行，直到接近目標點后重新開機制動；由推力曲線可以看出，推力幅值在每個階段都呈最大-最小-最大的趨勢。

圖5 平地飛躍分段凸優化的第一種情況Fig.5 The first case of piecewise convex optimization of flat hop

續圖5 Fig.5 Continued

第二種情況考慮了放寬飛躍段的末端速度約束，假設飛躍段末端飛躍器的質量仍為200 kg，如果采用最大推力著陸，容許速度為：

(40)

式中：vh為高度方向的容許速度；Tmax為最大推力；假設質量m恒定，因此加速度a恒定。

由式(40)，并考慮一定的工程余量，設置飛躍段末端垂向下降速率小于10 m/s，并將上升段末端速度v1設為自由，仿真結果如圖6所示。

圖6 平地飛躍分段凸優化的第二種情況Fig.6 The second case of piecewise convex optimization of flat hop

(2)全程凸優化

采用本文所提出的算法進行仿真，結果如圖7所示。表2列出了三種仿真條件下的燃料消耗，可以看出，放寬垂直上升段和飛躍段的末端約束可以減少燃料消耗，但總的來說，分段凸優化的燃料消耗始終高于全程凸優化的方法，這主要是因為全程凸優化是綜合考慮了所有約束，而分段凸優化只能保證在每個階段的燃耗是最優的。

圖7 平地飛躍全程凸優化仿真結果Fig.7 Full trajectory convex optimization simulation of flat hop

表2 平地飛躍不同條件下的燃料消耗

表3列出了三種仿真條件下的飛行時間。

表3 平地飛躍不同條件下的飛行時間

可以看出兩種分段凸優化的總飛行時間相差不大，而全程凸優化方法在飛躍段的飛行時間多于分段凸優化方法，這是因為更低的燃耗意味著用了更多時間的最小推力，因此增加了飛行時間。此外，注意到搜索得到的最優時間并不在給定搜索范圍的邊界處，因此在計算垂直上升、下降段時間的搜索邊界時，不考慮質量變化是合理的。

3.2 坑內飛躍

為進一步校驗算法的有效性，考慮如圖8所示工況，虛線為半多邊形約束，除垂直上升段和著陸段外，飛行軌跡不能與該半多邊形相交。

圖8 向坑內飛躍軌跡Fig.8 Hopping into the pit

首先直接用之前的算法進行仿真，得到如圖9所示結果?？梢宰⒁獾皆诔跏级芜`反了約束，因此要對約束條件進行修改。而由前面的仿真可知，飛躍器橫向運行1 km的時間為30 s，因此在全程凸優化框架內，對離散化后的前60個節點，增加約束使其高度大于4 030 m，仿真結果如圖10所示。可以看出，運動軌跡滿足約束，最優燃耗為25.069 2 kg，總飛躍時間171.244 6 s。由此工況可知，對本文算法的拓展應用只需根據實際情況修改式(31)～(33)參數即可。

圖9 坑內飛躍運動軌跡Fig.9 Trajectory of hopping into the pit

圖10 坑內飛躍仿真結果Fig.10 Simulation of hopping into the pit

續圖10 Fig.10 Continued

4 結束語

本文研究了飛躍器轉移過程中的燃料最優軌跡設計問題，通過對凸優化方法進行擴展，求解燃料最優問題從而得到最優軌跡。

1)算法通過假設垂直上升和下降段時間相同，利用黃金分割法搜索得到了上升、飛躍和下降段的時間，能將原問題中的非凸約束轉化為凸約束。

2)算法能在滿足約束的情況下實現精確飛躍轉移。對于平地飛躍情況，相同仿真參數下其燃耗為24.968 2 kg，優于兩種不同速度約束條件下分段凸優化的結果。

3)坑內飛躍的仿真結果表明，對于不同的工況，只需根據實際情況修改算法中的期望高度值，利用本文思想重新分段，便可得到燃耗最優軌跡。

由于考慮到計算效率的問題，本文算法近似了實際的上升、著陸段的飛行時間，因此，在不提升計算復雜性的情況下，精確考慮不同階段飛行時間的燃料最優解將是下一步的研究方向。