基于解析剖面的時間協同再入制導

2019-03-29 06:46:02王肖郭杰唐勝景祁帥

航空學報 2019年3期

王肖，郭杰，唐勝景，祁帥

北京理工大學宇航學院，北京 100081

高超聲速滑翔飛行器作為作戰武器使用時，具有速度快、射程遠、精度高、機動和突防能力強等特點，近年來受到世界各國的廣泛關注[1-2]。與此同時，針對高超聲速目標威脅，各國相繼研發了THAAD、宙斯盾、C-400等防空反導武器系統，使得單個高超聲速滑翔飛行器的突防能力和作戰效能大大降低[3]。在此背景下，發展多高超聲速滑翔飛行器協同打擊技術，可大大提高對防空反導武器系統的突防概率，實現對目標的飽和攻擊[4]。

以文獻[5-7]為代表的基于協調變量的雙層協同制導架構，以包含協調變量的協調策略為上層協調控制，以滿足飛行器特點的帶約束制導律為底層控制。以文獻[8-10]為代表的“領彈-從彈”協同制導架構，通過使從彈跟蹤領彈運動狀態以實現多導彈協同制導。文獻[11-12]提出了一種雙階段協同制導架構，第1階段通過一致性算法使過渡狀態一致，然后切換至比例導引實現終端協同，無需剩余時間估計。現有的協同制導大部分是基于以上3種協同架構。

然而，上述協同制導均研究的是末制導段。對于高超聲速滑翔飛行器而言，其再入飛行段占整個返回段飛行時間的90%以上[13]。若不考慮再入段的時間協同，則末制導段初始狀態可能相差很大，而末制導段的時間調節范圍有限，難以保證最終協同效果。因此，有必要研究時間協同的多飛行器再入制導，使多個飛行器在同一時刻到達指定位置，從而為協同末制導提供良好的交班條件[14-15]。

上述針對末制導段的協同方法也難以直接運用到協同再入制導中。這些方法一般研究的是定常速度下的飛行器運動學。然而對于再入飛行，由于整個過程中空氣密度變化很大，速度變化范圍大，飛行器動力學變化劇烈，熱流、過載、動壓等過程約束嚴苛，僅僅考慮定常速度下的飛行器運動學難以滿足再入飛行需求。此外，再入過程中“黑障”區的存在以及各飛行器間距離較大超過通訊距離，使得基于通訊的各類協同制導方法也不再適用[16]。因此，必須考慮再入飛行的特殊性，設計區別于傳統協同末制導的時間協同再入制導律。

目前關于時間協同再入制導的研究較少。文獻[17]基于模型預測靜態規劃設計了協同再入制導律，但該方法需預知各飛行器的飛行時間，從而相應地改變發射時間以實現終端時間一致，難以稱為協同制導。文獻[13]首次提出了一種時間可控再入制導律，并基于多飛行器時間協調信息設計了協同再入策略。該方法通過BP神經網絡在線預測剩余飛行時間，但需針對特定飛行器進行大量離線訓練。通過改變航向角走廊寬度增加側向機動以改變飛行時間，但由于再入飛行器的側向機動能力較弱，在不改變縱向動力學的情況下該方法調節時間能力較弱。該文獻顯示其時間調整范圍僅為整個再入時間4%～5%，這就大大降低了該方法的實際意義，且較寬的航向角走廊可能導致飛行器最終錯過目標。此外，其協調時間的計算還依賴于彈間通訊，不適用于再入飛行。

再入制導方法一般分為標準軌跡制導和預測校正制導兩類[18]。由于飛行時間是在線狀態量，可優先考慮使用預測校正制導，通過某種方法在線預測剩余飛行時間，校正飛行軌跡，從而實現時間約束以及協同飛行。

基于以上分析，本文提出一種基于高度-速度剖面的時間協同制導律。首先在高度-速度剖面內設計了由兩個軌跡參數確定的參考軌跡，則剩余航程和飛行時間可表示為兩個軌跡參數的函數。通過在線預測剩余飛行航程和時間并校正兩個軌跡參數，實現了時間約束再入制導。在此基礎上針對多飛行器協同再入任務設計了協同策略。該策略無需離線訓練和彈間通訊且時間可控范圍更大，更加適用于實際的再入過程。仿真結果說明了本文方法的有效性。

1 時間協同再入問題描述

1.1 運動方程

假設地球為均勻圓球，且不考慮自轉，則多個高超聲速滑翔飛行器的三自由度動力學方程描述為[19]

(1)

式中：i代表第i個飛行器;r為地心距;V為飛行速度;g為重力加速度;λ和φ分別為地球經度和緯度;θ為航跡角;ψ為航向角;σ為傾側角；L和D分別為升力加速度和阻力加速度，其計算公式可表示為

(2)

式中：m為質量;S為特征面積;ρ=ρ0e-h/7 110為大氣密度，ρ0=1.225 kg/m3為海平面大氣密度，h=r-R0為高,R0為地球半徑;CL和CD分別為升力和阻力系數，是迎角α和速度V的函數。文獻[20]給出了基于非線性最小二乘辨識的再入飛行器解析氣動系數模型，具有較高精度:

(3)

式中：CLi、CDi(i=0,1,2,3)為辨識出的氣動系數，具體數值參見文獻[20]。

1.2 再入約束

再入過程約束包括熱流約束、過載約束、動壓約束和準平衡滑翔條件[19]

(4)

(5)

(6)

(7)

終端約束包括終端高度約束、速度約束及經、緯度約束：

(8)

式中：tf,i為第i個飛行器的終端飛行時間；rf、Vf、λf和φf分別為相應的終端約束值。對于時間協同再入問題，還有終端飛行時間約束：

tf,1=tf,2=…=tf,n=tf

(9)

式中：tf為預設的協同飛行時間。

2 時間約束再入制導律

2.1 飛行剖面設計

再入過程中考慮到初始段熱保護要求和后續航程要求，通常采用分段函數形式的迎角方案：

(10)

式中：α1為一個較大的迎角；α2為最大升阻比迎角；Va和Vb為給定的速度值。

在給定迎角方案后，即可將熱流、過載、動壓約束和準平衡滑翔條件轉換為高度-速度平面內的再入飛行走廊：

(11)

初始下降段飛行器高度較高，氣動力作用很弱，難以對軌跡進行有效控制，通常采用常值傾側角開環制導，并當滿足一定條件后轉入滑翔段[21-22]。

針對滑翔段縱向制導，為充分利用走廊寬度和飛行器的再入能力，本文設計了如下的3段多項式形式的解析剖面：

(12)

式中：kij(i=1,2;j=0,1,2,3)分別為多項式系數；V1,V2∈(Vf,Vtran)為滑翔段內選定的兩個速度值；(Vtran、Htran)為初始下降段與滑翔段間的過渡點；H1、H2分別為V1、V2處的選定高度且滿足：

Hi=Hmax(Vi)-k1(Hmax(Vi)-Hmin(Vi))

i=1,2

(13)

其中：Hmax(Vi)、Hmin(Vi)分別為飛行走廊上對應速度Vi處的高度上邊界和下邊界；k1∈(0,0.98)為待設計的高度系數。當k1增大時，中段直線會更靠近走廊下界，整個H-V剖面也隨之下降。

k1確定后，中段直線段即確定。另外兩段均為三次多項式，確定三次多項式參考軌跡需要4組約束條件。當V1≤V≤Vtran時，考慮H-V剖面內滑翔段初始點(Vtran,Htran)及中段直線段連接點(V1,H1)處連續性和光滑性要求即可確定。

當Vf≤V

H3=Hmax(V3)-k2(Hmax(V3)-Hmin(V3))

(14)

式中：k2∈(0,0.98)為第2個待設計高度系數。

于是確定了整個H-V剖面，如圖1所示。該剖面分為3段。其中第1段多項式用于將飛行器快速拉升，避免超出走廊下邊界。中段直線段通過調整高度系數k1可充分利用走廊寬度，以調節飛行器航程和飛行時間，同時保證軌跡在走廊以內。第3段多項式一方面滿足終端高度、速度約束，另一方面通過高度系數k2調整飛行寬度，從而調節飛行器航程和飛行時間。

將再入段的H-V剖面設計為解析多項式是較為傳統的方法。但傳統方法均只有一個軌跡參數以滿足航程約束，相比之下本文設計的剖面有兩個待設計的軌跡參數，能夠更加充分地利用走廊高度和改變飛行器縱向動力學，從而同時滿足射程和飛行時間約束。

圖1 H-V剖面內的縱向軌跡Fig.1 Longitudinal trajectory in H-V profile

2.2 剩余飛行航程和時間預測

本文在每次軌跡預測中，同時預測剩余航程和剩余飛行時間。對于剩余航程sgo，在假設飛行軌跡近似于大圓弧下有

(15)

根據運動方程式(1)可得

(16)

則剩余航程對速度的導數為

(17)

考慮到dt=-dtgo，則由式(16)，剩余飛行時間tgo對速度的導數為

(18)

由式(2)知，阻力加速度D僅與高度、速度有關，在H-V剖面內僅與V有關。

由運動方程式(1)取高度對速度的微分有

(19)

(20)

(21)

進一步考慮在初始、終端速度確定的情況下，所設計的H-V剖面僅與兩個高度系數k1、k2有關，于是式(20)、式(21)可寫為

sgo=f(k1,k2)

(22)

tgo=g(k1,k2)

(23)

式中：f(·)、g(·)分別為剩余航程和飛行時間，為兩個高度系數的函數。

圖2和圖3分別展示了以CAV-H飛行器為例，在Vf=2 000 m/s、Hf=25 km約束和迎角方案式(10)下，數值仿真得到的航程、飛行時間與高度系數k1、k2間的關系。從圖中可知，隨著高度系數減小，航程、飛行時間均單調增大。由式(20)、式(21)分析知，高度系數減小，則飛行軌跡越接近走廊上邊界，飛行高度增大，阻力加速度減小，航程與飛行時間均增大。理論分析與仿真結果一致。

圖2 航程與高度系數的關系Fig.2 Relationship between range and altitude coefficients

圖3 飛行時間與高度系數的關系Fig.3 Relationship between flight time and altitude coefficients

2.3 校正算法

傳統預測校正算法往往采用牛頓迭代法或割線法以單航程約束校正傾側角。本文在每個制導周期內，以航程、時間雙約束，基于H-V剖面同時對兩個高度系數進行校正：

(24)

式中：s1為理想剩余航程，可由當前位置與終端位置求出；t1為理想剩余時間，可由終端時間減去已飛時間得到。上述兩個變量在每個制導周期內均為已知量，則方程組式(24)為關于高度系數k1、k2的二元非線性方程組

(25)

式中：F(k1,k2)=sgo-s1,G(k1,k2)=tgo-t1。利用二元非線性方程組求根的牛頓迭代法可快速求解方程組式(25)[23]:

(26)

(27)

式中：Δ為一正小量。利用式(26)計算得到下一次迭代的k1、k2值，重復迭代，直到F<ε1且G<ε2時停止迭代，得到滿足精度的k1、k2的解。ε1、ε2為允許的航程和時間誤差。

實際仿真中發現，當飛行過程中存在參數擾動且其影響較大時，方程組式(25)可能得不到滿足精度的解。此時，取目標函數

(28)

式中：s0為航程誤差的歸一化參數；t0為時間誤差的歸一化參數。于是將航程、時間雙約束的預測校正問題，轉化為使目標函數式(28)最小的參數優化問題。利用參數優化的牛頓迭代法可快速求解出使目標函數最小的高度系數k1、k2[24]:

(29)

高度系數確定之后則參考H-V剖面確定。根據運動學方程式(1)，選取控制輸入為cosσ，對高度求二階微分，得到動力學系統為

(30)

式中：a=-Dsinθ-g+V2cos2θ/r;b=Lcosθ。設計控制律跟蹤參考高度Hd:

(31)

式中：ξ、ω分別為阻尼比和自然頻率。則控制輸入cosσ為

(32)

于是得到了傾側角的幅值。

2.4 側向制導

側向制導采用經典的航向角走廊方法確定傾側角符號。當航向角誤差超過預設的誤差走廊時，改變傾側角符號使其重新回到走廊內；當航向角誤差未超過誤差走廊時，保持傾側角符號不變。航向角走廊方法相比于一次或兩次翻轉方法，傾側角符號翻轉次數較多，但對在線參數擾動魯棒性更好。

本文中的時間約束再入制導律的核心思路是將傳統再入制導律中航程約束的單軌跡參數搜索問題擴展為航程、時間雙約束的雙軌跡參數搜索問題。這種思路既適用于離線的軌跡規劃，也適用于在線預測校正制導；既可基于高度-速度剖面，也可基于阻力加速度-能量剖面乃至傾側角剖面。軌跡參數既可選為兩個高度系數，也可選擇任意兩個可同時影響軌跡高度的參數。

本文沿用了傳統預測校正制導的一般方法，通過在線數值積分計算式(20)、式(21)，但必須考慮計算實時性與精度之間的矛盾。

1) 注意到式(20)、式(21)以速度為積分變量，積分終止條件為終端速度約束，于是積分步長可取為(V-Vf)/N，V為當前速度，N為積分步數，可取為100。即采用變步長積分，而積分總步數一定，當距離目標較遠時，速度相差較大，對制導精度要求較低，積分步長較大；當距離目標較近時，積分步長較小，以提高預測精度。

2) 對于制導周期，即預測校正周期，在再入開始時刻，距離目標較遠，對制導精度要求較低，可將制導周期設置為多個仿真步長。隨著飛行器距離目標越來越近，每次預測校正的耗時逐漸減小，可逐漸減小制導周期至每一步預測校正以保證終端精度。

3) 傳統數值積分多采用四階龍格庫塔法，每次積分需要4次右端函數計算。本文采用Adams預估校正法，每次積分僅需兩次右端函數計算，在精度相當的情況下計算量為四階龍格庫塔法的一半，從而大大減少了計算時間。

3 多飛行器協同策略

3.1 飛行時間可調范圍

在多飛行器協同再入之前，需要確定協同飛行時間，為此先確定各飛行器的飛行時間可調范圍。圖3通過數值仿真得到了飛行時間與高度系數k1、k2間的關系曲面，但此曲面是在無航程約束下的。對于某個確定的再入任務，其航程約束一定，此時飛行時間與高度系數間的關系退化為三維空間中的一條曲線。

圖4為在Vf=2 000 m/s、Hf=25 km、航程sf=8 700 km約束下數值仿真得到的飛行時間與高度系數間的關系。從圖中可知，隨著高度系數k2減小、k1增大，飛行時間單調增大。

圖4 航程約束下的飛行時間與高度系數關系Fig.4 Relationship between flight time and altitude coefficients under a certain range constraint

將H-V剖面的第3段多項式軌跡稱為滑翔后半段，前兩段多項式軌跡稱為滑翔前半段。分析剩余航程與飛行時間的表達式(20)與式(21)知，兩者被積函數分母相同，但剩余航程的被積函數分子為Vcosθ，而剩余時間的分子為1。這意味著在相同分母下，剩余航程受速度的影響比剩余時間更大。由于滑翔前半段速度顯著大于滑翔后半段速度，可推測剩余航程受滑翔前半段影響較大，剩余時間則主要由滑翔后半段決定。因此，隨著滑翔后半段高度增大，高度系數k2減小，阻力加速度減小，飛行時間增大。同時，為使總航程保持一定，前半段高度系數k1增大。理論分析與圖4中的仿真結果一致。

于是，對于航程確定的再入任務，其最大飛行時間對應高度系數k2最小，即在滑翔后半段軌跡與走廊上邊界相切時得到；最小飛行時間對應高度系數k2最大，即在滑翔后半段軌跡與走廊下邊界相切時得到。據此思路，可在飛行之前求出兩條相切軌跡，求出對應飛行時間，作為飛行時間可調范圍。

圖5為在Vf=2 000 m/s、Hf=25 km、sf=8 700 km約束下的對應于最大、最小飛行時間的兩條軌跡，最大飛行時間與最小飛行時間相差227 s，時間可調范圍約占總飛行時間的15%，大大高于文獻[13]中的4%～5%。由于充分利用了縱向剖面，本文方法的時間調節范圍較大。

圖5 對應于最大、最小飛行時間的兩條軌跡Fig.5 Two trajectories corresponding to maximum and minimum flight times

3.2 協同策略

協同再入的目的是使多飛行器在同一時刻到達指定的末制導交班點，為此需要先確定協同飛行時間。為保證協同飛行時間存在，即各飛行器時間可調范圍有交集，應使再入段各飛行器的初始待飛航程大致相當。若各飛行器的初始待飛航程相差太大，無法實現協同再入。而高超聲速滑翔飛行器整個飛行過程分為助推段、再入段和末制導段[13]，再入段的初始條件是助推段的終端條件。于是在發射之前，通過設計助推段程序方案以保證再入段待飛航程大致相當，從而確保協同飛行時間存在以及協同再入的可行性。于是可求取如式(33)的協同飛行時間:

(33)

式中：tmax,i為第i個飛行器的最大飛行時間；tmin,i為第i個飛行器的最小飛行時間。

再入任務開始后，初始下降段各飛行器氣動力和機動能力較弱，采用常值傾側角開環制導，僅需將運動方向對準目標，不需考慮時間約束。在所有飛行器進入滑翔段后，氣動力和機動能力增強，于是各飛行器獨立地以tf為時間約束執行本地時間約束再入制導律，通過在線預測剩余航程和時間，實時校正飛行軌跡，最終實現時間協同再入任務。

考慮到再入過程中“黑障”區的存在以及各飛行器間距離較大，大大超過通訊距離，該協同策略避免了集中式或分布式通訊結構的使用，更加適用于實際的再入過程。

4 仿真分析

仿真以美國通用航空器CAV-H為對象[25]。該飛行器質量為907 kg，參考面積為0.483 9 m2。迎角方案參數:α1=20°,α2=10°,Va=6 500 m/s,Vb=5 000 m/s。高度-速度剖面參數:V1=7 000 m/s,V2=5 400 m/s,V3=3 000 m/s。歸一化參數s0=5 km,t0=5 s?？刂坡蓞郸?0.7,ω=0.1。航向角誤差門限選取為10°。過程約束選取為

本文分別在標稱條件和擾動條件下針對時間約束再入制導律進行仿真驗證，然后進行多飛行器協同再入仿真以驗證協同策略的有效性。

4.1 標稱條件下多任務仿真

首先在標稱條件下針對不同航程、時間約束約束下的再入任務驗證時間約束再入制導律的有效性。4個不同任務的初始、終端條件設置見表1。以終端速度為仿真截止條件計算終端誤差，終端誤差見表2。

圖6(a)為標稱條件下的高度-速度曲線。4個任務下的軌跡均比較平滑，無明顯振蕩，基本滿足準平衡滑翔條件，并能到達指定的終端速度高度。并且隨著航程的增大，H-V曲線逐漸抬高，始終保持在再入走廊內，滿足過程約束。圖6(b) 為標稱條件下的地面軌跡。表2為4個任務下的終端誤差?；诮馕鯤-V剖面的設計方法使得終端速度、高度、航程誤差均較小，終端時間誤差也保持在0.4 s以內，從而驗證了時間約束再入制導律對再入時間的可控性，為實現多飛行器協同再入奠定了基礎。

表1 再入任務Table 1 Entry mission

表2 終端誤差(標稱條件)Table 2 Terminal errors (standard conditions)

圖6 標稱條件下的高度-速度曲線和地面軌跡Fig.6 Height-velocity profile and ground tracks in standard conditions

4.2 擾動條件下時間約束仿真

為了驗證時間約束再入制導律在擾動條件下的精度和魯棒性，針對表2中任務2進行200次蒙特卡羅仿真。仿真中各初始狀態偏差和參數偏差假設符合正態分布，其偏差限如表3所示。

圖7(a)為擾動條件下的高度-速度曲線。擾動條件下的軌跡比較平滑，基本滿足準平衡滑翔條件，且均分布在再入走廊內，保證滿足過程約束。圖7(b)中的終端高度-速度誤差散布顯示終端高度誤差基本在400 m以內，速度誤差在0.3 m/s 以內，具有一定的精度。圖7(c)為擾動條件下的地面軌跡。圖7(d)為擾動條件下的落點散布。落點沿飛行方向散布，且散布誤差基本在10 km以內，滿足再入制導要求。圖8為終端飛行時間誤差直方圖。飛行時間誤差分布在3.5 s以內，且分布范圍比文獻[13]中小，從而驗證了時間約束再入制導律在擾動條件下對飛行時間的可控性。