RBCC高超聲速飛行器上升段軌跡快速優(yōu)化

閆循良,王舒眉,王培臣,劉海禮

(1.西北工業(yè)大學(xué) 航天學(xué)院,陜西 西安 710072; 2.中國商飛民用飛機試飛中心,上海 201323)

近年來,高超聲速飛行器的多樣化發(fā)展和復(fù)雜任務(wù)需求對助推上升段動力系統(tǒng)的設(shè)計提出了許多新的要求,如較強的短時加速、全程多次開關(guān)機、寬域工作能力等。與傳統(tǒng)火箭動力相比,火箭基組合循環(huán)(rocket based combined cycle,RBCC)動力系統(tǒng)將高推重比、低比沖的火箭發(fā)動機和低推重比、高比沖的沖壓發(fā)動機有機地組合在一起,具有低成本、技術(shù)先進、使用靈活等特征,可用于執(zhí)行廣空域、寬速域運載及新質(zhì)作戰(zhàn)等飛行任務(wù)[1-2]。

然而,RBCC動力系統(tǒng)的引入給高超聲速飛行器的上升段軌跡設(shè)計帶來了新的挑戰(zhàn)[3-4]。首先,RBCC動力系統(tǒng)工作模態(tài)多,各模態(tài)性能迥異、動態(tài)切換;其次,上升段狀態(tài)參數(shù)變化范圍大,動力性能與飛行狀態(tài)耦合程度高,且面臨多種復(fù)雜的過程及終端約束限制。上述因素使得上升段軌跡設(shè)計的可行域極小、設(shè)計難度極大,傳統(tǒng)火箭助推運載器的軌跡設(shè)計方法和經(jīng)驗規(guī)律已難以勝任,有必要針對其動力系統(tǒng)特點開展相應(yīng)的數(shù)值優(yōu)化求解算法及策略研究。

目前,典型的上升段軌跡優(yōu)化方法可分為間接法[5]和直接法[6-9]。其中,間接法由于模型精度低、通用性不高且收斂性較差,難以應(yīng)用于具有復(fù)雜、高保真模型的RBCC上升軌跡優(yōu)化問題。隨著計算機技術(shù)的發(fā)展,直接法中的配點法、偽譜法以及粒子群等數(shù)值優(yōu)化方法在上升段軌跡優(yōu)化設(shè)計中得到廣泛應(yīng)用[6-9],然而,上述方法均未能滿足上升段軌跡快速優(yōu)化的需求。因此,為提升軌跡優(yōu)化設(shè)計的計算效率,部分學(xué)者將凸優(yōu)化方法[10-14]引入上升段軌跡設(shè)計中。Szmuk等[12]以推力方向作為控制量,利用無損凸化技術(shù)求解運載火箭上升段軌跡優(yōu)化問題,但其忽略了運載器所受升力。Liu等[13]將推力方向和氣動系數(shù)作為控制量,但考慮的氣動力模型較為簡單。王嘉煒等[14]在對固體火箭助推飛行器軌跡優(yōu)化問題進行建模時,選擇攻角作為唯一控制量,有效處理了飛行器氣動力的計算。雖然基于凸優(yōu)化的方法已逐步用于解決上升段軌跡設(shè)計問題,然而現(xiàn)有公開文獻較少涉及RBCC動力上升段軌跡優(yōu)化問題。

因此,本文以RBCC動力高超聲速飛行器為研究對象,提出了一種考慮復(fù)雜約束、高非線性及強耦合模型限制的動力上升段軌跡快速優(yōu)化方法。針對攻角控制系統(tǒng)是否存在二階滯后情況,構(gòu)建了完整的RBCC動力上升段軌跡優(yōu)化模型;設(shè)計了基于序列凸優(yōu)化的上升段軌跡快速優(yōu)化策略及算法。仿真結(jié)果表明,該方法在保證可行性的同時,突破了傳統(tǒng)軌跡優(yōu)化方法求解該類復(fù)雜約束問題時效率低或無法找到可行解的局限,能夠提高RBCC動力飛行器的軌跡優(yōu)化設(shè)計效率。

1 上升段軌跡優(yōu)化問題描述

1.1 上升段運動模型

考慮飛行器的上升段運動主要保持在縱向平面內(nèi),故不考慮側(cè)向運動。假設(shè)地球為不旋轉(zhuǎn)圓球,建立縱向平面內(nèi)上升段質(zhì)心運動方程

(1)

式中:r,V,θ,m,α,g,ms,P,D,L分別為地心距、速度、當(dāng)?shù)貜椀纼A角、飛行器質(zhì)量、攻角、當(dāng)?shù)刂亓铀俣?、發(fā)動機燃料秒流量、推力、氣動阻力及升力。對于上述模型而言,其狀態(tài)量x=[r,V,θ,m]T,控制量u=[α,ms]T。

若考慮攻角控制存在二階滯后,則攻角實際值與指令值并不完全一致。此時,補充如下方程

(2)

氣動力D,L可由(3)式計算得到

(3)

式中:q=0.5ρV2為飛行動壓;ρ為大氣密度;Sref為飛行器氣動參考面積;升、阻力系數(shù)CL,CD均為高度、攻角和馬赫數(shù)的函數(shù),可通過插值得到。

本文研究的高超聲速飛行器采用RBCC動力系統(tǒng)助推,該系統(tǒng)由火箭發(fā)動機與沖壓發(fā)動機分系統(tǒng)組合而成,因此,其推力與燃料秒流量可表示為

(4)

式中:msr和msa分別為火箭發(fā)動機和沖壓發(fā)動機的燃料秒流量;PR和PA分別代表火箭發(fā)動機與沖壓發(fā)動機推力。

火箭發(fā)動機推力PR計算如(5)式所示

(5)

式中,Ispr為火箭發(fā)動機比沖,可表征為關(guān)于高度的函數(shù)fIspr(h)。

沖壓發(fā)動機推力計算為

(6)

式中,S,Er,φ,Ispa分別為沖壓發(fā)動機進氣道橫截面積、燃油流量比、流量系數(shù)及比沖?？梢?沖壓推力與高度、馬赫數(shù)、攻角及發(fā)動機參數(shù)均相關(guān),即動力性能與飛行狀態(tài)存在高度的非線性和強耦合關(guān)系。

1.2 約束條件描述

(7)

從實際可行性以及安全飛行角度考慮,狀態(tài)變量亦需要滿足一定的約束,即有

xmin≤x≤xmax

(8)

綜合考慮飛行性能、控制系統(tǒng)性能和發(fā)動機性能,可建立控制變量約束模型,有

umin≤u≤umax

(9)

此外,對于本文研究的問題而言,端點約束可分為初始約束與終端約束,即有

(10)

1.3 優(yōu)化問題描述

為挖掘飛行器運載潛力,將優(yōu)化目標(biāo)選取為上升段末端機械能最大,優(yōu)化目標(biāo)可以表示為

J=[0.5m·V2+m·g·h]|tf

(11)

綜上,RBCC上升段軌跡優(yōu)化問題可以描述為:尋找最佳狀態(tài)量x(t)及控制量u(t),使得目標(biāo)函數(shù)(11)式最大,同時滿足狀態(tài)方程約束(1)或(2)式,以及各類約束(7)～(10)式。

為了提升優(yōu)化求解的收斂性,往往需要對模型進行無量綱化處理,限于篇幅,此處不再贅述。

2 基于凸優(yōu)化的求解框架及策略設(shè)計

典型凸優(yōu)化問題即為尋找最優(yōu)控制量,使得

(12)

2.1 時間自由問題轉(zhuǎn)化處理

對于上升段總飛行時間未知的這類時間自由問題,首先定義新的自變量τ∈[0,1]和控制量ut=tf-t0,并且將原問題的時間區(qū)間映射到[0,1]上,得到

t=t0+(tf-t0)τ,τ∈[0,1]

(13)

同時,將原自變量時間作為新的狀態(tài)變量,則有

(14)

(15)

(16)

2.2 凸化處理

考慮到前述上升段軌跡優(yōu)化模型是非凸的,將其進行凸化處理是應(yīng)用凸優(yōu)化技術(shù)數(shù)值求解的關(guān)鍵。本節(jié)以姿態(tài)控制理想情況下上升段軌跡優(yōu)化問題為例采用逐次線性化等[11]技術(shù)進行模型凸化。

1) 動力學(xué)方程凸化

(17)

取其作為初始參數(shù)序列,將動力學(xué)方程(16)在該參考狀態(tài)序列下進行一階泰勒展開,得到線性化的動態(tài)方程約束,即

(18)

同時,為保證線性化的有效性和精度,需引入信賴域約束:

(19)

(20)

其動態(tài)方程約束可類比得到,此處不做贅述。

2) 過程約束處理

(21)

對其進行線性化可得

(22)

3) 指標(biāo)函數(shù)凸化處理

記機械能公式

J=[0.5m·V2+m·g·h]|tf=f10

(23)

為非凸形式。同理,采用逐次線性化方法,可將其轉(zhuǎn)化為

(24)

此外,注意到狀態(tài)量和控制量表達(8)～(10)式均為凸形式,故無需做凸化處理。至此,該上升段軌跡優(yōu)化即轉(zhuǎn)化為一典型的凸問題。

2.3 離散化處理

采用梯形法則[11]對(18)式進行離散,可得

(25)

上標(biāo)k表示第k次迭代,m表示離散點編號,m=0,1,…,M;Δτ=(τf-τ0)/M。

類似地,(22)式的離散化形式為

(26)

控制變量約束、狀態(tài)變量約束的離散形式為

(27)

狀態(tài)變量的端點條件可以進一步表述為

(28)

綜上所述,原軌跡優(yōu)化問題轉(zhuǎn)化為參數(shù)化凸問題,可表示為如下形式

(29)

2.4 凸優(yōu)化求解策略改進及流程

1) 終端約束處理

考慮終端高度等端點等式約束,若在優(yōu)化求解過程中直接限制為固定值,初期迭代過程中很難滿足該約束,甚至?xí)霈F(xiàn)不可行解進而導(dǎo)致優(yōu)化失敗。因此,可將這類約束轉(zhuǎn)化為性能指標(biāo)中的懲罰項[11]

(30)

式中:常系數(shù)ci>0,i為終端等式約束序號;I為約束總數(shù)量?？紤]該懲罰項為非凸形式,故通過引入松弛矢量R,將其轉(zhuǎn)化為(31)式的形式引入性能指標(biāo)中

p(R)=c·R

(31)

對應(yīng)的松弛約束為

(32)

其他終端等式約束即可采用該方法進行處理。

于是,上述參數(shù)化凸問題可松弛為

(33)

2) 優(yōu)化求解流程

考慮到上述凸化后的問題與原軌跡優(yōu)化問題之間存在偏差,通過序列凸優(yōu)化算法[15]迭代求解上述凸問題可逐步逼近原問題的解。

為保證線性化的有效性和解的收斂性及精度,迭代過程以當(dāng)前優(yōu)化軌跡不斷更新參考軌跡,并以相鄰2次迭代的凸優(yōu)化解對應(yīng)的狀態(tài)量最大偏差作為收斂準(zhǔn)則,即

(34)

式中,ε為收斂誤差限。

圖1 基于序列凸優(yōu)化的上升段軌跡優(yōu)化求解流程圖

3) 判斷是否滿足(34)式和回溯搜索條件,若滿足,則優(yōu)化結(jié)束;否則,令k=k+1,轉(zhuǎn)入步驟2);

3 算例仿真及分析

3.1 仿真條件

綜合考慮計算效率與解的精度,本文仿真過程中共取61個離散點,序列凸優(yōu)化算法的初始信賴域和收斂誤差限設(shè)置為

此外,回溯搜索算法中的2個參數(shù)均取0.8[15]。

所有仿真均在搭載Intel Core i7-8700 3.20 GHz Intel處理器的臺式機完成,仿真環(huán)境為MATLAB 2016b平臺,基于CVX工具包進行軌跡優(yōu)化算法開發(fā),并調(diào)用SDPT3求解器求解凸優(yōu)化子問題。

3.2 攻角控制系統(tǒng)理想情況下的上升段軌跡優(yōu)化

將控制系統(tǒng)理想情況下軌跡優(yōu)化模型簡記為模型一,針對該模型進行終端機械能最大的上升段軌跡優(yōu)化。凸優(yōu)化算法經(jīng)過12輪迭代后收斂,每輪迭代約耗時1.8 s,優(yōu)化過程共耗時19.88 s,優(yōu)化結(jié)果如圖2～3所示。

圖2 模型一對應(yīng)的部分狀態(tài)量優(yōu)化結(jié)果

圖3 模型一對應(yīng)的控制量優(yōu)化結(jié)果

分析仿真結(jié)果可知,上升段共飛行135.66 s,終端質(zhì)量為800 kg,終端高度為62.08 km,速度為4 866.99 m/s,滿足終端約束;動壓、過載、駐點熱流密度均滿足給定的過程約束限制。將所得控制量進行插值后代入運動方程進行積分,所得積分結(jié)果與優(yōu)化結(jié)果對應(yīng)的狀態(tài)相對偏差小于0.1%,驗證了優(yōu)化結(jié)果的可行性,表明本文所提序列凸優(yōu)化方法可逼近原軌跡優(yōu)化問題的解。圖2～3中虛線為凸優(yōu)化所使用的初始猜測軌跡,本文僅將初始猜測軌跡取為初、末端狀態(tài)猜測值的連線,這種情況下,凸優(yōu)化方法仍能有效、快速地收斂至優(yōu)化結(jié)果?？梢?本文所設(shè)計的優(yōu)化算法和求解策略可以快速、有效地求解RBCC高超聲速上升這一復(fù)雜、非線性、強耦合的非凸最優(yōu)控制問題。

此外,結(jié)合仿真結(jié)果可知,對于該問題而言,上升段飛行主要可分為3段。第一段為0～22 s,此時飛行馬赫數(shù)較低,因此沖壓發(fā)動機推力小、效率低,RBCC動力系統(tǒng)工作在沖壓+火箭混合模態(tài),且以火箭加速為主,此時馬赫數(shù)逐漸增加,但飛行高度變化較小。第二段為22～91 s,該段飛行高度、馬赫數(shù)都較適宜沖壓發(fā)動機工作,故此時火箭發(fā)動機關(guān)機,RBCC動力系統(tǒng)工作在沖壓模態(tài),飛行馬赫數(shù)和高度均緩慢變化,直至爬升至約30 km,飛行器在該段將達到動壓約束邊界并保持約17 s。第三段對應(yīng)為91 s至上升段結(jié)束,此時由于高度較高,沖壓發(fā)動機效率逐步降低,RBCC動力系統(tǒng)先工作在沖壓+火箭混合模態(tài),后逐步過渡到純火箭模態(tài)直至燃料耗盡,此時飛行高度和馬赫數(shù)快速增加。事實上,這種火箭動力峰-谷-峰的工作過程能夠最大限度地發(fā)揮火箭發(fā)動機工作范圍廣和沖壓發(fā)動機推進效率高的優(yōu)勢,從而使飛行器能夠在燃料一定的情況下獲得更大終端高度、馬赫數(shù),即實現(xiàn)終端機械能最大。

此外,為探究離散點數(shù)選取對凸優(yōu)化算法精度及計算效率的影響,取31和91個離散點分別進行軌跡優(yōu)化,并利用所得優(yōu)化控制量進行積分,所得優(yōu)化及積分結(jié)果對比如表1所示。

表1 不同離散點數(shù)量優(yōu)化及積分結(jié)果對比

由表1可知,3組優(yōu)化結(jié)果的迭代次數(shù)基本一致,而隨著離散點數(shù)的增加,計算耗時有所增加,優(yōu)化解與積分解之間的誤差顯著降低,表明優(yōu)化解的可行性逐漸提升。相較于離散點取61與91的情況,離散點取31時,終端狀態(tài)的優(yōu)化解與積分解誤差較大,表明該組優(yōu)化結(jié)果不可行;而相較于離散點91的情況,離散點取61的優(yōu)化結(jié)果在精度基本相當(dāng)情況下,計算耗時減少了30%。因此,綜合考慮計算效率和精度,取離散點為61較為合適。

3.3 考慮攻角控制系統(tǒng)二階滯后上升段軌跡優(yōu)化

將考慮控制系統(tǒng)二階滯后的軌跡優(yōu)化模型簡記為模型二,針對該模型進行上升段軌跡優(yōu)化。仿真中,取阻尼比ξ=0.7,固有頻率ωn=0.5。凸優(yōu)化算法經(jīng)過13輪迭代后收斂,每輪迭代約耗時1.8 s,優(yōu)化過程共耗時23.08 s,優(yōu)化結(jié)果如圖4～5所示。

圖4 考慮控制滯后的狀態(tài)量優(yōu)化結(jié)果

由結(jié)果可知,上升段共飛行136.20 s,剩余質(zhì)量為800 kg,終端高度為62.18 km、飛行馬赫數(shù)為15.59,滿足終端約束,動壓、過載、駐點熱流密度均滿足給定的過程約束限制。

由圖4可知,模型一、二的優(yōu)化結(jié)果較為接近,而由圖4d)的優(yōu)化曲線可知,模型二的彈道傾角變化存在一定滯后,尤其是在飛行前中期的幾個拐點位置,考慮是控制系統(tǒng)滯后導(dǎo)致實際攻角變化滯后于模型一,這一點亦可從圖4e)的攻角指令曲線得以印證。

由圖5a)可知,模型二的攻角變化率總是提前模型一變化,考慮是模型二控制系統(tǒng)存在二階滯后,需要提前改變控制指令,以補償滯后效應(yīng),使實際攻角曲線盡可能逼近模型一的最優(yōu)結(jié)果。此外,模型一終端機械能為952 380 kJ,模型二終端機械能為952 310 kJ,略低于模型一?？梢?模型二的控制滯后犧牲了一定的最優(yōu)性。

圖5 考慮控制滯后的控制量優(yōu)化結(jié)果

為進一步探究自然頻率ωn取值對軌跡的影響,分別取ωn=0.1,ωn=0.9進行仿真,結(jié)果如圖6及表2所示。由表2可知,相較于理想情況,考慮控制系統(tǒng)二階滯后特性軌跡優(yōu)化結(jié)果的最優(yōu)性有所降低,且最優(yōu)性隨著ωn的增大而提升。由圖6可知,當(dāng)ωn較小時,攻角曲線較為光滑,此時實際攻角對攻角指令的響應(yīng)速度較慢,攻角變化較指令值滯后情況明顯,因此需要提前給出指令信號以更好地操縱飛行。當(dāng)ωn較大時,實際攻角對攻角指令的跟蹤效果更好,此時飛行器具有較好的操縱性。通常來說,飛行器自然頻率由其靜穩(wěn)定性決定,但設(shè)計控制系統(tǒng)時,一般情況下要求自然頻率低于自動駕駛儀頻率以免出現(xiàn)共振。同時,由于控制系統(tǒng)的超調(diào)量及振蕩次數(shù)均只與阻尼比ξ相關(guān),而實際攻角的最大峰值并未超過指令攻角,因此設(shè)計姿態(tài)控制系統(tǒng)時應(yīng)依據(jù)實際需要進行控制參數(shù)選擇。

表2 不同ωn取值優(yōu)化結(jié)果

圖6 不同ωn取值優(yōu)化結(jié)果攻角對比

4 結(jié) 論

針對RBCC高超聲速上升段軌跡快速優(yōu)化問題,本文設(shè)計了一種基于序列凸優(yōu)化的軌跡優(yōu)化算法和策略,相關(guān)結(jié)論如下:

1) 所設(shè)計優(yōu)化算法可以快速、有效地求解復(fù)雜工作模態(tài)下RBCC高超聲速上升這一復(fù)雜、非線性、強耦合的非凸最優(yōu)控制問題;

2) RBCC發(fā)動機火箭模態(tài)峰-谷-峰的工作過程能夠最大限度地結(jié)合并發(fā)揮火箭發(fā)動機工作范圍廣和沖壓發(fā)動機經(jīng)濟性高的優(yōu)勢,從而使飛行器在燃料更省的情況下獲得更大的終端高度、馬赫數(shù)。

3) 相較于理想情況,考慮控制系統(tǒng)二階滯后特性軌跡優(yōu)化結(jié)果的最優(yōu)性有所降低,且降低程度與自然頻率ωn有關(guān),設(shè)計姿態(tài)控制系統(tǒng)時應(yīng)依據(jù)實際需要進行控制參數(shù)選擇。