吳嘉梁

(北京機(jī)電工程研究所，北京100074)

基于間接法的上升段軌跡優(yōu)化方法研究

吳嘉梁

(北京機(jī)電工程研究所，北京100074)

固體火箭上升段軌跡優(yōu)化設(shè)計具有重要的工程意義。針對此問題，提出了一種求解上升段最優(yōu)軌跡的可行方法。在零側(cè)滑角假設(shè)下構(gòu)建飛行器模型，以推力方向為最優(yōu)控制量，根據(jù)極小值原理推導(dǎo)一階最優(yōu)條件。采用間接法，將真空條件下上升段最優(yōu)軌跡的解作為初值，以狀態(tài)響應(yīng)方程構(gòu)造一種迭代的方法，最后在滿足攻角過程約束下，通過同倫算法獲得真實大氣環(huán)境下的最優(yōu)軌跡。仿真結(jié)果表明，該優(yōu)化算法能夠穩(wěn)定收斂，具有良好的魯棒性。

固體火箭；上升段；最優(yōu)控制；軌跡優(yōu)化；間接法

0 引言

固體火箭由于具有準(zhǔn)備時間短、機(jī)動性強(qiáng)的優(yōu)點，近年來應(yīng)用越來越廣泛，并且隨著其發(fā)射任務(wù)的多樣化和復(fù)雜化，在飛行中面臨著更強(qiáng)的過程約束與終點狀態(tài)約束，這就突出了上升段軌跡優(yōu)化設(shè)計的迫切需求。直接法和間接法是軌跡優(yōu)化最常用的兩種方法。間接法是基于龐德里亞金極小值原理推導(dǎo)最優(yōu)控制的一階必要條件，構(gòu)成了求解最優(yōu)軌跡的Hamiltonian邊值問題(HBVP)，由于不對性能指標(biāo)函數(shù)直接尋優(yōu)，因此稱為間接法；直接法采用參數(shù)化方法將連續(xù)空間的最優(yōu)控制問題求解轉(zhuǎn)化為一個非線性規(guī)劃(NLP)問題，通過數(shù)值求解非線性規(guī)劃問題來獲得最優(yōu)軌跡。

Lu提出了一種采用有限差分法求解最優(yōu)軌跡的方法，并對該方法進(jìn)行了進(jìn)一步的研究和發(fā)展[1-3]，宣穎等采用偽譜法進(jìn)行上升段飛行器軌跡優(yōu)化，能夠獲得正確有效的最優(yōu)軌跡[4-5]，Dukeman采用打靶法來求解軌跡優(yōu)化的兩點邊值變分問題，獲得上升段飛行器最優(yōu)軌跡[6]。綜合以上研究結(jié)果可以看出，直接法在求解優(yōu)化問題時，有著初值敏感度較低且收斂域更廣、不必推導(dǎo)一階最優(yōu)性條件的優(yōu)點，但是也有易收斂到局部最優(yōu)解、求解速度慢的缺點；間接法在求解最優(yōu)問題時，具有解的精度高、最優(yōu)解滿足一階最優(yōu)性必要條件、求解速度較快的優(yōu)點，但是相對于直接法，間接法對初值的選取較為敏感。因此，本文基于間接法，建立飛行器軌跡優(yōu)化模型，采用真空解為初始猜測，并通過狀態(tài)響應(yīng)方程迭代的思想求取最優(yōu)軌跡，通過仿真驗證算法的收斂性。

1 飛行器動力學(xué)模型

鑒于傳統(tǒng)基于地球固連坐標(biāo)系下的三自由度動力學(xué)方程在利用最優(yōu)控制理論推導(dǎo)協(xié)態(tài)方程時十分困難，為了提高算法效率，本文基于地心慣性坐標(biāo)系建立飛行器的三自由度動力學(xué)方程：

(1)

在式(1)中，R∈R3為地心慣性坐標(biāo)系下地心到質(zhì)心的向量，稱為矢徑向量，V∈R3為地心慣性坐標(biāo)系下速度向量，g(r)為重力加速度項，T為推力大小，m(t)為飛行器當(dāng)前時刻的質(zhì)量，Isp為比沖。lb為飛行器體軸方向單位向量，ln為飛行器縱對稱面內(nèi)法向方向單位向量，A和N分別為氣動力在lb與ln上的投影。

在本文上升段制導(dǎo)研究中，認(rèn)為飛行器的主對稱面和由單位矢量lb與飛行器空速向量Vr確定的平面始終保持重合，即側(cè)滑角始終保持為零。根據(jù)上述定義，已知彈體軸矢量lb和空速矢量Vr的條件下，定義彈體y軸空間矢量為ly

(2)

其中，‖lVr×lb‖=|sinα|，lVr=Vr/‖Vr‖是Vr的歸一化矢量，Vr=V-ωE×R-Vw，V為慣性速度矢量，ωE為地球自轉(zhuǎn)速度向量在地心慣性坐標(biāo)系中的投影，Vw為風(fēng)速矢量，由于風(fēng)速是不可知的，故在動力學(xué)模型中設(shè)Vw恒為0。

通過lb、ly和lVr可以得到彈體法向單位方向矢量ln和攻角α為：

ln=-lb×ly

(3)

(4)

通過式(2)、式(3)的定義可以確保飛行器在攻角正負(fù)變化時不會出現(xiàn)180°的滾動角跳變。之后將地心慣性系下的彈體坐標(biāo)系各軸空間矢量投影到地理系下，可以得到三個彈體姿態(tài)角：

lbt=Titlb,lnt=Titln

(5)

?=tan-1(lbty/lbtx),

(6)

其中，lbt、lnt分別代表地理系中l(wèi)b、ln的投影，Tit為地心慣性坐標(biāo)系到地理坐標(biāo)系的轉(zhuǎn)換矩陣，?、ψ、γ分別代表俯仰角、偏航角和滾動角。

(7)

2 最優(yōu)控制問題求解

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

飛行器軌跡優(yōu)化問題本質(zhì)上是一個最優(yōu)控制問題，即在多種約束的限制下，找到一個最優(yōu)推力方向，使得在此過程中軌跡的某個性能指標(biāo)最優(yōu)。本文選取終端速度最大，即終能量最大作為性能指標(biāo)

J=‖Vf‖

(8)

定義哈密頓函數(shù)為

(9)

其中，P∈R3為速度對應(yīng)的協(xié)態(tài)變量，Q∈R3為矢徑對應(yīng)的協(xié)態(tài)變量。選取彈體縱軸方向lb為控制量，根據(jù)龐德里亞金極小值原理，最優(yōu)駐點條件為

(10)

(11)

圖1 矢量示意圖Fig.1 Figure of vector

(12)

可以得到最優(yōu)攻角方程為

(13)

由式(13)通過一維數(shù)值搜索方法可以求出最優(yōu)攻角，之后帶入式(11)即可得到最優(yōu)彈體縱軸方向。

根據(jù)動力學(xué)方程列寫正則方程為

(14)

得到協(xié)態(tài)變量的微分方程：

(15)

(16)

2.2 過程及終端約束

為了實現(xiàn)飛行器終端的位置和姿態(tài)要求，本文建立末端時刻固定的終端約束方程

(17)

其中，rfs為終端高度約束，Rf為終端位置，Vf為終端速度，lf=Rf×Vf。并由最優(yōu)控制理論可得橫截條件

(18)

由式(18)易得如式(19)所示的包含狀態(tài)與協(xié)態(tài)變量終端值的等式

(19)

同時考慮飛行器實際結(jié)構(gòu)強(qiáng)度的要求，引入攻角過程約束：

C1=α-α0≤0

C2=-α-α0≤0

(20)

在式(20)中，顯含控制變量lb(攻角α為lb和lVr的函數(shù))，根據(jù)極小值原理，當(dāng)約束觸發(fā)時有：

(21)

其中，λα為攻角約束α對應(yīng)的增廣乘子，為得到式(21)的完整表達(dá)式，要求得λα。在滿足約束時，λα=0；當(dāng)不滿足約束時，最優(yōu)駐點條件變?yōu)?/p>

(22)

(23)

2.3 最優(yōu)軌跡求解

為了獲得需要的最優(yōu)軌跡，本文根據(jù)正則方程的形式，基于狀態(tài)響應(yīng)方程的思想設(shè)計了一種迭代方法來求取最優(yōu)軌跡。定義

(24)

則正則方程變?yōu)椋?/p>

(25)

(26)

其中：

N[cosφ-cosαcos(φ-α)]}

(27)

{-Acos(φ-α)+

N[cosφ-cosαcos(φ-α)]}

(28)

根據(jù)狀態(tài)響應(yīng)方程

(29)

通過選取合適大小的τ，使積分項用分段常值函數(shù)代替，即可以得到狀態(tài)和協(xié)態(tài)變量從ti時刻到ti+τ的狀態(tài)響應(yīng)表達(dá)式：

S(ti+τ)=Φ(τ)S(ti)+

(30)

(31)

其中：

(32)

(33)

首先，式(32)表示的狀態(tài)轉(zhuǎn)移矩陣只與時間有關(guān)，由式(3)、式(11)及式(30)的形式可知，協(xié)態(tài)變量初始值完全決定了協(xié)態(tài)變量的時間歷程，由式(31)可知，狀態(tài)的時間歷程由狀態(tài)初始值和積分項決定，而積分項又與協(xié)態(tài)變量的時間歷程相關(guān)，故而，協(xié)態(tài)變量的初始值完全決定了協(xié)態(tài)變量和狀態(tài)變量的時間歷程，所以將終端約束方程轉(zhuǎn)化為協(xié)態(tài)變量初值的方程組，利用具有6個未知數(shù)的高階方程來求解協(xié)態(tài)變量初值。鑒于現(xiàn)在數(shù)值算法已十分成熟，這種非線性求根問題具有很好的適應(yīng)性與魯棒性，能夠穩(wěn)定地計算出協(xié)態(tài)變量初值，進(jìn)而得到最優(yōu)軌跡。

為了開始迭代算法，需要給出協(xié)態(tài)變量初值的初始猜想，對于間接法而言初始猜測值直接影響到算法的收斂。火箭上升段軌跡優(yōu)化一般采用真空條件下的解作為初始猜想，具體的求解方法見文獻(xiàn)[7]，在得到初始猜想值后，通過同倫法逐步增加大氣密度直到獲得100%真實大氣環(huán)境下的最優(yōu)軌跡。

3 仿真分析

為了驗證軌跡優(yōu)化方法的實用性和魯棒性，取表1中三種不同上升情況進(jìn)行了仿真分析。仿真的初始條件為：高度1310m，速度84m/s，彈道傾角82°，并采用同一組起始值。

表1 仿真算例終端約束Tab.1 Terminal constraints for simulation

通過圖2～圖4的仿真結(jié)果可以看出，算法針對三種不同的終端約束，均可由同一組起始值得到最終結(jié)果，同時求得的最優(yōu)軌跡狀態(tài)量變化合理，能夠精確滿足不同的終端和過程約束，具有良好的適應(yīng)性和魯棒性。

圖2 時間-高度曲線Fig.2 Time-altitude curve

圖3 時間-攻角曲線Fig.3 Time-angle of attack curve

圖4 時間-彈道傾角曲線Fig.4 Time-flight path angle curve

4 結(jié)論

本文針對固體火箭飛行器上升段軌跡優(yōu)化問題設(shè)計了一種基于間接法的迭代求解方法，通過本文的研究可以看出：

1)對于不同的終端和過程約束，最優(yōu)算法均可以穩(wěn)定地得到最優(yōu)軌跡，具有較強(qiáng)的適應(yīng)性和魯棒性。

2)本文采用的基于狀態(tài)響應(yīng)方程的迭代方法，通過擬牛頓法進(jìn)行數(shù)值求解，可以較快地得出最優(yōu)軌跡，并且能夠保證收斂性，可以進(jìn)一步為在線最優(yōu)閉環(huán)制導(dǎo)方法提供基礎(chǔ)。

[1] Ping Lu，HongSheng Sun，Bruce Tsai.Closed-loop endo-atmospheric ascent guidance[J].AIAA 2002-4558:1-17.

[2] Ping Lu，Lijun Zhang，Hongsheng Sun.Ascent guidance for responsive launch a fixed point approach[J].AIAA 2005-6453，2005.

[3] PingLu，Binfeng Pan.Highly constrained optimal launch ascent guidance[J].Journal of Guidance，Control，and Dynamics. 2010，33(2):404-414.

[4] 宣穎，張為華，張育林.基于Legendre偽譜法的固體運(yùn)載火箭軌跡優(yōu)化研究[J].固體火箭技術(shù)，2008，31(5):425-429.

[5] 宗群，田伯苓，竇立謙.基于Gauss偽譜法的臨近空間飛行器上升段軌跡優(yōu)化[J].宇航學(xué)報，2010，31(7):1775-1781.

[6] Dukeman G A.Atmospheric ascent guidance for rochet-powered launch vehicles[R].AIAA 2002-4559.

[7] Calise A J，Nahum Melamed，Seungjae Lee.Design and evaluation of a three-dimensional optimal ascent guidance algorithm[J].Journal of Guidance，Control，and Dynamics，1998，21(6):867-875.

Indirect Method-based Ascent Trajectory Optimization

WU Jia-liang

(Beijing Electromechanical Engineering Research Institute，Beijing 100074，China)

Ascent trajectory optimization of solid rockets is important for engineering.Aiming at this problem，a feasible approach is proposed to obtain optimal ascent trajectory.Based on zero sideslip assumption to construct model，the direction of body longitudinal axis is chosen as control input.The minimum principle is used for optimization.Using indirect method，the vacuum solution is taken as initial guess，and based on state response to construct an iteration approach.With the path constraint of attack of angle，the optimal trajectory is gotten through homotopy algorithm.The results indicate that this method converges steadily and robustly.

Solid rockets；Ascent；Optimal control；Trajectory optimization；Indirect method

2015-03-26；

2015-05-18。

吳嘉梁(1989-)，男，工程師，主要從事導(dǎo)航制導(dǎo)與控制方面的研究工作。

V412.1

A

2095-8110(2016)02-0014-06