張合新，宮梓豐，蔡光斌，宋 睿

(1.火箭軍工程大學 導彈工程學院，西安 710025； 2.青州士官學院， 山東 青州 262500)

近年來，有關復雜約束的高超聲速滑翔飛行器軌跡優(yōu)化已經成了國內外研究的熱點[1]。高超聲速飛行器實際再入的過程中，不僅要考慮大氣壓強、熱流量、熱流密度等約束限制[2]，還要對再入環(huán)境狀態(tài)(自然條件和政治軍事因素)進行避障與突防[3-4]?，F(xiàn)有的制導控制方法可以解決再入過程中的問題，而對于高超聲速飛行器而言，要實現(xiàn)戰(zhàn)略突防[5]。現(xiàn)今研究對禁飛區(qū)與突防的研究成果相對較淺[6]，是軌跡優(yōu)化領域值得深入研究的課題。軌跡優(yōu)化問題具有較高的非線性和復雜度，常規(guī)的變分法和極大值原理難以在規(guī)定時間內高效求解出最優(yōu)解[7]，因而通常用數(shù)值法進行求解。偽譜法[8]是數(shù)值法的典型代表方法，大多數(shù)配點基于高斯積分規(guī)則[9]，采用契比雪夫或勒讓德多項式，而后對全局進行插值[10]。文獻[11]通過構建高超聲速飛行器物理模型，考慮飛行任務中的過程約束和終端約束，利用偽譜法解決了軌跡優(yōu)化問題；文獻[12]進一步討論了偽譜法再入軌跡制導的應用，但針對禁飛區(qū)約束條件沒有作出細化考慮。針對改進自適應偽譜法，對最優(yōu)控制問題進行轉化，即將整個區(qū)間進行網格細化[13]，避免出現(xiàn)極其精細的網格化和不合理的高逼近多項式，并綜合網格密度和插值多項式的階次，使系統(tǒng)滿足較快的收斂速度，提出最大軌跡空間概念，軌跡能夠完整繞開禁飛區(qū)，且提升了原有精度。通過仿真計算結果計算狀態(tài)量，評估軌跡狀態(tài)性能得出結論。

1 再入模型構建

復雜約束下的高超聲速滑翔飛行器再入軌跡優(yōu)化與飛行器物理結構、過程及終端約束和目標方程息息相關，下面就對模型及約束分別進行定義。

1.1 滑翔飛行器再入運動學模型

高超聲速滑翔飛行器再入過程受到很多約束的限制，信號通道間的耦合問題也較為突出，且非線性特征明顯，其動態(tài)變化可由一組非線性微分方程表示[9]。假定地球為均勻光滑的球體，則滑翔飛行器再入過程中的位置變化如下：

(1)

式(1)中，x、y、z分別是飛行器在航跡坐標系下三個軸所對應的位移分量；v是飛行器的線速度；ζ為航跡傾角；ξ為航跡方位角。飛行器的速度、航跡方位角和航跡傾角的一階微分滿足以下關系：

(2)

式(2)中，Q、S、g、m和σ分別代表飛行器的動壓、再入?yún)⒖济婷娣e、重力加速度、飛行器質量和傾斜角(bank angle)；KD、KL、KY分別是飛行器的阻力系數(shù)、升力系數(shù)和側向力系數(shù)。動壓Q滿足公式：

(3)

式(3)中，ρ表示大氣密度。其大小會隨飛行高度的變化而發(fā)生改變，一般用如下的公式進行擬合：

(4)

式(4)中:ρ0為海平面標準大氣密度；re為地球半徑；H為常值。

1.2 約束條件

為進一步模擬實際飛行條件，在軌跡優(yōu)化的過程中引入禁飛區(qū)約束。禁飛區(qū)即飛行器不允許飛入的空域，譬如防空識別區(qū)、導彈攔截區(qū)等，在設計軌跡的時候必須要繞過禁飛區(qū)。

根據(jù)需要，確定禁飛區(qū)路徑約束為：

C(x(t),u(t),t)≤0,t∈[t0,tf]

(5)

通常情況，熱流、動壓和過載被認為是路徑約束，以不等式的形式給出，比如說熱流約束表示為：

(6)

式(6)中，RN為滑翔飛行器頭部熱流駐點半徑；ρ0為海平面標準大氣密度；Vc=7.8×103m/s為地球第一宇宙速度；K1為常數(shù)。

過載約束的表達為：

(7)

動壓約束的表達式為：

(8)

1.3 禁飛區(qū)參數(shù)模型構建

禁飛區(qū)的形狀不盡相同，可以分為規(guī)則形狀禁飛區(qū)和不規(guī)則形狀禁飛區(qū)。常見的典型規(guī)則形狀禁飛區(qū)包括圓柱形、球形(半球形)、橢圓柱形和橢球形(半橢球形)，不規(guī)則形狀禁飛區(qū)可以通過多種、多個典型規(guī)則形狀的禁飛區(qū)形狀的組合來表征。下面通過具體分析，分別建立參數(shù)化數(shù)學模型來實現(xiàn)多種禁飛區(qū)覆蓋范圍的描述。

1.3.1 球體、橢球體禁飛區(qū)

在進行禁飛區(qū)坐標表示時，對于橢圓柱形禁飛區(qū)，首先確定禁飛區(qū)的中心點位置N的坐標，該點一般位于地球表面，一般以該點在地心赤道坐標系下的坐標(rN,λN,φN)為準。再以N點為原點建立禁飛區(qū)直角坐標系，其中Nx軸沿著橢圓柱面的短半軸方向，Ny軸垂直于水平面指向天，Nz軸沿著橢圓柱面的長半軸方向，Nx軸與正北方向的夾角為αN。參照文獻[13]，禁飛區(qū)中心點在地心赤道坐標系下的坐標可以表示為：

(9)

式(9)中，re為地球半徑；(r,λ,φ)為飛行器在地心求坐標系下的坐標。

由禁飛區(qū)坐標系到地心赤道坐標系的變換矩陣為：

F=F2[-(90°+αN)]F1[φN]F3[λN-90°]

(10)

其中

(11)

(12)

(13)

因此，飛行器在禁飛區(qū)坐標系下的坐標(xN,yN,zN)可以表示為：

(14)

由此，橢球體禁飛區(qū)參數(shù)坐標可以表示為：

(15)

式(14)中，a、b分別為橢球赤道面的長半軸和短半軸長度；c為橢球的極半徑長度。圖1為橢球體禁飛區(qū)示意圖。

圖1 橢球體禁飛區(qū)示意圖

當a=b=c時，橢球球心至橢球表面各點距離相等，此時的橢球體禁飛區(qū)退化成球體禁飛區(qū)，其半徑r=a=b=c。

1.3.2 圓柱、橢圓柱體禁飛區(qū)

令式(14)中的橢球體極半徑長度c=∞，此時的橢球體就變成了橢圓柱體，可以表達為：

(16)

橢圓柱體切面的長半軸、短半軸分別為a、b，分別與禁飛區(qū)坐標系Nxf、Nzf軸重合。若禁飛區(qū)有具體的高度界限，則可通過限定y的取值范圍來實現(xiàn)表示。特殊地，當橢圓柱體的長半軸和短半軸相等時，即a=b，此時的橢圓柱體即是圓柱體，其切面半徑r=a=b。圖2為橢圓柱體禁飛區(qū)示意圖。

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

禁飛區(qū)形狀具有多樣化的特點，而其他形狀的禁飛區(qū)可以類比以上兩種情況得出相關表達式，此處不再贅述。

飛行器再入時，其軌跡不得進入禁飛區(qū)內部，即飛行器與禁飛區(qū)的最短距離必須保持為正值，設L(x(t),t)為飛行器到禁飛區(qū)的最短距離，則禁飛區(qū)約束可以通過代數(shù)表達的方式計入計算，即使得L(x(t),t)≥εn。

2 自適應偽譜法解決最優(yōu)控制問題

2.1 自適應網格規(guī)劃策略

在判斷是否需要細化網格區(qū)間時，要確定誤差評估準則。在自適應偽譜法評估的基礎上，將動態(tài)約束和路徑約束方程在配點間的契合程度作為誤差評估標準[14]。

設[tk-1,tk]為第k個子區(qū)間，區(qū)間內共有Nk個配點數(shù)。取相鄰配點之間的等間距分布的s個點作為誤差評估采樣點，有：

(17)

式(17)中，i=1,…,N(k)；n是子區(qū)間內配點間評估采樣點的序數(shù)。評估采樣點的數(shù)量選取可以隨機確定，但不宜過多，一般取為2～4，因為過多的采樣點會造成計算代價的增加，而傳統(tǒng)的hp自適應偽譜法在采樣時僅僅取相鄰配點的中點作為采樣點，其誤差評判規(guī)則不夠精確。

(18)

2.2 考慮禁飛區(qū)的再入軌跡優(yōu)化求解策略

高超聲速飛行器再入軌跡的設置需要考慮復雜約束，而軌跡優(yōu)化問題經過自適應偽譜法的離散化后，最終軌跡的狀態(tài)變量和控制變量也變成了隨配點排布的離散量。而在復雜約束中，比如禁飛區(qū)約束，是針對連續(xù)的軌跡進行約束的。因此，經過自適應偽譜法解算出的最終結果，只能保證軌跡中的每個配點的狀態(tài)變量滿足復雜約束，而配點間的軌跡卻可能不滿足題設的約束條件。2.1節(jié)針對自適應偽譜法進行了改進，加入多誤差評估點插值和網格再細化的評判決策，最優(yōu)軌跡在精度上有了較好的進步，但改進算法仍然不能完全解決禁飛區(qū)連續(xù)性約束的問題。如圖3所示，誤差評估點之間的軌跡可能會越過禁飛區(qū)邊界，使得軌跡優(yōu)化問題不滿足禁飛區(qū)約束。

圖3 配點間不滿足禁飛區(qū)約束的情況

為了使得配點之間的軌跡滿足禁飛區(qū)約束，防止出現(xiàn)圖2的情況，就必須對禁飛區(qū)周圍的配點進行進一步優(yōu)化。可以利用2.1節(jié)，調節(jié)自適應網格細化策略的誤差門限值ε來實現(xiàn)，下面列出具體方案。

取L(x(t),t)最小的區(qū)間進行條件驗證。2.1節(jié)中，自適應偽譜法網格細化過程中，定義區(qū)間內的最大曲率向量為ρ，對應的，在L(x(t),t)上的最大曲率為ρl，則L(x(t),t)內，即由區(qū)間初始點Mk到區(qū)間終止點Mk+1軌跡的最大曲率半徑為：

(19)

由于飛行器再入過程中速度很快，每個子區(qū)間時間段很小，很難有大幅度的矢量角度變化，因此可將子區(qū)間內軌跡的曲率半徑看作恒定。以曲率半徑rl為半徑，作通過子區(qū)間軌跡初、末兩點Mk和Mk+1的劣弧Cl，定義Cl圍繞線段MkMk+1形成空間為最大軌跡空間。由于通過區(qū)間內軌跡的曲率必然比ρl小，可知最大軌跡空間可以包含飛行器在區(qū)間內所有可能的軌跡。

考慮Ξk為第k個區(qū)間內最大軌跡空間與禁飛區(qū)的距離，則若保證Ξk>ε，即可滿足L(x(t),t)>ε。

若取極端條件，即區(qū)間內的最大軌跡空間與禁飛區(qū)界面迫近相切時，此時軌跡滿足下式：

l(x)>L(x(t),t)>o(x)

(20)

式(19)中，O(X)=Ξk→0+。設區(qū)間中軌跡初、末兩點Mk和Mk+1的直線距離為sl，則sl≈l(x)，且有

(21)

至此，就可以保證再入軌跡全部滿足禁飛區(qū)約束。這一方法對解決離散問題路徑約束問題有較好的推廣性。

3 仿真算例與分析

本節(jié)的仿真實驗，使用CAV-H飛行器的數(shù)據(jù)進行模擬飛行，優(yōu)化仿真的硬件環(huán)境也是相同的?，F(xiàn)將軌跡優(yōu)化問題的條件描述如表1-3所示。

表1 狀態(tài)變量的端點狀態(tài)參數(shù)

表2 熱流率、動壓、過載約束

表3 控制量的取值邊界

除此之外，還要設立禁飛區(qū)約束。在得出上述條件的仿真軌跡后，在原始軌跡和優(yōu)化軌跡可能穿過的地域設立禁飛區(qū)約束，本節(jié)選取飛行實際情況中最為典型的禁飛區(qū)約束類型，針對橢圓柱體禁飛區(qū)約束和橢球體禁飛區(qū)約束模型進行仿真實驗。

建立2個區(qū)域的圓柱體禁飛區(qū)約束如下：

1) 以地球坐標系(30°E，5°N)為禁飛區(qū)水平切面中心，設立半徑r1=150 km的圓柱體禁飛區(qū)。

2) 以地球坐標系(25°E,5°N)為禁飛區(qū)水平切面中心，設立半徑r2=300 km的圓柱體禁飛區(qū)。

為了求解方便，兩個圓柱體禁飛區(qū)的高度可以定為飛行器最高的飛行高度，即70 000 km。

建立1個橢球體禁飛區(qū)約束如下：

1) 以地球坐標系(40°E,6°N)為禁飛區(qū)水平切面中心，橢球體禁飛區(qū)的長半軸a=440 km沿緯線方向，短半軸b=220 km沿經線方向，極半徑c=40 km垂直水平面指向天。

加入禁飛區(qū)約束后，得出原始軌跡和優(yōu)化軌跡的狀態(tài)變量隨時間變化如圖4所示。

由仿真實驗結果可以清晰看出，優(yōu)化后的軌跡完成避障任務時，彈道的跳躍性更強，每次跳躍幅度比原始軌跡的跳躍幅度高出50%以上，但最終在飛行精度上可以達到目標點。優(yōu)化前后的速度和彈道偏角變化區(qū)別不明顯，基本符合預期要求。

優(yōu)化前后的經緯度運行圖如圖5所示。

由圖4可知，相比于原始軌跡，優(yōu)化后的軌跡可以完全避開禁飛區(qū)區(qū)域，通過引入最大軌跡空間的概念，使得軌跡優(yōu)化配點間的軌跡完全處于禁飛區(qū)外部，通過設定誤差門限值，能夠在避障的同時完成最優(yōu)指標。從圖5可以看出，在距離禁飛區(qū)較近的軌跡區(qū)間內，其配點較為密集；在距離禁飛區(qū)較遠且軌跡平滑的階段，其配點較為稀疏。

優(yōu)化后軌跡的控制變量與時間的關系曲線如圖6所示。由圖6(a)可知，優(yōu)化后的軌跡在多數(shù)時刻內攻角控制量略微偏大，大攻角飛行可能會影響到升阻比，進而影響熱流率參數(shù)，但優(yōu)化方案的攻角仍在允許范圍內，飛行中25°以內的攻角不易產生交叉耦合現(xiàn)象。圖6(b)顯示優(yōu)化方案通過調節(jié)飛行器的傾側角實現(xiàn)程序轉彎，從而繞過禁飛區(qū)。

圖4 禁飛區(qū)約束下的軌跡狀態(tài)

圖5 飛行器航行軌跡

圖6 控制變量與時間的關系曲線

優(yōu)化后軌跡的熱流量、動壓值分別如圖7和圖8所示。

圖8 優(yōu)化軌跡的動壓

可以得出，優(yōu)化后的軌跡基本滿足熱流量和動壓約束。相比于原始軌跡，優(yōu)化軌跡在實現(xiàn)越過禁飛區(qū)轉彎時的熱流量會比較大，轉彎后則與原始方案基本持平，飛行的動壓約束在最終狀態(tài)達到峰值，主要是因為再入末端的空氣密度較大，使得動壓在低空飛行時猛增，總體來講，實驗數(shù)據(jù)在可接受范圍內，可以完成高超聲速飛行器再入快速軌跡優(yōu)化任務。

4 結論

通過對自適應偽譜法進行改進，綜合在配點間設立多個誤差采樣點策略，并提出最大軌跡空間的概念，解決了因問題離散化而導致配點間的軌跡不滿足禁飛區(qū)和其他路徑約束。通過仿真實驗可知，求解策略可以求得復雜約束條件下的軌跡最優(yōu)解。通過對比分析優(yōu)化前后方案，突出了采用最大軌跡空間法的優(yōu)化效果。