張科南，周 浩，陳萬春

（北京航空航天大學(xué) 宇航學(xué)院，北京100191）

近年來，隨著反導(dǎo)技術(shù)的不斷進步，機動突防設(shè)計成為現(xiàn)代導(dǎo)彈總體設(shè)計過程中必不可少的一環(huán).國內(nèi)外涌現(xiàn)出了許多關(guān)于彈道導(dǎo)彈、反艦導(dǎo)彈機動突防技術(shù)的研究成果[1～3].這些成果通常集中于研究導(dǎo)彈飛行過程中的某段（如彈道導(dǎo)彈的中段，反艦導(dǎo)彈末端），對該段導(dǎo)彈機動突防制導(dǎo)律進行設(shè)計，獲得制導(dǎo)律的最優(yōu)解.

高超聲速飛行器具有飛行速度快、射程遠和突防能力強等優(yōu)點.目前，關(guān)于高超聲速飛行器機動突防模式的成果還很少，本文旨在探討研究突防技術(shù)在這一全新領(lǐng)域的應(yīng)用.高超聲速飛行器因其自身的特點，采用彈道導(dǎo)彈、反艦導(dǎo)彈等傳統(tǒng)的機動突防技術(shù)效果未必會好.鑒于其速度快、射程遠的特點，本文提出了進行大空域機動甚至全程機動的策略.根據(jù)飛行任務(wù)中的航路點和禁飛區(qū)分布[4，5]，以及禁飛區(qū)類型不同，設(shè)計了多種機動突防模式：單S機動、半圓機動、橫向蛇形機動、滑翔跳躍機動、縱向蛇形機動.不同機動突防模式的設(shè)計問題可以看作是有特定約束要求的彈道優(yōu)化問題.本文受文獻[6]中關(guān)于再入機動彈道工程設(shè)計法的啟發(fā)，結(jié)合擬平衡滑翔飛行的特點，提出了一種簡便易行的直接法優(yōu)化策略[7]，用來求解彈道優(yōu)化問題.優(yōu)化過程中，不但考慮了終端條件約束和熱流、動壓、過載等過程約束，而且加入了航路點和禁飛區(qū)的限制.為便于工程實現(xiàn)，還限制了控制變量的變化幅度及變化率.在文獻[8]基礎(chǔ)上，提出了改進的擬平衡滑翔條件（Improved Quasi－Equilibrium Glide Condition，IQEGC），建立了攻角和傾側(cè)角間一種約束關(guān)系，減少了一個控制變量，從而提高了計算速度；而且能夠保證彈道的平緩以滿足熱流、動壓、過載等過程約束.

1 研究對象和動力學(xué)方程

本文采用國外公開的通用航空飛行器（Common Aero Vehicle，CAV）中的高升力體CAV－H為研究對象，完整的氣動參數(shù)和結(jié)構(gòu)參數(shù)參考文獻[9]的附錄.采用1976美國國家標準大氣模型.

本文設(shè)計的是CAV的滑翔段彈道.初始條件為主動段結(jié)束后，經(jīng)過一定變軌，滿足起滑點要求的條件.終端條件為滿足下壓段交班要求的條件.

假設(shè)地球為球體，不考慮地球自轉(zhuǎn)的影響，三自由度質(zhì)點動力學(xué)方程為

式中：r為地心到飛行器質(zhì)心的距離；θ，φ分別為經(jīng)度和緯度；v為對地速度；FD，F(xiàn)L為阻力和升力，即FD＝ρv2SrefCD/2，F(xiàn)L＝ρv2SrefCL/2，其中，Sref為飛行器氣動參考面積，CD，CL分別為阻力系數(shù)和升力系數(shù)；γ，σ，ψ分別為當?shù)睾桔E角、傾側(cè)角和航向角.

2 約束條件

本文考慮了復(fù)雜的約束條件，包括控制約束、邊界條件和各種過程約束：熱流、過載、動壓、航路點和禁飛區(qū)等.

2.1 控制約束

CAV的2個控制變量為傾側(cè)角σ和攻角α.為了保證飛行器的飛行穩(wěn)定，設(shè)定傾側(cè)角的變化極限是±50°，攻角范圍5°～20°，另外，對控制量變化率也進行約束，要求傾側(cè)角變化率不超過10（°）/s，攻角變化率不超過5（°）/s.

2.2 邊界條件

為滿足起滑點要求，初始條件設(shè)為高度h0＝60km，速度v0＝7km/s，彈道傾角γ0＝0°，彈道偏角ψ0＝90°.終端約束條件包括目標點的位置和指定的高度、速度限制等.為了滿足下壓段的交班條件，規(guī)定滑翔段末端的速度為2 km/s左右，高度為30～40km.

2.3 熱流約束

2.4 過載約束

式（8）左邊為飛行器受到的法向氣動過載，右邊為飛行器可以承受的最大過載，取nT，max＝2.5g，g為當?shù)氐闹亓铀俣?

2.5 動壓約束

式中，ρ為大氣密度；qmax為飛行器可以承受的最大動壓值，取qmax＝60kPa.

2.6 航路點約束

航路點是指飛行器在飛行途中必須經(jīng)過的一些位置點，進行導(dǎo)航信息確認、偵查拍照或者投遞物品等.設(shè)第i個航路點的位置（xi，yi），對飛行器要求在某一時刻ti，須有

2.7 禁飛區(qū)約束

禁飛區(qū)主要指預(yù)警雷達掃描區(qū)或其他防御系統(tǒng)的保護區(qū)域等.本文將禁飛區(qū)分為2類，可以避讓的禁飛區(qū)和無法避讓的禁飛區(qū).可避讓禁飛區(qū)可以簡化為無限高圓柱體型和有限高半球型.前者對應(yīng)預(yù)警雷達掃描范圍，要想避讓該區(qū)域，只能從掃描雷達的盲區(qū)飛行；后者對應(yīng)有射程和射高限制的反導(dǎo)系統(tǒng)作用區(qū)域，例如美國的愛國者系統(tǒng)，由于其反導(dǎo)導(dǎo)彈的射高受限制，也可以從該區(qū)域的上面飛過.對可避讓的禁飛區(qū)，要保證飛行過程中飛行器距離禁飛區(qū)中心在安全范圍以內(nèi)，即

第j個禁飛區(qū)的半徑為Rj，中心坐標（xj，yj）；飛行器到禁飛區(qū)中心的距離，x方向上Δxj＝xxj，y方向上Δyj＝y(tǒng)－yj.

第二類為無法避讓的禁飛區(qū)，如果目標位于禁飛區(qū)的中心，則無論如何都是沒辦法避讓的，這時需要通過一定的機動，改變彈道形狀，加大防御系統(tǒng)的攔截難度.

3 擬平衡滑翔條件

升力體飛行器的大部分再入彈道，彈道傾角γ通常很小，并且變化相對較慢.在式（5）中，令cosγ＝1，＝0，得：

這就是擬平衡滑翔條件QEGC.理論上，只要升力足夠，式（12）就可以保證彈道絕對平直.但實際上，隨著飛行時間變長，速度越來越小，由于攻角有上限制約，飛行器獲得的升力不能滿足式（12）時，＜0，彈道傾角γ會出現(xiàn)一個快速的下滑，不利于飛行器的穩(wěn)定控制.本文對式（12）進行了一定的修正，給出改進的擬平衡滑翔條件（Improved Quasi－Equilibrium Glide Condition，IQEGC）：

式中，ε是絕對值很小的負數(shù)，用來使彈道傾角γ很緩慢地變小，避免出現(xiàn)γ的快速下滑，保證飛行器的穩(wěn)定.另外，ε可以作為一個優(yōu)化變量參與彈道設(shè)計，通過設(shè)計ε，可以得到滿足指定終端條件要求的彈道.

由于FL＝ρv2SrefCL/2，而CL是攻角α的函數(shù)，可見，QEGC提供了2個控制量傾側(cè)角σ和攻角α的約束關(guān)系.若已知σ可以用σ來表示α，反之亦然.這個關(guān)系式用在優(yōu)化過程中可以減少優(yōu)化變量，提高優(yōu)化效率.

4 最優(yōu)控制問題的求解

根據(jù)參數(shù)化方法的不同，求解最優(yōu)控制問題的直接法分為2種基本類型[10]：①離散控制變量，將連續(xù)空間的控制變量參數(shù)化，狀態(tài)變量由數(shù)值積分獲得；②同時離散控制變量和狀態(tài)變量，通過若干數(shù)學(xué)變換將連續(xù)動力學(xué)微分方程轉(zhuǎn)換為代數(shù)方程.本文提出的方法屬于第一類直接法.由于飛行器按照QEGC飛行，選取傾側(cè)角作為獨立優(yōu)化變量，攻角可根據(jù)式（13）得到，反之亦然.按照航路點和禁飛區(qū)的分布，以及禁飛區(qū)類型的不同，將傾側(cè)角或攻角離散為分段常值函數(shù)，段數(shù)的多少根據(jù)機動突防模式?jīng)Q定.取每一段的幅值和持續(xù)時間作為離散后的優(yōu)化變量，然后利用序列二次規(guī)劃（Sequential Quadratic Programming，SQP）方法進行求解.例如，如果要躲避一個禁飛區(qū)并返回射面，只需要將傾側(cè)角分為兩段，兩段的傾側(cè)角符號相反.選擇每段的傾側(cè)角幅值和持續(xù)時間作為優(yōu)化變量即可.離散后的傾側(cè)角為

5 多機動突防模式彈道設(shè)計

不同的作戰(zhàn)任務(wù)，對應(yīng)不同的機動彈道.本文主要根據(jù)不同航路點和禁飛區(qū)的分布，以及禁飛區(qū)類型的不同，設(shè)計相應(yīng)的機動彈道.其中航路點是必須經(jīng)過的位置點.禁飛區(qū)分為可避讓和不可避讓2種情況.故文中將機動彈道分為橫向機動方式和縱向機動方式.

為了便于說明問題，不失一般性，假設(shè)CAV的初始位置經(jīng)緯度為（0°，0°），目標位于赤道平面內(nèi).CAV和目標之間存在若干禁飛區(qū)和航路點.

5.1 橫向機動方式

5.1.1 單S機動

假設(shè)在飛行器和目標點中間存在2個半徑為400km的禁飛區(qū)，中心分別位于經(jīng)緯度（40°，0°）和（70°，0°）處.其余約束條件如第3節(jié)所述.由于2個危險區(qū)間隔較遠，飛行器可以從2個危險區(qū)中間穿過擊中目標.假設(shè)目標位于赤道平面，經(jīng)度未指定.以縱程最大作為優(yōu)化目標函數(shù)，將傾側(cè)角σ離散為三段常值函數(shù)，相鄰兩段符號相反，以每段的幅值和持續(xù)時間為優(yōu)化變量.考慮傾側(cè)角變化率的約束，設(shè)定傾側(cè)角從一段常值轉(zhuǎn)換到另一段常值的過渡時間tf＝10s，這樣可以保證傾側(cè)角的變化率不超過10（°）/s.離散后的傾側(cè)角為

式中，σ1，σ2，σ3是常數(shù)，表示每一段上傾側(cè)角的幅值，相鄰2個符號相反.因為攻角的變化范圍比傾側(cè)角小得多，故只要傾側(cè)角滿足了變化率約束，攻角自然也滿足.

將σ1，σ2，σ3，t1，t2，tf作為優(yōu)化變量，利用 SQP方法求解轉(zhuǎn)換后的非線性規(guī)劃問題，結(jié)果如圖1～圖3所示.

圖1 單S機動軌跡地面投影

從圖1可以看出，優(yōu)化得到的軌跡成功地避開禁飛區(qū)擊中目標.從圖2中的速度和高度曲線可以看出，軌跡滿足終端約束條件.圖3中控制變量傾側(cè)角σ和攻角α的變化都在控制邊界內(nèi).σ變號過程，σ的絕對值先變小后變大，故根據(jù)QEGC決定的攻角α也伴隨出現(xiàn)先小后大的一種現(xiàn)象，反映在圖3中，每次傾側(cè)角σ變號，攻角α出現(xiàn)一個向下的小跳躍.同理可解釋圖2中彈道傾角γ小跳躍.傾側(cè)角σ變號以后，其幅值往往也發(fā)生變化，對應(yīng)的攻角α也會出現(xiàn)比較明顯的變化.圖4顯示機動彈道始終位于飛行走廊內(nèi)部，即駐點熱流、法向過載和動壓都滿足約束要求，這是按照QEGC飛行的優(yōu)點之一.

圖2 單S機動狀態(tài)變量

圖3 單S機動控制變量

圖4 單S機動飛行走廊

單S機動還適用于既要經(jīng)過航路點，又要避開禁飛區(qū)的情況，如圖5所示.

相關(guān)研究表明，影響數(shù)控機床熱變形誤差的主要原因是主軸部件熱變形誤差，由于實驗條件有限，并且這篇論文主要是驗證LWT-LSSVM建模預(yù)測方法的可行性與準確性，所以只對數(shù)控機床部分發(fā)熱部件進行實驗分析和研究。

圖5 單S機動過航路點軌跡地面投影

5.1.2 半圓機動

假設(shè)在5.1.1中的2個禁飛區(qū)中間（55°，0°）處增加一個半徑1 000km的大禁飛區(qū)，此時，從中間穿過禁飛區(qū)的方式難以實現(xiàn)，可以通過一次大空域半圓機動躲過所有禁飛區(qū).結(jié)果如圖6和圖7所示，圖中圖形的含義與5.1.1中相似.

圖6 半圓機動軌跡地面投影

圖7 半圓機動控制變量

5.1.3 橫向蛇形機動

前面幾種機動方式都是針對可以避讓的禁飛區(qū)給出的方案，如果目標位于防御區(qū)中心，此時僅靠避讓飛行是無法擊中目標的.為了以更高的概率突破敵方防御，命中目標，進行全程機動飛行，設(shè)計飛行器在整個滑翔飛行過程中進行蛇形機動，使彈道不規(guī)則，加大對方反導(dǎo)系統(tǒng)的攔截難度.將傾側(cè)角σ用N段常值函數(shù)進行離散，離散后每段的持續(xù)時間相等.其中，N是正整數(shù)，根據(jù)軌跡橫向機動的次數(shù)來確定，本文取N＝10.優(yōu)化變量選取每一段的傾側(cè)角大小及總的飛行時間.結(jié)果如圖8和圖9所示，圖中，R為縱程，X為橫程.

圖8 橫向蛇形機動軌跡地面投影

圖9 橫向蛇形機動控制變量

5.2 縱向機動方式

下面給出的2種機動方式均位于縱平面內(nèi)，飛行過程中傾側(cè)角為0°.

5.2.1 滑翔跳躍機動

針對2.7中介紹的第2類可避讓禁飛區(qū)，可以規(guī)劃一種更簡單的避讓機動方式.先讓飛行器按擬平衡滑翔飛行，當接近目標時通過一次跳躍，躲過目標防御系統(tǒng).擬平衡滑翔階段只進行彈道積分，不進行優(yōu)化計算.當接近目標時，對攻角進行離散處理，方法與橫向機動傾側(cè)角的離散方式相似，然后將離散后的分段攻角幅值和持續(xù)時間作為優(yōu)化變量，求解非線性規(guī)劃問題.

圖10為滑翔跳躍機動飛行軌跡，假設(shè)飛行器從（0°，0°）開始滑翔，禁飛區(qū)是位于（100°，0°）處半徑45km的半球，圖中用虛線表示；目標位于赤道面上禁飛區(qū)之后.圖10（b）為圖10（a）的局部放大圖.以縱程最大作為目標函數(shù)，進行彈道設(shè)計.

圖10 滑翔跳躍機動飛行軌跡

圖11 滑翔跳躍機動攻角

5.2.2 縱向蛇形機動

面對橫向機動方式5.1.3中的目標情況，還可以通過縱向全程跳躍機動，加大防御系統(tǒng)的攔截難度.有別于最大升阻比跳躍彈道，此處對彈道傾角進行了約束，最大幅值不超過1°，避免飛行器因跳躍幅度過大而失控.以擬平衡滑翔飛行彈道的攻角規(guī)律作為離散后攻角的初值，可以保證優(yōu)化出的彈道在擬平衡滑翔彈道附近做小幅振動，既便于工程實現(xiàn)，又能達到機動突防的效果.

此處，給定目標位置（104°，0°），以飛行時間最短作為目標函數(shù)，得到機動彈道如圖12～圖14所示.

圖12 縱向蛇形機動飛行軌跡

圖13 縱向蛇形機動攻角

圖14 縱向蛇形機動飛行走廊

6 結(jié)論

根據(jù)航路點、禁飛區(qū)和目標的分布情況，本文設(shè)計了5種機動突防模式彈道.其中，前3種屬于橫向機動方式，設(shè)計傾側(cè)角的變化規(guī)律，攻角可通過擬平衡滑翔條件得到.后2種屬于縱向機動方式，設(shè)計攻角規(guī)律，傾側(cè)角恒為0°.設(shè)計出的5種機動彈道的攻角規(guī)律、傾側(cè)角規(guī)律都比較簡單，而且得到的彈道滿足各種約束要求，驗證了本文設(shè)計方法的有效性.

[1]SHINAR J，STEINBERGT D.Analysis of optimal evasive maneuvers based on a linearized two－dimensional kinematic model[J].Journal of Aircraft，1977，14（8）：795－802.

[2]TROTTEMANT E J，SCHERER C W.Robust minimax strategies for missile guidance design[C].AIAA Guidance，Navigation and Control Conference and Exhibit.Honolulu，Hawaii：AIAA，2008.

[3]MAURO P，BRUCE A C.Optimal interception of evasive missile warheads：numerical solution of the differential game[J].Journal of Guidance Control and Dynamics，2008，31（4）：1 111－1 123.

[4]JORRIS T R，COBB R G.2－D trajectory optimization satisfying waypoints and no－fly zone constraints[J].Journal of Guidance，Control，and Dynamics，2008，31（3）：551－572.

[5]JORRIS T R，COBB R G.Three－dimensional trajectory optimization satisfying waypoints and no－fly zone constraints[J].Journal of Guidance，Control，and Dynamics，2009，32（2）：543－553.

[6]趙漢元.飛行器再入動力學(xué)和制導(dǎo)[M].長沙：國防科技大學(xué)出版社，1997.ZHAO Han－yuan.Vehicle reentry dynamics and guidance[M].Changsha：National University of Defense Technology Press，1997.（in Chinese）

[7]ZHANG K，CHEN W.Trajectory optimization for hypersonic vehicle satisfying maneuvering penetration[J].Applied Mechanics and Materials，2011，110：5 223－5 231.

[8]SHEN Z，LU P.On－board generation of three－dimensional constrained entry trajectories[J].Journal of Guidance，Control，and Dynamics，2003，26（1）：111－121.

[9]PHILIPS T H.A common aero vehicle（CAV）model，description，and employment guide[EB/OL].[2003].http：//www.dtic.millmatris/sbir041/srch/af031a.doc.

[10]BETTS J T.Practical methods for optimal control using nonlinear programming[M].Washington：Society for Industrial and Applied Mathematics，2001.