章吉力，佘智勇，樊雅卓，劉 凱，安帥斌

（1.大連理工大學(xué)航空航天學(xué)院，大連116024；2.北京空天技術(shù)研究所，北京100074）

1 引 言

隨著高超聲速航空航天技術(shù)的進(jìn)步與發(fā)展，新型的可重復(fù)使用空天飛行器為經(jīng)濟(jì)實(shí)用的空間往返提供了新的可行方案。空天飛行器具有速度快、可靠性高的特點(diǎn)，在民用和軍用領(lǐng)域均有重大應(yīng)用價(jià)值，因此在近幾年一直是世界各國(guó)研究的焦點(diǎn)。毫無(wú)疑問，空天飛行器將大大降低空間任務(wù)的運(yùn)營(yíng)成本，因?yàn)樗谌蝿?wù)之后可以進(jìn)行部分甚至全部的回收，將原有的一次性成本分?jǐn)偟蕉啻伟l(fā)射任務(wù)中去。

基于可重復(fù)使用與水平起降的要求，空天飛行器需要具有升力面，加上其本身飛行包線寬，環(huán)境不確定性大，導(dǎo)致其彈道設(shè)計(jì)十分困難，尤其是在空天飛行器的再入階段，歷經(jīng)真空和大氣環(huán)境，速度變化大，高度下降多，可調(diào)控制變量少，精度要求高。面對(duì)高精度再入制導(dǎo)的挑戰(zhàn)和要求，許多研究人員付出了巨大的努力，并提供了多種制導(dǎo)方法［1-7］。

數(shù)值預(yù)測(cè)校正制導(dǎo)算法是一種可以在線運(yùn)行制導(dǎo)方法。該算法不需要預(yù)先存儲(chǔ)參考軌跡，而是利用當(dāng)前狀態(tài)和最終目標(biāo)點(diǎn)信息給出制導(dǎo)指令。由于數(shù)值預(yù)測(cè)校正制導(dǎo)算法使用了這種制導(dǎo)邏輯，即使空天飛行器在再入階段遇到較大的擾動(dòng)，偏離了預(yù)先存儲(chǔ)的參考軌跡，制導(dǎo)算法也可以生成一條合適的傾側(cè)角指令，引導(dǎo)空天飛行器沿著新的軌跡飛向目標(biāo)點(diǎn)。

預(yù)測(cè)校正算法被應(yīng)用在許多方面，其可行性已在多個(gè)環(huán)境的仿真中得到了驗(yàn)證，主要包括可重復(fù)使用的運(yùn)載火箭再入，月球進(jìn)入及探月飛行器再入返回，火星進(jìn)入和精確著陸等。研究人員對(duì)預(yù)測(cè)校正制導(dǎo)算法進(jìn)行了有針對(duì)性的優(yōu)化?；A(chǔ)的預(yù)測(cè)校正制導(dǎo)算法由Xue S等在文獻(xiàn)［8］中提出，并通過X-33 測(cè)試驗(yàn)證了算法的有效性。Wang T 等［9］提出了一種基于模糊邏輯的預(yù)測(cè)校正制導(dǎo)方法，利用龍格-庫(kù)塔數(shù)值積分實(shí)現(xiàn)彈道預(yù)測(cè)。Xia Y 等［10］和Zheng Y 等［11］研究了預(yù)測(cè)校正算法在進(jìn)入火星過程中的應(yīng)用，驗(yàn)證了該算法在極低升阻比的火星著陸器上也是可行的。文獻(xiàn)［12］給出了一種基于降階運(yùn)動(dòng)方程的預(yù)測(cè)校正制導(dǎo)算法，降階后的方程利用較少的計(jì)算量，就可以快速生成可行的三維進(jìn)入軌跡。Wang T等［13］考慮了再入過程中的禁飛區(qū)約束，禁飛區(qū)和航路點(diǎn)被轉(zhuǎn)換成一系列的參考點(diǎn)。在參考點(diǎn)上都設(shè)計(jì)了一次傾側(cè)角反轉(zhuǎn)。王光綸［14］首先提出使用總航程而不是剩余航程來進(jìn)行傾側(cè)角指令的迭代，這樣可以避免校正算法的發(fā)散。

本文首先建立了空天飛行器在發(fā)射坐標(biāo)系下的三自由度動(dòng)力學(xué)模型，分析了再入階段的過程約束和終端約束。然后，基于當(dāng)前狀態(tài)和目標(biāo)點(diǎn)的終端狀態(tài)，使用在線的預(yù)測(cè)-校正再入制導(dǎo)算法對(duì)完整的再入軌跡進(jìn)行反復(fù)的數(shù)值迭代計(jì)算，得到所需的傾側(cè)角幅值指令，達(dá)到制導(dǎo)精度，同時(shí)，通過對(duì)航程的在線計(jì)算，給出了傾側(cè)角的反轉(zhuǎn)邏輯。此外，為了避免再入過程后期剩余航程接近零而導(dǎo)致的迭代發(fā)散，文中使用再入總航程代替剩余航程進(jìn)行迭代；為了避免迭代算法自身原因?qū)е碌牟皇諗浚闹羞€引入一個(gè)調(diào)節(jié)因子來自適應(yīng)地調(diào)整迭代步長(zhǎng)。最后，本文通過數(shù)值仿真驗(yàn)證了該算法的有效性，為了解決控制指令飽和問題、增強(qiáng)制導(dǎo)律氣動(dòng)不確定性下的魯棒性引入了分段方法，仿真結(jié)果表明，再入制導(dǎo)算法工作正常、性能良好、能夠?qū)崿F(xiàn)精度要求。

2 空天飛行器再入問題建模

2.1 發(fā)射坐標(biāo)系中的動(dòng)力學(xué)方程

典型的再入及返場(chǎng)過程如圖1所示。

圖1 典型的再入返場(chǎng)過程剖面Fig.1 Typical re-entry and back process

三維動(dòng)力學(xué)方程在發(fā)射系中建立，位置和速度信息均被分解成了三個(gè)方向上的分量來表示，發(fā)射坐標(biāo)系O-xyz定義如下：O為空天飛行器發(fā)射點(diǎn)，也是空天飛行器返回原著陸場(chǎng)的目標(biāo)點(diǎn)。Ox在水平面內(nèi)指向射向，Oy在豎直平面內(nèi)指向上方，與Ox垂直，Oz由右手定則確定。式（1）～（6）給出了空天飛行器的動(dòng)力學(xué)方程。

式中，x，y，z是位置坐標(biāo)，Vx，Vy，Vz是速度在發(fā)射坐標(biāo)系中的三個(gè)分量，X，Y，Z是氣動(dòng)力的分量，F(xiàn)x，F(xiàn)y，F(xiàn)z是離心力的分量，F(xiàn)gx，F(xiàn)gy，F(xiàn)gz是歌氏慣性力的分量，g是重力加速度，r是飛行器的地心距。

X，Y，Z的計(jì)算方法由式（7）～（9）給出。

式中，D是阻力，L是升力，θ是彈道傾角，ψ是彈道偏角。θ和ψ可以通過Vx，Vy，Vz求得，σ是傾側(cè)角。

一般來說，再入過程的控制變量是攻角α和傾側(cè)角σ，而攻角α的值一般由事先設(shè)定好的α-V剖面給出。預(yù)測(cè)校正制導(dǎo)算法能夠根據(jù)給定的傾側(cè)角預(yù)測(cè)后續(xù)航程，并基于此進(jìn)行迭代得到滿足航程要求的傾側(cè)角。在本文中，傾側(cè)角是唯一的控制變量。

2.2 再入的過程約束和終端約束

典型的再入過程約束由式（10）～（12）給出。

式中，˙是熱流率，q是動(dòng)壓，n是過載系數(shù)；、qmax、nmax即是對(duì)應(yīng)的熱流率、動(dòng)壓、和過載的約束。kQ是與飛行器本身相關(guān)的常數(shù)，Vc=是一個(gè)無(wú)量綱化參數(shù)，L是升力，D是阻力。

對(duì)于末端需要滿足的約束，終端的高度，速度以及航程需要滿足式（13）。

其中高度和速度又可以用能量來進(jìn)行統(tǒng)一表示：

因此式（13）表示的終端狀態(tài)約束還有一個(gè)等效的表達(dá)形式為：

式中，ef表示終端能量，sf是再入過程的總航程，eTAEM和sTAEM是相應(yīng)的能量和航程的約束。

3 預(yù)測(cè)校正制導(dǎo)算法

預(yù)測(cè)-校正制導(dǎo)方法（也稱預(yù)測(cè)制導(dǎo)法）是以消除實(shí)際飛行軌跡的預(yù)測(cè)落點(diǎn)和期望落點(diǎn)之間的偏差為目的的制導(dǎo)方法。該方法的基本思想是利用機(jī)載計(jì)算機(jī)在線預(yù)測(cè)飛行軌跡的終端點(diǎn)，并將求解出的終端點(diǎn)狀態(tài)與理想狀態(tài)比較得出預(yù)測(cè)終端誤差，制導(dǎo)系統(tǒng)根據(jù)預(yù)測(cè)終端誤差校正制導(dǎo)指令，使得飛行軌跡的預(yù)測(cè)終端誤差為零，其工作原理如圖2所示。預(yù)測(cè)-校正制導(dǎo)方法按照軌跡預(yù)測(cè)方式又可以分為解析法和數(shù)值法。本文重點(diǎn)討論數(shù)值預(yù)測(cè)校正制導(dǎo)方法。

圖2 預(yù)測(cè)-校正制導(dǎo)算法流程Fig.2 Predictor-corrector algorithm process

3.1 約束轉(zhuǎn)化與再入走廊

在再入過程中，飛行器始終要遵守約束條件。但是，在全程計(jì)算和判斷所有的約束是否都能滿足是難以實(shí)現(xiàn)的，這會(huì)極大地加大運(yùn)算量。事實(shí)上，由于各個(gè)約束中均含有傾側(cè)角，可以把過程中的約束都轉(zhuǎn)化為傾側(cè)角的約束，若設(shè)置平衡滑翔段，則還需要滿足平衡滑翔約束，利用準(zhǔn)平衡滑翔條件（QEGC）可以將高度和速度約束轉(zhuǎn)化為傾側(cè)角約束。

圖3 過程約束下的再入走廊Fig.3 Re-entry corridor under constraints

如圖3所示，在H-V剖面上，高度的上下界由各項(xiàng)約束決定，高度的上界由QEGC 條件決定，高度的下界由過載、動(dòng)壓、熱流約束的最大值決定，即：

再入過程中的約束決定了傾側(cè)角幅值的邊界，從圖4中可以看出，隨著傾側(cè)角幅值的增大，再入走廊也會(huì)收窄，為了不使平衡滑翔約束與熱流約束曲線相交并保留足夠的裕度，需要給傾側(cè)角設(shè)置合理的上界。

圖4 不同傾側(cè)角下QEGC邊界的變化情況Fig.4 QEGC’s changing under different bank angle

事實(shí)上，根據(jù)QEGC條件：

若已知r和V，就可以求出對(duì)應(yīng)的σ值。根據(jù)已經(jīng)得到的Hup和Hdown，可以求得對(duì)應(yīng)的地心距rup和rdown，對(duì)應(yīng)的傾側(cè)角邊界由式（18）給出：

所以，控制指令σ需要滿足：

其中，σi是由預(yù)測(cè)環(huán)節(jié)的迭代算法給出的傾側(cè)角幅值。需要注意，由于|σ|決定了飛行器在再入過程中的變向能力，|σ|min(V)的值不宜太小，通常取5°～15°。

3.2 縱向制導(dǎo)邏輯

如前文所述，數(shù)值預(yù)測(cè)校正制導(dǎo)算法的最終目的是生成傾側(cè)角指令。而縱向制導(dǎo)部分的任務(wù)是確定傾側(cè)角的幅值。

3.2.1 航程預(yù)測(cè)環(huán)節(jié)

傳統(tǒng)的再入制導(dǎo)研究往往在再入點(diǎn)建立坐標(biāo)系，使用剩余航程來迭代制導(dǎo)指令。然而，在再入過程接近終端時(shí)，剩余航程的值接近于0，這很可能會(huì)造成迭代發(fā)散，使算法無(wú)法獲得傾側(cè)角指令。同時(shí)，對(duì)剩余的導(dǎo)數(shù)進(jìn)行積分來計(jì)算剩余航程的過程會(huì)帶來巨大的計(jì)算量。為了解決這一問題，本文采用再入總航程來代替剩余航程。并且基于發(fā)射坐標(biāo)系到地心坐標(biāo)系的經(jīng)緯度位置信息，提出了一種不進(jìn)行積分直接計(jì)算總航程新方法。

式（1）～（6）給出了六個(gè)狀態(tài)變量的導(dǎo)數(shù)，通過對(duì)這六個(gè)方程的積分，我們可以得到空天飛行器在每個(gè)時(shí)刻的位置和速度信息。再入航程定義在飛行過程的縱平面內(nèi)。因此，通過發(fā)射坐標(biāo)系到地心坐標(biāo)系的轉(zhuǎn)換，可以實(shí)時(shí)得到飛機(jī)的經(jīng)度和緯度。通過給定一個(gè)初始的傾側(cè)角，預(yù)測(cè)校正算法的預(yù)測(cè)環(huán)節(jié)就可以預(yù)測(cè)飛行器的最終落點(diǎn)，再入航程可以用式（20）計(jì)算。

式中，λ0、φ0是再入點(diǎn)的經(jīng)緯度，λf、φf是由算法預(yù)測(cè)環(huán)節(jié)得到的目標(biāo)點(diǎn)的經(jīng)緯度，是對(duì)應(yīng)的地心角。

3.2.2 指令校正環(huán)節(jié)

在航程預(yù)測(cè)環(huán)節(jié)獲得了在給定傾側(cè)角下的航程之后，指令校正環(huán)節(jié)通過牛頓迭代來求解符合要求的傾側(cè)角指令，如式（21）所示。

在式（21）中，a是下山因子，一般情況下，為了保證計(jì)算效率，需要給它分配一個(gè)合適的值。然而，當(dāng)接近真正的解時(shí)，如果a的值太大，迭代可能陷入死循環(huán)，在真解附近振蕩。迭代中的另一個(gè)問題是目標(biāo)函數(shù)fi在局部的單調(diào)性可能與它的整體單調(diào)性不一致，由于局部的單調(diào)性異常，迭代可能會(huì)反向進(jìn)行從而偏離正確的解，這也會(huì)產(chǎn)生發(fā)散。為了解決這一問題，本文給出一種自適應(yīng)迭代算法，當(dāng)式（22）滿足時(shí)，a將縮減一定倍數(shù)。

式（21）中，sall是根據(jù)再入點(diǎn)和目標(biāo)落點(diǎn)的經(jīng)緯度求取的再入航程的參考值，具體求解方法由式（23）給出。

式中，λfinal、φfinal是目標(biāo)落點(diǎn)的經(jīng)緯度，Δall是sall對(duì)應(yīng)的地心角。

式（20）中，目標(biāo)函數(shù)fi的導(dǎo)數(shù)˙一般來說難以解析求取，因此在實(shí)際應(yīng)用中可以采用差分來代替。

迭代算法的流程如圖5所示。

3.3 側(cè)向制導(dǎo)邏輯

側(cè)向制導(dǎo)就是要通過合理的定義傾側(cè)角反轉(zhuǎn)邏輯來實(shí)現(xiàn)，在本節(jié)中反轉(zhuǎn)邏輯的設(shè)計(jì)是通過定義橫程和橫程邊界來實(shí)現(xiàn)。因此本文中的側(cè)向制導(dǎo)律的設(shè)計(jì)就是定義橫程并設(shè)計(jì)橫程邊界。橫程邊界的設(shè)計(jì)原則為使再入軌跡滿足終端位置約束，又不至于使側(cè)傾反轉(zhuǎn)過于頻繁。對(duì)于特定飛行器，可以通過多次仿真實(shí)驗(yàn)獲得適當(dāng)?shù)臋M程邊界參數(shù)。橫程有多種定義方法，對(duì)于不同定義的橫程，需要設(shè)計(jì)不同的橫程邊界。

為了定義橫程，首先在橫向定義橫向剩余航程，由式（24）給出。

橫程的定義如式（25）。

這種定義方法在傾側(cè)角變號(hào)時(shí)，橫程可以很快的響應(yīng)，有較好的控制效果，航向角誤差Δφ是以目標(biāo)落點(diǎn)為基準(zhǔn)，可以保證軌跡會(huì)逐漸趨向于目標(biāo)落點(diǎn)。

橫程的上下邊界設(shè)計(jì)為：

其中，k1、k2為可調(diào)參數(shù)。因此，傾側(cè)角符號(hào)翻轉(zhuǎn)邏輯為：

因此當(dāng)前使用傾側(cè)角指令變?yōu)椋?/p>

綜合縱向與側(cè)向制導(dǎo)，得到完整的預(yù)測(cè)校正制導(dǎo)邏輯如圖6所示。

3.4 分段制導(dǎo)邏輯

全程的預(yù)測(cè)-校正制導(dǎo)在標(biāo)稱情況下一般有很好的精度，但對(duì)于長(zhǎng)程的再入問題，氣動(dòng)不確定性在預(yù)測(cè)環(huán)節(jié)中的積累會(huì)影響傾側(cè)角幅值迭代的準(zhǔn)確性。為了進(jìn)一步提升算法應(yīng)對(duì)氣動(dòng)不確定性的性能，本文中的分段預(yù)測(cè)校正制導(dǎo)方法首先將再入過程分為兩大段，分別是初始下降段和預(yù)測(cè)校正段，由于再入段的仿真再入點(diǎn)位于再入走廊之外，需要經(jīng)歷一個(gè)高度迅速下降的階段才能進(jìn)入再入走廊，因此在再入段前期規(guī)劃初始下降段。

圖6 預(yù)測(cè)-校正制導(dǎo)邏輯流程Fig.6 Predictor-corrector guidance process

初始下降段下沉率大，且由于再入飛行器尚未進(jìn)入大氣環(huán)境，氣動(dòng)力的控制能力有限?；诔跏枷陆刀蔚倪@些物理特性，在該階段設(shè)計(jì)比較精確的傾側(cè)角方案，或通過改變這一過程中的傾側(cè)角來控制飛行軌跡，但效果并不明顯。所以，為了避免不必要的計(jì)算量，設(shè)置初始下降段的傾側(cè)角為一給定常值。該值由迭代算法確定，要求在使再入軌跡位于熱流約束曲線之上的前提下使傾側(cè)角盡量大，本文中設(shè)計(jì)對(duì)象以45°傾側(cè)角飛行直至滿足式（29）。

之后進(jìn)入預(yù)測(cè)校正段，在進(jìn)行傳統(tǒng)的再入段軌跡規(guī)劃時(shí)，往往根據(jù)物理過程把再入分為許多階段，但這種分段方式過于復(fù)雜，不利于快速實(shí)現(xiàn)。預(yù)測(cè)校正制導(dǎo)因其只需要初始和終端狀態(tài)即可進(jìn)行計(jì)算的特性，使其在標(biāo)稱情況下運(yùn)行時(shí)兼具了提供標(biāo)稱軌跡的功能，因此考慮在標(biāo)稱預(yù)測(cè)校正仿真的軌跡上選取若干航路點(diǎn)，以航路點(diǎn)為各段預(yù)測(cè)校正的終點(diǎn)來進(jìn)行新的分段預(yù)測(cè)校正制導(dǎo)仿真。

4 仿真分析與制導(dǎo)性能評(píng)估

本文研究空天飛行器的再入階段，初始條件由表1給出。表2給出了目標(biāo)落點(diǎn)的相關(guān)信息。

表1 再入點(diǎn)初始條件Table 1 The initial conditions of the re-entry point

表2 目標(biāo)落點(diǎn)相關(guān)信息Table 2 The conditions of the drop point

4.1 全程預(yù)測(cè)-校正制導(dǎo)仿真

由于每次迭代都要運(yùn)行預(yù)測(cè)環(huán)節(jié)，為了減少運(yùn)算時(shí)間，預(yù)測(cè)環(huán)節(jié)算法的步長(zhǎng)比狀態(tài)環(huán)節(jié)的步長(zhǎng)要長(zhǎng)，但是如果預(yù)測(cè)環(huán)節(jié)的狀態(tài)變化很大（狀態(tài)導(dǎo)數(shù)值偏大），算法會(huì)自動(dòng)減小步長(zhǎng)以規(guī)避誤差。

如圖7～8所示：制導(dǎo)算法能夠滿足速度、高度等終端約束。從圖9可以看出，在飛行的前期，傾側(cè)比較穩(wěn)定。在再入過程的后期，隨著橫程邊界的縮小和算法給出指導(dǎo)指令頻率的增加，傾側(cè)角開始反轉(zhuǎn)。在表3中，我們可以發(fā)現(xiàn)落點(diǎn)在誤差范圍內(nèi)。

圖7 再入軌跡高度-時(shí)間Fig.7 Altitude-time of the re-entry phase

4.2 分段預(yù)測(cè)-校正制導(dǎo)仿真

在上節(jié)的仿真結(jié)果中，飛行末端控制指令會(huì)出現(xiàn)飽和，且在氣動(dòng)不確定性存在的情況下制導(dǎo)精度下降明顯。

圖8 再入軌跡速度-時(shí)間Fig.8 Velocity-time of the re-entry phase

圖9 再入軌跡傾側(cè)角-時(shí)間Fig.9 Bank angle-time of the re-entry phase

圖10 再入軌跡地面軌跡-時(shí)間Fig.10 Subsatellite track-time of the re-entry phase

表3 落點(diǎn)誤差Table 3 Error of drop point condition

本節(jié)在標(biāo)稱預(yù)測(cè)校正仿真的軌跡上選取若干航路點(diǎn)，以航路點(diǎn)為各段預(yù)測(cè)校正的終點(diǎn)來進(jìn)行新的分段預(yù)測(cè)校正制導(dǎo)仿真，仿真彈道參數(shù)由圖11～14給出。

圖12 再入軌跡速度-時(shí)間Fig.12 Velocity-time of the re-entry phase

圖13 再入軌跡傾側(cè)角-時(shí)間Fig.13 Bank angle-time of the re-entry phase

圖14 再入軌跡攻角-時(shí)間Fig.14 Angle of attack-time of the re-entry phase

從圖15中可以看到，將再入過程分段后有效的緩解了末端控制指令飽和的問題，飽和時(shí)間明顯后移且不再有震蕩，落點(diǎn)誤差也被控制在3km以內(nèi)。

圖15 再入軌跡落點(diǎn)誤差Fig.15 Drop point error of the re-entry phase

5 結(jié)束語(yǔ)

通過對(duì)仿真結(jié)果的分析，可以發(fā)現(xiàn)該算法在空天飛行器再入階段工作良好。作為對(duì)數(shù)值預(yù)測(cè)校正方法的改進(jìn)，本文設(shè)置航路點(diǎn)對(duì)再入過程進(jìn)行分段，采用自適應(yīng)迭代算法來提高預(yù)測(cè)校正制導(dǎo)算法迭代環(huán)節(jié)的收斂性和精度，有效地抑制了迭代的發(fā)散。在未來的工作中，將考慮把深度學(xué)習(xí)應(yīng)用于航程預(yù)測(cè)環(huán)節(jié)，在使用本套算法獲取大量仿真數(shù)據(jù)的基礎(chǔ)上，利用AI 學(xué)習(xí)獲取傾側(cè)角與航程之間的輸入輸出規(guī)則，用基于數(shù)據(jù)模型的傾側(cè)角-剩余航程映射代替原有的預(yù)測(cè)環(huán)節(jié)，由此進(jìn)一步提高制導(dǎo)算法的計(jì)算效率。