(西南科技大學(xué) 計(jì)算機(jī)科學(xué)與技術(shù)學(xué)院,四川 綿陽(yáng) 621010)
自中國(guó)進(jìn)入“超級(jí)航空2018”,就預(yù)示著在未來(lái)幾十年中我國(guó)將在航空領(lǐng)域迎來(lái)更大的突破,而高超聲速技術(shù)則是突破航空領(lǐng)域的重要研究點(diǎn)。飛行器最優(yōu)軌跡規(guī)劃是高超聲速技術(shù)研究的關(guān)鍵點(diǎn)之一,其求解方法主要有兩大類:間接法和直接法[1]。偽譜法以較少的參數(shù)、較高的計(jì)算精度、指數(shù)級(jí)的收斂速度,同時(shí)兼具了直接法和間接法的優(yōu)勢(shì),被認(rèn)為是目前最具潛力的軌跡最優(yōu)控制求解方法[2]。
直接求解飛行器軌跡最優(yōu)控制問(wèn)題難度較大,利用偽譜法將軌跡最優(yōu)控制問(wèn)題進(jìn)行離散,對(duì)離散化得到的非線性規(guī)劃(Nonlinear Programming,NLP)問(wèn)題進(jìn)行插值獲得原問(wèn)題的解,可將求解最優(yōu)軌跡控制問(wèn)題上的難度轉(zhuǎn)移到NLP問(wèn)題的求解上,轉(zhuǎn)化得到的NLP問(wèn)題表現(xiàn)出優(yōu)化空間大、約束控制量多等特點(diǎn)[3-4]。近幾年,在選擇求解該類NLP問(wèn)題的策略上,很多研究者偏向使用序列二次規(guī)劃(Sequential Quadratic Programming,SQP)算法[5]來(lái)求解。由于SQP算法可以將一些約束添加到變量中去,對(duì)于強(qiáng)約束問(wèn)題的求解有很大優(yōu)勢(shì)。近兩年有研究者使用SQP算法對(duì)高超聲速飛行器周期巡航軌跡[6]、爬升軌跡[7]進(jìn)行求解,取得了較好的結(jié)果。
作為求解強(qiáng)約束NLP問(wèn)題的主流方法,SQP算法具有計(jì)算速度快,求解精度高等優(yōu)勢(shì),但當(dāng)NLP問(wèn)題規(guī)模過(guò)大時(shí),使用SQP算法求解存在初值敏感,收斂速度減慢的情況。因此求解較大規(guī)模的NLP問(wèn)題時(shí),選擇一個(gè)可靠的初值是十分必要的。K.Okuda等人[8]在求解帶翼火箭彈道優(yōu)化時(shí),提出使用具有全局最優(yōu)性的遺傳算法為SQP算法提供初始值,有效的解決初值選取問(wèn)題。但遺傳算法等全局優(yōu)化算法在求解約束性較強(qiáng)的NLP問(wèn)題時(shí),存在精度低,收斂速度慢,結(jié)果隨機(jī)抖動(dòng)等問(wèn)題。
本文結(jié)合偽譜法與SQP算法各自的優(yōu)勢(shì),提出基于偽譜法的滑翔軌跡多級(jí)迭代優(yōu)化策略。以X-51A[9]相似飛行器模型為研究對(duì)象,利用增量法與查表插值建立飛行器縱向氣動(dòng)力模型。使用偽譜法與SQP算法對(duì)滑翔段軌跡進(jìn)行求解,提出多級(jí)迭代優(yōu)化策略彌補(bǔ)軌跡規(guī)劃過(guò)程中,SQP算法求解大規(guī)模NLP問(wèn)題的不足之處,并與傳統(tǒng)方法求解出的狀態(tài)量與控制量仿真飛行狀態(tài)進(jìn)行對(duì)比。仿真結(jié)果表明本文提出的多級(jí)迭代優(yōu)化策略穩(wěn)定性好,計(jì)算精度高,收斂速度快,在實(shí)際工程應(yīng)用中起到了很好的支撐作用。
X-51A相似飛行器是氣動(dòng)推進(jìn)一體化飛行器,如圖1所示,模型總長(zhǎng)4.26 m,最大寬度0.58 m,最大高度0.566 m,翼展最大0.927 m,兩舵面與水平面成45°,內(nèi)流道寬0.225 m,進(jìn)氣道總收縮比為4.9,內(nèi)收縮比為2.4。
圖1 一體化飛行器數(shù)模
本文借鑒意大利航空航天中心(CIRA)[10]利用增量法實(shí)現(xiàn)PRORA-USV1高超聲速飛行器風(fēng)洞試驗(yàn)數(shù)據(jù)與數(shù)值計(jì)算數(shù)據(jù)的融合建模思路,使用數(shù)值方法模擬X-51A相似飛行器的數(shù)據(jù)集進(jìn)行建模,為仿真實(shí)驗(yàn)提供數(shù)據(jù)基礎(chǔ)。結(jié)合增量建模理論,在基準(zhǔn)數(shù)據(jù)集與舵效數(shù)據(jù)集的基礎(chǔ)上構(gòu)建增量模型,流程圖如圖2所示。
圖2 增量模型建模流程圖
(1)
(2)
式中,ρ∞為大氣密度,L為升力,D為阻力,m表示飛行器質(zhì)量,Sref=1m2,bref=1m2,V∞表示飛行速度。
滑翔段軌跡規(guī)劃是X-51A相似飛行器軌跡優(yōu)化的重要階段,滑翔飛行能延長(zhǎng)飛行器整體航程,實(shí)現(xiàn)遠(yuǎn)距離投送或打擊。由于高超聲速飛行器的動(dòng)力學(xué)系統(tǒng)較為復(fù)雜,其建模所涉及的因素繁多,為研究需要,將地球視作一個(gè)均勻球體,認(rèn)為滑翔飛行過(guò)程中飛行器是一個(gè)質(zhì)點(diǎn),質(zhì)量保持不變,并忽略外界環(huán)境的干擾?;瓒蔚倪\(yùn)動(dòng)學(xué)方程形式含多個(gè)狀態(tài)變量與非線性等式約束,且方程中隱含有控制變量。本文以最大滑翔距離作為優(yōu)化目標(biāo),其余指標(biāo)作為約束條件,構(gòu)建最優(yōu)軌跡控制問(wèn)題,以分析X-51A相似飛行器的滑翔彈道軌跡,運(yùn)動(dòng)模型表達(dá)式如下:
(3)
式中,x=[H,V,γ]為軌跡狀態(tài)變量,分別表示飛行高度、飛行速度、飛行航跡角;u=α為軌跡的控制變量,α為攻角;tf為終端時(shí)間;D為阻力方向的合力,L為升力,m為飛行器質(zhì)量,g為重力加速度,δ為最優(yōu)控制變量α對(duì)應(yīng)飛行姿態(tài)下保持縱向平衡的俯仰舵偏。J為目標(biāo)函數(shù),表示滑翔段的最大飛行距離。
飛行器滑翔段為無(wú)動(dòng)力飛行,為保證平衡飛行,在滑翔段軌跡設(shè)計(jì)過(guò)程中要考慮相應(yīng)的約束條件,本文的參數(shù)約束如表1所示。
表1 滑翔段約束條件
偽譜法源于求解微分方程的譜方法,其思想是在離散的狀態(tài)變量與控制變量基礎(chǔ)上,構(gòu)造逼近式,求得微分矩陣與積分算子,將最優(yōu)控制問(wèn)題轉(zhuǎn)換為NLP問(wèn)題,其核心技術(shù)包括:離散網(wǎng)格劃分策略的選擇;離散點(diǎn)插值函數(shù)對(duì)連續(xù)函數(shù)的逼近;插值函數(shù)的微積分對(duì)連續(xù)函數(shù)的微積分的近似。偽譜法參數(shù)化過(guò)程細(xì)化為以下幾步。
1)時(shí)域變換:
本文使用Legendre多項(xiàng)式作為插值基,其系數(shù){Ln(τ)}是τ∈[-1,1]上的正交多項(xiàng)式。當(dāng)n≥1時(shí)滿足以下遞推關(guān)系:
(4)
由于插值多項(xiàng)式的作用域?yàn)棣印蔥-1,1],所以需要將t∈[0,tf]時(shí)域變換到τ∈[-1,1]:
(5)
2)變量離散化:
通過(guò)離散網(wǎng)格劃分策略選取合適的配點(diǎn),能有效避免逼近可能出現(xiàn)的Runge現(xiàn)象,通過(guò)Gauss-Lobato求積,取n個(gè)節(jié)點(diǎn)將問(wèn)題離散,獲得離散時(shí)間:Tn=[τ1,τ2…τn+1]。離散時(shí)間對(duì)應(yīng)的狀態(tài)變量xn=[Hn,Vn,γn],Hn=[h1,h2…h(huán)n+1],Vn=[v1,v2…vn+1],γn=[γ1,γ2…γn+1],mn=[m1,m2…mn+1];控制變量un=[αn],αn=[α1,α2…αn+1]。
3)LGL配點(diǎn)近似逼近狀態(tài)與控制變量:
(6)
其中,Xi∈Rnx,Ui∈Rnu,φi(τ)為n階Lagrange插值基數(shù):
(7)
4)最優(yōu)控制轉(zhuǎn)化為NLP問(wèn)題:
通過(guò)前幾個(gè)步驟,可得到參數(shù)化后的NLP問(wèn)題,具體表達(dá)式如下所示:
(8)
其中: [Hn,Vn,γn]為軌跡狀態(tài)變量的離散向量,n為偽譜配點(diǎn)數(shù);D,L分別為各離散狀態(tài)量下的阻力與升力;D為微分矩陣;wi為積分算子在對(duì)應(yīng)第i個(gè)狀態(tài)變量下的積分權(quán)重:
(9)
初始值的選取對(duì)于NLP問(wèn)題的求解來(lái)說(shuō)具有很大影響。不可靠的初值,不僅會(huì)影響收斂速度,甚至有可能導(dǎo)致無(wú)法收斂,找不到可行解或最優(yōu)解的情況。并且隨著NLP問(wèn)題規(guī)模的增大,約束性的強(qiáng)增,初值對(duì)求解結(jié)果的影響也會(huì)不斷增大。偽譜法求解飛行器最優(yōu)軌跡的精度是由配點(diǎn)數(shù)量的多少所決定的,一般情況下配點(diǎn)數(shù)越多,求解精度越高。但是轉(zhuǎn)化得到的NLP問(wèn)題的規(guī)模也會(huì)變得更大,在使用SQP算法求解時(shí),會(huì)出現(xiàn)初值敏感、收斂速度減慢、陷入局部最優(yōu)的情況。而配點(diǎn)數(shù)過(guò)少,則會(huì)導(dǎo)致偽譜法求解精度無(wú)法滿足要求。
本文在滑翔軌跡問(wèn)題離散過(guò)程中取n個(gè)配點(diǎn)時(shí),會(huì)形成n+1個(gè)離散變量,對(duì)應(yīng)決策向量的維度為4n+5,其中狀態(tài)變量與控制變量構(gòu)成4(n+1)維決策變量,最后一維決策變量為終端時(shí)間tf,采用一個(gè)多級(jí)迭代優(yōu)化策略,對(duì)含高維度決策變量的NLP問(wèn)題求解初值進(jìn)行選取。
首先設(shè)置較少偽譜配點(diǎn),保證問(wèn)題規(guī)模不至于影響SQP算法的求解,對(duì)于得到的結(jié)果進(jìn)行插值,然后將插值后的結(jié)果作為增加配點(diǎn)后的大規(guī)模NLP問(wèn)題的初始值。詳細(xì)過(guò)程如下:
1)在變量上下邊界條件之間進(jìn)行初始猜測(cè)插值,獲得一個(gè)初始值x0;
2)設(shè)置少量偽譜配點(diǎn),通過(guò)偽譜法參數(shù)化構(gòu)建中小規(guī)模NLP問(wèn)題,以x0作為初始值,使用SQP算法求解NLP問(wèn)題得到粗略解x1;
3)增加配點(diǎn)數(shù),對(duì)獲得的粗略解x1進(jìn)行插值,將插值結(jié)果作為增大規(guī)模后的NLP問(wèn)題的初始值,求解NLP問(wèn)題,重復(fù),直到收斂。
多級(jí)迭代優(yōu)化策略流程圖如圖3所示。
圖3 多級(jí)迭代優(yōu)化策略流程圖
在初始猜測(cè)插值過(guò)程中,由于控制變量不光滑,如果與狀態(tài)變量一樣選擇使用Lagrange多項(xiàng)式插值,其變化過(guò)大時(shí)會(huì)出現(xiàn)Gipps現(xiàn)象,導(dǎo)致控制變量超出邊界。分段三次Hermite插值多項(xiàng)式不僅要求在節(jié)點(diǎn)上的值要等于目標(biāo)函數(shù)值,而且要求在節(jié)點(diǎn)上的導(dǎo)數(shù)值也要與目標(biāo)函數(shù)的導(dǎo)數(shù)值相等,光滑性更好。選擇使用分段三次Hermite插值對(duì)離散的控制變量進(jìn)行插值,可有效避免出現(xiàn)Gipps現(xiàn)象。
SQP算法是求解強(qiáng)約束NLP問(wèn)題最有效的方法之一,其思想是將原NLP問(wèn)題拆分為多個(gè)二次規(guī)劃子問(wèn)題,通過(guò)逐步求解二次規(guī)劃子問(wèn)題的解,為主迭代過(guò)程提供搜索方向,保證主迭代逐漸收斂于NLP問(wèn)題的最優(yōu)解。SQP算法的簡(jiǎn)要計(jì)算步驟如下:
1)初值設(shè)置:選擇初始點(diǎn)x0,初始矩陣H1是n×n的正定矩陣,令k=1,計(jì)算精度0<ε<Δ;
2)求解二次規(guī)劃子問(wèn)題,確定下次迭代搜索方向d:
(10)
得到滿足原問(wèn)題Kuhn-Tucker條件的近似解(d(k),λ(k));
3)若‖d(k)‖≤ε,則x(k+1)=x(k)+d(k)為求得的最優(yōu)解,停止計(jì)算;
4)如果P(x(k)+d(k),r)
5)一維搜索:對(duì)minP(x(k)+λd(k),r)進(jìn)行求解,得到最優(yōu)步長(zhǎng)λ(k),令x(k+1)=x(k)+λd(k);
6)更新矩陣H(k):根據(jù)擬牛頓公式更新Hessen矩陣H(k),確定新的正定矩陣H(k+1),令k=k+1,轉(zhuǎn)到第2)步。
通過(guò)MATLAB進(jìn)行仿真驗(yàn)證,仿真實(shí)驗(yàn)過(guò)程中設(shè)置X-51A相似飛行器質(zhì)量m=1 600 kg,考慮到滑翔飛行的終點(diǎn)與投射、打擊目標(biāo)之間還存在一定的距離,因此需要對(duì)高度與速度設(shè)置相應(yīng)的終端條件。設(shè)置初始速度V0=6 mach,初始高度H0=27 000 m,終端速度Vf=3 mach,終端高度Hf=1 000 m,為得到最遠(yuǎn)滑翔航程,不對(duì)終端時(shí)間做限制。
經(jīng)過(guò)反復(fù)計(jì)算,整理結(jié)果如表2~3所示。表2中n為配點(diǎn)數(shù),z為對(duì)應(yīng)的決策變量數(shù),t(s)為使用傳統(tǒng)方法求解對(duì)應(yīng)配點(diǎn)數(shù)規(guī)模的NLP問(wèn)題計(jì)算時(shí)間,titer(s)為使用多級(jí)迭代優(yōu)化策略下的計(jì)算時(shí)間。通過(guò)對(duì)比可以看出隨著偽譜法配點(diǎn)個(gè)數(shù)的增加,計(jì)算所需要的時(shí)間也在快速增加,多級(jí)迭代優(yōu)化策略可有效的減少計(jì)算時(shí)間。配點(diǎn)數(shù)越多,節(jié)省的計(jì)算時(shí)間越多,體現(xiàn)出的效果越好。
仿真實(shí)驗(yàn)過(guò)程中設(shè)置的最大迭代次數(shù)為1 000,表3中Kopt為最優(yōu)迭代次數(shù),L為計(jì)算求得的目標(biāo)值。通過(guò)分析比較,本文提出的多級(jí)迭代優(yōu)化策略能大幅度減少迭代次數(shù),加快收斂速度,獲得更精確的解。圖4、圖5為配點(diǎn)數(shù)n=40時(shí),傳統(tǒng)方法與使用多級(jí)迭代優(yōu)化策略計(jì)算下,目標(biāo)值的迭代情況??梢钥吹剑褂枚嗉?jí)迭代優(yōu)化策略時(shí),SQP算法能更快的收斂得到最優(yōu)值,減少計(jì)算時(shí)間。
表2 計(jì)算時(shí)間對(duì)比情況
表3 目標(biāo)值與迭代次數(shù)
圖4 傳統(tǒng)方法目標(biāo)值迭代情況
圖5 多級(jí)迭代優(yōu)化策略目標(biāo)值迭代情況
SQP算法對(duì)偽譜法的離散配點(diǎn)數(shù)量敏感,配點(diǎn)數(shù)越多決策變量維度越高,算法收斂時(shí)間越長(zhǎng)。但配點(diǎn)數(shù)過(guò)少,則會(huì)影響結(jié)果精度。對(duì)不同配點(diǎn)數(shù)求得的結(jié)果,進(jìn)行反解常微分方程,結(jié)果表明配點(diǎn)數(shù)過(guò)少時(shí),求得的結(jié)果會(huì)出現(xiàn)很大偏差,情況如圖6所示。
圖6 配點(diǎn)數(shù)對(duì)解的影響
由于控制變量的不光滑性,使用Lagrange多項(xiàng)式插值,在變量變化較大的情況下會(huì)出現(xiàn)Gipps現(xiàn)象,造成控制變量超出邊界,情況如圖7所示。采用分段三次Hermite插值方法對(duì)控制變量進(jìn)行插值,能保證光滑性,有效避免Gipps現(xiàn)象與Runge現(xiàn)象。
圖7 PCHIP插值克服Gipps現(xiàn)象
通過(guò)仿真實(shí)驗(yàn),對(duì)比本文方法與傳統(tǒng)方法求解的狀態(tài)變量與俯仰舵偏的飛行狀態(tài)。結(jié)果表明多級(jí)迭代優(yōu)化策略相比于傳統(tǒng)方法,在求解精度和收斂速度上有很大的優(yōu)勢(shì),圖8~11對(duì)比了兩種方法下配點(diǎn)數(shù)n=40時(shí),狀態(tài)變量和俯仰舵偏角的變化情況,虛線為使用傳統(tǒng)方法的變化情況。傳統(tǒng)方法求解下,飛行高度與終端條件設(shè)置的高度相差大約1818 m,速度相差大約0.3 mach。多級(jí)迭代優(yōu)化策略方法求解下,飛行高度與終端條件設(shè)置的高度相差大約86 m,速度相差大約0.01 mach。本文方法與傳統(tǒng)方法相比,求解偏差降低了很多,所得結(jié)果更加接近最優(yōu)軌跡。并且在多級(jí)迭代優(yōu)化策略求解下,飛行器滑翔時(shí)間更長(zhǎng),得到的滑翔水平距離比傳統(tǒng)方法遠(yuǎn)大約220 m,能獲得更遠(yuǎn)航程,實(shí)現(xiàn)更遠(yuǎn)距離的投送或打擊。通過(guò)以上各項(xiàng)指標(biāo),證明了本文提出的多級(jí)迭代優(yōu)化策略的有效性與高效性。
圖8 滑翔軌跡對(duì)比圖
圖9 速度變化對(duì)比圖
圖10 航跡角變化對(duì)比圖
圖11 俯仰舵偏變化對(duì)比圖
本文提出的多級(jí)迭代優(yōu)化策略有效的結(jié)合了偽譜法參數(shù)化最優(yōu)軌跡控制問(wèn)題與SQP算法求解強(qiáng)約束非線性規(guī)劃問(wèn)題的優(yōu)勢(shì)。偽譜法參數(shù)化后的NLP問(wèn)題具有大規(guī)模、強(qiáng)約束等特征,直接使用SQP算法求解時(shí)存在對(duì)初值敏感、收斂速度減慢等問(wèn)題。多級(jí)迭代優(yōu)化策略,為求解大規(guī)模NLP問(wèn)題提供一個(gè)較為有效的初始值,可有效解決SQP算法在求解過(guò)程中存在的不足。仿真驗(yàn)證表明,在多級(jí)迭代優(yōu)化策略下,NLP求解算法求解精度更高,計(jì)算時(shí)間更少。且該方法在實(shí)際工程應(yīng)用中取得了良好效果。
隨著高超聲速技術(shù)的發(fā)展,全軌跡優(yōu)化的研究也在火熱進(jìn)行中,所面對(duì)的問(wèn)題規(guī)模將更大,約束性將更強(qiáng),多級(jí)迭代策略體現(xiàn)的價(jià)值也將更大。