經濟型月球探測器精確定點軟著陸制導算法

高 峰，荊武興，高長生，李志剛，鐘 偉

(1. 哈爾濱工業大學航天工程系，哈爾濱150001；2. 上海宇航系統工程研究所，上海201109)

0 引 言

作為21世紀航天領域的研究熱點，深空探測吸引著越來越多的國家和科學技術人員進行深入的研究。而月球是人類開展外太空探測活動的第一顆地外天體，因此月球探測可以說是深空探測的開端。世界各航天大國對于月球探測的基本路線大體是“探”、“登”、“駐”[1]，所以建立月球基地、開發月球資源是實施月球探測計劃的最終目標。針對這一目標，未來月球探測器的研發需要考慮兩個方面的因素：一方面，需要運送大量的小型載荷到月球表面來建立月球通信網絡，以提高未來月球探測器的自主能力[2]；另一方面，對月球資源的開采和運輸需要探測器多次往返于地月之間。考慮到以上因素，若仍然采用傳統月球探測器，探月成本必然十分高昂，因此經濟型月球探測器將成為未來月球探測任務的不二之選，目前已知的相關項目包括NASA的“通用月球探測器”計劃(Common lunar lander,CLL)[3]和“月球網絡”計劃(International lunar network,ILN)等[4]。與傳統月球探測器不同，經濟型月球探測器采用固體火箭發動機作為制動減速的主發動機，相較于液體火箭發動機，固體火箭發動機結構簡單、加工和材料成本低廉，造價僅約為液體發動機的1/10～1/20，另外在主減速段結束之后，可以將其拋離，這樣也很大程度上減少了隨后過程的燃料消耗。但是采用固體發動機給制導律的設計帶來較大的困難，液體發動機可以通過開關機和變推力技術實現對探測器運動狀態的精確控制，而固體發動機不能做到隨時開關機，只能燃料耗盡關機，同時也很難實現推力大小的精確改變。綜上，需要對主發動機采用固體發動機的經濟型月球探測器的制導律進行深入研究和設計。

軟著陸任務對于制導算法有以下要求：燃耗最優或次優性、自主性、實時性、魯棒性以及避障性能[5]。目前，月球軟著陸制導方法主要有：月球垂線法、重力轉彎法、標稱軌跡法和顯式制導法，而目前的研究重點主要集中在標稱軌跡法和顯式制導法。標稱軌跡法包括標稱軌跡優化和軌跡跟蹤兩部分，孫軍偉等[6]通過將常推力月球軟著陸軌道離散化，將原軌跡優化問題轉化為非線性規劃問題，并基于SQP(Squential quadratic programmin)方法對該問題進行了求解。彭祺擘等[7]將控制變量和終端時間作為優化變量，同時離散控制變量與狀態變量，采用Gauss偽譜法對月球定點著陸最優軌跡進行快速優化設計。林曉輝等[8]在考慮敏感器視場約束及發動機推力大小約束的情況下，基于凸優化理論設計了月球軟著陸的最優標稱軌跡。對于軌跡跟蹤算法的研究，Liu等[9]基于H∞理論設計了可以魯棒追蹤標稱軌跡的最優反饋控制器。Wang等[10]提出了一種基于模糊神經網絡的非線性最優控制策略。相較于標稱軌跡法，顯式制導法表現出良好的自主性和實時性，同時具備抗干擾能力，能夠保證末端制導精度，具有一定的魯棒性，但由于在求解過程中使用了必要的假設，所以顯式制導算法的燃耗具有次優性。顯式制導法是根據探測器的實時運動狀態和終端狀態，按控制泛函的顯函數表達式進行實時解算的制導算法，該算法最早應用在日本SELENE-1探測器的設計方案中，隨后眾多學者對其進行了充分的研究。王大軼等[11]提出了一種燃料次優閉環制導方法，該方法以燃耗為最優性能指標，將制導指令表示為剩余時間的函數。王鵬基等[12]在顯式制導方法中引入前饋控制，使得制導律能夠消除初始偏差對末端狀態的影響。孫軍偉等[13]、Guo等[14]基于Pontryagin極大值原理設計了閉環多項式制導律。前期的顯式制導律均沒有考慮軟著陸終端水平位置的約束，李茂登[15]通過忽略動力學方程中的高階小量，采用變推力發動機，得到了能夠滿足末端全部位置約束的閉環顯式制導律。另外在再入飛行器領域，Lu[16]提出了一種單基線預測校正再入制導算法，月球軟著陸過程可借鑒該算法進行顯式制導律設計。盡管目前的顯式制導律具有自主性、實時性和易于工程實現的優點，但對于采用固體發動機作為主發動機進行減速制動的精確定點軟著陸任務，由于傳統顯式制導算法在假設條件及制導指令解算方法上的局限性，并不能夠滿足全部的終端狀態約束，因此有必要對其進行深入的研究。

基于上述分析，本文將對主發動機采用固體發動機的經濟型月球探測器的在線閉環制導算法進行研究。首先在北東地坐標系下建立軟著陸主減速段的質心動力學模型，然后依據動力學模型將探測器的運動分為縱向平面運動和橫向運動，分別設計在線閉環制導律，最后通過數值仿真驗證所設計算法的性能優劣。

1 月面著陸點全覆蓋的環月調相策略

在月球探測任務中，探測器經雙曲制動進入100×100 km環月軌道，然后經環月調相后最終機動到15×100 km環月軌道。在進行月面軟著陸制動減速之前，探測器進行環月調相的主要目的是使探測器到達近月點時，預定的著陸點隨月球自轉剛好進入環月軌道面內，若此時開啟主發動機進行主減速段制動，由于目標點已處于軌道面內，那么在探測器主減速段運動接近目標點的過程中，橫向偏差(軌道面法向)將會很小，這樣縱向平面(軌道面內)的運動通過固體主發動機產生的推力來控制，而橫向運動由于偏差很小，只需采用橫向小推力姿控發動機控制即可，從而為后續設計滿足全部終端狀態約束的在線閉環制導算法提供了基礎。另外，為了滿足未來建立月球通信網絡和月球資源的開采、運輸任務，需要探測器的著陸點能夠覆蓋月面包括極區在內的所有可著陸區域，因此選擇月球極軌道作為環月軌道是一個很好的選擇。通過環月調相等待，處于極軌道的探測器理論上可以軟著陸到月面任意一點。基于上述考慮，下面進行精確環月調相策略的設計。

如圖1所示，雙曲制動結束后探測器進入環月軌道A點，設此時探測器位于A點的角距為ω1, 15×100 km環月軌道近月點C點的角距為ω2(由預定著陸點和軟著陸過程的總航程決定)，根據幾何約束關系：

當ω2≥ω1時：

(1)

當ω2<ω1時：

(2)

則探測器由A點運動到15×100 km環月軌道遠月點B點所用的時間為：

(3)

同理，探測器在15×100 km環月軌道由B點運動到C點所用的時間為：

(4)

式中：a2為15×100 km環月軌道半長軸。

設預定著陸點D點的月心經度為λ1，探測器位于A點時環月軌道升交點(或降交點)經度為λ2，定義：

Δλ=λ2-λ1

(5)

則預定著陸點隨月球自轉運動到軌道面內的總時間：

(6)

式中：ωm為月球自轉角速度。

故可用于調相的時間為：

Δt=t-t1-t2

(7)

因此環月調相策略為：

Δt=ΣTi(i=1,2,3,…,)

(8)

式中：Ti(i=1,2,3,…,)為環月調相軌道周期。

圖1 環月調相策略示意圖Fig.1 Strategy of phasing

按照上述環月調相策略，只需設計合適的環月調相軌道，即可使探測器在到達15 km近月點時，剛好預定著陸點隨月球自轉進入軌道面內，此時開啟主發動機進行減速制動，同時啟動橫向小推力姿控發動機進行橫向控制。

2 經濟型月球探測器精確定點軟著陸問題建模

通過第1節的分析，探測器位于15×100 km環月極軌道近月點時開始主減速段制動。為了方便對月球探測器軟著陸制導方法進行研究，需要對動力學模型進行合理的簡化。考慮到月球的自轉角速度很小，而主減速段過程又很短，所以進行如下假設：1)忽略月球自轉；2)忽略月球非球形引力攝動；3)忽略其他天體的引力攝動。根據牛頓第二定律，在圖2所示的北東地坐標系下，探測器的質心動力學模型為：

(9)

式中：r=[0, 0, -r]T、v=[vx,vy,vz]T分別為探測器的位置和速度在北東地坐標系下的分量;φ為探測器所處的月心緯度;λ為探測器所處的月心經度;F=[Fx,Fy,Fz]T為發動機推力在北東地坐標系下的分量。

圖2 軟著陸主減速段坐標系Fig.2 Frames of main descent phase

考慮探測器質量的變化，其質量方程為：

(10)

為了方便后續制導算法的設計，這里定義等效引力加速度：

(11)

由于r相對于vx、vy、vz要大得多，所以探測器在主減速段運動過程中，g′x、g′y、g′z近似為常值。

3 精確定點軟著陸實時在線制導算法

對于配備固體發動機的經濟型月球著陸器來說，在燃耗最優著陸的情況下，只控制推力方向而不改變推力大小，傳統顯式制導算法不能同時滿足三個終端位置約束[11]。因此，根據式(9)以及第1節的分析，我們將探測器主減速段的運動分為縱向平面運動(x-z平面)和橫向運動(y方向)兩部分，然后分別設計在線閉環制導律。

如圖3所示，在北東地坐標系下，縱向平面的運動通過固體主發動機產生的推力來控制，控制量為俯仰角αbn，同時在探測器運動過程中，保持偏航角βbn=90°或270°。對于橫向運動控制，由第1節的分析可知，通過環月調相后，探測器在主減速段運動過程中橫向偏差很小，若采用大推力固體發動機進行控制，很容易造成控制過量和偏航角在零值附近“振蕩”的情況，所以橫向采用小推力姿控發動機進行控制，控制量為發動機開關機指令u，u的值域為{-1, 0, +1}。在此控制策略下，推力在北東地坐標系下的分量為：

(12)

式中：F為固體發動機推力大小，Tr為姿控發動機推力大小。

圖3 固體發動機推力在北東地坐標系下的分量Fig.3 Thrust of solid rocket motor in the NED frame

3.1 縱向平面制導算法

3.1.1改進顯式制導律

根據Pontryagin極大值原理，在固體發動機常值推力下，燃耗最優問題就等同于時間最優問題。令最優指標為：

(13)

引入Hamilton函數：

(14)

式中：λ=[λx,λvx,λz,λvz]T為縱向平面運動協態變量。

伴隨方程為：

(15)

在不考慮位置約束的情況下：

λx=λz=0

(16)

λvx=λvx0,λvz=λvz0

(17)

由極值條件：

(18)

對式(9)進行對時間積分，可以得到：

(19)

由上式可以得到不考慮位置約束下的俯仰角指令：

(20)

進一步，在考慮位置約束時，參考文獻[11]的線性近似方法，設計縱向平面的顯式制導律如下：

(21)

文獻[11]將(l1+l2t)視為小量并對控制角的正弦和余弦進行一階Taylor展開，然而在實際的仿真分析中發現(圖4、圖5)，在可達域內，除特定著陸點外，其余著陸點對應于位置約束的控制角常數項部分并不是小量，而只有一次項滿足小量假設，同時上述處理方式會在主減速段末段局部時間內造成姿態角突變，從而增大姿控系統的工作難度。通過上述分析，下面對傳統顯式制導律作式(22)的改進，僅把一次項l2t作為小量，而常數項l1不作為小量，然后進行Taylor展開，可以看到改進的顯式制導律滿足上述假設，同時解決了末段姿態角突變的問題。

(22)

圖4 l1變化曲線Fig.4 l1 curve

圖5 l2變化曲線Fig.5 l2 curve

根據探測器徑向zn的初始運動狀態和終端狀態約束，將式(22)代入式(9)并進行對時間一次積分和二次積分，可求得l1和l2的顯式表達式。

3.1.2最優點火角搜索

在改進的顯式制導算法中，控制角αbn的顯式解與當前時刻探測器的徑向位置坐標z0、終端徑向位置坐標zf、縱向平面內探測器當前時刻的速度v0和終端速度vf以及剩余時間tgo有關，而與航向位置x0、xf無關，因此所設計制導律缺乏對航向位置的控制能力，為了實現精確定點軟著陸，需要進一步設計探測器航向運動控制方法。

如圖6所示，當探測器在15×100 km環月極軌道近月點開啟點火制動時，若采用改進顯式制導律飛行，除特定著陸點外，終端實際著陸點將超前(或滯后)理想預定著陸點。為解決這一問題，參考有限推力制導思想[17]，對點火點位置進行修正，以達到控制主減速段航程的目的。

圖6 主減速段點火點修正示意圖Fig.6 Fire point correction of main descent phase

設點火點修正角為Δθ，通過仿真分析發現，在預定終端狀態一定時，隨著點火點位置的改變，主減速段采用改進顯式制導律飛行的探測器航程變化不大。

圖7 主減速段航程角隨點火點位置變化Fig.7 Flight angle with fire point in main descent phase

另外主減速段末端狀態誤差變化如圖8～圖11所示，可以看到存在最優點火角Δθopt，使得探測器終端位置誤差和速度誤差達到最小。

圖8 航向位置誤差隨點火點位置變化Fig.8 Course position error with fire point

圖9 徑向位置誤差隨點火點位置變化Fig.9 Radial position error with fire point

圖10 航向速度誤差隨點火點位置變化Fig.10 Course velocity error with fire point

圖11 徑向速度誤差隨點火點位置變化Fig.11 Radial velocity error with fire point

通過上述分析，根據主減速段改進顯式制導律的誤差傳播特性，設計如圖12所示的點火點修正角搜索算法(Fire angle search algorithm, FASA)，δ為定點軟著陸主減速段末端容許位置誤差。在實際任務中，為了與后續階段光滑過渡，一般要求在主減速段任務結束時，預定著陸點應位于探測器運動方向的前方，因此δ應大于0。另外，點火角搜索步長與環月軌道定軌精度有關，對于定軌精度位置偏差在300 m以內的月球探測器[18]，可設置點火角搜索步長為±0.01°(終端實際著陸點滯后理想預定著陸點時取+，終端實際著陸點超前理想預定著陸點時取-)。

圖12 點火角搜索算法Fig.12 Fire Angle Search Algorithm

3.2 基于ZEM/ZEV方法的橫向運動控制

軟著陸主減速段橫向運動控制采用ZEM/ZEV最優反饋制導律[19]，該方法具備燃耗最優、強魯棒性和易于工程實現的特點。定義橫向零控位置偏差yZEM和零控速度偏差vyZEV如下：

yZEM(t)=yf-y(tf)

(23)

vyZEV(t)=vyf-vy(tf)

(24)

式中：yf和vyf分別為預定著陸點的橫向位置及橫向速度，y(tf)和vy(tf)分別為在不施加控制情況下探測器在tf時刻的橫向位置和橫向速度。

對于式(9)和(12)，令ay=Tru/m，設計能量最優性能指標函數：

(25)

根據性能指標函數(25)得到Hamilton函數：

(26)

式中：λ=[λy,λvy]T為橫向運動協態變量。

由極值條件以及伴隨方程，可得最優控制量ay的表達式為：

ay=-λvy=-λy(tf)tgo-λvy(tf)

(27)

將式(27)代入式(9)并對時間積分，可以得到如下關系式：

(28)

(29)

式中：y0和vy0分別為探測器初始橫向位置和橫向速度。

由式(28)和(29)可得協態變量λy(tf)和λvy(tf)的解析表達式，將其代入式(27)可得ZEM/ZEV制導律為：

(30)

式(30)得到的制導律為連續控制變量，而橫向運動控制采用小推力液體發動機來實現，控制量為發動機開關機指令u(u={-1, 0, +1})，因此需要將其轉換為脈沖形式，下面基于PWPF對連續制導律進行調制[20]。

PWPF調制器工作原理如圖13所示，其由一階慣性環節、施密特觸發器及反饋回路組成，km和tm分別為一階慣性環節的增益和時間常數，Uon和Uoff分別為施密特觸發器的開、關閥值，ay為橫向運動連續制導指令，uy為調制脈沖輸出。發動機最小工作時間為：

(31)

式中：h=Uon-Uoff為調制器遲滯，Um為脈沖幅值。

圖13 PWPF調制器工作原理Fig.13 Structure of PWPF modulator

4 仿真算例

對于主發動機配備固體發動機的經濟型月球探測器，由于固體發動機啟動后必須持續工作到燃料耗盡，所以實際飛行中不能以運動狀態作為發動機關機條件，而必須以質量變化量或推進時間作為終止條件。而在實際任務中，需要提前設計固體發動機的燃料裝填量，因此仿真算例先以運動狀態作為積分終止條件，然后根據主減速段全程飛行時間給出固發燃料裝填質量參考值，最后以該參考值作為固體發動機終止條件來進一步驗證所設計制導算法的性能優劣。

主減速段仿真參數見表1：

表1 主減速段仿真參數Table 1 Simulation parameters of main descent phase

4.1 FASA算法優化性能

首先根據FASA算法優化搜索點火角修正值。為了驗證所設計點火角搜索算法的性能，在相同仿真條件下，選取終端航向位置誤差為代價函數，分別采用粒子群優化算法(Particle swarm optimization,PSO)[21]、樽海鞘群優化算法(Salp swarm algorithm,SSA)[22]和本文所設計的點火角搜索算法(FASA)進行50次獨立優化搜索，表2給出了三種算法代價函數值及平均搜索時間統計數據對比。

表2 點火角搜索結果統計Table 2 Fire angle search results statistics

由統計結果可以看出，三種算法的搜索結果都能夠滿足終端航向位置約束(容許位置誤差200 m)的要求，雖然FASA算法的誤差均值相比于PSO、SSA較大，但由于該算法的修正點火角初值是根據主減速段的誤差傳播特性得到，大大降低了初值猜測過程的隨機性和迭代搜索次數，同時搜索過程采用固定步長，因此相較于其它兩種算法多次搜索結果保持穩定。另外，FASA算法的搜索時間遠小于其他兩種算法，保證了在環月調相軌道周期內能夠快速完成多次優化搜索，即使在探測器環月飛行過程中臨時調整終端著陸點位置，也能夠實時搜索到滿足任務需求的點火角修正值，滿足制導算法實時性的要求。

4.2 主減速段全程仿真

下面采用所設計的制導算法，以運動狀態作為積分終止條件，進行主減速段全過程仿真，仿真中制導參數設置如下。縱向平面運動制導參數：制導周期T=0.08 s；點火角搜索結果：Δθcor=-0.08°。橫向運動制導參數：km=6，tm=0.48，Uon=3.6×10-5，Uoff=3×10-6。仿真結果見圖14～圖19。

圖14 徑向位移曲線Fig.14 Radial position curve

圖15 航向速度曲線Fig.15 Course velocity curve

圖16 橫向速度曲線Fig.16 Lateral velocity curve

圖17 徑向速度曲線Fig.17 Radial velocity curve

圖18 縱向平面運動控制俯仰角指令Fig.18 Thrust pitch angle in longitudinal plane motion

圖19 橫向運動控制發動機開關曲線Fig.19 Control instruction to lateral motion

為了驗證所設計制導算法的性能，采用相同的仿真條件，在縱向平面運動采用文獻[11]中的制導方法，橫向運動不控的情況下進行仿真，與本文所設計的制導律仿真結果進行對比，終端狀態誤差對比結果見表3：

表3 終端狀態誤差Table 3 Terminal state error

由仿真結果可以看出，軟著陸過程中采用本文設計的制導律飛行時，位置變化曲線和速度變化曲線平滑，三個方向上的終端著陸誤差能夠控制在很小的范圍內。同時主減速段結束時，探測器剩余速度大小為11.8 m/s，預定著陸點位于探測器運動方向的前方171.2 m處，能夠實現與后續階段的光滑銜接。另外，主減速段耗時259.04 s，固體發動機燃耗683.60 kg，將此數值作為發射前固發燃料裝填質量的參考值。

4.3 制導算法性能

為了進一步驗證所設計制導算法的魯棒性和適應性，根據上面的仿真結果，選取固發燃料裝填質量為683.60 kg，仿真過程終止條件設置為燃料耗盡關機，其它仿真條件不變，在目標著陸點緯度36 °N～37.6 °N的覆蓋范圍內進行蒙特卡羅仿真試驗，仿真結果如圖20、圖21所示。

圖20 終端位置誤差隨目標著陸點緯度變化曲線Fig.20 Terminal position error with latitude of target landing site

圖21 終端速度誤差隨著陸點緯度變化曲線Fig.21 Terminal velocity error with latitude of target landing site

從仿真結果可以看出，x方向的終端位置誤差散布在0～350 m的范圍內，且大部分結果滿足容許誤差200 m的要求；而終端速度誤差在目標著陸點緯度逐漸遠離原始預定著陸點緯度(36.8173 °N)的過程中，呈現出增大的趨勢。z方向的終端位置誤差受目標著陸點緯度變化的影響較小，能夠保持在很小的范圍內；而終端速度誤差的變化趨勢與x方向相同。對于y方向，由于橫向運動單獨采用小推力姿控發動機進行控制，因此在橫向偏差一定的情況下，y方向的終端位置誤差和終端速度誤差始終能夠控制在允許誤差范圍內。綜上所述，在一定的目標著陸區域內，終端位置誤差和終端速度誤差能夠控制在合理的范圍，所設計的制導算法表現出了很好的適應能力。

5 結 論

主發動機采用固體火箭發動機的月球探測器具有顯著的經濟優勢，本文針對這類經濟型月球探測器精確定點軟著陸主減速段的實時在線制導算法進行了詳細的研究，將探測器主減速段的運動分為縱向平面運動和橫向運動兩部分分別設計制導律。縱向平面運動采用改進顯式制導算法和點火角搜索技術，實現了對高程和航程的精確控制；對于橫向運動的控制，首先基于ZEM/ZEV最優反饋制導律給出橫向制導指令，然后根據發動機特性，采用PWPF將連續形式的制導指令調制為脈沖形式。仿真結果表明：

1)點火角搜索算法能夠在很短的時間內優化搜索得到滿足終端航向位置約束的點火角修正值，具有很好的實時性和穩定性；

2)飛行過程中運動狀態和制導指令變化曲線平滑，三個方向上的終端著陸誤差均控制在很小的范圍內，同時所設計的改進顯式制導算法解決了末段局部時間內姿態角突變的問題；

3)所設計制導算法具有自主性、實時性和一定的魯棒性，在一定的目標著陸區域內都具有很好的適應能力。