高天域，何睿智，溫長新，湯國建

(1.國防科技大學 空天科學學院空天工程系·長沙·410073；2.中國人民解放軍95899部隊·北京·100000)

0 引 言

近年來，隨著各國運載火箭技術(shù)不斷提升，運載能力不斷增強，對外層空間的開發(fā)與利用也在不斷深化，運載火箭需要更強的自主控制能力從而妥善處理突發(fā)故障[1]。發(fā)動機是運載火箭研制過程中的重中之重，其可靠性關(guān)乎整個飛行任務(wù)的成敗，然而新技術(shù)、新材料、新理論的應(yīng)用不可避免地造成其可靠性降低，據(jù)統(tǒng)計，運載火箭主動段動力系統(tǒng)故障占所發(fā)故障的60%[2]。例如日本H-2運載火箭由于一級發(fā)動機提前關(guān)機，升空后嚴重偏離預(yù)定軌道，導致發(fā)射任務(wù)失敗，造成巨大經(jīng)濟損失；印度GSLV-F06運載火箭發(fā)射后，2臺游動發(fā)動機點火失敗，導致任務(wù)失敗；同樣，我國也曾發(fā)生此類事故，長征五號遙二火箭因芯一級2臺發(fā)動機其中一臺故障導致發(fā)射任務(wù)失敗，對我國后續(xù)航天任務(wù)造成影響。隨著自適應(yīng)制導和軌跡在線規(guī)劃技術(shù)的不斷成熟，發(fā)動機非致命故障的情況下，仍能通過調(diào)整火箭飛行程序完成任務(wù)。SpaceX公司發(fā)射的Falcon 9運載火箭曾在一臺發(fā)動機因故障提前關(guān)機的情況下，通過增加一二級飛行時間，將主要載荷送入預(yù)定軌道，次要載荷進入較低軌道，基本完成任務(wù)；美國德爾塔4運載火箭曾利用制導控制系統(tǒng)重新生成故障后的飛行軌跡，對推力下降及時進行了補償，充分利用剩余燃料完成入軌。

針對火箭故障條件下的軌跡在線重規(guī)劃和制導問題，國內(nèi)外開展了一定的研究。韓業(yè)鵬[3]針對大氣層外動力故障，研究了入軌點更新迭代制導方法，提高了控制系統(tǒng)適應(yīng)故障的能力。王志祥等[4]以某型號捆綁火箭為研究對象，分析了單臺發(fā)動機推力下降對火箭飛行軌跡和姿態(tài)的影響，并利用Gauss偽譜法展開仿真研究。宋征宇等[5]提出了一種在線軌跡規(guī)劃策略，判斷火箭不可達時自主尋找最優(yōu)救援軌道。李源[6]針對二級運載火箭在不同飛行段發(fā)生不同程度的故障的軌跡規(guī)劃方法以及制導方法進行了研究，為運載火箭軌跡在線規(guī)劃和自適應(yīng)制導技術(shù)應(yīng)用提供了一定參考。H.Lee等[7]通過牛頓迭代法提升了火箭應(yīng)對機構(gòu)故障的能力。H.Sun[8]提出了一種基于有限微分的能夠適應(yīng)發(fā)動機推力故障的上升段制導方法。B.Beneditkter等[9]采用無損凸化技術(shù)和連續(xù)凸化技術(shù)相結(jié)合的方式，并在優(yōu)化過程中明確考慮自由終端時間的要求，提出了一種多級運載火箭主動段凸規(guī)劃方法。G.A.Dukema[10]和A.J.Calise等[11]提出了一種基于預(yù)測-校正的初值生成方法，能夠處理故障后返回地面的問題。除此之外，離線設(shè)計應(yīng)急軌道的方法也是提高運載火箭應(yīng)對故障能力的一種選擇，“挑戰(zhàn)者”號曾應(yīng)用離線設(shè)計好的應(yīng)急軌道，成功在一臺發(fā)動機故障后進入較低的安全軌道。A.L.Cowling[12]研究了不同時刻發(fā)動機停止工作后的中止軌跡。隨著箭載計算機的計算能力不斷提高，未來基于自主制導方法的新型閉環(huán)控制系統(tǒng)架構(gòu)以及多級分層軌跡優(yōu)化將成為主要的研究方向[13-15]。

針對運載火箭主動段非入軌飛行段出現(xiàn)推力下降的問題，本文提出了一種基于能量最優(yōu)的凸優(yōu)化軌跡規(guī)劃方法，將原非凸問題通過線性化估算和常值假設(shè)轉(zhuǎn)為凸優(yōu)化問題，并通過迭代得到原問題的最優(yōu)解，實現(xiàn)該問題下軌跡快速求解。

1 數(shù)學模型描述

1.1 坐標系定義

首先定義以下2個坐標系：1)發(fā)射慣性坐標系O-xyz，坐標原點為發(fā)射點，Ox軸在火箭發(fā)射點水平面內(nèi)指向發(fā)射瞄準方向，Oy軸垂直于發(fā)射點水平指向上方，Oz軸滿足右手坐標系；2)近焦點坐標系Op-xpypzp，坐標原點為地心，xp軸為地心和目標軌道近地點的連線，指向近地點，yp為xp在目標軌道平面內(nèi)沿軌道方向旋轉(zhuǎn)90°，zp滿足右手坐標系。近焦點坐標系如圖1所示。

圖1 近焦點坐標系Fig.1 Perifocal coordinate system

1.2 動力學模型

本文主要針對三級火箭的二級飛行段進行研究。假設(shè)故障發(fā)生于一二級分離時刻，此時火箭高度已達100km左右，大氣稀薄，可忽略氣動力的影響。假設(shè)地球為橢球，此時在發(fā)射慣性坐標系下建立運載火箭運動方程

(1)

式中，向量(r,v,a)分別為火箭在慣性系中的位置、速度和加速度；T為火箭在慣性系中3個方向的推力分量；Treal為運載火箭實際推力。一般情況下，火箭推力為常值，m為火箭質(zhì)量，隨時間變化，Isp為發(fā)動機比沖，g為引力加速度，g0為海平面重力加速度。

1.3 發(fā)動機推力故障模型

本文考慮的發(fā)動機推力下降模式為：由于秒耗量下降發(fā)動機推力迅速下降至某一固定值，且不考慮二次故障。表征故障的主要參數(shù)有故障發(fā)生時刻和發(fā)動機故障后輸出推力的大小，以某型號芯二級4臺發(fā)動機為例，各發(fā)動機發(fā)生推力下降后的數(shù)學模型為

Ti=kT(i=1,2,3,4)

(2)

其中，k為推力下降系數(shù)；Ti為故障后的推力；T為標準推力。

2 凸化算法

2.1 約束條件

2.1.1 初始狀態(tài)約束

假設(shè)發(fā)動機出現(xiàn)故障時刻為t0，后續(xù)軌跡優(yōu)化初始狀態(tài)即為t0時刻火箭運動狀態(tài)，X為火箭狀態(tài)變量，下標0表示優(yōu)化起始時刻，起始點等式約束即為

X0=[r0v0m0]

(3)

2.1.2 終端約束

由于研究對象為三級運載火箭，二級飛行段的主要任務(wù)為提升高度和速度，對終端約束并沒有嚴格的要求，為了保證后續(xù)飛行段的穩(wěn)定，增加兩個終端約束。

(1)軌道面約束

為防止后續(xù)飛行出現(xiàn)偏離目標軌道的情況，需要將運載火箭約束在目標軌道平面內(nèi)飛行。由于在發(fā)射坐標系下軌道平面相關(guān)參數(shù)計算復雜，因此將終端軌道面約束轉(zhuǎn)到近焦點坐標系下，可以有效提高計算效率，軌道面約束表現(xiàn)形式如下

rpfz=0
vpfz=0

(4)

式中，vpfz、rpfz分別為近焦點坐標系下終端速度和高度。

近焦點坐標系和發(fā)射慣性坐標系間的轉(zhuǎn)換矩陣為

(5)

式中，Ω為升交點經(jīng)度；ω為近地點幅角；i為軌道傾角；λ0為發(fā)射點經(jīng)度；B0為地理緯度；A0為射擊方位角。

(2)終端高度約束

由于發(fā)動機秒耗量下降，火箭在燃料冗余的情況下可以適當延長二級飛行時間，火箭在三級一次飛行段過后會進入一個過渡軌道，飛行時間的延長可能導致飛行高度超過或者接近這一過渡軌道，不利于后續(xù)三級一次飛行段的軌跡規(guī)劃和制導，因此在二級飛行段軌跡規(guī)劃時增加飛行高度約束。以標準軌跡二級飛行結(jié)束時刻地心距作為約束邊界條件。具體表現(xiàn)形式如下

(6)

其中，rf為故障條件下二級飛行結(jié)束時火箭的地心距；Rmax為標準軌跡二級飛行結(jié)束時火箭的地心距。

2.1.3 過程約束

當出現(xiàn)推力下降情況時，由于火箭推力下降，同一時刻的火箭速度、加速度和標準彈道相比也出現(xiàn)了下降的情況，且二級飛行段已飛出了稠密大氣層，大氣密度可以忽略不計。因此，與這三項相關(guān)的動壓、熱流、過載、彎矩等約束在同一時刻均能滿足約束要求，在優(yōu)化過程中不予考慮，主要考慮飛行過程中的推力約束與質(zhì)量約束。

(1)推力約束

該液體火箭發(fā)動機不具備推力調(diào)節(jié)能力，其推力固定，因此視為等式約束。

(7)

(2)質(zhì)量約束

為防止二級飛行段結(jié)束時出現(xiàn)計算錯誤，保證仿真符合物理規(guī)律，運載火箭在二級飛行時質(zhì)量必須大于其結(jié)構(gòu)質(zhì)量，具體形式如下

mdry≤m(t)

(8)

其中，mdry為火箭結(jié)構(gòu)質(zhì)量。

2.1.4 最優(yōu)化性能指標

本文所設(shè)計的最優(yōu)化性能指標為終端能量最大，即飛行過程中的能量損耗最小，具體指標形式如下

(9)

其中，μ為地球引力常數(shù)；vf為終端時刻的速度。

由于終端時刻在近焦點坐標系下z方向的速度、位置分量均為0，為了簡便計算，提高計算效率，將優(yōu)化命題在近焦點坐標系下表示，即為

(10)

其中，rpxf、rpyf、vpxf、vpyf分別為終端時刻運載火箭在近焦點坐標系下位置和速度在x軸和y軸的分量。

2.2 最優(yōu)控制問題的建立與求解

根據(jù)上述的性能指標及初始狀態(tài)、終端和過程約束，以及飛行過程中的動力學約束，同時定義狀態(tài)變量X=[rvm]和控制變量U=[TxTyTzTreal]，其中Tx、Ty、Tz分別為推力在近焦點坐標系下3個方向的分量。

建立如下最優(yōu)控制問題Problem 1：

(11)

式中的非凸項主要來源有：1)性能指標為非凸函數(shù)；2)動力學方程中的引力相關(guān)項為非線性。

2.2.1 非凸項凸化處理

由于性能指標函數(shù)為非凸，因此通過線性化手段進行凸化處理，將性能指標進行一階泰勒展開并做歸一化處理，省略掉無意義的常數(shù)項和高階項后，性能指標轉(zhuǎn)化為如下形式

(12)

其中，上標0表示標準軌跡下各個相關(guān)參數(shù)。

除此之外將控制量等式約束松弛為二階錐約束。即為

(13)

對于引力相關(guān)項，由于引力所引起的加速度與推力加速度相比為小量，同時軌跡規(guī)劃的高度變化一般在幾十千米以內(nèi)，遠小于地心距，所以對軌跡規(guī)劃的影響程度并不大，因此可以基于初始軌跡進行引力加速度的快速計算。根據(jù)初始軌跡將引力相關(guān)加速度項轉(zhuǎn)化為與地心距相關(guān)的線性函數(shù)，在軌跡規(guī)劃過程中根據(jù)各離散點運載火箭地心距插值得到引力加速度。

gk=f(r)

(14)

其中，gk表示在離散點處的引力加速度；f(·)為引力加速度和地心距之間的關(guān)系。

由于推進劑秒耗量為常值，因此各離散點的質(zhì)量mk可以直接計算

(15)

針對動力學方程，采用梯形離散法進行動力學方程的凸化，將飛行時間[t0tf]平均分為N段，離散間隔為Δt=(tf-t0)/N。具體表現(xiàn)形式如下

(16)

2.2.2 最優(yōu)控制模型

結(jié)合上述分析，建立如下二級飛行段軌跡規(guī)劃問題Problem 2：

(17)

2.2.3 求解流程

通過求解上述最優(yōu)控制問題，可以得到推力故障下的最優(yōu)軌跡，具體求解流程如圖2所示。

1)給定離散區(qū)間N以及終端精度要求εr、εv；

2)根據(jù)運載火箭故障時間及推力故障系數(shù)計算故障后飛行時間；

(18)

其中，twork為二級工作時長；t0為發(fā)生故障時刻；tf為二級飛行段結(jié)束時刻；k為推力下降系數(shù)。

3)根據(jù)所給參數(shù)積分求解初始軌跡；

4)根據(jù)初始軌跡進行引力項的凸化；

5)應(yīng)用原對偶內(nèi)點法，求解問題Problem 2，得到最優(yōu)解X*=[r*v*m*]，U*；

3 仿真驗證與結(jié)果分析

3.1 仿真條件設(shè)置

為了充分驗證本文推力故障下運載火箭非入軌段軌跡規(guī)劃算法的正確性和收斂性，重點針對運載火箭二級飛行段發(fā)生故障的情況進行仿真驗證。假設(shè)發(fā)生故障時刻的火箭位置速度參數(shù)均可根據(jù)導航系統(tǒng)或地面基站獲得。

圖2 求解流程Fig.2 Solving process

本文使用的仿真實驗軟硬件條件為Inter Corel i5-2300 CPU，3.10GHz，操作系統(tǒng)為Windows 7，求解軌跡規(guī)劃過程中的二階錐規(guī)劃問題。

仿真過程中的無故障下運載火箭參數(shù)設(shè)置如表1所示。

表1 初始任務(wù)參數(shù)

由于研究目標為三級運載火箭，在三級一次飛行段過后運載火箭會進入一個軌道滑行，三級一次飛行段進入的目標軌道參數(shù)如表2所示。

表2 目標軌道參數(shù)

發(fā)射場選為海南文昌衛(wèi)星發(fā)射中心，其相關(guān)參數(shù)如表3所示。

表3 發(fā)射場參數(shù)

運載火箭在整個三級飛行過程中引力和地心距的關(guān)系在發(fā)射坐標系下如圖3所示，在整個優(yōu)化過程中可以根據(jù)圖3的關(guān)系插值得到引力加速度，從而實現(xiàn)引力項的凸化。

圖3 引力變化曲線Fig.3 Gravitational force change curve

3.2 軌跡重規(guī)劃仿真驗證

針對二級點火時刻立即發(fā)生不同程度故障進行模擬，驗證基于能量最優(yōu)的軌跡在線重規(guī)劃算法。由于二次飛行段為運載火箭的非入軌段，其主要任務(wù)為爬高和加速，因此并無嚴格的終端精度要求，只需速度和高度達到要求即可。能量最優(yōu)性能指標可以使運載火箭充分利用所攜帶的燃料，以最小的能量損耗滿足運載火箭終端速度和高度的約束。取離散區(qū)間N=100，εr=10-5，εv=10-5。

針對二級點火時發(fā)動機發(fā)生故障，對推力、秒耗量同時下降10%的情況進行分析，二級飛行段在線軌跡重規(guī)劃仿真實驗結(jié)果如圖4所示。

(a)高度變化曲線

具體終端誤差如表4所示。

表4 終端誤差對比

為驗證推力經(jīng)松弛后其幅值不隨時間改變，推力幅值隨時間的變化曲線如圖5所示。

圖5 推力幅值變化曲線Fig.5 Thrust amplitude variation curve

仿真結(jié)果表明，應(yīng)用故障條件下基于能量最優(yōu)的軌跡在線重規(guī)劃算法得到的軌跡終端速度和高度誤差均控制在0.2%以內(nèi)，具有良好的準確度；且僅需迭代5次即可達到收斂條件，單次迭代時間在0.3s左右，具有良好的魯棒性和穩(wěn)定性。推力幅值經(jīng)松弛后不隨時間改變，驗證了其有效性。

4 結(jié) 論

本文研究了基于能量最優(yōu)的運載火箭非入軌段在發(fā)生推力故障后的軌跡重規(guī)劃問題。通過無損凸化、線性化等方法，將最優(yōu)控制問題轉(zhuǎn)化為二階錐規(guī)劃問題后，應(yīng)用原對偶內(nèi)點法精確求解。仿真結(jié)果表明，該方法能夠滿足終端速度和高度的要求，能夠完成運載火箭在該類故障下軌跡重規(guī)劃，計算效率高，能夠?qū)崿F(xiàn)軌跡的快速收斂，提高了運載火箭適應(yīng)故障的能力，具有在線應(yīng)用的潛力。