李海洋，寶音賀西

（清華大學(xué)航天航空學(xué)院，北京 100084）

引 言

電推進(jìn)與傳統(tǒng)的化學(xué)推進(jìn)相比，具有比沖高、燃料質(zhì)量消耗小等優(yōu)點(diǎn)，探測器可以攜帶更多的有效載荷。美國的“黎明號”探測器、日本的“隼鳥號”探測器等都采用了電推進(jìn)方式。在探測任務(wù)的初始設(shè)計階段，往往有著大量的可選軌道設(shè)計方案，每一屆的國際軌道競賽（Global Trajectory Optimization Com‐petition，GTOC）和國內(nèi)軌道競賽的問題，任務(wù)設(shè)計期初需要求解復(fù)雜的全局優(yōu)化問題。

小推力軌跡的設(shè)計與優(yōu)化問題，相對于脈沖軌跡的設(shè)計與優(yōu)化問題，求解較為困難，求解速度也較慢。由于計算能力與時間要求，不可能在全局優(yōu)化的過程中對每一個方案都進(jìn)行精確的小推力數(shù)值求解，這就需要對小推力轉(zhuǎn)移進(jìn)行快速準(zhǔn)確的估計。最重要的是快速準(zhǔn)確地估計小推力轉(zhuǎn)移的燃料消耗。

考慮交會問題的燃料最優(yōu)小推力轉(zhuǎn)移。目前最為常用的估計方法，是用雙脈沖轉(zhuǎn)移Lambert問題的速度增量來近似估計小推力轉(zhuǎn)移的速度增量，但是二者往往有著較大的誤差。除此之外，Gatto 等提出了基于Edelbaum 近似的估計方法[1]，李海洋等同樣基于Edelbaum 提出了考慮J2攝動下小推力多圈轉(zhuǎn)移的估計方法[2]。

本文提出一種基于軌道根數(shù)特征的機(jī)器學(xué)習(xí)的估計方法。機(jī)器學(xué)習(xí)作為人工智能的一個研究領(lǐng)域，在軌跡優(yōu)化領(lǐng)域也有應(yīng)用[3-4]。Hennes等[5]研究了機(jī)器學(xué)習(xí)方法估計小推力燃料消耗，比較了不同的機(jī)器學(xué)習(xí)方法對估計結(jié)果的影響，但是在機(jī)器學(xué)習(xí)方法中，特征的選擇對結(jié)果影響更為基礎(chǔ)與重要，主要創(chuàng)新點(diǎn)是研究不同的特征組合對于估計結(jié)果的影響。

本文使用第7 屆國際軌道競賽的小天體數(shù)據(jù)[6]，生成小推力轉(zhuǎn)移數(shù)據(jù)庫供機(jī)器學(xué)習(xí)使用，選擇不同的特征組合作為機(jī)器學(xué)習(xí)的特征，學(xué)習(xí)完成后即可對燃料最優(yōu)小推力轉(zhuǎn)移燃料消耗進(jìn)行準(zhǔn)確快遞估計。不同的特征組合結(jié)果不同，但都明顯優(yōu)于Lambert 估計方法。

1 小推力轉(zhuǎn)移數(shù)據(jù)庫生成

1.1 Lambert估計

使用第7屆國際軌道競賽的小天體數(shù)據(jù)，隨機(jī)選取出發(fā)天體A1和到達(dá)天體A2，在[0，3 000]天內(nèi)隨機(jī)確定出發(fā)時刻ts、在[100，300]天內(nèi)隨機(jī)確定轉(zhuǎn)移時間Δt，在[0.5，1]×minitial內(nèi)隨機(jī)確定初始質(zhì)量ms。燃料最優(yōu)小推力問題較難求解，不可能每次隨機(jī)生成的轉(zhuǎn)移都可以采用小推力實(shí)現(xiàn)。所以，首先求解雙脈沖Lambert 問題獲得轉(zhuǎn)移的速度增量，滿足一定條件的轉(zhuǎn)移再使用小推力實(shí)現(xiàn)。

若求解雙脈沖Lambert 問題所得的速度增量為ΔVL，小推力能產(chǎn)生的最大推力為Tmax，如果滿足下式

則認(rèn)為該段轉(zhuǎn)移可以轉(zhuǎn)為小推力求解。式中c的范圍一般取0.5～0.8之間，本文計算中取0.7。

Lambert估計的燃料消耗為

1.2 燃料最優(yōu)小推力轉(zhuǎn)移

經(jīng)過雙脈沖Lambert估計篩選后，可行的轉(zhuǎn)移即可進(jìn)行燃料最優(yōu)小推力求解[7]。

探測器初末狀態(tài)采用春分點(diǎn)軌道根數(shù)描述為

則小推力問題的動力學(xué)方程為

其中：M為 6×3 的矩陣，未列元素為 0；D為 6×1 的向量。

對于交會問題，約束條件為

對于燃料最優(yōu)問題，優(yōu)化指標(biāo)為

本文采用間接法求解燃料最優(yōu)小推力問題。為了更加快速有效地求解協(xié)態(tài)變量初值，采用對數(shù)同倫方法與協(xié)態(tài)變量歸一化技術(shù)[7]。

引入對數(shù)同倫參數(shù)ε和歸一化協(xié)態(tài)變量λ0，指標(biāo)函數(shù)變?yōu)?/p>

該最優(yōu)控制問題的哈密頓函數(shù)為

根據(jù)極大值原理，可以求得最優(yōu)控制率為

其中：開關(guān)函數(shù)ρ的表達(dá)式為

協(xié)態(tài)方程為

燃料最優(yōu)小推力轉(zhuǎn)移問題轉(zhuǎn)換為求解協(xié)態(tài)變量初值的兩點(diǎn)邊值問題，可以使用打靶法求解。打靶方程為

當(dāng)同倫參數(shù)ε從1 減到0，問題即變?yōu)槿剂献顑?yōu)問題。

2 機(jī)器學(xué)習(xí)方法

機(jī)器學(xué)習(xí)方法的本質(zhì)可以理解為函數(shù)的擬合方法。對于燃料最優(yōu)小推力轉(zhuǎn)移問題，其燃料消耗是初末軌道參數(shù)、初始質(zhì)量與轉(zhuǎn)移時間的函數(shù)，但是該函數(shù)目前無法使用數(shù)學(xué)推導(dǎo)解析地獲得準(zhǔn)確形式。機(jī)器學(xué)習(xí)的過程則是構(gòu)造擬合函數(shù)的過程，學(xué)習(xí)完成后構(gòu)造的擬合函數(shù)可以用來描述相應(yīng)的數(shù)學(xué)問題。

本文中的數(shù)學(xué)問題則是在二體太陽引力場中，給定小推力的最大推力和比沖，已知某組初末軌道參數(shù)、初始質(zhì)量與轉(zhuǎn)移時間，如何求解燃料消耗。只要滿足該數(shù)學(xué)問題的場景，均可以采用本文得到的機(jī)器學(xué)習(xí)擬合函數(shù)進(jìn)行燃料消耗估計，但是對于不同的場景，函數(shù)的擬合度可能不同，在實(shí)際使用過程中，還需要根據(jù)具體場景選擇合適的學(xué)習(xí)數(shù)據(jù)。使用機(jī)器學(xué)習(xí)方法進(jìn)行函數(shù)擬合的優(yōu)勢在于，機(jī)器學(xué)習(xí)方法對于復(fù)雜的函數(shù)關(guān)系可以有很好的描述，而解析公式必定會由于某些簡化從而失去復(fù)雜函數(shù)關(guān)系的某些特性。機(jī)器學(xué)習(xí)有著較為成熟的工具包，本文計算采用Py‐thon 語言的 scikit-learn 工具包[7]，參考為文獻(xiàn)[5]中得到的結(jié)論使用Gradient Boosting 回歸器，最大深度選擇為8。

根據(jù)第2 節(jié)的計算方法，一共生成6 萬組燃料最優(yōu)小推力轉(zhuǎn)移，每組轉(zhuǎn)移中包含了初末軌道參數(shù)、轉(zhuǎn)移時間、初始質(zhì)量和燃料消耗，其中，初末軌道參數(shù)、轉(zhuǎn)移時間、始質(zhì)量是基本變量，燃料消耗是希望通過機(jī)器學(xué)習(xí)獲得的擬合函數(shù)的函數(shù)值。其中5萬組用來進(jìn)行學(xué)習(xí)燃料最優(yōu)小推力轉(zhuǎn)移燃料消耗，其機(jī)器學(xué)習(xí)標(biāo)簽即為燃料消耗質(zhì)量；1萬組用來進(jìn)行檢驗(yàn)機(jī)器學(xué)習(xí)估計的燃料消耗與小推力精確計算得到的燃料消耗的差異。

機(jī)器學(xué)習(xí)最重要的是特征參數(shù)的選擇。燃料最優(yōu)的小推力轉(zhuǎn)移中，燃料消耗質(zhì)量由探測器初始質(zhì)量、轉(zhuǎn)移時間，初始位置速度和末端位置速度決定。表1參數(shù)可以作為機(jī)器學(xué)習(xí)的特征。

表1 可選特征Table 1 Possible attributes

本文選用表2中的7種特征組合進(jìn)行分析，組合1表示的是利用位置速度的軌道描述方法作為特征，組合2表示的是利用經(jīng)典軌道根數(shù)的軌道描述方法作為特征，組合3表示的是利用春分點(diǎn)軌道根數(shù)的軌道描述方法作為特征，組合4～6是在之前3組的基礎(chǔ)上，再加上Lambert估計的速度增量和燃料消耗作為特征。組合7是在組合6的基礎(chǔ)上考慮軌道積分中的能量與動量。

表2 6種特征組合Table 2 Six combinations of features

為了更為清楚地判斷估計結(jié)果是否準(zhǔn)確，采用估計燃料消耗與實(shí)際燃料消耗的均方差和均方根差作為估計精度的衡量標(biāo)準(zhǔn)。

3 算 例

本文使用第7屆國際軌道競賽的小天體數(shù)據(jù)，共16 256 個天體，軌道采用二體開普勒軌道，初始軌道根數(shù)的時刻為56 800MJD。探測器的初始質(zhì)量為2 000 kg，小推力產(chǎn)生的最大推力為0.3 N，比沖為3 000 s。每組燃料最優(yōu)小推力轉(zhuǎn)移生成中，從16 256個天體中隨機(jī)選取出發(fā)天體和到達(dá)天體，再隨機(jī)確定出發(fā)時刻，轉(zhuǎn)移時間和初始質(zhì)量，根據(jù)第2節(jié)中的方法判斷并生成小推力轉(zhuǎn)移。使用該方法共生成6萬組數(shù)據(jù)，其中1萬組用來檢驗(yàn)估計結(jié)果的小推力轉(zhuǎn)移軌跡信息如表3所示。

下面對Lambert 估計和6 種特征組合方式的機(jī)器學(xué)習(xí)結(jié)果進(jìn)行分析，結(jié)果如表4。

表3 小推力轉(zhuǎn)移軌跡信息Table 3 Information of low-thrust transfers

表4 估計結(jié)果比較Table 4 Estimation results

綜合來看，燃料最優(yōu)小推力轉(zhuǎn)移燃料消耗機(jī)器學(xué)習(xí)估計方法明顯優(yōu)于Lambert估計方法。機(jī)器學(xué)習(xí)估計方法中，帶有Lambert估計的速度增量和燃料消耗作為學(xué)習(xí)特征的特征組合明顯更優(yōu)。組合1、2、3的特征是構(gòu)建一段小推力轉(zhuǎn)移所需要的最基本的參數(shù)，但是作為機(jī)器學(xué)習(xí)的特征，為了得到更準(zhǔn)確的估計結(jié)果，這些基本參數(shù)還需要額外信息，比如Lambert估計的結(jié)果。

對于機(jī)器學(xué)習(xí)方法來說，Lambert 估計的結(jié)果提供了一個學(xué)習(xí)的方向，使得機(jī)器可以在Lambert估計結(jié)果的基礎(chǔ)上，更加快速準(zhǔn)確地逼近真實(shí)值；而對于沒有Lambert估計結(jié)果的特征組合，機(jī)器在學(xué)習(xí)過程中完全憑借轉(zhuǎn)移所需的最基本參數(shù)，可能需要耗費(fèi)更多的精力去學(xué)習(xí)，學(xué)習(xí)效果也相對較差。

對于不同的軌道描述方法，無Lambert估計特征的特征組合結(jié)果非常接近，有Lambert估計特征的特征組合中，春分點(diǎn)軌道根數(shù)的估計結(jié)果最優(yōu)，經(jīng)典軌道根數(shù)的估計結(jié)果次優(yōu)。機(jī)器學(xué)習(xí)最重要的是特征參數(shù)的選擇，特征的選擇決定了機(jī)器學(xué)習(xí)到的是問題的表面還是本質(zhì)。軌道根數(shù)的描述方法，相比于位置速度表示法，更能體現(xiàn)出軌道的本質(zhì)屬性，軌道根數(shù)表示法更能夠直接表現(xiàn)出2個軌道的差異，所以通過軌道根數(shù)學(xué)習(xí)得到的結(jié)果要明顯優(yōu)于通過位置速度學(xué)習(xí)得到的結(jié)果。春分點(diǎn)軌道根數(shù)相比于經(jīng)典軌道根數(shù)無奇異，更適合于小傾角近圓軌道，而本文使用的數(shù)據(jù)庫中小天體多為小傾角近圓軌道，故結(jié)果最優(yōu)。二體軌道中軌道積分應(yīng)該是最能夠反映軌道本質(zhì)的參數(shù)，但是在本文中考慮軌道積分后的結(jié)果并無明顯改觀，仍然需要進(jìn)一步進(jìn)行探究。

特征組合3和特征組合6的誤差分布圖如圖1～2所示，圖1中橫坐標(biāo)為機(jī)器學(xué)習(xí)估計值與實(shí)際值的誤差，縱坐標(biāo)為個數(shù)，藍(lán)色表示機(jī)器學(xué)習(xí)估計結(jié)果，橘紅色表示Lambert 估計結(jié)果。圖1中特征組合3 的機(jī)器學(xué)習(xí)結(jié)果誤差雖然平均值接近0，但是Lambert 估計的誤差分布更為集中，圖2中特征組合6的機(jī)器學(xué)習(xí)結(jié)果具有更優(yōu)的誤差平均值和集中程度。

圖1 特征組合3的誤差分布Fig.1 Error distribution of Features III

圖2 特征組合6的誤差分布Fig.2 Error distribution of features VI

4 結(jié) 論

本文采用機(jī)器學(xué)習(xí)的方法對燃料最優(yōu)小推力轉(zhuǎn)移的燃料消耗進(jìn)行估計，比目前最為常用的Lambert估計方法結(jié)果明顯更優(yōu)。不同的特征組合會影響機(jī)器學(xué)習(xí)的結(jié)果。相比于初末軌道、初始質(zhì)量與轉(zhuǎn)移時間的特征組合，增加了Lambert估計結(jié)果的特征組合估計結(jié)果更優(yōu)。在位置速度、經(jīng)典軌道根數(shù)、春分點(diǎn)軌道根數(shù)3種軌道描述方法中，春分點(diǎn)軌道根數(shù)描述方法作為特征的機(jī)器學(xué)習(xí)估計結(jié)果最為準(zhǔn)確。故帶有Lambert 估計特征的春分點(diǎn)軌道根數(shù)的特征組合為較好的機(jī)器學(xué)習(xí)特征組合。

本文只考慮了不同特征對估計結(jié)果的影響，不同的機(jī)器學(xué)習(xí)算法與軌道描述方法之間的關(guān)聯(lián)同樣也會對估計結(jié)果產(chǎn)生影響，未來可以考慮研究使用更為準(zhǔn)確的特征組合，以及不同的機(jī)器學(xué)習(xí)算法進(jìn)行估計。