高 峰， 陳少松， 牛津橋

(1.中國(guó)空間技術(shù)研究院， 北京 100094； 2.中國(guó)人民解放軍 63993 部隊(duì)， 北京 101149)

1 引言

隨著深空探測(cè)技術(shù)的不斷發(fā)展，人類(lèi)對(duì)月球探測(cè)的目的發(fā)生了本質(zhì)改變，由原來(lái)的飛行技術(shù)驗(yàn)證、地外科學(xué)實(shí)驗(yàn)向特定區(qū)域(月背/兩極/陰影區(qū))探測(cè)、表面大尺度巡視轉(zhuǎn)變，這對(duì)著陸系統(tǒng)的精確性、強(qiáng)魯棒、自適應(yīng)性提出了更高的要求。

近幾十年來(lái)，軟著陸制導(dǎo)技術(shù)獲得了長(zhǎng)足的發(fā)展。 從早期的蘇聯(lián)月球9 號(hào)采用的建立月球垂線法[1－2]到半開(kāi)環(huán)半閉環(huán)的重力轉(zhuǎn)彎制導(dǎo)方法[3－6]；從結(jié)合軌跡優(yōu)化[7]與軌跡跟蹤[8－11]的標(biāo)稱(chēng)軌跡制導(dǎo)方法到能夠在線應(yīng)用的顯示制導(dǎo)方法[12－14]以及最近興起的智能學(xué)習(xí)制導(dǎo)方法[15－18]。 從以上技術(shù)發(fā)展變化可以總結(jié)出軟著陸制導(dǎo)技術(shù)的發(fā)展趨勢(shì)呈現(xiàn)3 個(gè)階段[19]:基于人為智能的給定軌道直接著陸、基于機(jī)器智能的軌跡優(yōu)化以及基于人工智能的學(xué)習(xí)制導(dǎo)方法。

傳統(tǒng)的給定軌道及優(yōu)化思想實(shí)現(xiàn)的軟著陸制導(dǎo)很大程度上依賴于先驗(yàn)知識(shí)，比如模型、約束、環(huán)境等，同時(shí)容易受建模偏差、不確定擾動(dòng)等因素影響，且現(xiàn)有的障礙模型不具備較好的泛化能力。 因此在月球這樣多障礙、安全區(qū)域狹小和不確定性大的環(huán)境中，傳統(tǒng)方法難以滿足軟著陸所需的的精度要求和安全性要求。 近幾年，凸優(yōu)化技術(shù)逐漸成為研究的熱點(diǎn)，基于凸優(yōu)化的軌跡優(yōu)化方法成為解決上述問(wèn)題以及在線應(yīng)用 的 最 佳 方 案[20]。 Acikmese 等[21－22]提 出了將非凸優(yōu)化問(wèn)題進(jìn)行凸化轉(zhuǎn)化的思想，解決了以發(fā)動(dòng)機(jī)推力幅值的上下界作為約束條件的燃耗最優(yōu)問(wèn)題；Blackmore 等[23]的研究已經(jīng)擴(kuò)展到尋找最小的著陸誤差軌跡的燃耗最優(yōu)問(wèn)題；Bai 等[24]提出了利用非凸函數(shù)對(duì)障礙物進(jìn)行合并，將這些非凸障礙約束轉(zhuǎn)化為凸障礙約束的線性化方法。

嫦娥工程中面對(duì)復(fù)雜月背、極區(qū)、陰影區(qū)等條件，安全著陸區(qū)域很小，本文針對(duì)精確自主避障技術(shù)開(kāi)展深入研究，在凸錐模型的基礎(chǔ)上建立半球模型和通用凸包模型來(lái)描述月面障礙物，利用凸優(yōu)化技術(shù)，得到一種在多重障礙約束環(huán)境下的燃料最優(yōu)的制導(dǎo)算法，并進(jìn)行仿真實(shí)驗(yàn)。

2 軟著陸問(wèn)題描述

2.1 軟著陸數(shù)學(xué)模型

采用三自由度模型描述著陸器運(yùn)動(dòng)，考慮以月面固聯(lián)坐標(biāo)系為基礎(chǔ)進(jìn)行建模[21]，ox方向垂直于月面指向外空間，oy沿月面指向月球北極，oz與其他2 個(gè)坐標(biāo)軸構(gòu)成空間直角坐標(biāo)系，如圖1所示。

圖1 坐標(biāo)系定義Fig.1 Definition of the coordinate

基于此坐標(biāo)系，著陸器動(dòng)力學(xué)可以描述為式(1):

式中，r∈IR3是著陸器的位置向量，v∈IR3是著陸器的速度向量，g∈IR3是月球的引力加速度向量，Tc∈IR3是凈推力矢量，m是著陸器質(zhì)量，α是燃料消耗率。

推力器認(rèn)為是可變有限的，總推力受式(2)約束:

式中，T1、T2是推力的上下限。 在著陸器下降過(guò)程中，需要滿足式(3)的路徑約束:

考慮邊界條件、推力約束、著陸器傾角和表面高度，以燃料最優(yōu)為優(yōu)化指標(biāo)，軟著陸軌跡優(yōu)化問(wèn)題可描述為式(4)所示最優(yōu)控制問(wèn)題:

式中，t0為起始時(shí)刻，tf為終端時(shí)刻，m0為起始質(zhì)量，r0，rf分別為起始位置和終端位置，v0，vf分別為起始速度和終端速度。

2.2 障礙物數(shù)學(xué)模型

行星表面的障礙物具有不同的形狀和高度。Bai 等[24]將障礙物建模為高度可調(diào)、半錐角可調(diào)的凸錐，本文在此基礎(chǔ)上將其擴(kuò)展為半球和通用凸包模型。

在凸錐的情況下，著陸器不與障礙物碰撞的條件如式(5)所示:

式中，ne＝(1，0，0)T是單位向量，φ為半錐角，錐頂?shù)奈恢孟蛄繛镠＝(Hx，Hy，Hz)T，著陸器的位置向量為P＝(rx，ry，rz)T。

除了將障礙物描述為凸錐，還可以將其描述為包圍障礙物的以其幾何中心點(diǎn)S＝(Sx，Sy，Sz)為球心、RS為半徑最小外接半球，如圖2 所示。

圖2 半球障礙模型Fig.2 Hemisphere obstacle model

在這種情況下，只要著陸器在半球的半徑范圍外就是安全可行的，那么著陸器不與障礙物碰撞的約束條件如式(6)所示:

將障礙物建模為凸錐或者半球，雖然處理起來(lái)非常方便，然而缺少通用性，不能反映不同類(lèi)型障礙物的幾何特征，故需要一種統(tǒng)一的建模方式，即通過(guò)傳感器獲得的點(diǎn)云數(shù)據(jù)直接建立三維凸包。 本文使用隨機(jī)增量式算法來(lái)構(gòu)造點(diǎn)集T的凸包。 首先定義凸包可見(jiàn)區(qū)域和地平線，如圖3所示。

圖3 凸包可見(jiàn)區(qū)域和地平線Fig.3 Visible area and horizon of convex hull

假設(shè)Ti落在凸包lh(Ti－1) 之外，在點(diǎn)Ti的位置朝lh(Ti－1) 看去，在lh(Ti－1) 的表面上，那些可見(jiàn)的小平面連在一起組成的連通區(qū)域稱(chēng)之為點(diǎn)Ti在凸包lh(Ti－1) 上的可見(jiàn)區(qū)域。 而圍成這個(gè)區(qū)域的lh(Ti－1) 上的邊組成的折線稱(chēng)之為點(diǎn)Ti在凸包lh(Ti－1) 上的地平線[25]。 具體算法過(guò)程如圖4 所示。

圖4 凸包生成算法Fig.4 Generation algorithm of convex hull

容易證明在三維空間中，任意n個(gè)點(diǎn)的凸包的復(fù)雜度為o(n)[25]，得到凸包以及凸包所包含的點(diǎn)集后，給出建立障礙約束的方法。 假設(shè)最終得到的凸包由NT個(gè)凸頂點(diǎn)組成，定義函數(shù)f(P)＝∑‖P(t)－Ti‖/NT，其中Ti是凸頂點(diǎn)的位置向量。 該函數(shù)反映了點(diǎn)與凸包之間的相對(duì)距離關(guān)系，是IRn上的凸函數(shù)，也就是說(shuō)只要函數(shù)大于凸包邊界點(diǎn)函數(shù)的最大值，則點(diǎn)在凸包外。 為了更均勻地反映函數(shù)與點(diǎn)的位置關(guān)系，可將最大值用平均值替代，則著陸器不與障礙物相撞的條件如式(7)所示:

式中，NT為凸包障礙點(diǎn)集個(gè)數(shù)。 綜上，將一般性的月球軟著陸避障問(wèn)題描述為式(8)所示:

式中，H為障礙錐錐頂?shù)奈恢孟蛄?，S為半球障礙球心的位置向量，T為凸點(diǎn)集的位置向量，?是半錐角φ的補(bǔ)角。

3 基于凸優(yōu)化的軌跡優(yōu)化方法

3.1 凸化轉(zhuǎn)化

2.2 節(jié)式(8)的非凸性會(huì)造成計(jì)算求解困難，Acikmese 等[21]給出了將非凸約束轉(zhuǎn)化為凸約束的方法，令u＝Tc/m，σ＝Γ/m，z＝ln(m) ，其中Γ為推力松弛向量，則式(4)中的非凸約束將被轉(zhuǎn)化為式(9):

然而線性函數(shù)中省略了f(P) 的泰勒展開(kāi)式中的二階項(xiàng)，令向量(Pξ－Q) 和向量(P－P0)之間的夾角為l，容易求得f(P) 和其二階余項(xiàng)的相對(duì)誤差如式(12)所示:

通過(guò)在節(jié)點(diǎn)處連續(xù)線性化障礙約束，式(13)便可以表示為式(14)的形式:

式中，s為松弛因子。 根據(jù)式(12)可知f(P)的二階余項(xiàng)為正，這使得約束更加嚴(yán)格。 圖5 可以看出，著陸器的實(shí)際可飛行空間區(qū)域大于理論可飛行空間區(qū)域，故計(jì)算所得解是有效的。

圖5 障礙轉(zhuǎn)化前后關(guān)系Fig.5 Relationship before and after barrier transformation

線性化后的問(wèn)題由于縮小了解的范圍有可能導(dǎo)致問(wèn)題無(wú)解，故引入松弛因子來(lái)適當(dāng)放寬約束條件，以在避障和最優(yōu)解之間取得平衡。

3.2 二階錐優(yōu)化

在整個(gè)控制時(shí)域中，每個(gè)離散時(shí)間點(diǎn)的狀態(tài)滿足式(17):

4 實(shí)驗(yàn)仿真與結(jié)果

4.1 仿真條件

本文算法通過(guò)Matlab 和Yalmip 聯(lián)合求解，初始參數(shù)選擇如表1 所示。

表1 初始參數(shù)Table 1 Initial parameters

4.2 仿真驗(yàn)證

本文在障礙物分別為凸錐、半球和凸包的多障礙環(huán)境下對(duì)凸優(yōu)化算法進(jìn)行驗(yàn)證。

4.2.1 凸錐多障礙

選擇3 個(gè)障礙物作為約束條件，并且這3 個(gè)障礙物最好能夠排成1 個(gè)平面以凸顯避障的效果，其參數(shù)如表2 所示。

表2 凸錐障礙參數(shù)Table 2 Parameters of convex hull obstacle

仿真結(jié)果如圖6 所示，多重障礙物之間明顯存在可行路徑，但是安全區(qū)域較小。 從圖6 可以看出著陸器連續(xù)繞開(kāi)了3 個(gè)具有不同半錐角的障礙物，著陸位置為[0.0697，0.1070，0.9918]m，著陸速度為[－0.5184，－0.1682，－0.1468]m/s，消耗燃料384.4954 kg，得到了全局最優(yōu)解，說(shuō)明面對(duì)多重障礙時(shí)凸優(yōu)化算法有效。

4.2.2 半球多障礙

相對(duì)于凸錐障礙而言，多重半球障礙環(huán)境下著陸器的可行安全飛行空間更窄，通常情況下需要穿越2 個(gè)相鄰半球之間的縫隙，并且經(jīng)常需要采取爬升的策略進(jìn)行避障。 半球障礙參數(shù)如表3所示，其仿真結(jié)果如圖7 所示，可見(jiàn)飛行軌跡穿越黃色和橙色半球之間的縫隙，繞開(kāi)綠色半球從而到 達(dá) 目 標(biāo) 點(diǎn)； 著 陸 位 置 為[0.0697， 0.1071，0.9918] m， 著 陸 速 度 為[－0.4482， 0.1854，－0.121]m/s，消耗燃料286.8325 kg。

圖7 半球仿真結(jié)果Fig.7 Results of hemisphere simulation

4.2.3 凸包多障礙

點(diǎn)云圖是由計(jì)算機(jī)隨機(jī)生成的，以模擬傳感器所掃描的障礙的點(diǎn)云數(shù)據(jù)，如圖8(a)所示，由點(diǎn)云生成凸包的結(jié)果如圖8(b)所示。

圖8 凸包生成結(jié)果Fig.8 Results of convex hull generation

得到的凸包由2 個(gè)不規(guī)則多面體組成，相比前2 種模型而言，凸包更能反映障礙物的幾何特征，可飛行區(qū)域也是最精確的。 仿真結(jié)果如圖9所示，著陸器為了穿越狹窄的空間，首先采取了拉低的策略穿過(guò)低矮的峽谷，而后面對(duì)斜擋的障礙采取爬升的策略飛掠障礙表面，在極小的飛行空間內(nèi)找到了一條安全可行的軌跡；著陸位置為[0.0697， 0.107， 0.9918] m， 著 陸 速 度 為[－0.1762，－0.4675，－0.1171] m/s，消耗燃料361.8431 kg。

圖9 凸包仿真結(jié)果Fig.9 Results of convex hull simulation

凸包障礙模型的仿真結(jié)果表明對(duì)于任何形式的障礙物通過(guò)式(7)和式(14)建立的凸包約束均適用，不再受特定形式障礙約束的限制，大多數(shù)的物體都可以通過(guò)凸包來(lái)更好地近似，克服了凸錐和半球?qū)φ系K物近似效果不好的缺點(diǎn)。 盡管與凸錐和球體相比，凸包之間的計(jì)算更為復(fù)雜，但是相對(duì)于非凸的物體而言，還是要簡(jiǎn)單得多。

5 結(jié)論

1)在凸錐模型的基礎(chǔ)上建立了半球模型和通用凸包模型來(lái)描述月面障礙物，其中凸包模型能夠更準(zhǔn)確地對(duì)障礙物進(jìn)行近似；

2)利用凸優(yōu)化技術(shù)得到一種在多重障礙約束環(huán)境下的燃料最優(yōu)的制導(dǎo)算法，解決了月面多障礙環(huán)境下的精確自主避障問(wèn)題。