滑翔式高超声速目标可达区域计算方法*

张凯，熊家军

(空军预警学院 a. 研究生管理大队; b. 四系，湖北 武汉 430019)

滑翔式高超声速目标可达区域计算方法*

张凯a，熊家军b

(空军预警学院 a. 研究生管理大队; b. 四系，湖北 武汉 430019)

针对滑翔式高超声速目标能够利用强大的横向机动能力形成较大的可达区域的特点，着眼于系统阐述这类问题的解决方案，研究了几种可达区域计算方法。首先建立了滑翔式高超声速目标的三自由度动力学模型，并给出再入过程必须满足的相关约束。然后详细论述了轨迹优化、再入走廊规划、常倾侧角和近似椭圆拟合4种滑翔式高超声速目标可达区域计算方法，给出了几种方法具体的实现流程，并对方法的不足进行改进。通过算例仿真，对几种可达区域计算方法的优缺点与应用特点进行了比较分析。最后，以轨迹优化方法为例，分析了可达区域动态变化情况，验证方法的完备性和动态适应性。

高超声速飞行器；再入滑翔；可达区域；伪谱法；再入走廊；常值倾侧角；椭圆近似

0 引言

鉴于高超声速武器巨大的军事价值，美军在多个科技规划文件中均把“高超声速”列为未来改变游戏规则的技术方向[1]。X-51A与AHW的试验成功标志着美军高超声速飞行器相关技术的基本成熟，这对我国构建战略防御体系带来了巨大挑战。以HTV-2，AHW为代表的滑翔式高超声速目标因其高升阻比外形设计，相对于传统的惯性目标，此类飞行器具有较强的机动能力，能够利用强大的横向机动能力形成一个较大的打击区域，因此落点预报的可能性与弹道导弹相比大大降低。但是，由于飞行器再入过程中飞行状态跨度大、飞行环境恶劣、动力学特性复杂，必须满足各种复杂的约束条件，故再入轨迹并非无章可循。因此，防御方在获得空袭目标信息后，如能快速预报目标打击范围，这对尽早制定拦截打击方案意义重大。

针对计算滑翔式高超声速目标覆盖范围这一问题，早期Ngo将最优化方法应用于飞行器控制部件失效条件下的可达区域分析中，得到了一些有参考意义的结果[2]；Burkhardt以X-38为研究对象，利用一种最优化方法计算了飞行器可达区域，并展开了地面着陆场选择[3]。国内一些学者分别利用伪谱法对再入飞行器的可达区域进行了分析研究[4-7]。相对于利用优化方法计算覆盖范围，Saraf等人另辟途径，他们在阻力-能量(D-E)剖面中考虑了各种约束，研究了可达区域分析方法[8]。国内学者也在该方法的基础上改进了产生可达区域再入边界的插值方式和终端约束处理方法，使得计算结果更加精确[9-10]。近些年这些讨论再入覆盖范围的论文方兴未艾，但大多是站在攻击方角度讨论这类问题，同时，对再入覆盖问题均是基于单一方法的求解。本文着眼于系统阐述这类问题的解决方法，详细介绍了几类典型的再入可达区域计算方案，并对各自的优缺点与应用特点进行了比较分析，以期为滑翔式高超声速目标防御提供一定理论参考。

1 再入动力学建模

在半速度坐标系中建立飞行器三自由度动力学方程[11]，不考虑地球旋转，并令飞行器侧滑角为0。

(1)

式中：r,v,λ,φ,θ,σ为状态参数，分别为飞行器的地心距、速度、经度、纬度、速度倾角、航向角；D为阻力加速度；L为升力加速度；υ为倾侧角。

在再入过程中，考虑以下热流、法向过载、动压和平衡滑翔飞行约束条件：

(2)

(3)

(4)

(5)

2 可达区域计算方法

2.1 轨迹优化方法

针对飞行器再入可达区域问题，通常的做法是通过参数化，将优化问题转化为纵程不定，求横程最大的单参数优化问题。

伪谱法由于在计算效率上的巨大优势，在航天器轨迹优化领域得到了广泛应用[12-13]。hp-自适应伪谱法是一种融合了Radau伪谱法与hp-型有限元法优点的最优问题解决方法，一方面，其利用多项式在离散点(配点)上对状态变量和控制变量进行全局逼近；另一方面，配点处必须严格满足动力学方程中的路径约束和初始、终端条件，从而将微分方程转化为代数约束。同时，这种优化方法采用双层逻辑改善求解过程的效率——在迭代过程中自适应地改变单元长度h和插值多项式阶次p，相比于传统的伪谱法，具有更高的收敛效率和计算精度[14]。因此，本文选择利用hp-自适应伪谱法作为求解优化弹道的方法，其基本原理如图1所示，轨迹优化方法如下：

(1) 根据求解再入落点区域问题，其优化指标为不定纵程对应的最大横程。首先分别取优化指标为J=Lmax和J=Lmin，即计算落点区域的远边界和近边界。然后引入权重系数ω：

Li=Lmin+ω(Lmax-Lmin).

(6)

(2) 然后在确定的纵程Li情况下，取优化指标为J=Rmax，计算确定纵程对应的最大横程，将所有航程终点连接形成的区域即为打击可达区域。

图1 hp-自适应伪谱法原理图Fig.1 Schematic diagram of hp-adaptive pseudospectral method

2.2 再入走廊规划方法

工程上，通常使用再入走廊进行弹道设计，再入走廊规划方法的基本思想就是基于再入走廊的可行区域，在满足初始、终端约束的条件下，规划飞行器再入剖面从而计算落点区域。再入走廊一般在高度-速度(h-v)剖面、阻力-能量(D-E)剖面、阻力加速度-速度(D-v)剖面以及高度-能量(h-E)剖面内描述，不同剖面描述的物理意义根据相关公式可以互相转化。

图2 阻力加速度方案Fig.2 Resistance acceleration of reentry corridor method

如图2所示，本文以阻力加速度-速度(D-v)再入剖面为例，介绍再入走廊剖面规划方法的基本原理：

(1) 描述再入走廊。根据上文介绍的热流、动压、法向过载和平衡滑翔约束描述再入走廊D-v剖面。

(2) 改进再入走廊边界。首先计算初始阻力加速度D(v0)和终端阻力加速度D(vf)，确定再入剖面的起点和终点。再入走廊上边界为最小阻力边界，下边界为最大阻力边界。但阻力边界曲线存在不平滑的变化，可利用三阶多项式拟合改进最大和最小航程剖面。

(3) 设计再入剖面。对介于最大和最小航程剖面之间的航程，通过对两者D-v剖面进行线性插值，近似获得相应的再入剖面。设计一种非均匀线性插值阻力方法，避免等比例线性插值造成实际再入轨迹过于聚集或稀疏。

(7)

式中：Di(v)为第i条再入剖面对应的阻力加速度曲线；Dmax(v)为修正的最大阻力剖面；Dmin(v)为修正的最小阻力剖面；n为设计的再入剖面个数。

(4) 在线跟踪再入剖面。仿真时调整倾侧角大小控制飞行器的阻力加速度等于所选择的飞行剖面确定的阻力加速度，保持倾侧角的符号获得相应方向的最大横程，通过跟踪不同再入剖面得到再入落点区域。

一般而言，为使轨迹尽快偏离再入径向剖面，应当使倾侧角符号为正或为负保证飞行器始终朝一个方向转弯，然而一旦速度偏角超过了发射坐标系的正北或者正南方向，会造成航程回旋。为此文献[15]提出引入倾侧角υ反向策略，当速度方向与初始速度方向夹角大于设定值时改变倾侧角符号，增大横程。

2.3 常值倾侧角方法

一般而言，倾侧角为0时，飞行器航程最大，但横向机动能力最弱；倾侧角增大时，横向机动能力随之增强。因此，可以通过常倾侧角条件下的动力学仿真，分析不同倾侧角对应轨迹的纵程和横程，从而得到飞行器再入落点区域。在整个再入过程中，倾侧角的大小是固定不变的，因此计算速度快且过程简单。

当航向角超过了初始速度方向时，同样会造成航程回旋。为此，本文设计了倾侧角置零策略，一旦速度偏角指向初始速度方向，将倾侧角置零，此后飞行器将沿着垂直于初始速度的方向飞行，飞行器能充分利用空气动力达到増加横程的目的。

2.4 椭圆近似方法

椭圆近似方法的基本思想是利用飞行器再入落点区域近似为椭圆状的特点，建立关于纵程和横程近似符合椭圆分布的假设，通过优化算法获得最大、最小纵程及最大的横程轨迹，从而确定再入覆盖椭圆区域(图3)。该方法的基本方案如下：

(1) 计算椭圆边界点。采用优化方法计算飞行器最大纵程Smax、最小纵程Smin和最大横程Cmax对应的经纬度，其中O是再入点。

(2) 建立椭圆方程。以Cmax为短半轴，O’A为长半轴，O’为椭圆中心建立纵程与横程近似椭圆方程：

(8)

过最小纵程位置点D做直线OA的垂线，垂线截断的椭圆右侧部分即为椭圆近似方法计算的再入落点区域。

图3 可达区域近似椭圆Fig.3 Footprint of ellipse fitting

3 仿真分析

3.1 仿真参数设置

参考洛克希德·马丁公司的CAV-H为仿真实例模型[16]，具体参数如下：

初始参数：飞行器质量为1 500 kg，横截面积为0.605 m2。高度为50 km，速度为5 500 m/s，经纬度为(0°，0°)，速度倾角为0，航向角为-90°。

终端约束：取速度vf=1 500 m/s，高度为hf=20 km。

3.2 可达区域分析

图4为利用GPOPS优化软件包计算的弹道优化方法再入落点区域，GPOPS软件是佛罗里达大学Rao团队利用伪谱法分别开发了最优控制问题求解软件[17]，NASA在2015年上旬免费公开了该软件的第2代产品(GPOPSⅡ)，其简捷性和速度得到了大幅提升。根据图4可知，再入轨迹无法实现零射程，同时，靠近再入点的轨迹发生回旋，这是因为再入过程中受到热流、动压约束影响，飞行器无法实现高度和速度迅速减小。

图4 轨迹优化方法可达区域Fig.4 Footprint of trajectory optimization method

图5为再入走廊规划方法计算结果，其中方法1不考虑倾侧角反转策略，方法2考虑倾侧角反转。由图可知，方法2计算的落点区域比方法1略大，总体上2种方法计算结果基本相同。

图5 再入走廊规划方法可达区域Fig.5 Footprint of reentry corridor method

图6为常倾侧角方法计算结果，其中方法1不考虑倾侧角置零策略，方法2为考虑倾侧角置零策略。由图可知，方法2计算结果明显大于方法1，这是由于一旦速度方向为正北(南)，此时将倾侧角置零后，飞行器能充分利用空气动力达到増程的目的，从而增大横向距离。该方法的主要缺点是，首先未考虑过程约束条件，图7中倾侧角取值范围为[0,70]，当常倾侧角超过50°，再入状态会超过再入走廊下边界，即热流、动压约束不满足；其次，未考虑终端约束，导致再入过程终端不满足高度约束条件。

图6 常值倾侧角方法可达区域Fig.6 Footprint of constant bank angle method

图7 常值倾侧角方法高度-速度曲线Fig.7 Footprint of constant bank angle method constant bank angle

如图8所示为4种可达范围计算方法的计算结果对比，其中方法1为轨迹优化方法，方法2为修正的再入走廊规划方法，方法3为修正的常倾侧角方法，方法4为椭圆近似方法。

由图8可知，方法1与方法4结算结果基本一致，这是因为方法4的基本参数(长半轴和短半轴)均根据方法1中的hp-自适应伪谱法计算得到。由于实际进行落点预测时，很难确切了解敌方再入飞行器的全部参数，因此，相对方法1，利用方法4可在保证一定精度的情况下可以实现快速预测可达区域。

图8 可达区域对比图Fig.8 Comparison between the four methods

方法2与方法3的落点区域范围明显小于方法1，4，这是因为方法2，3 均在各自的变化量取值范围内线性插值，并没有对可行域内的取值进行最优化求解。其次，方法2，3的最小航程小于方法1，4，这是因为方法1，4考虑了航迹回旋的轨迹，方法2，3经过修正后均未产生回旋轨迹。

结合表1中几种方法的计算时间(Matlab仿真时间)可知，再入走廊规划方法求解耗时较短，但求解的可达区域最小，此方法的优势是工程上一般采用该方法规划再入轨迹。轨迹优化方法耗时较长，但求解的可达区域最精确，能够避免可达区域的遗漏。

表1 方法计算时间

为了进一步研究飞行器再入过程中可达区域动态变化情况，以轨迹优化方法为例，对再入过程动态可达区域进行了仿真。图9为再入过程中不同状态节点对应的可达区域，状态节点参数如表2所示。可见，随着再入仿真过程的推进，可达区域不断减小，这是因为速度减小导致机动能力降低；同时，实时可达区域均在前一个状态节点覆盖范围内，说明方法计算的可达区域具有动态适应性。

表2 状态节点参数

图9 再入过程动态可达区域Fig.9 Dynamic footprint of reentry process

4 结束语

本文首先系统论述了现有的几种滑翔式高超声速目标可达区域计算方法，同时，针对方法的不足提出了一些改进。通过对相关方法的仿真，分析了各自的优缺点和应用特点，并对再入过程中可达区域动态变化情况进行了研究。论文验证了几种方法计算的可达区域的完备性和动态适应性，为滑翔式高超声速目标防御提供一定理论参考。

[1] 刘杨, 李继勇, 赵明. 国外高超声速武器技术路线分析及启示[J]. 战术导弹技术, 2015(5):47-50. LIU Yang, LI Ji-yong, ZHAO Ming. Technical Route Analysis and Inspiration of Foreign Hypersonic Weapons [J]. Tactical Missile Technology, 2015(5):47-50.

[2] NGO A, DOMAN D. Footprint Determination for Reusable Launch Vehicles Experiencing Control Effector Failures [C]∥AIAA Guidance, Navigation, and Control Conference and Exhibit，2002.

[3] BURKHARDT J. Impact of Mission Constraints on Optimal Flight Trajectories for the Lifting Body X-38[C]∥Atmospheric Flight Mechanics Conference, 1999.

[4] 李瑜, 崔乃刚. 洲际助推-滑翔导弹可达区域优化[J]. 战术导弹技术,2010(1):40-44. LI Yu, CUI Nai-gang. Reachable Domain Optimization for Intercontinental Boost-Glide Missile [J]. Tactical Missile Technology, 2010 (1):40-44.

[5] 肖翔, 杨业, 郭涛, 等. 升力式再入飞行器两种可达区域计算方法的探讨[J]. 航天控制, 2015, 33(2):39-43. XIAO Xiang, YANG Ye, GUO Tao, et al. Discussion of Two Algorithms for Generating Landing Footprint for Lifting Reentry Vehicles [J]. Aerospace Control, 2015, 33(2),39-43.

[6] 王路, 邢清华. 再入高超声速滑翔飞行器可达区域快速预测方法[J]. 军事运筹与系统工程, 2015, 29(1):40-44. WANG Lu, XING Qing-hua. [J]. Military Operations Research and Systems Engineering, 2015, 29(1):40-44.

[7] 张凯, 郭振云. 基于优化弹道的BPNN预测制导方法[J]. 导弹与航天运载技术,2015, 339(3):5-9. ZHANG Kai, GUO Zhen-yun. BPNN Predictor-Corrector Guidance Based on Optimized Trajectory [J]. Missiles and Space Vehicles,2015, 339(3):5-9.

[8] SARAF A. Landing Footprint Computation for Entry Vehicles[C]∥AIAA Guidance,Navigation,and Control Conference and Exhibit, 2004.

[9] 常松涛, 杨业, 王永骥, 等. 计算再入飞行器可达区域的快速算法[J]. 华中科技大学学报:自然科学版, 2012, 40(7):1-5. CHANG Song-tao, YANG Ye, WANG Yong-ji, et al. A Rapid Algorithm for Generating Landing Footprint for Entry Vehicles [J]. Huazhong University of Science & Technology: Natural Science ed, 2012, 40(7):1-5.

[10] 孙尚, 杨博, 李大伟, 等. 高超声速飞行器动力学规划的落点区域计算[J]. 航天控制, 2011, 29(1):48-53. SUN Shang, YANG Bo, LI Da-wei, et al. The Calculation for Landing Footprint of Hypersonic Vehicle Based on Dynamics Programming [J]. Aerospace Control, 2011, 29(1):48-58.

[11] 赵汉元. 飞行器再入动力学与制导[M]. 长沙：国防科技大学出版社,1997. ZHAO Han-yuan. Vehicle Reentry Dynamic and Guidance[M]. Changsha: National University of Defense Technology Press, 1997.

[12] FAHROO F, ROSS I M. Direct Trajectory Optimization by a Chebyshev Pseudospectral Method [J]. Journal of Guidance, Control and Dynamics, 2002, 25(1):160-166.

[13] FAHROO F, DOMAN D B, NGO A D. Footprint Generation for Reusable Launch Vehicles Using a Direct Pseudospectral Method [C]∥American Control Conference, 2003.

[14] DARBY C L, HAGER W W, RAO A V. An hp-Adaptive Pseudospectral Method for Solving Optimal Control Problems[J]. Optimal Control Applications and Methods, 2011, 32(4):476-502.

[15] 冯必鸣, 聂万胜, 李柯. 再入飞行器可达区域近似算法及地面覆盖研究[J]. 航天控制,2012, 30(6):43-49. FENG Bi-ming, NIE Wan-sheng, LI Ke. Research on Closest-Aproach of Footprint and Coverage for Reentry Vehicle [J]. Aerospace Control, 2012, 30(6):43-49.

[16] PHILLIPS T H. A Common Aero Vehicle (CAV) Model, Description, and Employment Guide [R]. Schafer Corporation for AFRL and AF-SPC, 2003.

[17] RAO A V, BENSON D A, DARBY C, et al. GPOPS, A Matlab Software for Solving Multiple-Phase Optimal Control Problems Using the Gauss Pseudospectral Method [J]. Acm Transactions on Mathematical Software, 2011, 38(1):163-172.

Algorithm of Footprint Area for Gliding Hypersonic Target

ZHANG Kaia, XIONG Jia-junb

(Air Force Early Warning Academy, a. Graduation Management Department; b. No.4 Department, Hubei Wuhan 430019, China)

In view of the characteristics of the large footprint area formed by gliding hypersonic target, four algorithms are researched to describe the solutions of this problem systematically. Firstly, the three degree of freedom dynamics of glide hypersonic target is established, and multiple constraints of reentry process are given. Secondly, four footprint calculation algorithms (trajectory optimization method, reentry corridor planning method, constant bank angle method and ellipse fitting method) of gliding hypersonic target are discussed, the processes of algorithms are given in detail and the bugs are fixed. Then the advantages/disadvantages and characteristics of the algorithms are compared and analyzed by simulation. Finally, taking the trajectory optimization algorithm as an example, the dynamic variation of the footprint area is analyzed, and the completeness and dynamic adaptability of the algorithm are verified.

hypersonic vehicle; reentry glide; footprint area; pseudospectral method; reentry corridor; constant bank angle; ellipse fitting

2016-04-13；

2016-07-21 作者简介：张凯(1990-)，男，湖北红安人。博士生，主要从事预警情报分析。

10.3969/j.issn.1009-086x.2017.02.0010

V412.4；O242；TP391.9

A

1009-086X(2017)-02-0067-07

通信地址：430019 湖北省武汉市江岸区黄浦大街288号研究生管理大队 E-mail：christophzhang@163.com