熊少鋒，魏明英，趙明元，熊 華，王衛(wèi)紅，周本春

(1. 北京電子工程總體研究所，北京 100854；2. 中國航天科工集團防御技術(shù)研究院,北京 100854；3. 北京航空航天大學(xué)自動化科學(xué)與電氣工程學(xué)院,北京 100191)

0 引 言

攔截戰(zhàn)術(shù)彈道導(dǎo)彈(Tactical ballistic missile,TBM)的制導(dǎo)律設(shè)計是一個富有挑戰(zhàn)性的課題[1]。TBM為代表的高速機動目標，其再入大氣層時的速度可以達到3 km/s，遠大于攔截彈的速度，傳統(tǒng)的尾追制導(dǎo)方式難以讓導(dǎo)彈擊中目標。為了提高攔截概率，需采用迎面逆軌攔截的打擊方式。理想的逆軌攔截有兩個優(yōu)點：1)導(dǎo)彈速度和目標速度均與彈目視線重合，故無論是導(dǎo)彈速度還是目標速度大小的變化，都無需導(dǎo)彈付出額外的法向過載，只是攔截時間變長或變短而已;2)逆軌攔截時，導(dǎo)彈速度與彈目視線重合，導(dǎo)彈的法向過載與彈目視線垂直，其將全部被用于改變視線轉(zhuǎn)率，效率最高。在制導(dǎo)末端，為使導(dǎo)彈具有較充裕的機動能力，要求在制導(dǎo)末端需用過載盡可能小。

綜合而言，逆軌攔截高速機動目標對制導(dǎo)律提出了攻擊角度和過載兩方面的約束[2-3]。常見的制導(dǎo)律設(shè)計方法有比例制導(dǎo)、滑模制導(dǎo)和最優(yōu)制導(dǎo)，下面分別對其進行介紹。

針對固定目標，文獻[4]和[5]分別設(shè)計了虛擬目標比例導(dǎo)引律和偏置比例導(dǎo)引律。其中，文獻[4]有效地避免了末端控制量發(fā)散的問題，而文獻[5]通過設(shè)計一種盲區(qū)控制方案以減少終點處的法向過載。針對高速運動目標，文獻[6]將碰撞角約束作為反饋指令，使用正比例系數(shù)和負比例系數(shù)分別設(shè)計了逆軌和順軌攔截目標的角度約束制導(dǎo)律。

滑?？刂疲惴ê唵吻覍Ω蓴_具有魯棒性，因此被廣泛應(yīng)用于制導(dǎo)律設(shè)計。文獻[7]和[8]分別采用擴張狀態(tài)觀測器和基于Kriging的無模型有限長單位沖激響應(yīng)濾波器估計目標機動的加速度，設(shè)計了帶角度約束的滑模制導(dǎo)律。考慮導(dǎo)彈自動駕駛儀的動態(tài)特性，文獻[9-10]基于滑?？刂圃O(shè)計了角度約束制導(dǎo)律。

最優(yōu)控制因為對各種約束條件的考慮，在制導(dǎo)律設(shè)計方面有著較大的應(yīng)用潛力[11-12]。針對固定目標，文獻[13]和[14]分別應(yīng)用Schwarz不等式和間接Gauss偽譜法設(shè)計了帶角度約束的最優(yōu)制導(dǎo)律。文獻[15]考慮中制導(dǎo)末端的視場角約束，提出一種新的虛擬目標設(shè)置方法，并基于此設(shè)計了滿足視場角約束的最優(yōu)中制導(dǎo)律。文獻[16]應(yīng)用終端投影變換法對高階制導(dǎo)模型進行降階處理，根據(jù)線性二次型控制和微分對策理論設(shè)計了具有角度約束的制導(dǎo)律。在文獻[16]的基礎(chǔ)上，文獻[17]考慮時變的加速度界限，設(shè)計了一種打擊機動目標的角度約束制導(dǎo)律?？紤]導(dǎo)彈自動駕駛儀的二階動態(tài)特性，文獻[1]應(yīng)用高斯偽譜法設(shè)計了多約束條件下的三維最優(yōu)中制導(dǎo)律，且將剩余飛行時間n次冪的倒數(shù)作為控制權(quán)重系數(shù)，使得過載在制導(dǎo)末端趨于一個小值。該制導(dǎo)律需要使用算法軟件包SNOPT求解非線性規(guī)劃問題。

綜上所述文獻在設(shè)計制導(dǎo)律時，一般假設(shè)導(dǎo)彈速度恒定不變。但在實際工程中，導(dǎo)彈由于受到發(fā)動機推力、空氣動力和重力的作用，其速度是時變的，而且變化范圍較大，從Ma1以下上升到Ma5～6，故需要研究考慮導(dǎo)彈速度時變的制導(dǎo)律設(shè)計。

針對固定目標，文獻[18]將剩余飛行距離引入到指標函數(shù)，設(shè)計了帶導(dǎo)引頭視場角約束的制導(dǎo)律，使得過載在制導(dǎo)末端趨于小值，在仿真時考慮了導(dǎo)彈速度是時變的，效果良好。針對勻速運動目標，文獻[19]應(yīng)用間接高斯偽譜法設(shè)計了考慮導(dǎo)彈速度時變的角度約束制導(dǎo)律。對于低速運動的坦克目標，文獻[20]將導(dǎo)彈速度的變化率和目標機動視為總擾動，利用干擾觀測器對其進行估計，基于滑?？刂圃O(shè)計了角度約束制導(dǎo)律。針對機動目標，文獻[21]應(yīng)用伴隨變換法設(shè)計了考慮導(dǎo)彈速度時變的最優(yōu)制導(dǎo)律。文獻[22]在文獻[16]的基礎(chǔ)上，應(yīng)用線性時變系統(tǒng)理論設(shè)計了考慮導(dǎo)彈速度時變的角度約束最優(yōu)制導(dǎo)律。文獻[21-22]的制導(dǎo)律均包含了導(dǎo)彈速度從當前時刻到終端時刻的積分這一項。為解決導(dǎo)彈未來速度未知的問題，文獻[21]忽略誘導(dǎo)阻力，給出了預(yù)測導(dǎo)彈未來速度的一次性方法；文獻[22]則給出了預(yù)測導(dǎo)彈未來速度的多步迭代方法。文獻[21-22]預(yù)測導(dǎo)彈未來速度的方法均需要耗費大量的計算時間。

在傳統(tǒng)的最優(yōu)制導(dǎo)律中[16,22]，脫靶量和角度約束的權(quán)重系數(shù)被設(shè)為一個較大的常數(shù)，這會使得導(dǎo)彈在初始階段的過載指令較大。但很多時候，導(dǎo)彈在初始階段速度較低，還處于加速過程中，難以提供較大的可用過載。為解決這個問題，文獻[23]對雙曲正切函數(shù)進行改造，設(shè)計了雙曲正切函數(shù)的變種，它隨時間從零開始增大，且增大的速度越來越快，并將其設(shè)為脫靶量和角度約束的權(quán)重系數(shù)。這么做的思想是，在制導(dǎo)的初始階段，導(dǎo)彈離目標較遠，不需要對脫靶量和攻擊角度施加太強的約束，導(dǎo)彈可以“自由”朝目標飛去；隨著時間的推移，導(dǎo)彈逐漸接近目標，這時就得逐步強化對脫靶量和攻擊角度的約束。雙曲正切函數(shù)的變種符合這一要求特性，降低導(dǎo)彈在初始階段的過載需求。

上述文獻里的制導(dǎo)律基本都是二維的，但實際中導(dǎo)彈和目標是在三維空間里飛行的，需要研究三維制導(dǎo)律設(shè)計。文獻[24]基于李亞普諾夫穩(wěn)定性理論設(shè)計了控制撞擊角度和時間的三維導(dǎo)引律，其打擊對象為固定或者慢速移動目標。

針對導(dǎo)彈逆軌攔截來襲高速機動目標，本文將[23]的制導(dǎo)律從二維擴展到三維。攻擊角度定義為導(dǎo)彈速度矢量和目標速度矢量之間的夾角，在三維情況下，難以像文獻[25]直接得到期望視線角和攻擊角度之間的解析關(guān)系式，增加了設(shè)計難度。從便于工程應(yīng)用的角度出發(fā)，對三維制導(dǎo)方程進行簡化處理，將三維制導(dǎo)模型投影到垂直和橫向兩個二維平面。在投影過程中，導(dǎo)彈在三維空間里的彈道傾角與其在二維垂直平面內(nèi)的彈道傾角并不相等，兩者之間存在一個等效關(guān)系，因此使得垂直平面內(nèi)制導(dǎo)模型的建立較文獻[23]中單純的二維制導(dǎo)變得復(fù)雜，增加了三維制導(dǎo)律的設(shè)計難度。截至目前，作者未曾在公開發(fā)表的文獻中看到類似的逆軌攔截機動目標的三維制導(dǎo)律推導(dǎo)過程，本文中的垂直平面和偏航平面內(nèi)的制導(dǎo)方程都是作者原創(chuàng)推導(dǎo)的。本文設(shè)計的三維制導(dǎo)律具有解析表達式，可實時在線運行，滿足工程中導(dǎo)彈攔截TBM目標時的逆軌要求，且導(dǎo)彈在剛開始加速階段，過載指令不太大。

1 彈目攔截方程

先定義一下慣性坐標系，慣性坐標系以導(dǎo)彈發(fā)射車為原點；x軸正向為發(fā)射車的車頭方向；y軸在當?shù)劂U垂面內(nèi)，向上為正；z軸在當?shù)厮矫鎯?nèi)，且與x軸和y軸滿足右手定則。

在慣性坐標系中，彈目攔截幾何如圖1所示：

圖1 彈目攔截幾何Fig.1 Missile target engagement geometry

在圖1中，M和T分別表示導(dǎo)彈和目標，r表示彈目視線，qε和qβ分別表示視線傾角和視線偏角。飛行器(包括導(dǎo)彈和目標)在慣性坐標系中的速度如圖2所示。

圖2 飛行器速度矢量Fig.2 Velocity vector of aircraft

其中，v,γ和ψ分別表示飛行器的速度、彈道傾角和彈道偏角。導(dǎo)彈速度vm和目標速度vt在三維空間里的投影分別為

(1)

將三維空間分解成鉛垂和水平兩個二維平面，進而將三維制導(dǎo)律簡化為兩個二維制導(dǎo)律來設(shè)計。在鉛垂面內(nèi)設(shè)計制導(dǎo)律時，假設(shè)導(dǎo)彈的彈道偏角ψm(t)在t時刻保持瞬時恒定；在水平面內(nèi)設(shè)計制導(dǎo)律時，假設(shè)導(dǎo)彈的彈道傾角γm(t)在t時刻保持瞬時恒定。這里需要強調(diào)的是，假設(shè)導(dǎo)彈的彈道傾角γm和彈道偏角ψm保持瞬時恒定，這樣做僅僅只是為了便于將三維制導(dǎo)分解成兩個二維制導(dǎo)，從而簡化三維制導(dǎo)律的推導(dǎo)過程，但是在執(zhí)行制導(dǎo)律時，依然還是用γm和ψm隨時間變化的實時值。另外，一般來說，彈道傾角和彈道偏角相對于制導(dǎo)指令而言是慢時變量，所以推導(dǎo)過程中假設(shè)它們保持瞬時不變。穩(wěn)定性方面，制導(dǎo)系統(tǒng)的狀態(tài)變量收斂很快，將三維制導(dǎo)問題解耦成兩個二維制導(dǎo)問題不會影響制導(dǎo)系統(tǒng)的穩(wěn)定性，后面的仿真結(jié)果也表明解耦設(shè)計的制導(dǎo)律是穩(wěn)定的。另外，還可以參考文獻[26]，介紹了三維制導(dǎo)律可以分解成兩個二維制導(dǎo)律來設(shè)計，而不影響制導(dǎo)系統(tǒng)的穩(wěn)定性。

1.1 鉛垂面內(nèi)的制導(dǎo)方程

在鉛垂平面內(nèi)，彈目相對運動如圖3所示。

圖3 oxy平面內(nèi)的彈目攔截幾何Fig.3 Missile target engagement geometry in oxy plane

彈目距離沿垂直于初始視線方向的分量y的動力學(xué)方程為：

(2)

對式(2)求導(dǎo)可得

(3)

由式(1)可得導(dǎo)彈速度vm和目標速度vt在oxy平面上的投影分別為

(4)

對式(4)求導(dǎo)可得

(5)

根據(jù)式(1)可得導(dǎo)彈和目標在oxy平面內(nèi)的等效彈道傾角分別為

(6)

對式(6)求導(dǎo)可得

(7)

定義如下的兩個變量

(8)

則式(7)可以表示為

(9)

將式(4)和式(5)及式(9)代入式(3)可得

(10)

定義如下兩個變量為

(11)

將式(11)代入式(10)可得

sin(γmxy-q10)-ηmxyay

(12)

(13)

根據(jù)式(12)可建立如下的狀態(tài)空間方程

(14)

1.2 水平面內(nèi)的制導(dǎo)方程

在水平面內(nèi)，彈目相對運動如圖4所示。

圖4 oxz平面內(nèi)的彈目攔截幾何Fig.4 Missile target engagement geometry in oxz plane

彈目距離沿垂直于初始視線方向的分量z的動力學(xué)方程為：

(15)

對式(15)求導(dǎo)可得

(16)

根據(jù)式(1)可得導(dǎo)彈速度vm和目標速度vt在oxz平面上的投影分別為

vmxz=vmcosγm,vtxz=vtcosγt

(17)

對式(17)求導(dǎo)可得

(18)

根據(jù)式(1)可得導(dǎo)彈和目標在水平面內(nèi)的等效彈道偏角分別為

ψmxz=ψm,ψtxz=-ψt±π

(19)

對式(19)求導(dǎo)可得

(20)

將式(18)和(20)代入式(16)可得

(21)

定義一個新的變量ftxz為

(22)

將式(22)代入式(21)可得

cos(ψm-q20)az

(23)

(24)

根據(jù)式(23)可建立如下的狀態(tài)空間方程

(25)

2 最優(yōu)制導(dǎo)律設(shè)計

將導(dǎo)彈逆軌攔截機動目標的要求轉(zhuǎn)化為對彈目速度交會角的約束，基于最優(yōu)控制理論設(shè)計制導(dǎo)律。

鉛垂面內(nèi)的最優(yōu)制導(dǎo)律，設(shè)計性能指標函數(shù)為

(26)

式中:tf表示終端時刻，tgo表示剩余飛行時間，a為脫靶量的權(quán)重系數(shù)，b為攻擊角度的權(quán)重系數(shù)。針對式(26)，構(gòu)造Hamiltonian函數(shù)為

(27)

根據(jù)極小值原理可得最優(yōu)控制解為

(28)

協(xié)態(tài)方程為

(29)

對式(29)求解微分方程可得

(30)

(31)

從式(31)可以看到，最優(yōu)控制解里面包含了終端時刻的狀態(tài)量x1tf和x4tf，接下來給出這兩個終端狀態(tài)量的求解過程。根據(jù)線性系統(tǒng)理論可得

(32)

(33)

求解式(33)可得

(34)

其中，

(35)

(36)

定義一個新的狀態(tài)向量為

(37)

ZEM(Zero effort miss)和ZEAE(Zero effort angle error)分別為零控脫靶量和零控角誤差，它們的數(shù)學(xué)表達式為

(38)

式中：

(39)

將式(34)代入式(32)右邊的第二個部分，并從t到tf積分可得

(40)

式中：

(41)

將式(37)和(40)代入式(32)可得

(42)

求解方程(42)可得

(43)

(44)

將式(43)和式(44)代入式(31)可得最優(yōu)制導(dǎo)律為

(45)

式中：

(46)

(47)

式(45)中的最優(yōu)制導(dǎo)律包含了導(dǎo)彈速度從t時刻到tf時刻的積分，但導(dǎo)彈在未來時刻的速度是未知的，這給制導(dǎo)律的執(zhí)行帶來了困難，從工程易于執(zhí)行的角度出發(fā)，簡單地假設(shè)導(dǎo)彈速度在t時刻到tf時刻是不變的?；谶@一假設(shè)，可直接得到下面的關(guān)系式

(48)

假設(shè)導(dǎo)彈和目標的相對位置偏離碰撞三角形較小，也即(q1-q10)<<1，則下列近似關(guān)系式成立

y=rxy(q1-q10)

(49)

對式(49)求導(dǎo)可得

(50)

一般來說，目標未來的加速度是難以預(yù)知的，因此假設(shè)在時間域[t,tf]內(nèi)目標加速度是不變的，故可得到下列關(guān)系式

(51)

將式(48)，式(50)～(51)代入式(38)可得零控脫靶量和零控角誤差的簡化表達式為

(52)

采用同樣的方法，可推導(dǎo)設(shè)計水平面內(nèi)的最優(yōu)制導(dǎo)律，為節(jié)省篇幅，此處略去具體的推導(dǎo)過程。

考慮到前面所做的易于工程實施的簡化處理，現(xiàn)將y/z向過載指令的計算公式演繹、綜合如式(53)和式(54)所示。

(53)

(54)

3 雙曲正切函數(shù)

雙曲正切函數(shù)的表達式為

(55)

其隨時間變化的曲線如圖5所示。

圖5 雙曲正切函數(shù)Fig.5 Hyperbolic tangent function

設(shè)計雙曲正切函數(shù)的變種函數(shù)為

(56)

令n=2，ζ=2，τ=3，雙曲正切函數(shù)的變種函數(shù)隨時間變化的曲線如圖6所示。

采用雙曲正切函數(shù)的變種設(shè)計指標函數(shù)如下

(57)

4 仿真校驗

通過數(shù)字仿真以說明所設(shè)計最優(yōu)制導(dǎo)律的性能。

假設(shè)導(dǎo)彈是軸對稱的，其動力學(xué)方程為式(58)。

(58)

導(dǎo)彈動力學(xué)方程(58)闡明了導(dǎo)彈所受到的力，導(dǎo)彈速度的變化關(guān)系式，在后續(xù)仿真中將由式(58)實時計算導(dǎo)彈的速度，升力、阻力和側(cè)向力。在式(58)中，T為發(fā)動機推力，X，Y和Z分別為阻力、升力和側(cè)向力，α和β為攻角和側(cè)滑角，ρ和S為空氣密度和特征面積，K，CX0，CYα和CZβ為空氣動力系數(shù)。

仿真參數(shù)：導(dǎo)彈的初始位置為(2000 m, 1000 m, -1000 m)，目標的初始位置為(80000 m，75000 m, -10000 m)，導(dǎo)彈的初始速度為269.2582 m/s，目標的速度、彈道傾角和彈道偏角的初值分別為3471.3 m/s，-43.7396°和-175.4261°，目標沿x，y和z三軸的加速度分別是-3g，-2g和-2g，導(dǎo)彈的初始質(zhì)量為734 kg，燃料的燃燒時間和消耗速率分別為10 s和-24.4 kg/s，在燃料燃燒期間發(fā)動機的推力T為110 kN，燃料消耗完畢后發(fā)動機推力T為0 N，導(dǎo)彈的特征面積S為0.146 m2，重力加速度g=9.8 m/s2，攻角和側(cè)滑角的最大值為35°，仿真步長為0.01 s，空氣動力系數(shù)K，CX0，CYα和CZβ由表1查表獲得。

表1 空氣動力系數(shù)Table 1 Aerodynamic coefficient

為便于性能對比驗證，將脫靶量和角度約束權(quán)重系數(shù)設(shè)為常值的傳統(tǒng)最優(yōu)制導(dǎo)律(Traditional optimal guidance law,TOGL)[16]也一并仿真。本文制導(dǎo)律和文獻[16]的制導(dǎo)律，它們的仿真條件是相同的，導(dǎo)彈速度是時變的，其變化關(guān)系由式(58)確定。對于本文提出的制導(dǎo)律，控制參數(shù)(俯仰通道和偏航通道均采用相同的控制參數(shù))為：N=1，n=2，μ=4，ζ=1。

(59)

TOGL的控制參數(shù)為：N=1，a=103，b=107。

導(dǎo)彈的初始彈道傾角和彈道偏角分別為55°和15°，仿真結(jié)果如圖7～12所示。脫靶量為0.5633 m，導(dǎo)彈的末速度為1274.9 m/s。彈目速度縱向和橫向交會角分別為-1.4799°和-179.9103°，與逆軌攔截要求的交會角0°和-180°相差-1.4799°和0.0897°，滿足逆軌攔截對彈目交會角的約束。從圖10～11看出，初始階段導(dǎo)彈沿Y軸和Z軸的過載較小，分別為-4.1463gm/s2和2.9191gm/s2，在制導(dǎo)末段，Y向和Z向過載指令均可以趨于零附近的小值，滿足制導(dǎo)律對過載指令在初始段和末段不能太大的約束。Y/Z向過載指令在第10 s時存在尖點的現(xiàn)象，這是因為第10 s的時候，導(dǎo)彈燃料消耗完畢，發(fā)動機推力由110 kN突降到0 N所致。第18 s時的尖點現(xiàn)象，則是由控制參數(shù)τ在第18 s時連續(xù)但不可導(dǎo)(不光滑)造成的。在圖12中，由于發(fā)動機推力的作用，導(dǎo)彈前10 s處于加速過程，第10 s時，導(dǎo)彈速度達到峰值，后面由于發(fā)動機推力為零，導(dǎo)彈只受到空氣阻力和重力的作用，其速度逐漸減小。

圖7 三維空間內(nèi)的彈目運動軌跡Fig.7 Trajectories of missile and target in threedimensional space

圖8 X-Y平面內(nèi)的彈目運動軌跡Fig.8 Trajectories of missile and target in X-Y plane

圖9 X-Z平面內(nèi)的彈目運動軌跡Fig.9 Trajectories of missile and target in X-Z plane

圖10 Y向過載Fig.10 Guidance command along Y axis

圖11 Z向過載Fig.11 Guidance command along Z axis

圖12 導(dǎo)彈速度Fig.12 Velocity of missile

5 結(jié) 論

針對大氣層內(nèi)導(dǎo)彈逆軌攔截高速運動目標的三維制導(dǎo)律設(shè)計問題，本文考慮導(dǎo)彈速度時變的情況，運用最優(yōu)控制理論和雙曲正切函數(shù)設(shè)計了帶角度約束的最優(yōu)制導(dǎo)律。通過合理的簡化處理，本文設(shè)計的制導(dǎo)律具有解析形式，可實時在線運行。仿真結(jié)果表明該最優(yōu)制導(dǎo)律可以提升導(dǎo)彈的末速度，增強導(dǎo)彈的殺傷動能。