尤志鵬 邵 干 李 洋 楊 勇

中國(guó)運(yùn)載火箭技術(shù)研究院，北京 100076

0 引言

空天飛行器再入過(guò)程具有約束多、大氣及氣動(dòng)不確定性高等特點(diǎn)，設(shè)計(jì)精度高、自主性好、自適應(yīng)性強(qiáng)的再入制導(dǎo)算法仍然充滿挑戰(zhàn)。航天飛機(jī)再入制導(dǎo)律是最早得到應(yīng)用的升力式再入制導(dǎo)律，經(jīng)過(guò)了多次飛行考驗(yàn)，從未發(fā)生致命性故障，被后續(xù)可重復(fù)使用飛行器(RLV)所繼承并得到進(jìn)一步發(fā)展[1-2]。此后，更加先進(jìn)的再入制導(dǎo)算法得到廣泛研究，自主性、快速性、適應(yīng)性進(jìn)一步增強(qiáng)[3-5]。

現(xiàn)有再入制導(dǎo)算法主要分為2類：標(biāo)準(zhǔn)軌跡制導(dǎo)和預(yù)測(cè)校正制導(dǎo)。目前，標(biāo)準(zhǔn)軌跡制導(dǎo)應(yīng)用最廣泛，Mease等[6-7]提出衍化加速度再入制導(dǎo)律(EAGLE)，比傳統(tǒng)航天飛機(jī)再入制導(dǎo)律具有更強(qiáng)的自適應(yīng)能力，可以應(yīng)用于大橫程再入。王鵬等[8]進(jìn)一步考慮了基于阻力加速度-能量剖面再入制導(dǎo)算法的速度傾角控制問(wèn)題，提升了制導(dǎo)精度和適應(yīng)性。趙頔等[9]通過(guò)數(shù)值迭代確定滿足終端約束的速度-高度剖面，制導(dǎo)精度較高。但是，標(biāo)準(zhǔn)軌跡跟蹤制導(dǎo)算法需要離線設(shè)計(jì)標(biāo)準(zhǔn)剖面，自主性較差，不利于飛行器性能發(fā)揮。

預(yù)測(cè)校正制導(dǎo)通過(guò)比較預(yù)測(cè)終端條件與期望終端條件偏差形成制導(dǎo)指令，不需要事先規(guī)劃飛行剖面。早期預(yù)測(cè)校正制導(dǎo)迭代參數(shù)較多，F(xiàn)uhry等[10]通過(guò)解算傾斜角指令和傾斜角變號(hào)時(shí)機(jī)消除預(yù)測(cè)落點(diǎn)與期望落點(diǎn)偏差，但該算法計(jì)算復(fù)雜，計(jì)算資源消耗多。為降低迭代復(fù)雜度，Lu等[11]將攻角剖面固定，每制導(dǎo)周期僅迭代獲取滿足航程要求的傾側(cè)角剖面，并通過(guò)傾側(cè)角反號(hào)控制橫程，得到一種通用預(yù)測(cè)校正制導(dǎo)算法[12-13]。劉剛等[14]進(jìn)一步考慮了預(yù)測(cè)校正制導(dǎo)律設(shè)計(jì)過(guò)程中的攻角調(diào)節(jié)問(wèn)題，增強(qiáng)了算法適應(yīng)性。但是，上述制導(dǎo)算法在每個(gè)制導(dǎo)周期均需要在縱向運(yùn)動(dòng)平面內(nèi)執(zhí)行多次彈道積分，算法效率仍存在較大提升空間。

為降低預(yù)測(cè)校正制導(dǎo)計(jì)算復(fù)雜度，本文提出一種基于速度-高度剖面的預(yù)測(cè)校正再入制導(dǎo)算法。在每個(gè)制導(dǎo)周期通過(guò)卡爾曼濾波估計(jì)速度中點(diǎn)對(duì)應(yīng)的無(wú)量綱高度，預(yù)測(cè)滿足航程要求的速度-高度剖面，避免傳統(tǒng)數(shù)值預(yù)測(cè)校正的多次數(shù)值迭代。同時(shí)，針對(duì)算法在再入飛行末端擬合系數(shù)奇異問(wèn)題，通過(guò)切換為標(biāo)準(zhǔn)軌跡制導(dǎo)予以避免。

1 運(yùn)動(dòng)方程及約束

1.1 運(yùn)動(dòng)方程

假設(shè)地球是均質(zhì)圓球，三維質(zhì)點(diǎn)再入運(yùn)動(dòng)無(wú)量綱方程為：

(1)

L=ρ(VcV)2SrefCL/(2mg0)

(2)

D=ρ(VcV)2SrefCD/(2mg0)

(3)

式中:ρ表示大氣密度；Sref,m分別表示參考面積和飛行器質(zhì)量；CL,CD分別表示升力系數(shù)和阻力系數(shù)，與攻角α有關(guān)，而攻角通常按速度進(jìn)行裝訂。

1.2 再入約束

再入過(guò)程約束主要包括動(dòng)壓約束、熱流約束、過(guò)載約束、平衡滑翔約束，分別通過(guò)下式計(jì)算：

(4)

(5)

n=|Lcosα+Dsinα|≤nmax

(6)

(7)

(8)

(9)

(10)

通過(guò)牛頓迭代求解式(7)，得到該約束下r和V的對(duì)應(yīng)關(guān)系，至此，給定攻角剖面后，即可以得到再入走廊。

末端約束主要包含末端高度、末端速度、末端經(jīng)緯度約束，對(duì)于基于速度-高度設(shè)計(jì)的飛行剖面，速度是自變量，因此高度約束和經(jīng)緯度約束是末端主要約束。表達(dá)如下

r(Vf)=rf

(11)

θ(Vf)=θf(wàn)

(12)

φ(Vf)=φf(shuō)

(13)

式中：rf、θf(wàn)和φf(shuō)分別表示末端飛行器質(zhì)心至地心距離、末端經(jīng)度和末端緯度。

2 制導(dǎo)算法設(shè)計(jì)

飛行器再入段包括初始再入段和擬平衡滑翔再入階。初始再入段氣動(dòng)力作用較弱，通常在滿足航程、熱流等約束作用下，以固定傾角飛行。當(dāng)飛行器再入軌跡進(jìn)入再入走廊并滿足擬平衡滑翔條件時(shí)，飛行器進(jìn)入擬平衡滑翔階段，該階段飛行時(shí)間長(zhǎng)、飛行距離遠(yuǎn)、狀態(tài)變化大，是再入制導(dǎo)律主要起作用階段，本文從該段開(kāi)始設(shè)計(jì)再入制導(dǎo)律。由于進(jìn)入該階段后，飛行速度隨飛行時(shí)間逐漸降低，具有良好的單調(diào)性，因此可將速度作為運(yùn)動(dòng)方程的自變量，實(shí)現(xiàn)運(yùn)動(dòng)方程降階，提高預(yù)測(cè)校正制導(dǎo)計(jì)算效率。

2.1 速度-高度飛行剖面擬合

再入飛行過(guò)程中，飛行器除需要滿足再入走廊約束外，縱向軌跡還需要滿足航程、終端速度和終端高度約束，除此之外，為保證飛行軌跡平滑，可進(jìn)一步引入終端高度對(duì)終端速度的導(dǎo)數(shù)約束。飛行器當(dāng)前縱向軌跡狀態(tài)包含當(dāng)前速度、當(dāng)前高度、當(dāng)前航跡傾角。由于

(14)

式中：h表示無(wú)量綱飛行高度，因此

(15)

即航跡傾角約束可轉(zhuǎn)化為飛行高度對(duì)飛行速度的導(dǎo)數(shù)約束。因此可利用三次分段多項(xiàng)式形式的速度-高度飛行剖面擬合縱向飛行剖面，即

(16)

式中：aj,bj,j=0,…,3是擬合系數(shù)，Vm表示當(dāng)前速度至交班點(diǎn)速度的速度中點(diǎn)，即Vm=(Vi+Vf)/2，在第i個(gè)制導(dǎo)周期，已知狀態(tài)主要包括當(dāng)前點(diǎn)的速度Vi、高度hi、航跡傾角γi，再入段交班點(diǎn)的速度Vf、高度hf。初始及終端狀態(tài)構(gòu)成四組約束，即(Vi,hi)、(Vi,(dh/dV)i)、(Vf,hf)、(Vf,(dh/dV)f)。此外為滿足連續(xù)性要求，需要兩個(gè)分段在Vm處飛行高度相等且高度相對(duì)于速度的導(dǎo)數(shù)及二階導(dǎo)數(shù)相等，形成三組約束。同時(shí)，航程要求可通過(guò)調(diào)整Vm處飛行高度hm實(shí)現(xiàn)，至此形成了8組約束，可對(duì)擬合系數(shù)進(jìn)行求解。在速度-高度剖面內(nèi)航程隨速度變化表達(dá)為：

(17)

將式(15)代入，可得飛行航程與期望待飛航程之差可以表示為

(18)

式中：Rexp表示期望待飛航程。

選擇合適的hm，使得S(hm)=0，從而滿足航程約束。

2.2 縱向制導(dǎo)指令

為得到期望的縱向航程，需要在每個(gè)制導(dǎo)周期獲得合適的hm，從而確定飛行剖面并得到傾側(cè)角幅值。

通過(guò)卡爾曼濾波對(duì)每個(gè)制導(dǎo)周期內(nèi)滿足式(18)的hm進(jìn)行辨識(shí)，利用辨識(shí)結(jié)果求解擬合剖面的擬合系數(shù)，進(jìn)而產(chǎn)生該制導(dǎo)周期所需要的制導(dǎo)指令，避免了傳統(tǒng)數(shù)值預(yù)測(cè)算法需要多次積分縱向運(yùn)動(dòng)方程的不足，增強(qiáng)了算法實(shí)時(shí)性。

建模過(guò)程如下：

狀態(tài)方程：

hm_i+1=hm_i+εi

(19)

觀測(cè)方程：

Zi=S(hm_i)+νi

(20)

式中：i表示第i個(gè)制導(dǎo)周期，hm_i+1,hm_i表示第i+1和第i個(gè)制導(dǎo)周期速度中點(diǎn)對(duì)應(yīng)的無(wú)量綱飛行高度；εi及νi為互不相關(guān)零均值高斯白噪聲，其協(xié)方差分別記為Qi,Wi。

為獲得使再入剖面滿足航程約束的hm_i，可將觀測(cè)值設(shè)置為Zi+1≡0。hm_i必須滿足再入走廊約束并且具有一定的裕度，將速度中點(diǎn)對(duì)應(yīng)飛行高度距離再入走廊邊界的無(wú)量綱距離取為不小于κ=5×10-4，修正后的辨識(shí)結(jié)果見(jiàn)式(21)。

(21)

觀測(cè)方程是非線性的，濾波過(guò)程中需要通過(guò)式(22)所示的差分計(jì)算觀測(cè)矩陣Hi+1。

(22)

在每個(gè)制導(dǎo)周期，當(dāng)hm_i確定后，可通過(guò)式(23)求解該制導(dǎo)周期對(duì)應(yīng)的速度-高度剖面擬合系數(shù)。

(23)

為獲得傾側(cè)角指令，將飛行高度對(duì)速度求二階導(dǎo)數(shù)，得

(24)

其中

(25)

(26)

(27)

當(dāng)前飛行狀態(tài)接近終端狀態(tài)時(shí)，式(23)求逆過(guò)程會(huì)出現(xiàn)奇異。為解決該問(wèn)題，當(dāng)前飛行速度與終端速度之差小于100m/s時(shí)，停止預(yù)測(cè)校正制導(dǎo)，不再更新飛行剖面，改為跟蹤最后一次產(chǎn)生的飛行剖面，即式(28)。

cosσ=

(28)

式中：ξ,ω分別是阻尼比和自然頻率，可取為固定值；href為參考軌跡對(duì)應(yīng)的高度。為進(jìn)一步增強(qiáng)跟蹤能力，考慮引入傾側(cè)角反饋，即

Lcosσ′=Lcosσ-K(γ-γref)

(29)

式中：K是反饋系數(shù)；γref為參考剖面對(duì)應(yīng)的航跡傾角。

2.3 橫向制導(dǎo)

橫向制導(dǎo)采用基于航跡偏角偏差走廊的傾側(cè)角反轉(zhuǎn)實(shí)現(xiàn)，表達(dá)如下

(30)

式中：Δψu(yù)p和Δψdown表示航跡偏角走廊上邊界和下邊界；Δψ表示航跡偏角偏差，計(jì)算如下

Δψ=ψ-ψLOS

(31)

ψLOS即當(dāng)前位置至目標(biāo)點(diǎn)的理想視線角，計(jì)算如下

(32)

至此完成算法設(shè)計(jì)，流程見(jiàn)圖1。

圖1 算法流程

3 仿真校驗(yàn)

仿真模型參考文獻(xiàn)[12]，攻角剖面取為

(33)

3.1 標(biāo)稱狀態(tài)仿真

對(duì)4種不同航程下的算例進(jìn)行仿真驗(yàn)證。它們初始再入位置不同，交班點(diǎn)高度為25km，初始狀態(tài)如表1所示。

表1 初始條件

仿真結(jié)果如圖2～3所示?？梢?jiàn)在標(biāo)稱狀態(tài)下，基于卡爾曼濾波的速度-高度剖面預(yù)測(cè)校正制導(dǎo)律能夠滿足飛行約束條件，軌跡變化平緩，精度較高。

圖2 不同算例下標(biāo)稱再入軌跡

相比于文獻(xiàn)[12]所展示的數(shù)值預(yù)測(cè)校正制導(dǎo)算法(NPC)，本文制導(dǎo)指令生成速度更快。NPC制導(dǎo)復(fù)現(xiàn)時(shí)，制導(dǎo)周期選擇為1s，預(yù)測(cè)步長(zhǎng)選擇為5s，積分方式采用歐拉積分，本文算法制導(dǎo)指令生成速度與之對(duì)比如下。

表2 全飛行段生成一個(gè)制導(dǎo)指令所需時(shí)長(zhǎng)

圖3 無(wú)量綱高度濾波結(jié)果及制導(dǎo)指令

3.2 參數(shù)擾動(dòng)狀態(tài)仿真

考慮參數(shù)初始狀態(tài)偏差及不確定性參數(shù)辨識(shí)初值偏差，針對(duì)算例4，利用蒙特卡洛仿真校驗(yàn)算法魯棒性，仿真次數(shù)設(shè)定為300次。初值偏差均取為隨機(jī)偏差(表3)，升力系數(shù)、阻力系數(shù)、大氣密度、攻角偏差取為隨機(jī)偏差和固定偏差的組合(表4)。

表3 初值參數(shù)偏差

表4 過(guò)程參數(shù)偏差

仿真結(jié)果如圖4～6所示，圖4是偏差條件下飛行軌跡，可見(jiàn)制導(dǎo)律能夠適應(yīng)偏差工況。圖5是對(duì)應(yīng)的辨識(shí)結(jié)果和制導(dǎo)指令，算法制導(dǎo)指令變化較為平滑，無(wú)發(fā)散現(xiàn)象。圖6是再入交班點(diǎn)分布，可見(jiàn)大部分工況下交班點(diǎn)在理論交班點(diǎn)3km以內(nèi)，制導(dǎo)精度較高。

圖4 蒙特卡洛仿真再入軌跡

圖5 蒙特卡洛仿真辨識(shí)結(jié)果及制導(dǎo)指令

圖6 蒙特卡洛仿真交班點(diǎn)分布

4 結(jié)論

基于卡爾曼濾波的速度-高度剖面預(yù)測(cè)校正再入制導(dǎo)算法，每次制導(dǎo)指令的生成僅需要2次單變量數(shù)值積分，相比于傳統(tǒng)數(shù)值預(yù)測(cè)校正制導(dǎo)，避免了多次積分縱向運(yùn)動(dòng)方程，提高了實(shí)時(shí)性。針對(duì)算法末端擬合精度下降的問(wèn)題，設(shè)置了更新終止條件并跟蹤最后一次規(guī)劃產(chǎn)生的標(biāo)準(zhǔn)飛行剖面，同時(shí)引入航跡傾角反饋，提升了算法在飛行末端的收斂性。

仿真表明，算法制導(dǎo)精度高，計(jì)算速度快，具有良好的適應(yīng)性和魯棒性。