運載火箭推力故障下軌跡在線重規劃方法*

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

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

0 引 言

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

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

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

1 數學模型描述

1.1 坐標系定義

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

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

1.2 動力學模型

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

(1)

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

1.3 發動機推力故障模型

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

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

(2)

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

2 凸化算法

2.1 約束條件

2.1.1 初始狀態約束

假設發動機出現故障時刻為t0，后續軌跡優化初始狀態即為t0時刻火箭運動狀態，X為火箭狀態變量，下標0表示優化起始時刻，起始點等式約束即為

X0=[r0v0m0]

(3)

2.1.2 終端約束

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

(1)軌道面約束

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

rpfz=0
vpfz=0

(4)

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

近焦點坐標系和發射慣性坐標系間的轉換矩陣為

(5)

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

(2)終端高度約束

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

(6)

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

2.1.3 過程約束

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

(1)推力約束

該液體火箭發動機不具備推力調節能力，其推力固定，因此視為等式約束。

(7)

(2)質量約束

為防止二級飛行段結束時出現計算錯誤，保證仿真符合物理規律，運載火箭在二級飛行時質量必須大于其結構質量，具體形式如下

mdry≤m(t)

(8)

其中，mdry為火箭結構質量。

2.1.4 最優化性能指標

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

(9)

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

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

(10)

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

2.2 最優控制問題的建立與求解

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

建立如下最優控制問題Problem 1：

(11)

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

2.2.1 非凸項凸化處理

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

(12)

其中，上標0表示標準軌跡下各個相關參數。

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

(13)

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

gk=f(r)

(14)

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

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

(15)

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

(16)

2.2.2 最優控制模型

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

(17)

2.2.3 求解流程

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

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

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

(18)

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

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

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

5)應用原對偶內點法，求解問題Problem 2，得到最優解X*=[r*v*m*]，U*；

3 仿真驗證與結果分析

3.1 仿真條件設置

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

圖2 求解流程Fig.2 Solving process

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

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

表1 初始任務參數

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

表2 目標軌道參數

發射場選為海南文昌衛星發射中心，其相關參數如表3所示。

表3 發射場參數

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

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

3.2 軌跡重規劃仿真驗證

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

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

(a)高度變化曲線

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

表4 終端誤差對比

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

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

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

4 結 論

本文研究了基于能量最優的運載火箭非入軌段在發生推力故障后的軌跡重規劃問題。通過無損凸化、線性化等方法，將最優控制問題轉化為二階錐規劃問題后，應用原對偶內點法精確求解。仿真結果表明，該方法能夠滿足終端速度和高度的要求，能夠完成運載火箭在該類故障下軌跡重規劃，計算效率高，能夠實現軌跡的快速收斂，提高了運載火箭適應故障的能力，具有在線應用的潛力。