卢宝刚，傅 瑜，崔乃刚，单文昭

(1.哈尔滨工业大学 航天学院，150001哈尔滨;2.北京宇航系统工程研究所，100076北京)

高超声速飞行器是指飞行马赫数大于或等于5的飞行器，能够依靠其较大的升阻比在大气层内做较长时间的无动力滑翔飞行，具有飞行速度快、突防能力强、飞行距离远等优势.目前，主要有两大类高超声速飞行器，即可重复使用运载器(RLV)与通用航空飞行器(CAV).航天飞机是经过实际工程验证的最典型的RLV，在此之后，美国开展了X-33、X-34、X-37、X-40和X-43等一系列新一代RLV验证机的研究工作，其中X-37B于2010年4月第一次获得了成功飞行试验，2011年3月，第二架X-37B飞行器升空.为了实现低成本的远程快速打击或者物资投送，美国空军司令部对军用太空飞机的概念进行了论证，提出了CAV概念飞行器，并于2010年进行了HTV-2的飞行试验，该飞行器可看作是CAV的延续.

飞行器在再入过程中要受到动压、过载与气动热的限制，与此同时，需保证飞行器到达指定的位置(经度、纬度和高度)和速度(速度的大小和方向)，满足以上约束条件的再入制导问题一直为近年来航空航天领域的研究热点.从20世纪70年代起，以航天飞机为背景的大升阻比再入制导问题研究开始兴起［1］，20世纪90年代以来，新的制导方法层出不穷［2-5］，主要分为两大类:标称轨迹制导和预测制导.标称轨迹制导包括轨迹规划与轨迹跟踪两部分，标称轨迹可以基于不同的剖面进行规划，如阻力加速度剖面［1］、高度剖面［6］，动压剖面［7］等.轨迹跟踪可以分为剖面跟踪［1］及全状态跟踪［8］两种类型.预测制导是根据终端状态的期望值与预测值之差修正控制量，按照对终端状态预测方法的不同，分为解析法和数值法.解析法一般需作大量简化，当某些关键参数被简化处理后，会对再入轨迹影响较大，预测的精度比较低.数值法预测的精度高，在计算机高度发展的今天，数值预测计算量较大的问题也得到了有效解决，所以，数值预测制导受到广泛关注［9-10］.在众多再入制导方法中，美国马歇尔空间飞行中心对它们作了测评［11］，其中，能够在线生成标称轨迹的拟平衡滑翔制导是顺利通过全部测试的再入制导方法.

本文结合拟平衡滑翔条件与数值预测的优势，提出一种高超声速飞行器再入轨迹规划方法，该方法可离线设计标称轨迹，也可在线规划轨迹(受弹载计算水平限制).考虑到攻角对大升阻比飞行器飞行轨迹的影响较大，单纯的依靠调节倾侧角无法满足需求，因此，通过对攻角剖面的设计满足终端速度约束，通过对飞行路径角角剖面的设计满足航程约束，然后利用拟平衡滑翔条件计算倾侧角进行轨迹保持.此时，三维约束再入轨迹规划问题被简化为2个参数的搜索问题，轨迹规划效率高，适用于大升阻比高超声速飞行器的平滑轨迹规划.

1 运动数学模型

考虑地球自转且假设地球为均质的圆球，建立如下再入飞行器三自由度运动方程［12］:

其中:r表示地心距离;γ表示飞行路径角;V表示再入飞行器相对于地球的速度;ψ表示速度方位角，从北向顺时针量测为正;θ和φ分别为经度和纬度，以上为再入飞行器的状态变量，且以时间t为微分自变量．S表示特征参考面积;m表示质量，以上为与飞行器相关的总体参数．ω表示地球的自转角速度;μ表示地球万有引力常数;ρ=ρ0e-h/hs表示当前高度的大气密度，以上为地球参数及大气密度模型．CD与CL为再入飞行器的阻力系数和升力系数，是攻角α与飞行马赫数Ma的函数;σ为倾侧角，其正向是指从后方看飞行器顺时针转动的角度．

地球自转会产生附加力，但是与气动力相比可视为小量［13］，所以，在轨迹规划的预测阶段，忽略地球自转，再入飞行器简化运动方程为

其中:LR为再入航程;h为飞行高度;地心距r=h+Re，Re为地球半径;g=μ/r2为引力加速度．

2 轨迹规划问题描述

为保证再入飞行器安全，需满足动压、过载、热流密度等过程约束:

其中:Qs为热流密度;ks为与飞行器外形相关的常值系数;q为动压;N为法向过载;D和L分别为阻力加速度和升力加速度，即

轨迹规划问题可描述为:给定飞行器再入初始条件，在复杂的过程约束下，确定攻角及倾侧角的变化过程，使飞行器能够无动力滑翔至末端交班点位置，同时满足终端速度大小约束并消除航向误差，即

其中:tf为再入飞行时间，不加限定;Δψ表示航向误差;Δψ*f为其允许值;X(tf)表示X在tf时刻的值;X*

f表示X的终端约束条件(X表示式(7)中的高度h、经度θ、纬度φ及速度V)．

3 拟平衡滑翔轨迹规划算法

3.1 拟平衡滑翔条件

拟平衡滑翔条件存在两种定义方式:令再入飞行过程中的飞行路径角保持常数，即=0，这是普遍采用的拟平衡滑翔条件;令飞行过程中的高度变化率为常数，即=0，该条件应用于火星着陆解析预测制导中．本文根据第一种定义，在运动方程(2)中令=0，得到拟平衡滑翔条件(quasi equilibrium glide condition，QEGC):

其中L为飞行器的升力加速度，由式(6)得到．

3.2 初始下降段

通常情况下，飞行器的再入初始条件不满足拟平衡滑翔条件.为了使再入轨迹平滑切换到拟平衡滑翔状态，避免产生跳跃轨迹，需要引入初始下降段.平滑切换的评判准则是在高度-速度剖面内，寻找一个满足QEGC的转移点，使再入弹道的斜率dh/dV与QEGC的斜率(dh/dV)QEGC保持一致，并且平滑切换的成功表示了初始下降段的结束.其平滑切换的评判准则为［6］

其中ε为事先确定的小量．根据简化运动方程(1)与(3)可以得到

其中D为阻力加速度，由式(5)得到．(dh/dV)QEGC为当前点(h，V)对应的QEGC的斜率，根据式(17)对高度h求速度V的导数，可以得到

在采用常值倾侧角σ0积分运动方程过程中，当飞行高度小于一定高度而切换条件(9)仍不满足时，需要重新迭代计算σ0．然而σ0越大，轨迹下降越快，同时热流密度越容易超出约束范围．当超出约束的σ0仍不能满足切换条件时，由再入飞行器的气动特性可知，增大攻角能够提升气动力的减速效果，从而降低再入过程中的热流密度约束．所以，初始下降段一般选择常值大攻角飞行．

3.3 拟平衡滑翔段

传统的基于拟平衡滑翔条件的再入轨迹规划方法，首先给定攻角-速度剖面，然后只对倾侧角进行调整［6］.但是，对于大升阻比的高超声速再入飞行器，攻角-速度剖面对其气动参数影响较大，单纯的调节倾侧角有时无法完成再入需求.因此，本文选取攻角与倾侧角作为控制量，进而实现再入轨迹规划任务.

首先，基于飞行路径角与航程间的数值关系，确定满足航程要求的飞行路径角-高度剖面;然后，基于攻角与终端速度间的数值关系，确定满足终端速度约束的攻角-速度剖面;最后，对具有该飞行路径角的拟平衡滑翔轨迹由倾侧角进行保持控制.另一方面，为消除航向误差，倾侧角的符号由通过设计航向误差走廊控制倾侧角反转确定.此外，过程约束借助于拟平衡滑翔条件(QEGC)转化为对倾侧角的约束，从而确定出攻角和倾侧角，它们共同维持拟平衡滑翔轨迹飞行，最终实现再入飞行器终端位置和速度同时满足要求.在推导过程中，忽略了地球自转的影响，所以在一定时间周期内需要重新规划参考轨迹.拟平衡滑翔再入轨迹规划原理如图1所示.

图1 拟平衡滑翔再入轨迹规划原理

3.3.1 航程控制

航程约束是再入飞行器主要约束之一，因此，本节基于运动数学模型，推导了航程与飞行路径角的关系，由运动方程(4)与(3)可以得到

由于拟平衡滑翔轨迹的高度呈现近似线性变化趋势，所以选取高度为微分自变量是合理的.根据式(10)可以得知，再入航程是飞行路径角γ的函数．当假设再入过程中γ始终是小量且·γ=0的情况下，可以得到再入航程与飞行路径角的解析关系．但是整条再入轨迹不可能一直保持该假设条件，会造成一定的航程预测误差，故本文采用数值积分的方式预测再入航程．

定义了飞行路径角-高度剖面，为保证再入飞行路径角的平滑及连续性，该剖面选取二次解析多项式形式:

由于式(11)存在3个未知量，需要3组参数进行求解.3组参数分别为:当前点(hcur，γcur)，终端点，另外定义一个中间点(htr，γtr)且满足

其中:hct、kct为人为选取的常值参数;γdes为剖面设计参数，当γdes确定后，飞行路径角-高度剖面即可确定，再入航程LR也随之确定．

采用牛顿迭代方法搜寻参数γdes，使再入航程LR与待飞航程LRtogo相等．LRtogo根据当前点的经纬度(θcur，φcur)及目标点的经纬度，采用球面三角原理确定:

此时，再入航程控制问题转化为单参数的搜索问题，用于数值预测再入航程的微分方程(10)，涉及参数少，计算量小，效率高.

3.3.2 速度控制

终端速度约束是再入飞行器另一个主要约束，由运动方程(1)与(3)可以得到

由于飞行路径角-高度剖面已经通过航程控制确定，根据式(12)，终端速度可以通过调节隐含在阻力加速度中的攻角进行调节.为简化设计思路，定义攻角-速度剖面:

其中:α1、V1、V2为选定的常值参数(通常α1取最大攻角αmax);α2为调节终端速度的设计参数;kα、bα是与α1、V1、V2、α2相关的参数．

当α2确定后，攻角剖面即确定，数值积分式(12)可以预测终端速度．为使终端速度满足约束条件，采用牛顿迭代的方法搜寻参数α2．当搜索到的α2小于最小攻角αmin(或大于αmax)时，可以进一步调节参数V2，确保攻角α∈［αmin，αmax］．攻角剖面调节原理如图2所示．此时，再入速度控制问题也同样转化为相应的单参数的搜索问题．

图2 攻角剖面调节原理

3.3.3 拟平衡滑翔轨迹保持

当满足航程约束的飞行路径角剖面与满足终端速度约束的攻角剖面给定后，根据式(8)可以得到保持拟平衡滑翔轨迹飞行的倾侧角，即

其中升力加速度Lpf由攻角剖面αpf插值气动参数后计算得到．

再入飞行的过程约束，根据拟平衡滑翔条件可以将它们转化为对倾侧角的约束，即

取上述约束最小值作为倾侧角的约束边界，即

所以，能够保持拟平衡滑翔轨迹飞行的并考虑过程约束的倾侧角:

3.3.4 航向控制

航向控制是通过调节倾侧角符号进而消除航向误差，航向误差表示为

其中:ψcur为当前速度方位角;ψLOS为当前目标视线角，根据球面三角原理可以得到

消除航向误差的主要方法有:倾侧角单次反转模式、倾侧角两次反转模式以及航向误差走廊反转模式.前两种方法需要迭代搜索反转时机，能够提高精度，但是同时也增加了计算量．基于航向误差走廊的反转模式为:当航向误差Δψ大于航向误差边界Δψd时，倾斜角要反号;在边界区域内，倾侧角保持原来的符号．倾侧角符号表达式如下:

其中sign(σi-1)为前一时刻倾侧角的符号．

为降低倾侧角反转次数并提高精度，航向误差边界Δψd一般设计为“漏斗”式，本文采用该方法消除航向误差．航向误差边界如下:

4 仿真算例分析

4.1 仿真条件

本文采用美国的CAV-H作为研究对象［14］，该飞行器在再入飞行过程中无动力、无侧滑，转弯通过倾侧实现，气动参考面积为0.483 9 m2，质量为907 kg，最大升阻比为(L/D)max=3.5.

再入初始条件:初始高度h0=70 km，速度V0=6 500 m/s，飞行路径角γ0=0°，经度θ0=0°，纬度φ0=0°，速度方位角ψ0=80°;过程约束:热流密度约束Qmax=2 MW/m2，动压约束qmax=200 kPa，过载约束Nmax=3;控制量攻角和倾侧角的约束:最大攻角αmax=20°，最小攻角αmin=5°，最大倾侧角|σmax|=90°，在弹道计算过程中由过程约束进一步加以限制.

4.2 仿真验证

为验证本文提出的轨迹规划算法的有效性，设置终端约束条件为:高度20 km，速度V=1 500 m/s，经度=70°，纬度φ=0°，终端航向误差|Δψ(tf)|≤5°.为保证终端航向误差约束，设定航向误差走廊的相关数据为:Vd=2 000 m/s，Δψd1=10°，Δψd2=5°.

图3显示了攻角变化曲线，在再入初期以最大攻角进行飞行，但是为了达到期望的终端速度，在速度小于5 020 m/s时，以常值攻角6.2°飞行.图4显示了倾侧角的变化规律，为消除航向误差，倾侧角反转3次.此外，由于本文研究的再入飞行器具有大升阻比的气动特性，在保持拟平衡滑翔轨迹飞行阶段，为平衡重力需要较大的倾侧角.终端速度为1 496 m/s，偏差4 m/s，该约束条件的精度可以增加，轨迹规划次数进一步提高.从图5可以看出，再入轨迹非常平滑，由于本文将终端高度作为弹道积分终止条件，所以高度约束可以精确满足.图6显示了三维再入轨迹，终端经纬度分别为(69.97°，-0.007°)，该约束条件的精度可以通过调整航向误差走廊进一步提高.从图7可以看出，在初始下降段之后，飞行路径角在拟平衡滑翔飞行过程中变化缓慢，基本保持在0°附近.图8显示了过程约束的变化规律，由于本文再入轨迹规划算法是以QEGC为核心的，所以，热流密度、动压、过载等约束都能够比较容易满足约束条件.

图3 速度-攻角曲线

图4 速度-倾侧角曲线

图5 航程-高度曲线

图6 三维再入轨迹曲线

图7 航程-飞行路径角曲线

图8 速度-路径约束曲线

为验证本文提出的轨迹规划算法具有较强的适应能力，在初始条件不变的情况下，选取不同的终端纬度=(0°，5°，10°，15°，20°)进行仿真分析.图9显示了不同目标点的三维再入轨迹，再入轨迹平滑，且满足终端约束与路径约束条件，限于篇幅，没有给出相应的仿真结果曲线.表1给出了不同目标点的轨迹规划结果，能够达到比较理想的交班条件，而且在计算机(运行平台为Lenovo Y430，操作系统为Windows XP，编译环境为Matlab 6.5，CPU T5800/2 GHz，内存为DDR2 2 G)完成上述仿真的时间消耗均小于1 s，故采用此方法能快速求解再入参考轨迹.

图9 不同目标点的三维再入轨迹

表1 不同目标点的轨迹规划结果

5 结 语

本文根据大升阻比高超声速再入飞行器的飞行特点，将攻角和倾侧角作为控制量，充分利用QEGC，提出了一种全新的再入轨迹规划算法.该算法是数值预测与QEGC的完美结合，一方面获得了再入航程及速度的精确数值预测，进而借助QEGC进行轨迹保持控制，使再入滑翔轨迹具备平滑弹道的优良特性;另一方面基于QEGC实现了飞行过程约束向倾侧角约束的转化，成功解决了再入过程约束问题.仿真结果表明，本文提出的拟平衡滑翔再入轨迹规划算法，能够满足复杂的多约束再入条件，并具有较强的适应能力，轨迹规划算法简单，物理意义清晰，易于工程实现，能够为再入飞行器的轨迹规划问题提供参考.

［1］HARPOLD J C，GRAVES C A.Shuttle entry guidance［J］.Journal of Astronautical Sciences，1979，27(3):239-268.

［2］XUE S B，LU P.Constrained predictor-corrector entry guidance［J］.Journal of Guidance Control and Dynamics，2010，33(4):1273-1281.

［3］LU P，XUE S B.Rapid generation of accurate entry landing footprints［J］.Journal of Guidance，Control，and Dynamics，2010，33(3):756-767.

［4］BRUNNER C M，LU P.Comparison of numerical predict or corrector and Apollo skip entry guidance algorithms［R］.Toronto:AIAA Guidance，Navigation，and Control Conference，2010:AIAA-2010-8307.

［5］沈振，胡钰，任章，等.一种新型RLV再入轨迹在线规划算法［J］.宇航学报，2011，32(8):1670-1675.

［6］SHEN Z J，LU P.Onboard generation of three-dimensional constrained entry trajectories［J］.Journal of Guidance，Control，and Dynamics，2003，26(1):111-121.

［7］吴了泥，黄一敏，贺成龙.基于动压剖面的再入弹道解析解［J］.弹箭与制导学报，2009，29(6):173-176.

［8］DUKEMAN G A.Profile-following entry guidance using linear quadratic regulator theory［R］.Monterey:AIAA Guidance，Navigation，and Control Conference and Exhibit，2002:AIAA-2002-4457.

［9］YOUSSEF H，CHOWDHRY R S，LEE H，et al.Predictorcorrector entry guidance for Reusable launch vehicles［R］.Montreal:AIAA Guidance，Navigation，and Control Conference，2001:AIAA-2001-4043.

［10］POWELL R W.Numerical roll reversal predictor-corrector aerocapture and precision landing guidance algorithms for the mars surveyor program 2001 missions［R］.Boston:AIAA Atmospheric Flight Mechanics Conference，1998:AIAA-98-4574.

［11］HANSON J M，JONES R E.Test results for entry guidance methods for Reusable launch vehicles［R］.Reno:42nd AIAA Aerospace Sciences Meeting and Exhibit，2004:AIAA-2004-701.

［12］赵汉元.飞行器再入动力学和制导［M］.长沙:国防科技大学出版社，1997:74-81.

［13］阮春荣.大气中飞行的最优轨迹［M］.茅振东，译.北京:宇航出版社，1987:30-35.

［14］李瑜.助推-滑翔导弹弹道优化与制导方法研究［D］.哈尔滨:哈尔滨工业大学，2009.