劉 飛,王松艷,楊 明,晁 濤
(哈爾濱工業(yè)大學(xué)航天學(xué)院控制與仿真中心,哈爾濱 150001)
近些年來(lái),固體運(yùn)載火箭(Solid launch vehicle,SLV)上升段制導(dǎo)受到了國(guó)內(nèi)外學(xué)者的廣泛關(guān)注[1]?;趦?yōu)化理論,文獻(xiàn)[2]利用間接法對(duì)固體火箭的上升段軌跡進(jìn)行了優(yōu)化,將軌跡優(yōu)化問(wèn)題轉(zhuǎn)化為Hamilton兩點(diǎn)邊值問(wèn)題,并用數(shù)值方法求解。然而,間接法存在兩個(gè)明顯的缺點(diǎn):由最優(yōu)條件組成的約束方程非常復(fù)雜;數(shù)值方法對(duì)初值猜想的敏感度較高。相比于間接法,直接法(如高斯偽譜法[3])對(duì)于具有快時(shí)變和非線性特性的SLV上升段具有更大的工程實(shí)踐意義。自適應(yīng)制導(dǎo)算法是近年來(lái)的研究熱點(diǎn)[4],但在動(dòng)態(tài)參數(shù)變化較大的極端情況下,制導(dǎo)算法的可靠性無(wú)法得到保證。為了提高軌跡生成效率,許多學(xué)者對(duì)無(wú)損凸優(yōu)化技術(shù)進(jìn)行了研究,已成功應(yīng)用到運(yùn)載火箭[5]、航天器[6]、高超聲速滑翔飛行器[7]的軌跡優(yōu)化中。然而,大多數(shù)過(guò)程約束和非線性動(dòng)態(tài)都不能無(wú)損凸化。近些年來(lái),計(jì)算制導(dǎo)引起了學(xué)者們的廣泛關(guān)注。預(yù)測(cè)校正方法是一種典型的計(jì)算制導(dǎo)方法[8],但其對(duì)算法的實(shí)時(shí)性要求很高。
SLV的一個(gè)顯著特點(diǎn)是發(fā)動(dòng)機(jī)推力和關(guān)機(jī)時(shí)間不可控。對(duì)此,文獻(xiàn)[9]采用了耗盡關(guān)機(jī)技術(shù),通過(guò)改變飛行姿態(tài)角來(lái)調(diào)節(jié)所需剩余速度。在過(guò)去的幾十年里,學(xué)者們致力于SLV能量管理算法的研究[10]。Yao等[11]提出了一種單向姿態(tài)能量管理方法。文獻(xiàn)[12]研究了基于零射程線的耗盡關(guān)機(jī)能量管理算法,其優(yōu)點(diǎn)是降低了算法對(duì)關(guān)機(jī)點(diǎn)的靈敏度,但并沒(méi)有考慮過(guò)程約束。Patha等[13]研究了基于交替俯仰角變化的姿態(tài)能量管理(Attitude energy management,AEM)方法,但姿態(tài)角變化曲線需要離線計(jì)算。Wang等[14]提出了一種改進(jìn)的姿態(tài)能量管理,減小由于忽略角速度而產(chǎn)生的跟蹤誤差。上述方法均未考慮不確定性,在有參數(shù)攝動(dòng)的情況下魯棒性較差。Zarchan[15]提出的通用能量管理(General energy management,GEM)方法是一種閉環(huán)制導(dǎo)算法,其速度控制精度較高,但起始點(diǎn)和結(jié)束點(diǎn)的攻角較大,可能會(huì)違反過(guò)程約束。文獻(xiàn)[16]采用樣條能量管理方法構(gòu)建約束方程,利用數(shù)值方法在線求解。文獻(xiàn)[17]提出了一種基于動(dòng)態(tài)逆的能量管理方法。這兩種方法都能實(shí)現(xiàn)對(duì)終端速度的精確控制,但不能調(diào)整終端高度。文獻(xiàn)[18]提出了基于推力矢量控制的能量管理的完整方案。上述算法在縱向平面進(jìn)行能量管理,避免了與側(cè)向運(yùn)動(dòng)的耦合,但減弱了速度調(diào)節(jié)能力。文獻(xiàn)[19]引入側(cè)向能量管理消耗剩余能量,從而精確控制終端速度。然而由于縱、側(cè)向機(jī)動(dòng)的耦合,側(cè)向制導(dǎo)指令可能會(huì)對(duì)終端高度產(chǎn)生影響。文獻(xiàn)[20]針對(duì)滿足多約束條件的標(biāo)稱軌跡,提出了一種離線軌跡優(yōu)化方法。雖然該算法的精度較高,但不具備在線自適應(yīng)能力?,F(xiàn)有的能量管理算法僅在大氣層外進(jìn)行[21],這導(dǎo)致終端速度調(diào)節(jié)能力較差。因此,需要開(kāi)發(fā)一種大氣內(nèi)能量管理算法以延長(zhǎng)能量管理時(shí)間,提高速度調(diào)節(jié)能力。
結(jié)合上述國(guó)內(nèi)外研究現(xiàn)狀分析,本文研究的大氣層內(nèi)三維能量管理存在如下難點(diǎn):
1)由于氣動(dòng)力的存在,難以準(zhǔn)確計(jì)算剩余速度能力。
2)由于發(fā)動(dòng)機(jī)工作環(huán)境的差異,導(dǎo)致發(fā)動(dòng)機(jī)的比沖和秒流量與地面試驗(yàn)存在差異。理論計(jì)算的局限性以及風(fēng)洞實(shí)驗(yàn)與實(shí)際飛行條件不同,使得計(jì)算時(shí)用到的氣動(dòng)系數(shù)與真實(shí)值存在偏差。這些參數(shù)不確定性將影響制導(dǎo)精度。
3)與大氣層外能量管理相比,大氣層內(nèi)能量管理需要考慮動(dòng)壓約束。對(duì)于在真空中飛行的SLV,過(guò)載約束只與當(dāng)前速度有關(guān)。而動(dòng)壓約束不僅與速度有關(guān),還與高度有關(guān)。動(dòng)壓約束下的控制量可行域的計(jì)算比由過(guò)載約束構(gòu)成的可行域[22]更為復(fù)雜。
4)當(dāng)發(fā)動(dòng)機(jī)關(guān)機(jī)時(shí),SLV需要以零側(cè)向速度回到發(fā)射平面。三維制導(dǎo)帶來(lái)了終端側(cè)向位移和側(cè)向速度約束,這使得終端約束的數(shù)量增多,從而增加了制導(dǎo)的復(fù)雜性。
基于上述問(wèn)題,本文提出一種改進(jìn)的上升段制導(dǎo)算法——基于能量管理的魯棒三維制導(dǎo)算法(Robust three-dimensional energy management, R3DEM),用以解決考慮模型參數(shù)不確定的上升段在線軌跡規(guī)劃制導(dǎo)問(wèn)題。引入側(cè)向能量管理對(duì)發(fā)動(dòng)機(jī)剩余速度能力進(jìn)行耗散,提高速度調(diào)節(jié)能力。將閉環(huán)制導(dǎo)轉(zhuǎn)化為求解關(guān)于攻角和側(cè)向速度能力曲線參數(shù)的非線性方程組問(wèn)題。此外,還考慮了控制量變化率約束、過(guò)載約束和動(dòng)壓約束。針對(duì)氣動(dòng)系數(shù)、發(fā)動(dòng)機(jī)比沖和秒流量存在的攝動(dòng),利用容積卡爾曼濾波進(jìn)行模型參數(shù)辨識(shí)提高了制導(dǎo)算法的魯棒性。
在發(fā)射坐標(biāo)下建立SLV上升段運(yùn)動(dòng)模型,將三維空間運(yùn)動(dòng)分解為縱向平面運(yùn)動(dòng)和側(cè)向平面運(yùn)動(dòng),如圖1所示。v為SLV的速度;vxy和vxz分別為縱向平面和側(cè)向平面上的速度投影;Pe為發(fā)動(dòng)機(jī)推力,其方向沿SLV體軸;Pey和Pez分別為推力在縱向平面和側(cè)向平面上的投影;Θ為當(dāng)?shù)貜椀纼A角;Φ為SLV彈體軸線與當(dāng)?shù)厮矫娴膴A角,稱為當(dāng)?shù)馗┭鼋?;α為攻角;σ為?dāng)?shù)貜椀榔?;Ψ為體軸與當(dāng)?shù)乜v向平面的夾角,稱為當(dāng)?shù)仄浇牵沪聻閭?cè)滑角;D,L和R分別為阻力、升力和側(cè)力;G為重力。
圖1 上升段受力分析示意圖Fig.1 Diagram of force analysis in ascending phase
為了簡(jiǎn)化計(jì)算,對(duì)模型進(jìn)行了若干假設(shè)。第一,將SLV視為質(zhì)點(diǎn),忽略旋轉(zhuǎn)。第二,不考慮擺動(dòng)噴管的影響。第三,地球被視為一個(gè)球體。第四,忽略科氏力和離心力。第五,重力加速度為常值?;谏鲜黾僭O(shè),大氣層內(nèi)SLV上升段三維運(yùn)動(dòng)方程為
(1)
式中:h為高度;z為側(cè)向位移;Pe和m分別表示發(fā)動(dòng)機(jī)推力大小和SLV質(zhì)量;g為重力加速度。
氣動(dòng)力的計(jì)算方式如下:
(2)
式中:cD,cL和cR分別為升力系數(shù)、阻力系數(shù)和側(cè)向力系數(shù);ρ為大氣密度;Sm為氣動(dòng)有效面積。
初始條件:根據(jù)上升段三維動(dòng)力學(xué)模型(1),初始條件可表示為
(3)
過(guò)程約束:對(duì)于大氣層內(nèi)的上升段運(yùn)動(dòng),SLV受過(guò)載、攻角和側(cè)滑角變化率約束,以及動(dòng)壓約束。
其中,縱向平面和側(cè)向平面的法向過(guò)載約束分別為
(4)
式中:nv和nh分別為縱向平面和側(cè)向平面的法向過(guò)載;nvmax和nhmax分別為nv和nh允許最大值絕對(duì)值。
攻角側(cè)滑角的變化率約束分別為
(5)
動(dòng)壓約束為
q=0.5ρv2≤qmax
(6)
式中:qmax為允許的最大動(dòng)壓值。
終端約束:包含終端高度、速度、終端當(dāng)?shù)貜椀纼A角、側(cè)向位移約束以及當(dāng)?shù)貜椀榔羌s束
(7)
不確定性:SLV飛行過(guò)程中面臨的主要不確定性來(lái)自于發(fā)動(dòng)機(jī)參數(shù),即發(fā)動(dòng)機(jī)比沖和質(zhì)量秒流量
(8)
對(duì)于大氣層內(nèi)的飛行,除了式(8)中的不確定性外,還考慮氣動(dòng)系數(shù)不確定性。為了簡(jiǎn)化問(wèn)題,本文只考慮cL和cD的不確定性:
(9)
R3DEM算法的制導(dǎo)方案示意圖如圖2所示,圖中物理量的意義將在后文進(jìn)行說(shuō)明。R3DEM各部分功能如下:
在線參數(shù)辨識(shí):根據(jù)加速度計(jì)的測(cè)量值,將ΔcD, ΔcL和Pe/m視作擴(kuò)展?fàn)顟B(tài)量。通過(guò)對(duì)模型的合理簡(jiǎn)化,建立上升段參數(shù)辨識(shí)的近似模型。利用容積卡爾曼濾波對(duì)模型中的不確定性參數(shù)進(jìn)行辨識(shí),并實(shí)時(shí)修正模型參數(shù)。詳細(xì)過(guò)程見(jiàn)3.1。
在線軌跡規(guī)劃:實(shí)現(xiàn)軌跡的在線更新。構(gòu)建攻角和剩余速度能力曲線的數(shù)學(xué)形式,將終端高度、速度、當(dāng)彈道傾角、側(cè)向位移和當(dāng)?shù)貜椀榔潜硎緸殛P(guān)于曲線參數(shù)的代數(shù)方程。提出一種改進(jìn)的側(cè)向能量管理,在傳統(tǒng)真空能量管理的基礎(chǔ)上,考慮氣動(dòng)力的影響,分別給出發(fā)動(dòng)機(jī)和氣動(dòng)力產(chǎn)生的剩余速度能力的計(jì)算方法。詳細(xì)過(guò)程見(jiàn)3.3。
軌跡解算:采用牛頓法實(shí)時(shí)求解攻角參數(shù)和速度能力曲線,實(shí)現(xiàn)在線軌跡規(guī)劃。詳細(xì)過(guò)程見(jiàn)3.3。
過(guò)程約束:將動(dòng)壓、過(guò)載和控制量幅值及變化率約束轉(zhuǎn)化為攻角曲線和側(cè)向速度能力曲線的可行域。若在線軌跡規(guī)劃得到的攻角曲線和側(cè)向速度能力曲線在可行界內(nèi),則按照規(guī)劃的路徑飛行;否則,沿著可行域邊界進(jìn)行飛行,保證軌跡實(shí)時(shí)滿足路徑約束。詳細(xì)過(guò)程見(jiàn)3.4。
圖2 制導(dǎo)方案示意圖Fig.2 Diagram of the guidance scheme
由于大氣層內(nèi)上升段面臨發(fā)動(dòng)機(jī)參數(shù)和氣動(dòng)系數(shù)等多種不確定性,常用的擴(kuò)張狀態(tài)觀測(cè)器[22]對(duì)所有不確定性的總體影響進(jìn)行觀測(cè),無(wú)法將每種不確定性分開(kāi)觀測(cè)。因此,本文采用濾波方法對(duì)各參數(shù)不確定性進(jìn)行辨識(shí)。參數(shù)辨識(shí)中常用的非線性濾波方法,如擴(kuò)展卡爾曼濾波(Extended Kalman filter,EKF),是將非線性方程線性化處理。EKF算法需要計(jì)算擴(kuò)展?fàn)顟B(tài)方程的雅可比矩陣,不但增加了計(jì)算量,且容易引起濾波發(fā)散[23]。相比之下,容積卡爾曼濾波(Cubature Kalman filter,CKF)通過(guò)三階容積積分規(guī)則,利用容積點(diǎn)的加權(quán)對(duì)高斯積分進(jìn)行逼近,可達(dá)到三階Taylor展開(kāi)的精度,具有較高的濾波穩(wěn)定性[24]。
(10)
由于SLV在第二級(jí)時(shí)質(zhì)量較大,有如下近似關(guān)系:
(11)
在實(shí)際飛行中,上升段動(dòng)力學(xué)參數(shù)的攝動(dòng)具有隨機(jī)性,難以用精確的模型描述。為了簡(jiǎn)化問(wèn)題,本節(jié)將參數(shù)攝動(dòng)模型近似如下:
(12)
式中:δD,δL>0,為常數(shù)。
(13)
將待辨識(shí)的模型參數(shù)與原系統(tǒng)狀態(tài)量合并,擴(kuò)展成新的狀態(tài)量。將式(1)和式(12)組合后,采用四階龍格-庫(kù)塔法則離散化,得到擴(kuò)展后的狀態(tài)方程:
xk=F(xk-1,uk-1)
(14)
式中:xk=[v,Θ,σ,h,z,ΔcD,ΔcL,ΔPe/m]T|t=k是系統(tǒng)在k時(shí)刻的狀態(tài)向量;uk=[α,β]T|t=k是系統(tǒng)在k時(shí)刻的輸入向量。
采用四階龍格-庫(kù)塔法則將式(13)離散化,得到觀測(cè)方程:
zk=H(xk,uk)+wk
(15)
利用文獻(xiàn)[24]中的容積卡爾曼濾波,對(duì)由式(14)和(15)構(gòu)成的離散系統(tǒng)進(jìn)行辨識(shí),可實(shí)現(xiàn)對(duì)ΔPe/m, ΔcD和ΔcL的在線估計(jì)。
2.2.1大氣層內(nèi)側(cè)向能量管理的基本原理
大氣層內(nèi)能量管理的基本原理如圖3所示。橫軸和縱軸分別代表當(dāng)?shù)厮胶彤?dāng)?shù)貍?cè)向速度。vz和vf分別為當(dāng)前時(shí)刻側(cè)向速度矢量和終端時(shí)刻期望側(cè)向速度矢量;Δvf為期望速度增量;DP和RP分別是阻力和側(cè)向力在側(cè)向平面推力方向上的分量;對(duì)于耗盡式關(guān)機(jī)的SLV,側(cè)向平面內(nèi)能量管理將發(fā)動(dòng)機(jī)推力和氣動(dòng)力在側(cè)向平面所提供的剩余速度能力進(jìn)行耗散。VCAPz是側(cè)向平面剩余速度能力,可以分解為VTCAPE(由發(fā)動(dòng)機(jī)推力提供的剩余速度能力)、VACAPP和VACAPL(由氣動(dòng)力提供的剩余速度能力)。根據(jù)能量管理基本原理[22],有如下關(guān)系式:
(16)
式中:VTCAPz為側(cè)向平面內(nèi),由發(fā)動(dòng)機(jī)推力和側(cè)向力共同產(chǎn)生的負(fù)荷剩余速度能力曲線。
圖3 側(cè)向能量管理原理Fig.3 The principle of lateral energy management
注1.根據(jù)SLV在上升過(guò)程中的飛行特點(diǎn),v在推力作用下增大。SLV在發(fā)射的過(guò)程中,姿態(tài)由垂直到水平,當(dāng)?shù)貜椀纼A角從90°降至0。因此,vx=vcosΘ是單調(diào)遞增的,可以用來(lái)作為飛行過(guò)程的自變量。
2.2.2氣動(dòng)力的分解
傳統(tǒng)的縱向平面大氣層外能量管理分別計(jì)算推力和重力產(chǎn)生的剩余速度能力[16]。基于該思想,分別將D和R沿著側(cè)向平面內(nèi)發(fā)動(dòng)機(jī)推力方向和當(dāng)?shù)貍?cè)向方向分解,如圖4所示。DP,DL,RP和RL分別是D和R沿著這兩個(gè)方向的分量。Dxoz表示阻力在當(dāng)?shù)厮矫娴姆至俊?/p>
圖4 側(cè)向平面氣動(dòng)力分解Fig.4 Aerodynamic force decompositionin lateral plane
通過(guò)上述分解,可以得到R和D在側(cè)向平面內(nèi),沿推力方向和當(dāng)?shù)貍?cè)向方向上的分量為
(17)
2.2.3剩余速度能力分量的計(jì)算
將積分變量由t替換為vx,則的表達(dá)式如下:
(18)
式中:M(vx)為質(zhì)量關(guān)于vx的函數(shù)。
由于氣動(dòng)力遠(yuǎn)小于飛行第二階段的推力且固體SLV采用耗盡式關(guān)機(jī),可以得到以下近似關(guān)系:
(19)
由式(18)可知,在發(fā)動(dòng)機(jī)剩余速度能力計(jì)算時(shí),將Pe/M(vx)作為一個(gè)整體??赏ㄟ^(guò)2.1節(jié)對(duì)Pe/m進(jìn)行更新,獲得較為精確的發(fā)動(dòng)機(jī)剩余速度能力。
參考VTCAPE的計(jì)算方法,有:
(20)
式中:
(21)
2.2.4VTCAPz的起始點(diǎn)與終點(diǎn)
設(shè)VTCAPz的表達(dá)式為S(vx),則S(vx)當(dāng)前時(shí)刻下的起始點(diǎn)為(vtcosΘtcosσt,vtcosΘtsinσt)
S(vtcosΘtcosσt)=vtcosΘtsinσt
(22)
(23)
為了保證終端當(dāng)?shù)仄浇铅?0,S(vx)在終點(diǎn)(vf,0)的導(dǎo)數(shù)為:
S′(vf)=0
(24)
當(dāng)速度到達(dá)vf時(shí),為了使推力和氣動(dòng)力在側(cè)向平面的分量產(chǎn)生的剩余速度能力完全耗盡,有:
(25)
2.2.5側(cè)向位移增量Δz計(jì)算
設(shè)A(vxA,S(vxA))和B(vxB,S(vxB))分別為VTCAPz上的兩點(diǎn)。當(dāng)速度沿著VTCAPz從A至B時(shí),側(cè)向位移ΔzA→B可以按照如下計(jì)算:
(26)
由于直接求取?ΔzA→B/?vx較為困難,可以將其寫(xiě)成導(dǎo)數(shù)乘積的形式:
(27)
在圖1中,有如下近似關(guān)系:
β=Ψ-σ
(28)
(29)
將式(27)和(28)代入到式(1),有
(30)
因此,當(dāng)速度向量沿著VTCAPz由A點(diǎn)運(yùn)動(dòng)到B點(diǎn)時(shí),產(chǎn)生的側(cè)向位移ΔzA→B為
ΔzA→B=
(31)
則從當(dāng)前時(shí)刻到終端時(shí)刻,由VTCAPz和VACAPP共同引起的側(cè)向位移Δz為
(32)
2.2.6高度增量和當(dāng)?shù)貜椀纼A角增量的計(jì)算
(33)
(34)
因此,Θ和h均可以表示成關(guān)于vx的方程:
(35)
終端高度和當(dāng)?shù)貜椀纼A角約束可表示為
(36)
由于終端高度和當(dāng)?shù)貜椀纼A角通過(guò)攻角進(jìn)行控制,為了保證攻角解的唯一性,構(gòu)造了一條包含兩個(gè)待定參數(shù)的攻角曲線。假設(shè)SLV將在如下攻角的數(shù)學(xué)形式下進(jìn)行剩余飛行:
α(vx)=γ1vx+γ2
(37)
式中:γ1和γ2是制導(dǎo)周期內(nèi)需要求解的參數(shù)。
2.3.1約束代數(shù)方程的建立
由上述設(shè)計(jì)過(guò)程可知,VTCAPz=VTCAPE+VACAPP。根據(jù)式(36),(22)-(25)和(32),三維軌跡應(yīng)滿足:
(38)
由于SLV在終端時(shí)刻需要返回發(fā)射面,因此終端側(cè)向位移為零,即Δzf=0。
注3.為了保證式(38)有唯一解,α(vx)和S(vx)中未知數(shù)的總數(shù)應(yīng)該等于式(37)中的方程個(gè)數(shù)。由于α(vx)包含了2個(gè)未知數(shù)[γ1,γ2]T,S(vx)應(yīng)包含5個(gè)未知數(shù)。有無(wú)窮多數(shù)學(xué)形式的S(vx)可以滿足式(38)。對(duì)于該問(wèn)題,多項(xiàng)式曲線有兩個(gè)優(yōu)點(diǎn):首先,多項(xiàng)式函數(shù)的特點(diǎn)是它可以任意增加未知數(shù)的數(shù)目來(lái)匹配方程組的數(shù)目。其次,它的導(dǎo)數(shù)形式比較簡(jiǎn)單。因此,可用5階多項(xiàng)式曲線來(lái)擬合S(vx):
S(vx)=A(vx-vxf)4+B(vx-vxf)3+C(vx-vxf)2+
D(vx-vxf)+E
(39)
通過(guò)將攻角和剩余速度能力曲線設(shè)計(jì)為特定的數(shù)學(xué)形式,將上升段能量管理制導(dǎo)問(wèn)題轉(zhuǎn)化為代數(shù)方程求解問(wèn)題,進(jìn)而在線獲得制導(dǎo)指令。
2.3.2S(vx)和α(vx)參數(shù)計(jì)算
假設(shè)所有狀態(tài)量都是可測(cè)的,則式(38)是關(guān)于[γ1,γ2,A,B,C,D,E]T的方程組。在式(38)中:
S(vf)=VACAPL,S′(vf)=0
(40)
顯然有:
(41)
因此,式(38)可以化簡(jiǎn)為
(42)
設(shè)K=[K1,K2,K3,K4,K5]T,Π=[γ1,γ2,A,B,C]T。為了便于計(jì)算,采用Simpson法對(duì)積分進(jìn)行近似。
注4.一些改進(jìn)的數(shù)值迭代算法同樣可以有效地解決本問(wèn)題,兼顧計(jì)算效率和收斂性能,本文選擇牛頓法進(jìn)行迭代計(jì)算:
(43)
本文考慮過(guò)載、和控制量變化率、控制量幅值以及動(dòng)壓約束約束。
注6.側(cè)向平面法向過(guò)載約束和側(cè)滑角的變化率約束可行域求解方法與文獻(xiàn)[22]的相同,這里不再贅述。
2.4.1攻角變化率約束
根據(jù)導(dǎo)數(shù)分解原則,有
(44)
則
(45)
式中:
(46)
2.4.2縱向平面法向過(guò)載約束
縱向平面內(nèi)的法向過(guò)載約束如下:
(47)
其邊界情況分為:(a)α>0和(b)α<0??紤]到篇幅的限制,這里僅以α>0的情況為例進(jìn)行說(shuō)明。
根據(jù)式(47),當(dāng)觸發(fā)法向過(guò)載約束時(shí)滿足:
cosΨmax1cosΦmax1sinΘt+D
(48)
設(shè)υ=sinΦmax,則式(48)可以寫(xiě)作:
S2(vxt))cLSm
(49)
對(duì)式(49)進(jìn)行牛頓迭代,可以獲得υ1:
υ(k-1)
(50)
式中:υ(k-1)和υk分別表示相鄰兩次牛頓迭代的υ值。
根據(jù)式(44),有:
|Φ|=|γ1vx+γ2+Θt|<|γ1vx+Θt|+|γ2|=
γ2 (51) 因此,有: |γ2| (52) 由此可以得出,縱向平面內(nèi)的法向過(guò)載約束對(duì)|γ2|施加了限制。 2.4.3動(dòng)壓約束 大氣層內(nèi)能量管理不但要考慮過(guò)載約束、攻角和側(cè)滑角的幅值、變化率約束,還要滿足動(dòng)壓約束。由于側(cè)向位移相對(duì)高度而言較小,近似認(rèn)為動(dòng)壓只與高度和速度有關(guān)。假設(shè)大氣密度隨高度變化如下: ρ(h)=ρ0e-ζh (53) 式中:ρ0是零海拔處的大氣密度;ζ>0為常數(shù)。 根據(jù)式(53),需要實(shí)時(shí)滿足: (54) 定義動(dòng)壓約束觸發(fā)判斷方程(55)。第一項(xiàng)表示在某一速度vxt+τ時(shí),觸發(fā)最大動(dòng)壓約束qmax的最低臨界高度。第二項(xiàng)表示速度為vxt+τ時(shí)的預(yù)測(cè)高度。這里,τ=(vfcos Θf-vxt)/100是關(guān)于vx的積分步長(zhǎng)。當(dāng)預(yù)測(cè)高度高于臨界高度時(shí),不違反動(dòng)壓約束。 h(vxt+τ) (55) 若Jq>0,則說(shuō)明當(dāng)vx=vxt+τ時(shí),動(dòng)壓超出約束值qmax。此時(shí),對(duì)式(42)中的K5(vx)進(jìn)行替換,使得vx=vxt+τ時(shí)滿足動(dòng)壓約束: (56) 式中: (57) 另一方面,若當(dāng)前時(shí)刻Jq<0,則K5(vx)保持不變。 本節(jié)通過(guò)仿真校驗(yàn)來(lái)驗(yàn)證R3DEM算法的有效性。首先,給出仿真的初始條件。然后,驗(yàn)證該算法在不同任務(wù)中的適用性,以及在含有模型參數(shù)不確定性情況下的魯棒性。最后,通過(guò)與模型預(yù)測(cè)靜態(tài)規(guī)劃算法和改進(jìn)粒子群優(yōu)化算法的比較,驗(yàn)證算法的優(yōu)越性。 本文的研究背景是三級(jí)SLV的第二級(jí),表1給出了制導(dǎo)任務(wù)滿足的部分約束條件。 表1 初始條件、終端約束和路徑約束Table 1 Initial, terminal and path constraints 為了驗(yàn)證R3DEM算法的適應(yīng)性,對(duì)終端高度63000~72000 m,終端速度2550~3000 m/s下的四種制導(dǎo)任務(wù)進(jìn)行測(cè)試。終端當(dāng)?shù)貜椀纼A角設(shè)置為6.5°。 4種終端高度和終端速度下的仿真結(jié)果如圖5所示。圖5(a)和(c)為h-v曲線和當(dāng)?shù)貜椀纼A角隨時(shí)間變化曲線??梢缘贸?,R3DEM算法可以滿足4種不同任務(wù)下的終端狀態(tài)約束。圖5(e)為攻角隨時(shí)間變化曲線。由圖5(d)中當(dāng)?shù)貜椀榔请S時(shí)間變化曲線可知,當(dāng)速度達(dá)到終端約束值時(shí),當(dāng)?shù)貜椀榔鞘諗恐?。由圖5(b)中的z-v曲線可知,側(cè)向位移z在終端時(shí)刻收斂至0,這表明SLV最終回到了發(fā)射平面。此外,圖5(b)中側(cè)向位移變化率收斂至0,圖5(d)中終端當(dāng)?shù)貜椀榔鞘諗恐?,二者均表明終端時(shí)刻側(cè)向速度收斂至0。終端側(cè)滑角收斂至0,表示速度方向與彈體體軸方向重合,體軸調(diào)整到發(fā)射平面,滿足了第三級(jí)助推的起始條件。圖5(f)表明,為了達(dá)到精確控制速度的目的,在側(cè)向平面內(nèi)機(jī)動(dòng)消耗了多余的能量,且終端高度和速度越小,需要消耗的發(fā)動(dòng)機(jī)能量越多。圖5(g)和(h)分別為縱向平面法向過(guò)載絕對(duì)值和側(cè)向平面法向過(guò)載絕對(duì)值隨時(shí)間變化曲線。圖5(i)為動(dòng)壓隨時(shí)間變化曲線。第二級(jí)飛行過(guò)程始終滿足所有的過(guò)程約束,且觸發(fā)了過(guò)程約束的邊界條件。盡管如此,由于速度能力曲線始終在其可行域內(nèi),因此約束始終沒(méi)有超出邊界值。 圖5 不同高度、速度任務(wù)下的仿真結(jié)果Fig.5 Simulation results for various combinations of altitude and velocity 圖6(a)和(b)分別為蒙特卡洛仿真試驗(yàn)下,h-v和z-v的包絡(luò)線。所有軌跡均滿足終端高度、速度和當(dāng)?shù)貜椀纼A角約束,且攝動(dòng)軌跡均在標(biāo)稱情況附近。終端高度、速度、當(dāng)?shù)貜椀纼A角最大誤差值分別為13.3 m, 5.44 m/s, 0.013°,滿足精度指標(biāo)要求。 圖6 含有不確定性的蒙特卡洛仿真結(jié)果Fig.6 Monte Carlo simulation results with uncertainty 為了檢驗(yàn)R3DEM的速度調(diào)節(jié)能力,分別給出了R3DEM和MPSP算法的終端可達(dá)區(qū)域。由圖7(a)可知,R3DEM的可達(dá)區(qū)域完全覆蓋了MPSP。因此,與MPSP相比,R3DEM具有更強(qiáng)的終端適應(yīng)能力。IPSO適應(yīng)度隨迭代代數(shù)的變化曲線如圖7(b)所示。IPSO經(jīng)過(guò)120次迭代后收斂,得到最優(yōu)解,平均收斂時(shí)間超過(guò)5 min,這顯然不適用于在線制導(dǎo)。相比之下,R3DEM僅需10 s左右。圖7(c)和(d)分別為高度和速度隨時(shí)間變化曲線,三種算法都可以達(dá)到期望的終端速度和高度。由圖7(e)可知,MPSP產(chǎn)生的攻角大于R3DEM和IPSO。R3DEM和IPSO算法在縱、側(cè)向兩個(gè)平面對(duì)速度進(jìn)行耗散,因此產(chǎn)生的攻角和側(cè)滑角較小。由于MPSP的速度耗散完全依賴于縱向平面,導(dǎo)致了大范圍的攻角變化,進(jìn)而產(chǎn)生圖7(f)中較大的法向過(guò)載。過(guò)載約束是制約速度調(diào)節(jié)能力的重要因素,MPSP算法觸發(fā)過(guò)載約束的時(shí)間較長(zhǎng),導(dǎo)致其速度調(diào)節(jié)能力較低。R3DEM和IPSO算法的過(guò)載在整個(gè)過(guò)程中都處于約束值之下。 圖7 R3DEM,MPSP,IPSO的仿真結(jié)果對(duì)比Fig.7 Simulation results comparison among R3DEM, MPSP and IPSO 本文提出了一種基于能量管理的三維在線軌跡規(guī)劃制導(dǎo)算法,仿真結(jié)果表明,該算法適用于耗盡關(guān)機(jī)的SLV大氣層內(nèi)上升段。本文的主要結(jié)論如下: 1) 提出了一種SLV上升段制導(dǎo)算法,將求解空間軌跡曲線的問(wèn)題轉(zhuǎn)化為求解攻角和速度能力曲線參數(shù)的問(wèn)題,能以較高的精度滿足終端約束。 2) 本文提出的制導(dǎo)算法將傳統(tǒng)的縱向平面制導(dǎo)擴(kuò)展到三維空間,擴(kuò)大了終端速度的可調(diào)范圍。當(dāng)終端在縱向平面上達(dá)到規(guī)定的高度和速度時(shí),其側(cè)向位移和速度收斂至零。與模型預(yù)測(cè)靜態(tài)規(guī)劃方法相比,拓寬了終端速度調(diào)整范圍。與改進(jìn)的粒子群算法相比,計(jì)算效率得到了提高。 3) 可在過(guò)載、動(dòng)壓和控制量變化率約束,以及模型參數(shù)不確定性情況下,能以較高精度完成制導(dǎo)任務(wù)。3 仿真校驗(yàn)
3.1 仿真條件
3.2 適應(yīng)性驗(yàn)證
3.3 魯棒性驗(yàn)證
3.4 與現(xiàn)有方法比較
4 結(jié) 論