多禁飛區(qū)高超聲速滑翔飛行器再入機動制導

高 楊,蔡光斌,2,張勝修,徐 慧

(1.火箭軍工程大學 導彈工程學院, 西安 710025; 2.西北工業(yè)大學 航天學院, 西安 710072)

近年來，高超聲速滑翔飛行器的再入機動制導研究逐漸成為世界各大國航天工業(yè)發(fā)展的一個熱點[1-5]。高超聲速滑翔飛行器再入飛行必須滿足復雜的約束條件，其中包括傳統(tǒng)的熱流率、動壓、過載約束和終端約束，也包括自然、軍事等因素造成的復雜禁飛區(qū)約束[4]。飛行器的禁飛區(qū)是一種路徑約束，在飛行器最優(yōu)軌跡的解算中，如果禁飛區(qū)的數(shù)量越多或者模型比較復雜，那么飛行器在軌跡優(yōu)化的路徑約束就越多，優(yōu)化軌跡的解算難度就越大?，F(xiàn)今研究禁飛區(qū)突防的成果較少，且主要集中在離線求取最優(yōu)解，其中，偽譜法是數(shù)值求解方法的代表，采用不同的配點和多項式插值方法達到全局優(yōu)化的目的，這些方法無法應用到在線的機動制導中，尤其是在禁飛區(qū)信息較復雜或者信息缺失的情況下[6]。因而，預測反饋的制導機制在飛行中作用明顯[7]，一些基于誤差校正的神經網絡與系統(tǒng)辨識方法被廣泛應用[8-11]，然而這些研究中對于未知形狀的禁飛區(qū)約束條件研究并不成熟，當預測數(shù)據(jù)因為禁飛區(qū)等因素急劇變化時，會導致反饋誤差變大而導致制導錯誤。部分學者在研究昆蟲生物學觸角反饋機制的基礎上，使用較少的“觸角”探測前方未知區(qū)域并進行反饋，從而制定制導方案[12-14]。文獻[12]在移動機器人路徑規(guī)劃研究中，首先提出了基于觸角的方法，其中觸角代表了機器人可能的運動方向。文獻[14]使飛行器連續(xù)產生左右雙觸角對前方進行探測，通過跟蹤標準軌跡并使用側滑角正負反轉策略規(guī)避不同形狀的禁飛區(qū)。

鑒于現(xiàn)有研究結果中存在問題，本文改進了軌跡預測時觸角的設計與側滑角的瞬變策略，引入了航向角誤差限制條件和延時計數(shù)開關，計算了規(guī)避策略的優(yōu)先級，軌跡能夠順利規(guī)避禁飛區(qū)并到達終點。通過魯棒性仿真驗證了該策略預測過程計算時間足夠小，滿足機動制導策略的要求。

1 再入模型構建

高超聲速滑翔飛行器的動力學模型以及再入飛行過程中的約束條件模型是高超聲速滑翔飛行器機動制導過程的基礎，下面對這些模型與條件分別進行定義。

1.1 高超聲速滑翔飛行器再入運動學建模

針對高超聲速滑翔飛行器在再入飛行階段的機動制導過程，如圖1所示，在地心赤道旋轉坐標系內，忽略地球自轉角的影響，建立高超聲速滑翔飛行器的動力學方程[20]：

其中，r是高超聲速滑翔飛行器地心距，V是飛行器的地球相對速度，ψ與γ是飛行器的航向角與航跡角，θ和φ分別是飛行器所處的經度緯度，m是飛行器的質量，g是重力加速度，σ是飛行器的側滑角。D=ρV2SrefCD/2與L=ρV2SrefCL/2是飛行器在飛行過程中氣動阻力與升力，其中ρ是空氣密度，Sref是飛行器的參考橫截面積，CL與CD分別是與飛行器攻角有關的空氣動力學參數(shù)。不同的飛行器由于設計的不同，導致升力阻力的空氣動力學參數(shù)不同以及升阻比的不同。

圖1 地心赤道旋轉坐標系下飛行器動力學方程建模

1.2 常規(guī)約束條件建模

在臨近空間的長距離、無動力飛行過程中，高超聲速滑翔飛行器會受到多種約束的影響，這些因素都會影響制導的結果。本文橫向制導策略在預測的過程中，每條預測“觸角”上預測時的每一個積分步長端點為一個飛行器約束檢測點，預測計算時每到一個檢測點都要檢測是否遵循常規(guī)約束條件：

一是起始狀態(tài)與終止狀態(tài)的端點約束，它是對飛行器起始與終止位置、速度和姿態(tài)的約束，也是飛行器不同狀態(tài)轉換的交接班條件，起始條件通常為給定的任務起始狀態(tài)[19]，令X=[r,θ,φ,V,γ,ψ]為飛行器狀態(tài)矩陣，εX為較小的常值矩陣，終端約束模型建立如下：

|X(tf)-Xf|≤εX

二是飛行器的常規(guī)路徑約束，分別是飛行器的熱流率約束、過載約束和動壓約束[10]，建立模型分別為

Q=KQρ0.5V3.15≤Qmax

三是航向角約束。在高超聲速滑翔飛行器因禁飛區(qū)等限制原因機動飛行的過程中，為了保證飛行器轉彎規(guī)避禁飛區(qū)后仍能到達目的地，需要對飛行過程中的航向角進行約束。這種約束使得飛行器能夠有足夠的距離通過對側滑角和攻角的控制，把飛行器制導到終點位置。首先將航向角約束定義為在軌跡檢測點上的飛行器實時航向與目標地點在經緯度平面的夾角。圖2中先求出要研究的飛行器的最大轉彎軌跡，在該軌跡上求出各檢測點Pi(i=1,2,…,m)的速度方向與終點連線在經緯度平面上的夾角，為速度航向角誤差限制Ψi(i=1,2,…,m)，即

Ψi=ψi-arctan((θi-θf)/(φi-ff))

式中，[θi,φi,ψi]為檢測點Pi位置的經緯度坐標和航向角，[θf,φf]是終點的經緯度坐標。

圖3中是兩種飛行器Shuttle與Common Aero Vehicle-HPMARV(CAV_H)根據(jù)最大轉彎軌跡計算的航向角約束以及在文獻[4，14]中按照經驗設定的對CAV-H飛行器的航向角約束。

圖2 最大轉彎軌跡與檢測點上的航向角

圖3 速度-航向角約束變化

1.3 多禁飛區(qū)約束建模

在飛行器實際飛行中禁飛區(qū)的形狀是多樣的，用解析式來表示很難，而用多個圓柱形禁飛區(qū)組合表示是比較方便的做法。假設這些禁飛區(qū)都是如圖4所示的圓柱形的區(qū)域，在俯視視角下如圖5所示，高超聲速滑翔飛行器面對的多個小型禁飛區(qū)的復雜約束情形以及各禁飛區(qū)與飛行器的距離。

圖4 圓柱體禁飛區(qū)示意圖

圖5 高超聲速滑翔飛行器多禁飛區(qū)約束情形

將圓柱形禁飛區(qū)約束簡化為在經緯度上的二維約束，針對n個圓柱形禁飛區(qū)中的第i(i=1,2,…,n)個禁飛區(qū)，令禁飛區(qū)的中心經緯度坐標為(θi,φi)，半徑為Ri，其在經緯度上的禁飛區(qū)邊界表示為

飛行器與各禁飛區(qū)中心之間的距離定義為{d1,d2,…,dn}，當di≤Ri(i=1,2,…,n)時，飛行器被判定進入第i個禁飛區(qū)，飛行任務失敗。在未知禁飛區(qū)具體形狀的情況下，飛行器只能探測到禁飛區(qū)的邊界。

2 基于“觸角”預測的側滑角瞬變機動制導策略

2.1 基于線性二次調節(jié)器的縱向制導方法

在機動制導過程中，制導算法分為縱向和橫向兩個制導方法，縱向平面跟蹤方法通過線性二次調節(jié)器(Linear Quadratic Regulator,LQR)實現(xiàn)[17-19]，即

其中，Δα和Δσ為攻角和側滑角的調整量，[Δr,ΔV,Δγ]為實際軌跡與標準軌跡的地心距、速度和航跡角差值，因為是縱向平面的軌跡跟蹤，所以狀態(tài)量上只選擇了這3個關系密切的變量。

橫向制導方法是一種側滑角符號轉換方法，它與縱向制導策略同時進行，橫向制導策略的周期一般大于縱向制導策略周期，使側滑角在飛行器縱向制導的跟蹤過程中保持不變，即

ΔTLat=kΔT

其中，ΔTLat與ΔT分別是橫向和縱向制導策略的計算周期，k取整數(shù)，k值越大，橫向機動策略計算的次數(shù)就越少，總機動策略計算的平均時間就越小，k值越小，對禁飛區(qū)的探測效果越好。本文中兩個計算周期都取1 s。

2.2 基于動力學方程的三“觸角”探測模型

在對高超聲速滑翔飛行器機動制導的過程中，只有最開始解算出的標準軌跡是已知的，而前方飛行通道的實時禁飛狀況由于不能提前探知，需要飛行器采取預測校正的制導方法進行探測。本文采取的基于“觸手”的方法需要對前方發(fā)出3個觸手來預測軌跡的飛行狀態(tài)。這種預測軌跡是由計算機根據(jù)飛行器的動力學方程積分產生的，其中積分的步長越小，積分產生的預測路徑越符合實際軌跡。而高超聲速滑翔飛行器動力學方程較為復雜，預測軌跡的積分速度較慢，因而產生的預測軌跡越少越好。在文獻[14]中，每一個預測周期產生兩條“觸角”軌跡線，極大地縮減了計算量，但是缺少了直線前進的“觸角”預測軌跡線，導致在可以無機動地直線飛行時，采取了左右擺動的復雜制導策略，側滑角正負瞬變對控制系統(tǒng)造成的壓力也較大。

本文采取的是前向“觸角”與左右兩側“觸角”結合的三“觸角”預測方式，如圖6，當飛行器在前向觸角遇到禁飛區(qū)時，即飛行段OA段，采取三觸角與雙觸角的反饋結果是基本相同的，飛行器只能在兩側的觸角中選擇其中的較優(yōu)路線進行機動。當飛行器在前向觸角沒有遇到禁飛區(qū)而停止時，即飛行段AB段，飛行器需要判斷前向觸角和另外兩側觸角的反饋信息，文獻[14]中提到的方法，在類似A點的飛行狀態(tài)時，由于左右觸角給出的反饋信息的實時變化，飛行器會進行左右搖擺的機動，直到因為航向角約束的逐漸減小而最后使軌跡逼近終點；而采取類似于B點的三觸角探測策略，飛行器的前向觸角在AB段沒有禁飛區(qū)約束，因而觸角的終點更貼近終點，飛行器在經過短暫的左右機動調整后可以沿直線飛行到終點。左右兩側觸角預測軌跡的側滑角符號一正一負，數(shù)值上使用飛行器的最大側滑角度數(shù)。如果使用度數(shù)較小的側滑角進行兩側軌跡的預測，大多時候會因為機動角度過小而無法使飛行器飛躍禁飛區(qū)；如果在原有三條“觸角”的基礎上，再加兩條小角度的預測線，會增加計算的時間，為橫向制導策略的選擇增加困難，同時因為高超聲速滑翔飛行器高靈敏的特性而增加控制的難度。

圖6 三“觸角”預測示意圖

飛行器每一個預測周期的三個“觸角”的起點是相同的，它們都以飛行器當前時刻仿真積分一個時間步長的狀態(tài)為預測起點，當預測結束并且獲得相應的控制指令后，飛行器也剛剛飛至下一時刻，由于積分時間步長固定且較小，預測起點與實際狀態(tài)誤差較小。三個“觸角”的終止條件如下：

C1：積分時間超過了理論上飛行的最大任務時間，本文采取的最大任務時間是3 000 s；

C2：“觸角”的末端觸碰到了禁飛區(qū)的邊緣或者因為積分步長的原因進入了禁飛區(qū)；

C3：預測軌跡線上的航向角超出了該速度所對應的最大航向角約束；

C4：預測軌跡順利到達終點位置。

其中，前3個條件屬于禁忌類條件，飛行器應避免飛行到這三類“觸角”的終點，第4個終止條件屬于允許類條件，飛行器應優(yōu)先選擇這樣終止的“觸手”制導飛行。其中的禁飛區(qū)終止條件是禁忌條件中禁忌程度最高的，而時間終止條件的成立幾率較小，由此得到終止條件所對應的“觸角”選擇優(yōu)先級Ka，如表1所示。

表1 終止條件優(yōu)先級

3條軌跡預測只能以這4種情況結束，并將終止時的狀態(tài)和終止原因反饋給橫向制導方法，由橫向制導方法選擇三“觸角”的其中一條進行機動制導。

2.3 基于側滑角瞬變的橫向制導方法

文獻[4]和文獻[14]等提到的側滑角逆變策略是一種提供側滑角符號的方法，側滑角的絕對值是基本不變的，其中文獻[14]在這一基礎上根據(jù)LQR方法其側滑角的數(shù)值進行了微調，但都屬于側滑角逆變策略，即在機動時將飛行器的側滑角變號，改變飛行器的空氣動力方向。本文將側滑角歸零這一選擇加入到側滑角的瞬變策略中，使得飛行器能夠在不需要機動的時候選擇前向飛行，或者在側滑角進行正負轉換時，作為穩(wěn)定控制系統(tǒng)的過渡狀態(tài)，因此飛行器的側滑角度數(shù)有3種可能，即

σ=sgn(σ)·σmax,sgn(σ)∈{-1,0,1}

其中，σmax是側滑角最大值，左右觸角分別對應符號函數(shù)sgn(σ)的-1與1。符號函數(shù)sgn(σ)的賦值由觸角反饋的總優(yōu)先級K決定?？們?yōu)先級K由終止條件優(yōu)先級Ka和觸角終止點與終點的距離優(yōu)先級Kb共同決定。因此，每一個觸角的總優(yōu)先級如下式所示：

Ki=-Ka·Kb,i=Left,Mid,Right

其中，Kb根據(jù)觸角終止位置的經緯度坐標與飛行路徑終點經緯度坐標的距離排序，距離越小，該觸角的終止位置越靠近終點，該觸角越應該被選擇，即根據(jù)這一距離的遠中近3個等級，Kb分別賦值為{1.1,1.0,0.9}。式中Ki的優(yōu)先級越大，代表該觸角更應該被選擇作為下一個時刻的制導指令，即：

為了避免飛行器側滑角瞬變角度過大或者因為優(yōu)先級的振蕩變化造成的側滑角振蕩變化而控制系統(tǒng)不能完成這樣的變化，加入一種延時鎖定開關。如圖 7所示，由側滑角符號函數(shù)計算出來的結果進入延時鎖定開關后，如果開關關閉，延時計數(shù)開關開始計數(shù)，控制系統(tǒng)中的側滑角不變；如果開關打開，則重新判斷。

圖7 側滑角計算流程框圖

當前時刻為j時刻，預測周期為ΔTLat，則控制系統(tǒng)獲得的側滑角σ′由下式計算得到：

其中，σj、σj-ΔTLat與σj+ΔTLat代表j、j-ΔTLat與j+ΔTLat三個時刻的側滑角符號函數(shù)輸入值，延時開關的延時時間最小值根據(jù)飛行器自身特性來決定，本文中的延時時間統(tǒng)一取10 s。在延時過程中，統(tǒng)計式(17)中出現(xiàn)的側滑角同類振蕩變化次數(shù)，振蕩次數(shù)超過一定限制，則打開延時開關，并將σj+ΔTLat輸入控制系統(tǒng)。

3 仿真算例與分析

本節(jié)采用CAV-H和Shuttle的高超聲速滑翔飛行器模型在多禁飛區(qū)條件下驗證基于三觸手的側滑角瞬變機動制導策略。飛行器Shuttle和CAV-H的模型數(shù)據(jù)分別來源于文獻[19]和文獻[20]，表2是兩個飛行器標準軌跡的飛行數(shù)據(jù)[14,19]：

根據(jù)表2初始數(shù)據(jù)獲得兩飛行器的初始標準軌跡高度曲線，如圖8所示。本文兩飛行器的標準軌跡中側滑角設為0，攻角的變化幅度較小，導致飛行時間較長。

3.1 未知多禁飛區(qū)情況下的機動制導仿真

為檢測機動制導算法的有效性，設置了11個半徑為1°的小型圓柱禁飛區(qū)，它們的中心經緯度坐標如表3所示，組成了復雜的多禁飛區(qū)模型。

表2 飛行器飛行數(shù)據(jù)

圖8 標準軌跡高度曲線

序號經度緯度序號經度緯度130-3765-1.523038651.5350-3980345001080055031180-36400

根據(jù)表3的禁飛區(qū)信息，按照圖9的算法流程采取兩種方法求解并比較。一種是利用GPOPS II工具包求解，該求解方法是在離線狀態(tài)下進行的，因此使用已知禁飛區(qū)具體信息求解，同時側滑角在最小值與最大值構成的區(qū)間內自由取值。另一方面采取本文的機動制導策略進行在線求解，求解過程只能通過觸角的探測，反饋給制導策略禁飛區(qū)部分輪廓信息，側滑角的取值只有最大值、最小值和零。如圖10和圖11所示，基于“觸角”的方法與基于GPOPS II的方法都可以使兩種飛行器順利規(guī)避禁飛區(qū)到達終點，同時，如圖12和圖13所示，由于基于“觸角”的方法的側滑角的選取區(qū)間相對簡單，路徑圖中的軌跡相對平滑。

圖9 仿真算法流程框圖

圖10 CAV-H飛行器規(guī)避禁飛區(qū)路徑

圖11 Shuttle飛行器規(guī)避禁飛區(qū)路徑

圖12 CAV-H控制變量變化

圖14和圖15顯示，基于“觸角”的機動制導策略在規(guī)避禁飛區(qū)的同時，由于速度和高度嚴格跟蹤標準軌跡，雖然在算法中觸角終止條件沒有判別是否符合3個經典路徑約束條件，機動制導的結果依舊滿足約束條件。

圖14 CAV-H熱流率、動壓、過載變化

圖15 Shuttle熱流率、動壓、過載變化

由圖16中可以看出：飛行過程中的單次預測周期計算時間是呈下降趨勢的，其中的大小變化主要是由預測時飛行器位置決定的，飛行器距離禁飛區(qū)、航向角約束或終點越近，計算的時間就越短。

圖16 飛行過程中單次預測周期計算時間

3.2 模型的魯棒性與計算時間

為了驗證本文基于三“觸角”的側滑角瞬變機動制導模型，對兩個飛行器進行魯棒性仿真實驗。采用蒙特卡洛方法，對初始狀態(tài)中的緯度、航跡角、航向角以及飛行器的質量、面積實施3°、3°、10°和5%、5%的離散，經仿真獲得了每個飛行器的1 000次結果，從圖17和圖18中可以看出雖然有一些仿真路徑的機動幅度較大，但是所有仿真軌跡均能順利規(guī)避本文所設計的復雜禁飛區(qū)到達各自的標準軌跡終點。

圖17 CAV-H蒙特卡洛仿真結果

圖18 Shuttle蒙特卡洛仿真結果

本文進行實驗時使用的計算機CPU主頻為2.6 GHz，軟件為Matlab。在仿真實驗中，統(tǒng)計每一個仿真實驗的單次預測周期計算時間的最大值、最小值、平均值，得到圖19和圖20，可以看出，本文提出的機動制導策略在預測時最大預測時間不超過0.25 s，而且平均預測時間集中在0.01 s，完全符合機動制導的時間要求。同時，圖21與圖22展示了仿真中計算得到的熱流率、動壓和過載最大值均滿足表2的約束條件。

圖20 Shuttle單預測周期計算時間統(tǒng)計圖

圖21 CAV-H熱流率、動壓、過載最大值點集

圖22 Shuttle熱流率、動壓、過載最大值點集

4 結論

本文通過對基于“觸角”的預測校正制導方法的改進，綜合航向角約束和側滑角延時計數(shù)開關設計，提出了基于三“觸角”的側滑角瞬變機動制導策略，解決了高超聲速滑翔飛行器在多禁飛區(qū)信息未知情況下再入機動制導的問題。通過仿真驗證了該方法能夠有效解算出多禁飛區(qū)情況下高超聲速滑翔飛行器機動制導的在線飛行方案，通過蒙特卡洛仿真實驗驗證，突出體現(xiàn)了該方法計算時間上能夠滿足機動制導的要求。