呂名添，周 祥，張洪波

(國(guó)防科技大學(xué)空天科學(xué)學(xué)院空天工程系，長(zhǎng)沙 410073)

0 引言

小型滑翔飛行器是一種由飛機(jī)或助推器加速，在多種高度投放、進(jìn)行亞/跨聲速飛行的無(wú)動(dòng)力飛行器，具有體積小、發(fā)射時(shí)間短、使用靈活等特點(diǎn)。為保證小型滑翔飛行器能夠完成預(yù)定任務(wù)，快速生成一條滿足各類約束的參考軌跡至關(guān)重要。

根據(jù)是否對(duì)目標(biāo)函數(shù)直接尋優(yōu)，軌跡規(guī)劃方法分為直接法和間接法。由于直接法約束處理更靈活、問(wèn)題適應(yīng)性更強(qiáng)等優(yōu)勢(shì)，近年來(lái)得到了廣泛研究和應(yīng)用，其代表方法有偽譜法、凸優(yōu)化等，本文主要關(guān)注偽譜法。偽譜法采用正交多項(xiàng)式的根作為離散點(diǎn)對(duì)狀態(tài)變量和控制變量同時(shí)離散，將軌跡優(yōu)化問(wèn)題轉(zhuǎn)化非線性規(guī)劃問(wèn)題，再利用序列二次規(guī)劃(Sequential Quadratic Programming，SQP)等底層求解算法求解。SQP通過(guò)迭代求解一個(gè)二次規(guī)劃問(wèn)題序列，使得迭代解逐步收斂至原問(wèn)題的最優(yōu)解。在首次構(gòu)造迭代子問(wèn)題時(shí)，必須借助初始軌跡近似處理約束和性能指標(biāo)，以后的迭代則采用上一次的解構(gòu)造新的子問(wèn)題，直至滿足收斂條件。因此，初始軌跡會(huì)影響序列二次規(guī)劃算法的迭代子問(wèn)題構(gòu)造，進(jìn)而影響偽譜法的收斂性和求解效率。

文獻(xiàn)中提出的初始軌跡生成方法包括線性假設(shè)法、智能算法等。周祥等提出利用線性假設(shè)法生成初始軌跡，對(duì)非凸的三維剖面規(guī)劃問(wèn)題進(jìn)行凸化處理，實(shí)現(xiàn)了序列凸化算法迭代求解原始非凸問(wèn)題。孫志遠(yuǎn)等針對(duì)火星進(jìn)入軌跡規(guī)劃問(wèn)題采用線性假設(shè)法生成初始軌跡，將其作為自適應(yīng)偽譜法的輸入，實(shí)現(xiàn)了最優(yōu)軌跡的快速生成。Kim等提出引導(dǎo)策略搜索方法，每次迭代采用序列凸規(guī)劃算法生成軌跡，再通過(guò)反饋控制律生成其他軌跡，所有軌跡共同用于更新神經(jīng)網(wǎng)格策略，從而提高初始軌跡精度。Zhang等應(yīng)用隨機(jī)快速搜索樹算法生成初始軌跡，減小隨機(jī)選擇點(diǎn)與預(yù)期點(diǎn)的距離，并將其作為高斯偽譜法的輸入從而得到最優(yōu)軌跡。馬林設(shè)計(jì)了啟發(fā)式初始化策略，以基于低密度網(wǎng)格的軌跡優(yōu)化問(wèn)題最優(yōu)解作為高密度網(wǎng)格軌跡規(guī)劃的初始軌跡，改進(jìn)了非線性規(guī)劃問(wèn)題的收斂性。Wu等提出了一種基于混合粒子群算法的初始軌跡生成方法，通過(guò)該算法在優(yōu)化外循環(huán)中作用，生成有利于偽譜法收斂的初始軌跡。從已有結(jié)果來(lái)看，雖然線性假設(shè)法生成初始軌跡的效率較高，但存在精度較低、收斂性較差的問(wèn)題；智能算法可以獲得接近最優(yōu)解的初始軌跡，精度較高、收斂性好，但初始解生成的效率較低，從而降低了偽譜法的整體效率。因此，需要研究能在精度與效率間取得平衡、具備工程應(yīng)用潛力的初始軌跡生成方法。

在一定的假設(shè)條件下，滑翔飛行器的運(yùn)動(dòng)方程存在解析解，能夠近似描述飛行器的運(yùn)動(dòng)特性。王肖等基于準(zhǔn)平衡滑翔條件建立了航程與能量、傾側(cè)角的解析關(guān)系，從而得到航跡角和攻角解析解。曾亮基于一定假設(shè)，推導(dǎo)了速度、速度傾角、傾側(cè)角關(guān)于高度的解析解，有效提高了終端狀態(tài)的預(yù)測(cè)精度和效率。本文將解析解的研究思路拓展到軌跡規(guī)劃問(wèn)題中，研究以飛行軌跡的近似解析解作為偽譜法迭代初值的方法，從而提高偽譜法的收斂速度和計(jì)算效率。

1 滑翔飛行器軌跡規(guī)劃問(wèn)題建模

1.1 運(yùn)動(dòng)方程

考慮到小型滑翔飛行器的航程較短、高度跨度較小、速度較低等特點(diǎn)，提出以下幾點(diǎn)假設(shè)，作為運(yùn)動(dòng)學(xué)建模基礎(chǔ)：

1)假設(shè)地球?yàn)橐粓A球，半徑為其平均半徑；

2)略去地球自轉(zhuǎn)的影響，=0；

3)假設(shè)引力加速度受高度的影響可以忽略，設(shè)為常值=9.8 m/s；

4)飛行器以純空氣舵控制，認(rèn)為飛行器的縱軸始終處于起滑點(diǎn)速度矢量與地心矢所在的平面內(nèi)，側(cè)滑角=0。

結(jié)合上述假設(shè)，小型滑翔飛行器在半速度坐標(biāo)系中的微分方程為

(1)

式中，為飛行器相對(duì)地球速度，為當(dāng)?shù)厮俣葍A角，為航跡偏航角，規(guī)定順時(shí)針為正，為地心距，為經(jīng)度，為緯度。為阻力系數(shù)，為升力系數(shù)，為氣動(dòng)參考面積，為飛行器質(zhì)量。

1.2 軌跡規(guī)劃問(wèn)題描述

飛行器在半速度系的射程變化率為

(2)

式中，為飛行過(guò)程中的速度傾角。

本文研究縱程最大，故構(gòu)造目標(biāo)函數(shù)為

(3)

為滿足飛行器控制系統(tǒng)要求，考慮攻角、傾側(cè)角變化率約束，將原控制量擴(kuò)維至狀態(tài)量中，在方程(1)的基礎(chǔ)上增加兩個(gè)微分方程

(4)

飛行器滑翔過(guò)程滿足動(dòng)力學(xué)方程、邊界約束、過(guò)程約束等條件，以控制量、終端時(shí)間為優(yōu)化變量，可將軌跡優(yōu)化問(wèn)題寫為以下形式

(5)

2 初始軌跡在偽譜法中的影響分析

偽譜法求解軌跡規(guī)劃問(wèn)題的主要思路是將原問(wèn)題轉(zhuǎn)化為非線性規(guī)劃問(wèn)題，并通過(guò)序列二次規(guī)劃算法進(jìn)行求解。初始軌跡主要用于首次構(gòu)造二次規(guī)劃問(wèn)題。

令={,}，={,，…,()}，，(=1,2，…,)分別為離散點(diǎn)處的狀態(tài)量和控制量，則非線性規(guī)劃問(wèn)題可描述為

(6)

式中，為等式約束的個(gè)數(shù)，為所有約束個(gè)數(shù)。偽譜法處理原問(wèn)題的過(guò)程可參考文獻(xiàn)[13-14]。

SQP算法是擬牛頓法解無(wú)約束問(wèn)題的推廣，可用來(lái)求解上述非線性規(guī)劃問(wèn)題，其基本思路是把式(6)所示的非線性問(wèn)題轉(zhuǎn)為二次規(guī)劃(Quadratic Programming，QP)子問(wèn)題，通過(guò)依次求解每個(gè)子問(wèn)題得到下降方向，重復(fù)以上步驟直至滿足誤差容限，最終獲得規(guī)劃問(wèn)題的最優(yōu)解。

在每次迭代過(guò)程構(gòu)造二次規(guī)劃問(wèn)題

(7)

上述構(gòu)造主要是將原始非線性動(dòng)力學(xué)約束和過(guò)程約束在初始軌跡附近進(jìn)行線性展開，得到式(7)所示的線性化矩陣。在完成首次QP問(wèn)題構(gòu)造后，求解該QP問(wèn)題會(huì)得到一組新的迭代解，將該迭代解作為下一次迭代時(shí)求解矩陣的參考軌跡，從而啟動(dòng)SQP算法。因此，初始軌跡對(duì)于后續(xù)序列二次規(guī)劃的收斂進(jìn)程將產(chǎn)生影響，有必要研究一種能夠平衡精度與效率的初始軌跡生成方法。

初始軌跡在偽譜法求解軌跡規(guī)劃問(wèn)題中的作用如圖1所示。

圖1 初始軌跡在軌跡規(guī)劃問(wèn)題中的影響示意圖Fig.1 Schematic diagram of the influence of the initial trajectory in the trajectory planning problem

3 初始軌跡生成的解析解

小型滑翔飛行器動(dòng)力學(xué)方程具有強(qiáng)非線性，控制量與狀態(tài)量相互耦合，直接推導(dǎo)解析解是非常困難的。本節(jié)在第1部分的微分方程基礎(chǔ)上做如下假設(shè)，推導(dǎo)用于生成初始軌跡的解析解：

1)準(zhǔn)等溫大氣模型與實(shí)際大氣模型計(jì)算氣體參數(shù)比較接近，假設(shè)大氣密度模型為指數(shù)模型，即

()=-

(8)

式中，為標(biāo)準(zhǔn)大氣密度，為大氣密度參數(shù),=1720 01 m；

2)采用常值攻角、傾側(cè)角剖面；

3)假設(shè)氣動(dòng)模型在整個(gè)飛行過(guò)程為常值，不隨速度、攻角變化。

由圓球假設(shè)，地心距與高度的關(guān)系為

=+

(9)

由式(1)和式(9)得

(10)

對(duì)當(dāng)?shù)厮俣葍A角微分方程進(jìn)行簡(jiǎn)化

(11)

令=2，且設(shè)函數(shù)(,)為

(12)

則式(11)等價(jià)于

(13)

若(,)始終小于0，與呈正比，且釋放點(diǎn)處的當(dāng)?shù)厮俣葍A角小于0，則整個(gè)滑翔過(guò)程的當(dāng)?shù)厮俣葍A角都將小于0，從而保證飛行器高度單調(diào)下降。

方程(11)右端3項(xiàng)經(jīng)過(guò)數(shù)值計(jì)算，發(fā)現(xiàn)氣動(dòng)升力項(xiàng)遠(yuǎn)大于引力項(xiàng)與附加力項(xiàng)，根據(jù)式(12)中氣動(dòng)升力項(xiàng)的表達(dá)式，通過(guò)設(shè)計(jì)傾側(cè)角幅值，可構(gòu)造函數(shù)(,)為

(,)=--

(14)

在上式構(gòu)造的函數(shù)中，由于大氣密度參數(shù)與指數(shù)函數(shù)恒正，則參數(shù)的符號(hào)決定了函數(shù)(,)的正負(fù)。若參數(shù)取正值，則飛行器的高度在滑行飛行段將單調(diào)下降。

結(jié)合式(13)和(14)，可得到當(dāng)?shù)厮俣葍A角的解析表達(dá)式為

(15)

式中，為滿足初始狀態(tài)的積分常數(shù)，0和0分別為釋放點(diǎn)處的當(dāng)?shù)厮俣葍A角和高度。

同理，推導(dǎo)速度關(guān)于高度的解析解，結(jié)合式(1)和式(10)，對(duì)速度微分方程簡(jiǎn)化得到

(16)

式(16)中由于引力項(xiàng)與氣動(dòng)升力項(xiàng)量級(jí)相當(dāng)，難以略去，因此速度大小不能獲得解析解，只能得到上述半解析解。

至此，利用上述解析解實(shí)現(xiàn)初始軌跡生成的流程如下：

1)根據(jù)初始高度和終端高度確定離散點(diǎn)個(gè)數(shù)，給出高度下降曲線；

2)利用式(15)計(jì)算速度傾角關(guān)于高度的解析變化曲線；

3)利用Runge-Kutta積分方法計(jì)算方程(16)，得到速度關(guān)于高度的變化規(guī)律；

4)假設(shè)其他狀態(tài)量及控制量為常值，得到完整的初始軌跡。

4 仿真分析

4.1 仿真條件

本節(jié)利用所提解析初值進(jìn)行小型滑翔飛行器軌跡快速規(guī)劃，對(duì)解析方法的有效性進(jìn)行仿真驗(yàn)證。調(diào)用MATLAB平臺(tái)下的IPOPT求解器對(duì)軌跡規(guī)劃問(wèn)題進(jìn)行求解，相對(duì)誤差精度設(shè)為1×10。飛行器質(zhì)量為200 kg，參考面積為0.1 m，升阻比為2左右。初始任務(wù)參數(shù)為：=3，=-1°，=90°，=25 km，=0°，=0°。

在初始高度與終端高度之間均勻確定250個(gè)離散點(diǎn)，根據(jù)所提解析方法獲得速度與速度傾角在每個(gè)高度處的大小，航跡偏航角、經(jīng)緯度、傾側(cè)角均設(shè)為常值0°，攻角設(shè)為常值15°，至此可得到完整的初始軌跡，序列二次規(guī)劃算法可以啟動(dòng)。

4.2 仿真結(jié)果

基于偽譜法軌跡規(guī)劃得到的控制量，以四階Runge-Kutta積分結(jié)果作為最優(yōu)軌跡精確值，驗(yàn)證所提解析方法的可行性。圖2～5給出了解析解作初始軌跡對(duì)應(yīng)的狀態(tài)量變化曲線，藍(lán)色實(shí)線是偽譜法規(guī)劃結(jié)果，紅色短線是Runge-Kutta積分精確結(jié)果，綠色虛線為解析方法計(jì)算的初始軌跡。

圖2 速度變化曲線Fig.2 Speed curve

圖3 速度傾角變化曲線Fig.3 Velocity inclination curve

圖4 規(guī)劃與積分軌跡三維圖Fig.4 3D plot of planning and integrating trajectories

圖5 航跡偏航角變化曲線Fig.5 Track angle curve

由圖2～4中的規(guī)劃解與積分解可以看出，速度傾角全程小于0°，終端高度滿足邊界約束，速度、高度呈三段式持續(xù)減小，說(shuō)明飛行前期、后期升力較小，不滿足平衡滑翔條件，而飛行中期，大氣密度、速度在一個(gè)合適范圍內(nèi)，升力大于重力。圖4、圖5中，規(guī)劃解的航跡偏航角變化小于0.5°，飛行器側(cè)向機(jī)動(dòng)范圍較小，軌跡整體沿正東方向從而得到初始飛行能力下的最大縱程。規(guī)劃結(jié)果與積分結(jié)果航程為166.48 km，終端狀態(tài)量偏差較小，說(shuō)明優(yōu)化解的狀態(tài)量與控制量基本滿足動(dòng)力學(xué)約束，軌跡規(guī)劃方法是正確的。對(duì)比初始軌跡與規(guī)劃結(jié)果，雖然解析解與規(guī)劃軌跡相差較遠(yuǎn)，但偽譜法仍能找到滿足迭代精度的可行解，說(shuō)明偽譜法對(duì)初始軌跡具有較好的適應(yīng)性，能夠確保方法自身收斂。

從圖6、圖7可以看出，控制量全程滿足幅值約束與角速率約束，對(duì)控制執(zhí)行機(jī)構(gòu)產(chǎn)生的壓力較小。

圖6 攻角和傾側(cè)角變化曲線Fig.6 Variation curve of attack angle and bank angle

圖7 攻角和傾側(cè)角變化率曲線Fig.7 Variation rate curve of attack angle and bank angle

在偽譜法初始軌跡生成中，常用線性假設(shè)和簡(jiǎn)化數(shù)值積分兩種方法，其中線性假設(shè)法是以初始狀態(tài)、終端狀態(tài)為端點(diǎn)的線性變化曲線作初始軌跡，簡(jiǎn)化數(shù)值積分則是假定控制量為常值，通過(guò)數(shù)值積分得到初始軌跡。

表1給出了3種初始軌跡進(jìn)行軌跡規(guī)劃所得收斂情況的對(duì)比結(jié)果。CPU時(shí)間是指調(diào)用IPOPT求解器所用時(shí)間，總規(guī)劃時(shí)間為初始軌跡生成時(shí)間與CPU時(shí)間之和，迭代次數(shù)是指相對(duì)誤差滿足網(wǎng)格細(xì)化精度的迭代次數(shù)，間隔數(shù)是指網(wǎng)格細(xì)化分段內(nèi)的間隔數(shù)目。

由表1可以看出，3種方法實(shí)際計(jì)算精度均小于誤差容限1×10，得到精度較高的收斂解；解析方法與簡(jiǎn)化數(shù)值積分法需要大約1.68 s，迭代21次得到最優(yōu)軌跡；解析方法總規(guī)劃時(shí)間最短。

表1 不同初始軌跡生成方法對(duì)應(yīng)的收斂結(jié)果Tab.1 Convergence results corresponding to different initial trajectory generation methods

分析表1所示結(jié)果，主要原因在于初始軌跡的精度和與最優(yōu)解的距離。線性假設(shè)生成的初始軌跡不滿足動(dòng)力學(xué)約束，精度較低，且與最優(yōu)解距離較大，所以將其作為初始軌跡使得偽譜法迭代次數(shù)較多。簡(jiǎn)化數(shù)值積分生成的初始軌跡嚴(yán)格滿足動(dòng)力學(xué)約束，精度較高，解析方法生成的初始軌跡在簡(jiǎn)化條件下滿足動(dòng)力學(xué)約束，精度略低，因此這兩種方法的收斂迭代次數(shù)與CPU時(shí)間基本一致，計(jì)算精度、求解效率均優(yōu)于線性假設(shè)法。

總規(guī)劃時(shí)間不同是因?yàn)樯沙跏架壽E的時(shí)間不同，簡(jiǎn)化數(shù)值積分法需對(duì)8個(gè)狀態(tài)量數(shù)值積分才能獲得初始軌跡，而本文所提解析方法僅需對(duì)速度一項(xiàng)進(jìn)行積分，速度傾角通過(guò)解析公式計(jì)算，其他量通過(guò)常值確定，所以計(jì)算量小于數(shù)值積分法。因此，可認(rèn)為本文所提解析法是一種平衡了初始軌跡生成效率與求解精度的較好方法。

5 結(jié)論

本文針對(duì)傳統(tǒng)初始軌跡生成方法難以兼顧效率與精度的問(wèn)題，以小型滑翔飛行器為對(duì)象，提出了一種基于解析初值的軌跡快速規(guī)劃方法。首先基于一定假設(shè)建立了簡(jiǎn)化的運(yùn)動(dòng)學(xué)微分方程，其次分析了初始軌跡對(duì)偽譜法收斂性的影響過(guò)程，最后推導(dǎo)了速度傾角、速度關(guān)于高度的解析/半解析解。仿真中，相比線性假設(shè)法，解析方法能夠提升偽譜法的收斂迭代精度和計(jì)算效率；相比常值積分法，解析方法生成初始軌跡的計(jì)算量小且求解精度、效率與其相當(dāng)。因此，解析初值平衡了初始軌跡生成效率與求解精度，基于該解析初值的軌跡規(guī)劃方法具備較強(qiáng)的工程適用潛力。