遲哲敏，李俊峰，蔣方華，寶音賀西

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

0 引 言

目前，各國正在積極開展深空探測任務，對具有科學價值的天體的探測已經成為深空探測領域的研究熱題。我國計劃在2020年發(fā)射第一個火星探測器“天問”一號，同時也籌劃著未來的小行星探測任務。遠距離探測活動的開展離不開推進系統(tǒng)的選擇。1998年NASA成功發(fā)射的“深空”一號首次驗證了電推進技術在深空探測領域的應用。由于電推進的比沖是化學推進的10倍左右，所以電推進在減少燃料消耗、增加有效載荷質量方面更具有優(yōu)勢。電推進相較于化學推進可以提供“連續(xù)”的小推力，因此發(fā)動機的壽命較長，可以工作15000小時以上，適合完成長距離、長時間的飛行任務。

常見的采用電推進方式的發(fā)動機有三種，首先是美國艾德阿斯特拉火箭公司在研的磁等離子體火箭(VASIMR)[1]，它在恒定的大功率下工作，通過控制比沖來調節(jié)推力等參量的大小，計劃采用核能作為其動力來源，因此有利用VASIMR完成30天到達火星的構想。另外兩種電推進發(fā)動機——霍爾效應推進器和離子電推進發(fā)動機均已經應用于深空探測任務中?；魻栃七M器利用基于電磁角力的霍爾效應原理來提供飛行動力[2]。2003年歐空局發(fā)射的第一個月球探測器“智慧”一號(SMART-1)就是將霍爾效應推進器作為主推進器，成功完成了月球探測任務[3]。而離子電推進發(fā)動機是在高溫下將工質電離，再令離子在電場中加速噴出從而產生推力，目前普遍選擇氙氣作為工質。近幾年，離子電推進發(fā)動機越來越廣泛地應用于深空探測任務中，NASA先后研制了NSTAR(NASA Solar Technology Application Readiness)和NEXT(NASA Evolutionary Xenon Thruster)離子發(fā)動機，NEXT發(fā)動機在推進效能和燃料節(jié)省方面均優(yōu)于NSTAR發(fā)動機，表1展示的是NSTAR和NEXT發(fā)動機的部分性能參數[4]?！吧羁铡币惶栄b載一個NSTAR 發(fā)動機[5]，發(fā)動機工作的功率由太陽能帆板提供。“黎明”號攜帶3個NSTAR發(fā)動機先后探測了灶神星和谷神星[6]，成為了第一個環(huán)繞2顆地球外天體的探測器。NEXT發(fā)動機也在2010年測試成功，原定被應用于2024年探測彗星67P的“凱撒”號任務中，然而任務競標失敗。它的商業(yè)化版NEXT-C計劃于2021年作為“雙小行星重定向測試”任務(Double Asteroid Redirection Test, DART)[7]的主推進器。JAXA的小行星采樣任務Hayabusa和Hayabusa2均采用μ10離子電推進發(fā)動機[8]，Hayabusa在2005年11月完成了對小行星Itokawa的采樣任務，而Hayabusa2也于2019年11月在小行星Ryugu上進行采樣，Hayabusa2所采用的μ10推進器是Hayabusa任務推進器的升級版[9]，它們均是變推力變比沖發(fā)動機。英國的T6離子發(fā)動機[10]正應用于在飛的Bepicolombo任務中，Bepicolombo借助一次地球、兩次金星和六次水星引力輔助計劃于2025年到達水星。

表1 NSTAR推進器和NEXT推進器的性能特征對比

上述電推進發(fā)動機的比沖均是變化的，因此變比沖電推進模型是更符合工程實際的推進系統(tǒng)模型。而對于變比沖連續(xù)小推力軌跡優(yōu)化的研究是具有工程應用意義的。由于連續(xù)小推力自身的特點和太陽能提供功率的限制，以及考慮多目標或多次引力輔助探測問題而引起的多內點約束，變比沖電推進探測任務的設計仍面臨著困難和挑戰(zhàn)。傳統(tǒng)的直接法和間接法也可以適用于解決變比沖電推進軌跡優(yōu)化問題，同時也需要針對具體模型和具體問題進行方法上的改進。本文對目前常見的變比沖電推進模型以及常用的變比沖連續(xù)小推力軌跡優(yōu)化方法的研究現狀進行了相關介紹。

1 變比沖電推進模型

目前，在深空軌跡計算中，比較常用的變比沖電推進簡化模型有三種，分別對應著不同真實發(fā)動機系統(tǒng)的工作原理和方式。

首先是比沖在區(qū)間內變化的推進模型。在這一模型下，將比沖作為控制量，控制其在區(qū)間范圍內變化，優(yōu)點是若在需要減小燃料消耗的時候，令比沖取范圍中的較大或是最大值，若需要在一個位置加速，那么又會取比沖范圍中相應的較小值，根據不同情況可以靈活調節(jié)[11]。最符合這一模型的是磁等離子體電推進發(fā)動機，因為考慮采用核能供電，所以發(fā)動機輸入功率是恒定的，比沖是可控的。Seywald等[12]就是構建了VASIMR發(fā)動機模型來求解共面圓軌道之間的燃料最優(yōu)軌跡，Park等[13]也是應用VASIMR模型求解到達越地小行星的交會軌跡，為地球防御任務的實施提供參考；Kechichian[14]也應用了恒定功率變比沖模型，在他的研究中，將春分點軌道根數代替位置速度來表示用于優(yōu)化的狀態(tài)量。然而也存在另一種功率模型假設：在太陽能電推進模式下，功率是星日距離的函數，但比沖仍是在區(qū)間內變化并可被當作控制量[15-17]。Casalino等[18]在此變比沖模型的基礎上構建了雙比沖模型，即比沖只能取范圍的最大和最小值，再對比了兩種模型下的燃料消耗；他們還將變比沖模型用于求解從近地軌道向地球同步軌道轉移的多圈問題[19]。Mengali等[20]在變功率變比沖模型的基礎上還考慮了推進效率的影響，推進效率是比沖的多項式函數，根據擬合程度有線性、二次、三次等函數關系，陳楊[21]在博士論文中將其中的線性關系用于行星際間的軌跡優(yōu)化。Taheri等[22]在研究中考慮了多模式的推進系統(tǒng)模型，在他們的模型構建中，分別將發(fā)動機輸入功率和比沖作為控制量來求解燃料最優(yōu)軌跡。

另一種是比沖隨發(fā)動機輸入功率變化的數學模型。這種情況下，發(fā)動機是由太陽能供電，因此比沖和推力也是星日距離的函數。此模型是最貼合光柵離子電推進發(fā)動機工作原理的，發(fā)動機的輸入功率作用于兩部分：等離子體的產生和離子的加速。等離子體生成器是在等離子體生成功率的作用下工作的，在燃燒室中，燃料中的中性原子在電子轟擊下電離，生成的離子在一組光柵中靜電加速，這組光柵的控制電壓為高電壓，相對應的功率取決于高電壓和靜電加速的電流。一旦給定一個高電壓，性能參數：推力和比沖，確定為是發(fā)動機輸入功率的函數。圖1是NSTAR光柵離子發(fā)動機的原理圖[23]。

圖1 NSTAR離子發(fā)動機[23]

發(fā)動機試驗得出的結果是一系列相關的推力隨功率以及比沖隨功率變化的離散點(參見圖2)，因此相應的函數關系是通過曲線擬合得到，最常見的是多項式擬合。需要解釋的是，一般通過多項式擬合的是推力和質量流率對輸入功率的函數關系，而比沖是通過推力和質量流率計算得到。William等[24]在1997年已經對將應用于“深空”一號探測任務的NSTAR離子電推進發(fā)動機的參數進行線性擬合和評估，Petropoulos等[25-26]采用此線性擬合模型進行小推力引力輔助軌跡設計。William等[27]在火星探測軌跡的求解中應用了NSTAR的質量流率是功率的二次多項式擬合函數模型。更精確一些的NSTAR和NEXT發(fā)動機的四次多項式擬合曲線[28]也用于求解行星際轉移的燃料最優(yōu)軌跡。David等[29]分別根據用于NASA的五個深空探測任務的NEXT發(fā)動機的參數離散點給出相應的“大推力”和“高比沖”四次多項式擬合模型。Saripalli等[30]分別擬合了不同電壓下的、“大推力”和“高比沖”情況下的推力-功率以及比沖-功率四次曲線，并且將它們用于任務計算來對比各種模型得到的軌跡的最優(yōu)性。一系列霍爾效應發(fā)動機的參量擬合關系也用于軌跡優(yōu)化的計算中[31-34]。Casalino等[35]在電推進三次曲線擬合模型下研究了功率作為控制量的多模式推進系統(tǒng)的軌跡優(yōu)化問題。國內的孫沖等[36]應用了太陽帆電推進的二次多項式擬合模型進行虛擬引力場的軌道設計。除此之外，Karthik等[37]對微小衛(wèi)星上的電推進系統(tǒng)在考慮了電壓電流等各種因素的影響下進行了精確建模。

(a)推力變化圖

最后是功率離散點模型或稱為功率分檔模型。這種模型對應兩種情況：首先是一個發(fā)動機具有多個工況點，推力和比沖隨著功率分檔變化；第二種情況是一個推進系統(tǒng)由多個恒定推力恒定比沖的發(fā)動機組成，發(fā)動機工作的功率由太陽能帆板提供，所以隨著功率減小，發(fā)動機依次關閉，也就導致推力隨著功率分檔變化[38]。Brophy等[39]給出了NSTAR發(fā)動機所有工作周期內的6個工況點數據，Chi等[38]應用此6個離散點的發(fā)動機模型求解火星探測的時間最優(yōu)與燃料最優(yōu)軌跡。Polk等[40]提供了用于“深空”一號任務的NSTAR的16個工況點試驗數據，并且將它們用于任務軌跡的驗證和評估。Quarta等[41-42]應用NEXT發(fā)動機的41個工況離散點[43]求解行星際轉移的燃料最優(yōu)軌跡，他們對代入離散點數據(功率、推力和比沖)的哈密頓函數進行排序來優(yōu)化軌跡。在不考慮功率和比沖的情況下，Axelrod等[44]研究了帶有推力離散檔位的發(fā)動機模型下的軌跡優(yōu)化方法，該方法既適用于多推力檔位的推進系統(tǒng)又適用于多個單向發(fā)動機組成的推進系統(tǒng)。

2 變比沖小推力軌跡優(yōu)化方法

2.1 軌跡優(yōu)化方法

變比沖電推進軌跡優(yōu)化問題的本質是連續(xù)小推力軌跡優(yōu)化問題。這里的“連續(xù)”和“小”是相對于脈沖化學推進而言的。常見的求解連續(xù)小推力最優(yōu)軌跡的數值方法分為兩種：一是將問題進行轉換，比如直接法，將原先連續(xù)控制問題轉化為參數優(yōu)化問題；第二種是針對問題的特點給出方法，連續(xù)小推力軌跡優(yōu)化問題等效于泛函極值問題，所以在此可應用基于變分法和最優(yōu)控制原理的間接法。下面分別介紹一下直接法和間接法的基本原理和特點。

2.1.1 直接法

直接法實際上是將連續(xù)控制問題轉化為參數優(yōu)化問題。對于包含狀態(tài)量、控制量、過程約束和性能指標的小推力軌跡優(yōu)化問題，直接法是將狀態(tài)量和(或)控制量離散，并將狀態(tài)約束條件和性能指標寫為數值微分或是數值積分的代數形式，所以原問題轉化為求解多個離散后序列的非線性規(guī)劃問題(Nonlinear Programming, NLP)[45]。常見的求解非線性規(guī)劃問題的方法是序列二次規(guī)劃算法(Sequence Quadratic Program, SQP)，常用的SNOPT軟件就是基于此算法來求解大規(guī)模的非線性規(guī)劃問題[46]。

通過上述描述，直接法包含一個離散的過程，這就涉及到兩個問題，首先是離散的過程中并不能確保狀態(tài)量和控制量的匹配，也不能保證得到的解可以滿足最優(yōu)控制的一階必要條件，所以直接法不能保證解的最優(yōu)性。再者為了提高求解的精度，需要增加離散點數量，會導致計算量增大，相應地也就會降低計算效率。但直接法也有很多優(yōu)點，比如離散過程簡單，解的收斂性高等。

直接法包括幾種常用的方法，比如直接配點法和直接打靶法。它們的區(qū)別體現在前者將動態(tài)約束條件寫為數值微分(有限差分或拉格朗日插值)的形式，后者是通過數值積分的形式來實現。偽譜法就是直接配點法的一種[47-49]，它是利用基于全局近似的拉格朗日插值的數值微分方法重構動力學方程。偽譜法可以改善上一段所敘述的應用直接法存在的問題，首先偽譜法自身構造的一階必要條件與最優(yōu)控制的一階必要條件等價[50]，所以可以提高解的最優(yōu)性；其次是在相同數量離散點的情況下，偽譜法相比其他直接法的求解精度更高[51]。Rao等[52-53]基于偽譜法和自適應離散網格法創(chuàng)建了MATLAB的GPOPS(Gauss Pseudospectral Optimization Software)軟件包，目前該軟件在軌跡優(yōu)化中廣泛應用。偽譜法還有一個優(yōu)勢在于，由于一階必要條件的等價性，所以偽譜法的解可以為間接法求解提供協(xié)態(tài)初值猜測[54]。但是利用偽譜法求解得到的Bang-Bang控制曲線存在數值振蕩等問題，所以結果的精度還是比間接法弱一些，Bai等[55]和郭鐵丁等[50,56]對此提出了相應的改進方法。

另外，基于直接法原理還開發(fā)了其他常用的程序包，比如MALTO、EMTG等，在后文的敘述中會對它們做簡單介紹。

2.1.2 間接法

應用間接法的求解思路是構造對應于狀態(tài)變量的協(xié)態(tài)變量，在性能指標的基礎上建立哈密頓函數，再推導歐拉-拉格朗日方程(協(xié)態(tài)微分方程)。同時積分狀態(tài)和協(xié)態(tài)微分方程，在邊界點處打靶，打靶方程包括狀態(tài)約束以及與協(xié)態(tài)變量相關的一階必要條件。通過龐得里亞金極大或是極小值原理推導得到最優(yōu)控制律。因此，“間接”是通過構造輔助優(yōu)化的協(xié)態(tài)變量來體現的。

間接法是將最優(yōu)控制問題轉化為兩點邊值問題(Two-Point Boundary-Value Problem, TPBVP)或多點邊值問題(Multi-Point Boundary-Value Problem, MPBVP)，MPBVP對應于多內點約束問題，比如多目標交會或是多次引力輔助問題。其本質都是求解非線性方程組，在此可以應用牛頓迭代法[57-58]等。Jiang等[59]采用Fortran里面的Minpack-1程序包求解TPBVP和MPBVP，Minpack-1是基于結合了牛頓法和梯度法的Powell混合法的非線性求解器[60]。Zhang等[61]通過構造解析的雅各比矩陣來估算打靶過程的梯度信息，以此提高牛頓迭代的收斂效率。

間接法可以有效克服應用直接法產生的困難，比如采用間接法得到的解可以滿足一階必要條件，至少能夠保證是局部最優(yōu)的，再者因為無需進行多參數離散求解，計算效率會相應地提高。但是因為引入協(xié)態(tài)變量且協(xié)態(tài)變量自身物理意義不明確，所以無法解析給定，如果對其猜測不準確，那么將一直無法收斂到最優(yōu)解。對于多內點約束問題，隨著協(xié)態(tài)初值數量的增多，收斂難度增大。因此應用間接法的最大困難是由于引入協(xié)態(tài)變量導致的較低的收斂效率。目前采用協(xié)態(tài)初值估計[59,62-63]、同倫方法[15,59,64-66]以及開關檢測[59,61,67]等一系列策略來提高間接法求解的計算精度和收斂效率。

最后簡單介紹一下混合法，它是直接法和間接法的結合，用直接法離散協(xié)態(tài)變量，通過間接法推導最優(yōu)控制律，最后利用非線性規(guī)劃算法對問題進行求解?；旌戏朔酥苯臃ê烷g接法的一些弊端，比如解的最優(yōu)性難以保證，協(xié)態(tài)初值猜測困難等，但是因為計算參數較多，求解效率仍不夠高。Kluever和Gao等[68-71]提出并改進了混合法用于求解行星際間轉移的軌跡優(yōu)化問題。

2.2 比沖區(qū)間變化的軌跡優(yōu)化方法

對于比沖區(qū)間變化的模型下的軌跡優(yōu)化問題，在所查找到的文獻中基本都應用間接法來處理[12-15,17-20,72]。應用間接法在處理區(qū)間變比沖與恒定比沖軌跡優(yōu)化問題的不同之處在于，在此將比沖作為優(yōu)化過程中的控制量，因此需要用龐得里亞金極大值原理推導關于比沖的最優(yōu)控制律。在使用間接法的基礎上輔助有其他策略，可以提高求解效率。比如對于VASIMR恒定功率變比沖模型，Seywald等[12]在求解共面圓軌道之間的最優(yōu)轉移軌跡時，先假設飛行器的最優(yōu)推力方向始終與速度方向同向或反向；Kechichian等[14]在春分點軌道根數下推導了采用間接法的最優(yōu)控制律；Casalino等[19]將由經典軌道根數表示的Edelbaum方法擴展到變比沖最優(yōu)軌跡的求解中。對于太陽能電推進下的變功率變比沖的情況，Mengali等[20]和陳楊[21]在推導最優(yōu)控制律的判斷條件上考慮了推進效率的影響，推進效率是比沖的函數；Chi等[15-16]結合了間接法、同倫法和協(xié)態(tài)歸一化方法，提高了對于TPBVP和MPBVP的求解效率，并且他們改進了對數同倫函數以適用于比沖區(qū)間變化的情況。對于這種比沖區(qū)間變化的情況，將比沖作為控制量，所以軌跡優(yōu)化過程中至少有三個控制量存在，Taheri等[17,22,72]提出了適用于多控制量優(yōu)化的復合平滑控制技術，它是以雙曲正切函數近似的平滑思想[66]為基礎，應用該策略可以輔助間接法更有效地求解多控制量以及多模式情況下的優(yōu)化問題。

對于多次引力輔助軌跡優(yōu)化問題，Chi等[16]用間接法推導了多次引力輔助的一階必要條件，通過隨機協(xié)態(tài)初值猜測進行軌跡優(yōu)化。然而對于三次及以上引力輔助問題，隨機猜測難以收斂，因此采用單段拼接的方法[73]獲得全局可行解。并且他們在單次引力輔助最優(yōu)軌跡的計算中，利用前文提到的基于偽譜法的GPOPS軟件包[52-53]進行了算例對比。Chi等[74]之后采用協(xié)態(tài)變量轉換的方法[65]重新推導了一階必要條件，因此單段的最優(yōu)解可以為多段協(xié)態(tài)初值提供參考，有效地提高了收斂效率。

2.3 比沖隨距離變化的軌跡優(yōu)化方法

直接法和間接法均曾用于求解比沖隨距離變化的發(fā)動機模型下的軌跡優(yōu)化問題。其中的一些研究是應用直接法軟件MALTO(Mission Analysis Low Thrust Optimization)[75]來計算變比沖最優(yōu)軌跡[28,34,76]，MALTO是根據Sims等[77]提出的將軌跡按控制點劃分成多段，每一段都假設在中途位置施加一個脈沖的離散方法而創(chuàng)建的，之后再通過SNOPT軟件包求解此離散變量的非線性規(guī)劃問題，此外，MALTO非常適用于求解小推力多目標交會以及多次引力輔助問題，控制點可以根據交會天體和甩擺天體來確定。Conway等[27]提出了用解析梯度代替直接配點法中采用的有限差分近似，在求解連續(xù)小推力的軌跡優(yōu)化問題中，無論在收斂率還是求解效率方面，解析梯度法更具有優(yōu)勢。Saripalli等[30]應用Englander等[32,78]開發(fā)的EMTG(Evolutionary Mission Trajectory Generator)軟件重新計算NEXT發(fā)動機多種擬合曲線下的Dawn任務和NEARER任務軌跡，EMTG是基于直接法分為外環(huán)和內環(huán)優(yōu)化程序集合，外環(huán)采用非支配排序遺傳算法離散變量，內環(huán)通過單調域跳躍(MBH)和非線性規(guī)劃算法得到收斂解。Petropoulos等[25-26]提出了用基于形狀的標稱軌跡算法篩選出合適的軌跡，為直接法的精確求解提供初值。Karthik等[37]在對電推進發(fā)動機精確建模的基礎上，用MATLAB的工具包DIRETTO對參量進行離散，再用工具包IPOPT求解此非線性規(guī)劃問題。

Casalino等[35]和Mengali等[33]應用間接法求解比沖隨距離變化模型下的燃料最優(yōu)軌跡。Williams等[27]采用混合法，即直接的遺傳算法優(yōu)化為間接法求解提供初值參考，在此應用軟件包SETOP(Solar Electric Propulsion Trajectory Optimization Program)進行求解，SETOP是間接法軟件VARITOP[79]的改進版。Taheri等[80]將他提出的復合平滑控制技術[22]用于求解光柵離子發(fā)動機模型下的軌跡優(yōu)化問題。另外，孫沖等[36]在虛擬引力場下通過平均算法推導了軌道設計的解析解，該解析方法可以為直接法或間接法提供初值。

對于復雜的多次引力輔助軌跡優(yōu)化問題, 目前的很多研究應用形狀曲線近似比沖隨距離變化的連續(xù)小推力多次引力輔助軌跡的方法。Petropoulos等[26,81]應用基于形狀的指數正弦曲線逼近最優(yōu)軌跡，該方法不僅可以對軌跡進行快速的全局搜索，也可以為直接法或間接法軌跡優(yōu)化提供初值。普渡大學開發(fā)的一款軌跡優(yōu)化軟件GALLOP(Gravity-Assist Low-thrust Local Optimization Program)[82]就是根據基于形狀的全局搜索和局部直接法優(yōu)化相結合的思想而編制的。

2.4 功率分檔軌跡優(yōu)化方法

最后對于功率分檔發(fā)動機模型，很難找到相關的應用直接法的研究。Quarta等[41-42]在NEXT發(fā)動機41個工況離散點模型下采用間接法求解時間最優(yōu)軌跡，41個工況點對應各自的功率、比沖和推力值，因此軌跡優(yōu)化的控制量為41個工況點以及推力方向，在滿足功率約束的前提下，對哈密頓函數按照離散點進行排序并尋找到使哈密頓函數最大的工況點。Chi等[38]應用間接法結合同倫方法以及功率檔位檢測方法優(yōu)化行星際轉移的時間最優(yōu)和燃料最優(yōu)軌跡，這套方法既適用于單個發(fā)動機具有多個工況檔位的情況，也適合于多個恒定比沖恒定推力發(fā)動機帶有功率約束的情況。Axelrod等[44]針對離散推力模型提出了基于間接法求解的閉環(huán)控制技術，并分別討論了應用離散控制律和次優(yōu)控制律的求解效果。

3 結 論

本文對常用的變比沖電推進簡化模型和變比沖連續(xù)小推力軌跡優(yōu)化方法進行了簡介和研究現狀綜述。變比沖電推進簡化模型有比沖在區(qū)間內變化、比沖隨距離變化和功率分檔變化三種，國內外的很多學者針對此三種模型進行了軌跡優(yōu)化研究。常見的軌跡優(yōu)化方法分為直接法和間接法，直接法收斂效率相對較高，但是求得的解無法保證局部最優(yōu)性；應用間接法可以滿足局部最優(yōu)的一階必要條件，但是引入協(xié)態(tài)變量會導致收斂效率較低。由于方法的不足和模型的特殊性，很多研究針對具體問題提出了相應的處理辦法。例如，對于比沖在區(qū)間內變化的模型，采用同倫法以克服初值猜測的敏感性，但由于考慮比沖作為控制量，所以需要對常用的同倫函數進行改進；對于功率分檔發(fā)動機模型，提出功率檔位檢測技術來解決被積項離散不連續(xù)的問題，從而可以提高積分精度和收斂效率。此外，對于任務設計中的多次引力輔助問題，采用相應的直接法程序包進行求解，或在應用間接法的過程中增加提高計算效率的策略。

目前沒有針對變比沖電推進模型的初始軌道設計方法，所以仍采用適用于連續(xù)小推力的脈沖搜索和形狀曲線近似等方法。并且，對于考慮更為實際的變比沖發(fā)動機模型和更為復雜的多目標探測任務，連續(xù)小推力軌跡優(yōu)化問題的求解仍是一個難題，高效的求解方法和策略仍是需要研究的目標。