李天任，雷建長，王宇航，周 華，黃 佩

0 引 言

近年來，國內(nèi)外學(xué)者對高超聲速飛行器再入軌跡優(yōu)化與制導(dǎo)開展了大量研究：一類是傳統(tǒng)的返回艙式低升阻比再入飛行器，常采用傾側(cè)角反轉(zhuǎn)的軌跡優(yōu)化和制導(dǎo)方法；另一類是近年來發(fā)展的高升阻比滑翔飛行器，一般采用攻角、傾側(cè)角作為控制變量，通過設(shè)計兩者在再入過程中隨時間、速度、高度的變化規(guī)律來滿足過程、終端狀態(tài)約束，同時達(dá)到某一性能指標(biāo)最優(yōu)。

目前，離線進行優(yōu)化的偽譜法是較為成熟的再入體軌跡設(shè)計方法。偽譜法通過在整條軌跡的配點處離散，將再入動力學(xué)微分方程轉(zhuǎn)化為代數(shù)方程，最后采用序列二次規(guī)劃（Sequential Quadratic Programming，SQP）等方法求解非線性規(guī)劃問題。其在單機上進行一次軌跡設(shè)計的計算時間較長，一般用于離線生成標(biāo)稱軌跡。

為了發(fā)揮高升阻比飛行器的在線機動能力，水尊師[1]和肖紅[2]等人發(fā)展了基于航路點的局部預(yù)測校正制導(dǎo)方法，但其離線設(shè)計的航路點是固定的，削弱了飛行器軌跡的變更、再生成能力?；贖arpold[3]的航天飛機再入軌跡設(shè)計思路，周浩[4]、李強[5]、劉運鵬[6]等人采用常值、分段線性等攻角設(shè)計方案，用于軌跡的初步設(shè)計和在線預(yù)測。另一方面，胡錦川[7]的平穩(wěn)滑翔軌跡以及張冉[8]均要求攻角曲線連續(xù)且可導(dǎo)，在使軌跡平穩(wěn)的同時降低飛行器姿態(tài)控制器的設(shè)計要求。

下降段是再入飛行器進入大氣層后的初始飛行階段，也被稱為引入段。為了在下降段軌跡的設(shè)計中簡化控制量曲線的設(shè)計，同時增強控制量剖面的光滑度，易于與平衡滑翔段軌跡光滑對接，除了現(xiàn)有的常值、分段線性以及復(fù)合指數(shù)函數(shù)的剖面設(shè)計方案外，有必要引入新的函數(shù)形式。本文采用改進的復(fù)合三角函數(shù)，其具有幅值可調(diào)、波峰的位置和峰形可變、導(dǎo)數(shù)連續(xù)且光滑等優(yōu)點，適合在控制量剖面設(shè)計中采用。

1 再入動力學(xué)模型

針對無動力再入飛行器，在位置系下建立動力學(xué)模型，O-XYZ為地心固連坐標(biāo)系，將其繞Z軸正向旋轉(zhuǎn)經(jīng)度θ，再繞新的 Y軸旋轉(zhuǎn)緯度φ得到位置坐標(biāo)系O-xyz，如圖1所示。

圖1 三自由度運動參數(shù)示意Fig.1 3D Motion Parameters Diagram

采用圓球無旋地球假設(shè)，飛行器動力學(xué)方程[9]為

式中 r，θ，φ，V，γ，ψ分別為地心距、經(jīng)度、緯度、速度、彈道傾角、航向角；g為重力加速度；m為飛行器質(zhì)量；σ為傾側(cè)角；L，D分別為升力和阻力，且，

式中 ρ為大氣密度；S為參考面積； CL， CD分別為升力系數(shù)、阻力系數(shù)。

m，S，CL, CD均參考國外滑翔飛行器CAV-H，m=907 kg，S=0.4839 m2。對升力系數(shù) CL和阻力系數(shù) CD進行擬合，得到兩者對攻角α的關(guān)系式[10]：

在再入飛行器初始下降段，常采用零傾側(cè)角、大攻角拉起的方案，在飛行器狀態(tài)參數(shù)滿足一定條件時，切換到平衡滑翔階段，隨后采用平衡滑翔段的攻角、傾側(cè)角剖面。本文采用的切換條件為

式中1K為彈道傾角的權(quán)系數(shù)；ε為一小量。

采用下降段零傾側(cè)角、單波峰攻角剖面，平衡滑翔段的攻角剖面由擬平衡滑翔條件（Qusi Equilibrium Glide Condition，QEGC）推導(dǎo)得到。根據(jù)擬平衡滑翔條件：

根據(jù)當(dāng)前狀態(tài)已知的 V ,r,g，聯(lián)立式（2）、式（3）、式（6），得到擬平衡滑翔攻角 αQEGC為

式中 A，B為式（3）中攻角的一階項系數(shù)和常數(shù)項。

考慮到氣動系數(shù)誤差，在 αQEGC上反饋彈道傾角的誤差項 Δ α= K2(γ -γ*)，其中，γ*為平衡滑翔段跟蹤的彈道傾角， K2為反饋增益系數(shù)。由于本文研究的重點不是平衡滑翔段，可取γ*=-0.05°，從而有：

仿真結(jié)果表明，該方法能夠在考慮各項偏差的情況下生成擬平衡滑翔軌跡。

2 下降段控制量剖面設(shè)計方法

2.1 下降段控制量設(shè)計思路

下降段攻角剖面常采用最大攻角、分段線性等攻角方案[3～5]，在與滑翔段的對接處難以連續(xù)和可導(dǎo)，使軌跡平穩(wěn)性降低、同時對姿態(tài)控制器的設(shè)計提出了很高要求。此外，大量的彈道仿真表明不連續(xù)或不可導(dǎo)的控制量剖面常導(dǎo)致對接點過載和熱流峰值較大，難以滿足過程約束。

考慮在再入下降段采用改進的復(fù)合三角函數(shù)，以增強對接點光滑性、降低最大過載和熱流峰值。

2.2 單波峰攻角剖面生成方法

再入飛行器下降到一定高度，采用如下攻角時序：

式中mα為飛行器最大可用攻角。根據(jù)CAV-H的氣動系數(shù)表[10]，mα取20°，()Ft采用時間t的二次多項式函數(shù)：

通過調(diào)節(jié)多項式系數(shù) ai(i = 0 ,1,2)來調(diào)整波長、波形和波峰位置，生成廣義正弦函數(shù)簇，見圖2。

圖2 不同波長、波峰位置的三角函數(shù)簇Fig.2 Trigonometric Functions with Diffierent Wavelength and Peak Position

2.3 搜索過渡點光滑的攻角剖面

為了增強下降段與平衡滑翔段過渡點處的攻角剖面光滑性，需要優(yōu)化控制下降段攻角剖面形狀的ai(i = 0 ,1,2)3個參數(shù)。本文采用改進的差分進化算法（Differential Evolution，DE）[11]進行參數(shù)尋優(yōu)。目標(biāo)函數(shù)見式（11）。

式中 j表示第j個積分區(qū)間，積分區(qū)間個數(shù) N = tf/h ，tf，h為積分結(jié)束時間和積分步長。本文在Storn[11]的交叉操作之后增加拉回步驟，拉回不在解空間的新個體，以增強算法的魯棒性。

3 數(shù)值仿真及分析

3.1 攻角剖面光滑性驗證

CAV-H 再入點參數(shù)取 H=70 km，V=7.5 km/s，γ=-2°。下降段采用復(fù)合三角函數(shù)的攻角剖面，平衡滑翔段采用擬平衡滑翔攻角，以控制量光滑性為目標(biāo)函數(shù)進行參數(shù)優(yōu)化，得到如圖3所示的攻角曲線。

圖3 分段的攻角曲線Fig.3 Segment Curve of Angle of Attack

檢查攻角剖面光滑性，見圖4、圖5。結(jié)果表明，攻角剖面在下降段前期與傳統(tǒng)最大正攻角方案一致，而在過渡點附近與平衡滑翔段光滑銜接，攻角、攻角一階導(dǎo)數(shù)均連續(xù)。攻角剖面一階導(dǎo)數(shù)大小在過渡點附近不超過1 (°)/s，二階導(dǎo)數(shù)大小在過渡點附近不超過2 (°)/s2，基本滿足張冉[8]對攻角光滑性的要求。

圖4 攻角一階導(dǎo)數(shù)Fig.4 First Derivative of Angle of Attack

圖5 攻角二階導(dǎo)數(shù)Fig.5 Second Derivative of Angle of Attack

采用優(yōu)化后的攻角剖面對CAV-H飛行器的再入下降段和平衡滑翔段進行數(shù)值仿真，見圖6、圖7。

圖6 高度-時間曲線Fig.6 Hight-Time Curve

圖7 彈道傾角及其導(dǎo)數(shù)曲線Fig.7 Curves ofγand its First Derivative

結(jié)果表明，再入飛行器能夠平滑地進入擬平衡滑翔狀態(tài)，彈道傾角由初始值過渡到零附近的過程連續(xù)且二階可導(dǎo)，一階導(dǎo)數(shù)絕對值不超過0.05 (°)/s。

3.2 不同攻角剖面方案的光滑性對比

采用式（11）中的目標(biāo)函數(shù)f作為控制量曲線的光滑性能指標(biāo)，對單波峰、分段線性、分段常值3種攻角剖面進行光滑性對比，見圖8和表1。

圖8 攻角曲線對比Fig.8 Comparision of Diffierent AOA Curves

表1 不同攻角剖面方案的光滑性對比Tab.1 Comparison of Smoothness with Diffierent AOA Profiles

結(jié)果表明，復(fù)合三角函數(shù)提高了攻角剖面的光滑性，使得曲線在分段過渡點處保持了較好的連續(xù)性。即使不對參數(shù)進行優(yōu)化，單波峰攻角剖面相比另外兩者有先天性的光滑性優(yōu)勢。此外，以上 3種剖面參數(shù)化方法的參數(shù)數(shù)量是相同的（均為3個），因此在軌跡優(yōu)化時的復(fù)雜度基本一致。

3.3 不同攻角剖面方案的軌跡特性對比

為了檢驗單波峰攻角方案的優(yōu)越性，將仿真結(jié)果與 Harpold[3]分段線性攻角方案、周浩[4]分段常值攻角方案得到的軌跡進行對比。采用胡錦川[7]提出的性能指標(biāo)來評價軌跡平滑性：

式中 aε為縱向加速度，aε=V γ˙；tf為積分時間，tf=160 s。檢查最大過載 nmax、峰值熱流 Q˙max、過渡點動壓q等軌跡特性，見表2。

表2 3種攻角剖面的軌跡特性對比Tab.2 Comparison of Trajectory Characteristics with Three Diffierent AOA Profiles

結(jié)果表明，經(jīng)過優(yōu)化的單波峰攻角剖面得到的軌跡最平滑，其他軌跡性能指標(biāo)與Harpold[3]攻角方案差別不大，而分段常值攻角剖面的軌跡特性最差。分析可知，分段線性攻角方案在α-t和α-V剖面均連續(xù)但不可導(dǎo)，導(dǎo)致軌跡的平滑性在過渡點處下降，而分段常值攻角方案難以保證軌跡平滑，無法滿足再入軌跡精細(xì)化設(shè)計的任務(wù)需求。

4 結(jié)束語

本文通過引入新的基函數(shù)，對原有控制量剖面進行了改進。新方法在保證參數(shù)數(shù)量不增加的情況下改善了控制量剖面的光滑性，基本滿足胡錦川[7]、張冉[8]對控制量曲線的光滑性要求。通過差分進化算法對基函數(shù)進行參數(shù)優(yōu)化，使得過渡點處的攻角曲線光滑性較好，提高了軌跡平滑性，下降段的熱流、過載等環(huán)境也有一定程度改善。

復(fù)合三角函數(shù)具有在前、后端平緩變化的特點，適合在高速、大動壓情況下控制飛行器的姿態(tài)角，可以在其他場景中加以進一步的應(yīng)用。

[1] 水尊師, 周軍, 葛致磊. 基于高斯偽譜方法的再入飛行器預(yù)測校正制導(dǎo)方法研究[J]. 宇航學(xué)報, 2011, 32(6)： 1249-1255.

Shui Zunshi, Zhou Jun, Ge Zhi-lei. On-line predictor-corrector reentry guidance law based on gauss pseudospectral method[J]. Journal of Astronautics, 2011, 32(6)： 1249-1255.

[2] 肖紅, 王勇, 國海峰, 王征, 臧旭. 基于航路點的高超聲速飛行器預(yù)測校正制導(dǎo)律設(shè)計[J]. 電光與控制, 2013, 20(10)： 33-37.

Xiao Hong, Wang Yong, Guo Haifeng, Wang Zheng, Zang Xu. Design of waypoint based predictor-corrector guidance law for hypersonic vehicles[J]. Electronics Optics＆Control, 2013, 20(10)： 33-37.

[3] Harpold J C, Graves C A. Shuttle entry guidance[J]. The Journal of the Astronautical Sciences, 1979, 37(3)： 239-268.

[4] 周浩, 陳萬春, 殷興良. 高超聲速飛行器多段彈道優(yōu)化設(shè)計[J]. 宇航學(xué)報, 2008, 29(4)： 1205-1209.

Zhou Hao, Chen Wanchun, Yin Xingliang. Multi-phase trajectory optimization for a hypersonic vehicle[J]. Journal of Astronautics, 2008,29(4)： 1205-1209.

[5] 李強. 高超聲速滑翔飛行器再入制導(dǎo)控制技術(shù)研究[D]. 北京： 北京理工大學(xué), 2015.

Li Qiang. Study on reentry guidance and control method for hypersonic glide vehicle[D]. Beijing： Beijing Institute of Technology, 2015.

[6] 劉運鵬, 李伶. 高超聲速導(dǎo)彈高空再入減速段彈道優(yōu)化設(shè)計[J]. 航天控制, 2010.11, 28(6)： 3-6.

Liu Yunpeng, Li Ling. The optimization design for high-altitude reentry deceleration trajectory of hypersonic missile[J]. Aerospace Control, 2010,28(6)： 3-6.

[7] 胡錦川, 陳萬春. 高超聲速飛行器平穩(wěn)滑翔彈道設(shè)計方法[J]. 北京航空航天大學(xué)學(xué)報, 2015, 41(8)： 1464-1475.

Hu Jinchuan, Chen Wanchun. Steady glide trajectory planning method for hypersonic reentry vehicle[J]. Journal of Beijing University of Aeronautics and Astronautics, 2015, 41(8)： 1464-1475.

[8] 張冉, 李惠峰. 基于平衡滑翔的升力再入攻角設(shè)計分析及應(yīng)用[J].宇航學(xué)報, 2012, 33(11)： 1557-1563.

Zhang Ran, Li Huifeng. Angle of attack operation zone analysis and application based on equilibrium glide condition for lifting reentry[J].Journal of Astronautics, 2012, 33(11)： 1557-1563.

[9] Vinh N X, Busemann A, et al. Hypersonic and planetary entry flight mechanics[M]. Ann Arbor MI： Univ of Michigan Press, 1980.

[10] Zhang K, Chen W. Reentry vehicle constrained trajectory optimization[R]. AIAA-2011-2231, 2011.

[11] Storn R, Price K. Differential evolution-a simple and efficient heuristic for global optimization over continuous spaces[J]. Journal of Global Optimization, 1997, 11 (4)： 341-359.