基于高斯偽譜法的火星表面上升燃耗最優(yōu)軌跡設(shè)計

柯森锎，李爽，肖東東，王衛(wèi)華，聶欽博

（1. 南京航空航天大學(xué) 航天學(xué)院，南京 210016；2. 南京航空航天大學(xué) 航天新技術(shù)實驗室，南京 210016；3. 上海航天控制技術(shù)研究所，上海 201109）

0 引 言

火星是離地球最近的類地行星，對火星的探測有助于人類進(jìn)一步理解認(rèn)識太陽系的起源和演變，將宇宙科學(xué)推上新臺階[1]。在人類成功實施火星著陸探測之后，采樣返回將是火星探測的下一個里程碑[2-3]?；鹦巧仙鳎∕ars Ascent Vehicle，MAV）發(fā)射總質(zhì)量是關(guān)系到火星采樣返回任務(wù)可行性的關(guān)鍵因素，不僅關(guān)系到所能攜帶的火星樣本的質(zhì)量，而且直接影響火星表面上升軌跡設(shè)計。同時，它也關(guān)系到火星發(fā)射支持系統(tǒng)的質(zhì)量，進(jìn)一步影響到火星登陸器的質(zhì)量，最后傳導(dǎo)到地球上火箭的運載能力要求上。MAV每增加1 kg，則地球發(fā)射質(zhì)量將相應(yīng)增加約8 kg[4]。而在當(dāng)前技術(shù)條件下每將1 kg的有效設(shè)備送上火星需要花費約100萬美元的成本[5]。所以MAV的設(shè)計不僅需要考慮火箭的運載能力、火星進(jìn)入下降著陸（Entry，Descent and Landing，EDL）能力的可行性約束[6]，還面臨著經(jīng)濟(jì)性的考慮。

為了實現(xiàn)最輕的火星上升器設(shè)計，必須對MAV分級參數(shù)和上升軌跡進(jìn)行優(yōu)化。上升器參數(shù)包括各級火箭的干重、燃料質(zhì)量、推力、比沖等。軌跡優(yōu)化的經(jīng)典應(yīng)用領(lǐng)域為火箭的最優(yōu)控制問題[7-10]，在該問題中給定火箭參數(shù)，要求以燃耗最優(yōu)為目標(biāo)函數(shù)進(jìn)行軌跡優(yōu)化以得到最優(yōu)的控制量。Betts對軌跡優(yōu)化問題的求解方法進(jìn)行了概括分析，將其分為解析法和數(shù)值法，其中數(shù)值法又可以分為直接法和間接法[11]。高斯偽譜法是直接法的一種，Patterson提出了一種自適應(yīng)的高斯偽譜優(yōu)化方法并開發(fā)了相應(yīng)的Matlab軟件GPOPS-II[12]。國內(nèi)的楊希祥基于高斯偽譜法做了固體火箭快速軌跡優(yōu)化研究[13]。如果火箭分級參數(shù)設(shè)計得不合理，那么軌跡優(yōu)化問題就無法收斂得到理想的優(yōu)化結(jié)果，針對這一問題Coskun提出了一種基于高斯偽譜法的優(yōu)化方法，能夠?qū)⒒鸺旨墐?yōu)化和軌跡優(yōu)化耦合[14]。Benito也將高斯偽譜法引入到MAV的上升軌跡優(yōu)化問題中來[15]。

本文以火星上升器總發(fā)射質(zhì)量最小為目標(biāo)函數(shù)，基于高斯偽譜法設(shè)計了一種針對兩級MAV的分級與軌跡耦合多階段優(yōu)化算法，最終得到MAV總發(fā)射質(zhì)量最小的火星上升器分級參數(shù)以及一條燃耗最優(yōu)上升軌跡。

1 問題描述

1.1 飛行方案

在火星表面發(fā)射入軌所需的速度增量與地球相比大大減小，適合采用兩級MAV設(shè)計。本論文中上升器從火星表面起飛入軌的飛行方案如圖1所示，第1級點火開始是一個很短的大推力的動力上升段；然后是長時間的無動力滑行段一直到上升軌道的遠(yuǎn)拱點處；進(jìn)一步在遠(yuǎn)拱點處第2級火箭點火進(jìn)行軌道圓化進(jìn)入目標(biāo)圓軌道，然后釋放火星樣本艙，最終由火星軌道器實現(xiàn)對火星樣本艙的交匯捕獲以及返回地球[16-18]。

圖1 火星上升器飛行方案Fig. 1 MAV fly profile

1.2 坐標(biāo)系定義

1.3 兩級MAV質(zhì)量模型

一個一般化的兩級MAV的質(zhì)量模型如圖2所示。每一級的質(zhì)量由干重和燃料質(zhì)量組成[19]。而每級干重又由固定質(zhì)量和與燃料質(zhì)量相關(guān)的可變質(zhì)量組成。

圖2 一般化的兩級MAV質(zhì)量模型Fig. 2 The universal model for 2-stage MAVs

燃料質(zhì)量由每級的干重、載荷質(zhì)量、速度增量以及燃料比沖決定，燃料質(zhì)量為

第1級的載荷包括第2級總質(zhì)量加上有效載荷的質(zhì)量。

1.4 動力學(xué)模型

假設(shè)火星為旋轉(zhuǎn)的標(biāo)準(zhǔn)圓球，只考慮空氣阻力作用，在火心慣性直角坐標(biāo)系下建立MAV上升段的動力學(xué)方程為

2 兩級MAV的分級與軌跡耦合多階段優(yōu)化算法

將兩級MAV的總發(fā)射質(zhì)量最小化問題轉(zhuǎn)化為兩級MAV的分級與軌跡耦合多階段優(yōu)化問題。根據(jù)兩級MAV的飛行方案將兩級MAV的分級與軌跡耦合優(yōu)化問題分為3個階段：動力段、滑行段、圓化段。

2.1 目標(biāo)函數(shù)與優(yōu)化變量

對于MAV的總發(fā)射質(zhì)量最小化問題，其目標(biāo)是使得總發(fā)射質(zhì)量最小，等效于最小化第一階段的初始質(zhì)量，其目標(biāo)函數(shù)為

2.2 動力學(xué)約束

動力學(xué)約束由第1節(jié)中的動力學(xué)模型給出。由于MAV總發(fā)射質(zhì)量最小化問題需要對各級燃料質(zhì)量進(jìn)行優(yōu)化，所以增加了消耗的燃料質(zhì)量作為狀態(tài)量，需要增加式（10）所示的微分方程作為動力學(xué)約束。

其中：p為消耗的燃料質(zhì)量。

2.3 邊界約束

本節(jié)將給出MAV總發(fā)射質(zhì)量最小化問題各階段的邊界約束條件。

1）第1級推力段

2）無動力滑行段

無動力滑行段的初始條件由第1級推力段的終端狀態(tài)給出，即由2個階段的連接點約束決定。

為了增加MAV進(jìn)入火星穩(wěn)定軌道的概率，MAV將一直滑行到滑行軌道的遠(yuǎn)拱點再進(jìn)行軌道圓化，即該階段的終端約束為徑向速度為零，即

3）軌道圓化

該初始條件由滑行段的終端狀態(tài)給出，由兩個階段的連接點約束決定。在此階段的開端MAV進(jìn)行一、二級火箭分離，拋棄第1級的同時第2級火箭點火進(jìn)行軌道圓化，所以第2級軌道圓化段與無動力滑行段間的連接點約束要考慮第1級干重質(zhì)量。連接點約束為

將第2級軌道圓化段的初始燃料消耗質(zhì)量設(shè)置為零，即

該階段的終端條件為進(jìn)入目標(biāo)軌道，即滿足給定的軌道要素，與軌跡優(yōu)化問題類似，但是需要增加終端質(zhì)量約束為

升交點赤經(jīng)可以通過改變發(fā)射時間來調(diào)整，因此本文中不考慮升交點赤經(jīng)約束，目標(biāo)軌道為圓軌道所以為目標(biāo)軌道高度。雖然圓軌道是數(shù)學(xué)上最簡單的軌道，但是在現(xiàn)實中由于各種擾動很難實現(xiàn)精確的圓軌道，如果將離心率設(shè)為零將會帶來一些數(shù)值計算上的困難，比如偽震蕩和不收斂。對于引入誤差范圍如下式所示。

2.4 路徑約束

1）第1級推力段

飛行高度約束

動壓約束

推力角約束如式（31），推力角為推力方向與相對火星速度的夾角，其計算方程如式（32）所示，由于使用3DoF模型無法表示攻角，用推力角近似攻角。

2）無動力滑行段

該階段的飛行高度約束以及動壓約束與第1級推力段相同，因為整個過程沒有推力，所以并沒有推力方向和推力角約束。

3）軌道圓化段

該階段在大氣層外所以無需考慮動壓和推力角約束，但推力方向和飛行高度約束與第1級推力段相同。

3 高斯偽譜法

高斯偽譜法屬于直接法中的配點法，在處理分級與軌跡耦合優(yōu)化問題時，通過在一系列拉格朗日高斯（Legendre Gauss，LG）點上構(gòu)建拉格朗日多項式來逼近狀態(tài)變量和控制變量[13]，進(jìn)一步對得到的多項式求導(dǎo)來近似狀態(tài)量的導(dǎo)數(shù)，將微分方程轉(zhuǎn)換為代數(shù)方程。對于終端條件和目標(biāo)函數(shù)中的積分項使用高斯積分來替換，將積分方程也轉(zhuǎn)換為代數(shù)方程，把分級與軌跡耦合優(yōu)化問題離散后轉(zhuǎn)換為非線性規(guī)劃問題[進(jìn)行]求解。

對（34）求導(dǎo)可以得到狀態(tài)量關(guān)于時間的導(dǎo)數(shù)為

基于式（36）可以將LG點處的動力學(xué)方程約束轉(zhuǎn)換為代數(shù)方程形式為

終端狀態(tài)也需要滿足動力學(xué)方程約束，將動力學(xué)方程從初始時刻積分至終端時刻可以得到終端狀態(tài)量，使用高斯積分近似積分項，可以得到終端狀態(tài)量為

加上各離散點處的邊界條件和路徑約束就可以將兩級MAV分級與軌跡耦合優(yōu)化問題離散化為非線性規(guī)劃問題，求得滿足約束條件的狀態(tài)量和控制量使得目標(biāo)函數(shù)最小。

4 仿真實驗與分析

對提出的基于高斯偽譜法的兩級MAV的分級與軌跡耦合多階段優(yōu)化算法仿真，目標(biāo)軌道參數(shù)、發(fā)射點位置以及兩級MAV質(zhì)量模型參數(shù)的詳細(xì)輸入值都由表1給出[19]。仿真程序在Matlab環(huán)境下運行，計算機CPU為2.67 GHz Core i5，仿真計算得到兩級MAV的分級優(yōu)化結(jié)果如表2所示，同時得到一條標(biāo)稱上升軌跡如圖3～9所示。

表1 兩級MAV分級與軌跡耦合優(yōu)化算法仿真輸入?yún)?shù)Table 1 Input parameters of two-stage MAV coupled stagingtrajectory optimization algorithm simulation

表2 兩級MAV分級優(yōu)化結(jié)果Table 2 The result of two-stage MAV staging optimization

圖3 上升高度隨時間變化曲線Fig. 3 Altitude vs. time

圖4 飛行速度隨時間變化曲線Fig. 4 Speed vs. time

圖5 軌道傾角隨時間變化曲線Fig. 5 Inclination vs. time

由圖3～5可知滑行段時間約為565 s，燃耗最優(yōu)的軌跡并不是在第1級耗盡時達(dá)到目標(biāo)軌道傾角，而是留下約2.5°的軌道傾角誤差留給第2級進(jìn)行調(diào)整；由圖6～7第1級控制量可得發(fā)射方位角為42.5°，仰角為40.5°，第2級進(jìn)行軌道圓化的推力方向接近于固定方向；由于MAV的氣動阻力系數(shù)Cd在這里被設(shè)定為簡單的恒定值0.2[21]，所以得到的動壓值只能作為一個概略的定性值，由圖8知由于燃燒時間短，第1級耗盡點高度較低，在第1級耗盡點處達(dá)到了動壓的峰值，可以考慮增加燃燒時間以提高耗盡點高度，以降低耗盡點的動壓使其為在峰值之后的一個較小值，可以減小滑行段姿態(tài)保持的難度；圖9中第1級的推力角起始值約為60°，之后快速減小并保持在5°以下，這是因為推力角是推力與MAV相對火星速度的夾角，而在起始時刻MAV相對火星速度為零，這里的起始推力角是由數(shù)值計算誤差得到的。

圖6 第一階段的控制量隨時間變化曲線Fig. 6 Control vs. time for phase 1

圖7 第三階段的控制量隨時間變化曲線Fig. 7 Control vs. time for phase 2

圖8 動壓隨時間變化曲線Fig. 8 Dynamic pressure vs. time

圖9 第1級推力角隨時間變化曲線Fig. 9 Thrust angle vs. time for phase 1

通過仿真計算可知，該優(yōu)化算法的初值在給定范圍內(nèi)隨機選取都能在約7 min內(nèi)收斂得到最優(yōu)結(jié)果，計算精度達(dá)到10-6量級。

5 結(jié) 論

以火星采樣返回任務(wù)中的兩級火星上升器為研究對象，利用高斯偽譜法設(shè)計了在復(fù)雜多約束條件下兩級MAV分級與軌跡耦合多階段優(yōu)化算法，以求得總發(fā)射質(zhì)量最小的兩級MAV分級設(shè)計以及燃耗最優(yōu)的上升軌跡。該算法將MAV的分級優(yōu)化與軌跡優(yōu)化耦合在一起，解決了由于不合理的分級參數(shù)導(dǎo)致軌跡優(yōu)化無法收斂的問題。通過仿真驗證可以得到如下結(jié)論。

1）為達(dá)到燃耗最優(yōu)的目的，兩級MAV的第1級將留下一定的軌道傾角偏轉(zhuǎn)量由第2級完成。

2）為減小兩級MAV在滑行段的姿控難度可以適當(dāng)減小第1級火箭的推力進(jìn)而減小第1級耗盡點的動壓。

3）提出的優(yōu)化算法對于初值選擇的敏感度較小，具有較強的魯棒性，且收斂速度較快，方便于進(jìn)行快速的兩級MAV設(shè)計優(yōu)化。

本文僅考慮了總發(fā)射質(zhì)量最小和燃耗最優(yōu)作為性能指標(biāo)的火星表面上升段軌跡設(shè)計問題，而在實際應(yīng)用中還需要考慮飛行時間、火星軌道器觀測跟蹤范圍等其他因素，未來需要針對MAV上升軌跡的多性能指標(biāo)優(yōu)化開展進(jìn)一步的研究。