考虑机动系数的标准轨迹再入制导方法

乔 浩，李新国，常武权

(1. 西北工业大学航天学院，西安 710072；2. 陕西省空天飞行器设计重点实验室，西安 710072；3. 中国运载火箭技术研究院，北京 100076)

0 引 言

再入飞行是沟通太空与陆地的重要飞行阶段。自上个世纪的阿波罗计划开始，各国的科技工作者对再入过程中的轨迹设计、制导与控制等技术进行了持续而深入的研究[1-2]。按照是否跟踪标准轨迹，制导方法主要可分为标准轨迹制导与预测校正制导两大类。后者具有较强的任务适应性，但对快速计算提出了较高要求[3]，对于常规的再入问题，前者以其简洁与稳定性具有较高的应用价值[4]。在标准轨迹再入制导问题中，如何设计所需的标准轨迹是较为基础的研究内容，对后续制导、控制方法的选取以及性能有着深刻的影响[5]。

再入飞行需要考虑飞行器的防热以及结构安全问题。反映到轨迹设计过程中，即需要满足热流、过载以及动压等约束。再加上对再入末端状态的要求，轨迹设计问题转变为一个有限可行域的多约束求解问题。按照求解的出发点可分为最优解与可行解两类方法。一些学者选择采用最优化方法寻找最优解，取得了大量的研究成果[6]。最优解的质量一般比较好，但这类方法有一个较为明显的缺陷，即寻优过程不可控。大量的迭代寻优往往对应着较长的求解时间，在某些情况下，计算还会陷入局部寻优而导致求解失败，因此一般用在飞行前的评估与地面准备阶段[7]。

另一求解思路即寻求满足任务需求的可行解。这类方法通常将三维轨迹设计分解为纵向与侧向轨迹设计两部分。首先在纵向设计满足过程约束以及终端约束的标准轨迹，之后加入侧向制导逻辑控制倾侧角反转位置，使再入终端满足任务需求。这类方法最早在航天飞机的再入制导中获得应用，即经典的阻力加速度制导方法。具体分为轨迹设计与轨迹跟踪两部分。前者创造性地引入了航程与阻力加速度之间的关系，将航程转化成再入走廊内分段阻力加速度曲线，大大简化了轨迹的设计流程[8]。在完成航程预测之后，侧向通过设置航向偏差走廊控制倾侧角的反转。X-33、X-34以及X-37B等的再入制导，都是以航天飞机的制导方法为基础设计的。航天飞机制导方法最大的问题在于采用了大圆弧假设，没有考虑轨迹的曲率问题，因此末端航程有较大的偏差。

在此背景下，Saraf等[9]改进了轨迹设计流程，提出了衍化的阻力加速方法(Evolved acceleration guidance logic for entry，EAGLE)。EAGLE将轨迹设计分为轨迹长度与轨迹曲率两个子问题，相比于航天飞机的制导方法，主动考虑了轨迹长度与曲率的耦合关系，飞行末端的航程误差大幅减小，所设计的轨迹更接近于真实轨迹。EAGLE在AG&C项目中获得了很高的评价。但从本质上讲，在初始化与轨迹长度更新过程中该方法仍未完全脱离大圆弧假设，导致每次更新后的轨迹长度和真实长度都有一定偏差，需要通过多次迭代进行更新。

除EAGLE方法之外，其他一些学者在轨迹设计过程中也考虑了轨迹的曲率问题。Shen等[10]在利用准平衡滑翔条件进行再入轨迹的设计过程中，得出在航程预测公式中保留航向偏差对大横程飞行有益的结论，并根据初始与终端的航向偏差假设了线性近似的航向偏差剖面，改善了剩余航程的预测效果。Liang等[11]详细讨论了航向偏差对于轨迹设计的影响，并从宏观的角度出发首次提出机动系数的概念，给出了一种新的再入轨迹设计流程。沈振等[12]改进了传统阻力加速度剖面，采用两段折线与走廊边界快速设计标准轨迹，提高了标准轨迹制导方法的实时性。王涛等[13]将阻力加速度剖面与预测校正方法结合，由走廊上下边界插值获得标准轨迹，插值参数由预测校正方法确定，并在三维空间对标准轨迹进行校正，提升了再入飞行能力。

本文提出一种改进的三维再入轨迹设计方法。与EAGLE类似，该方法将轨迹设计分为轨迹长度、轨迹曲率两个子问题。首先在轨迹长度子问题中对于轨迹长度相关的定义进行梳理，对相应计算公式进行分析归纳；给出改进的机动系数定义形式；之后给出完整的考虑机动系数的三维轨迹设计方法；最后设计了用于制导的轨迹跟踪控制器。

1 再入动力学模型

为方便后续讨论，在开展再入问题的研究之前，首先给出需要用到的动力学模型以及相关约束。假设飞行器为质点，地球为均匀球体，忽略地球自转。再入飞行器的总能量E定义为：

(1)

式中：V为速度大小，r为地心距，μ为地球引力场量，R0为均匀球体假设下的地球半径。在实际求解中将能量归一化，定义初始能量为0，终端能量为1，再入能量表示为(E-E0)/(Ef-E0)。以总能量E为自变量，则对应的动力学模型为[14]：

(2)

式中：θ，φ分别为经、纬度；γ，ψ分别为航迹倾角与航向角；g为地球引力加速度；σ为倾侧角；L与D分别为升力、阻力加速度，其表达式如下：

(3)

式中：m为飞行器质量，ρ为大气密度，SA为参考面积，CL，CD分别为升力系数与阻力系数。

需要考虑的主要路径约束为：

(4)

对于升阻比较大的再入飞行器，通常为改善轨迹质量，还需要考虑平衡滑翔约束，取平衡滑翔约束中的倾侧角σ=0°，有平衡滑翔条件为：

L-g+V2/r=0

(5)

此外，还要考虑终端高度hf、速度Vf、距目标点距离Stogo,f以及航向偏差Δψf等约束。

2 轨迹长度相关定义

再入问题可分为轨迹长度问题与轨迹曲率问题两部分，其中前者是再入制导的核心。对于依赖标准轨迹的制导方法，轨迹长度是规划参考阻力加速度剖面的依据，标准轨迹的轨迹长度和实际轨迹长度越接近，后续跟踪过程中由末端航程偏差导致的更新与迭代计算量越小；对于预测校正等不依赖标准轨迹的制导方法，每次预测阶段对剩余航程的预测越准确，指令校正收敛越快，也就降低了弹载计算机的计算量。

轨迹长度常见的描述方法有射程、航程、剩余航程、纵程、横程等。其中，射程一般用于描述弹道导弹等从初始星下点到当前星下点的大圆弧长度；航程在大多数研究中指星下点轨迹的真实长度；剩余航程是指飞行器当前位置与目标点连线对应的星下点轨迹长度，通常被用于描述待飞距离。纵程与横程是将整个再入飞行分解为沿着纵向目标平面与垂直纵向目标平面两部分，前者对应纵程，后者对应横程。在航天飞机的制导问题中，假设航向偏差较小，飞行器基本沿着纵向平面飞行，此时航程与纵程采用了相同的定义方式。常见的三种轨迹长度描述方式为：

(6)

式中：R0为地球半径，Δψ为航向偏差：

Δψ=ψ-ψLOS

(7)

式中：ψLOS为飞行器当前位置与终端位置连线对应的视线角。由球面三角公式可得：

(8)

ΔΨ定义为关于初始纵向平面的航向偏差：

ΔΨ=ψ-ψLOS,0

(9)

式中：ψLOS,0为初始视线角。

Stogo=R0cos-1[cosφcosφfcos(θ0-θf)+sinφ0sinφf]

(10)

值得注意的是，该定义方式下Stogo对应了每一瞬时飞行器与目标点连线的大圆弧，并不是真实的剩余航程，这种差异在具有横向大机动的场合中尤为明显，导致制导过程需要不断修正对剩余航程的估计。

(11)

图1 S1～S3几何关系示意图

其中,δS1为直角三角形的斜边，对应的轨迹长度最长，描述了真实的航程情况。δSm为δS1在当前视线方向的垂直分量。δS2与δS3分别对应了δS1在视线方向与初始纵向平面内的投影。二者之间的差异是由于Δψ与ΔΨ的不一致造成的，当且仅当飞行器位于星下点轨迹与初始视线在地面投影的交点位置(包括初始与终端点)时δS2=δS3。

(12)

式(10)中Stogo的计算方式不变，将其定义修改为剩余距离或剩余视线长度，而剩余航程的定义改为真实的剩余轨迹长度。对于定点再入问题，成功的判定标准为：在满足高度、速度约束下，末端剩余视线长度为零；或再入末端横程为零，纵程等于初始大圆弧长度。

3 再入机动系数

文献[11]在定义机动系数时，将δS2作为瞬时纵程δStogo，将δSm作为瞬时横程，此时：

(13)

式中：δSm/δStogo对应了瞬时横向运动与纵向运动的比例关系，用于描述当前瞬时状态下的机动水平。

之后，对式(14)进行积分，获得总航程Rflight的表达式为：

(14)

进一步根据积分中值定理，将上式转化为：

(15)

(16)

保留中间项是因为在一些诸如再入导弹的任务分析初期，虽然已知初始位置与目标位置，但不确定给定的制导指令是否能导引导弹到达目标位置，此时可采用中间项衡量机动系数，考察机动系数是否超过限制。当给定的制导指令使Rfligth≥Rl，采用右端项简单的形式求解机动系数。

4 三维轨迹设计

4.1 轨迹长度子问题

轨迹长度是跟踪标准轨迹制导方法的前提，不同于EAGLE求解轨迹长度时的迭代方法，本文在已知Rl并给出合理Pm之后，即可确所需的真实轨迹长度Rflight。因此，本节主要研究Pm的确定方法。一个合理的飞行任务，首先必须保证目标点在可达域之内，此时Pm是一非空集合；当Pm=φ时，表明没有到达目标的可行轨迹，任务不可行，需要更新初始条件或重新选择目标点。本文后续研究假设飞行任务可行，此时Rfligth≥Rl，对应机动系数Pm≥1。但在实际飞行过程中，Pm=1表示无侧向运动，飞行器始终在初始目标平面内飞行，在给定攻角剖面的前提下，对应了可达域的最远边界，此时机动系数区间只有唯一解。因此对于一般的飞行任务，需要同时确定Pm的上、下界。

文献[11]采用枚举方法确定机动系数的可行集。在分析初期给定一系列机动系数，之后对每个机动系数对应的轨迹长度规划阻力加速度剖面，求出对应的倾侧角指令并加上侧向制导逻辑生成完整轨迹，判断是否满足终端误差要求。这种方法计算量较大，需要进行大量的前期准备工作，导致由机动系数确定轨迹长度的算法相比EAGLE等方法优势不明显。而Pm的边界问题相当于求解给定纵程Rl情况下，航程Rfligth的最大值、最小值，因此可以转化为简单的边界优化问题。采用提前给定的攻角剖面αref，控制量为σ，Rfligth，max对应最优化问题：

minJ(σ)=-Rflight

(17)

类似的，Rfligth,min对应的最优化问题为：

minJ(σ)=Rflight

(18)

式中：sT为剩余距离约束，Stogo,f由末端经、纬度与目标点经、纬度确定。值得注意的是，给定初始状态与终端状态后，Pm的取值范围随之确定，即：Pm∈[Rflight,min/Rl,Rflight,max/Rl]。后续只要初始状态与终端约束不改变，则可直接根据机动性需求选取合适的Pm，确定出实际需要的轨迹长度，对于轨迹的快速生成有较好的应用潜力。

4.2 阻力加速剖面生成

确定轨迹长度之后，即可在阻力加速度与能量空间按照纵向轨迹的规划方法，规划满足实际轨迹长度的阻力加速度剖面[15]。对E关于时间求导，有：

(19)

(20)

式中：E0与Ef分别为初始与终端能量。由于r-R0≼r，且γ较小，因此引入假设：R0/r≈1，cosγ≈1，Rflight关于D的关系转化为：

(21)

热流率、过载、动压对应的阻力加速度边界为：

(22)

对应的准平衡滑翔边界为：

(23)

参考阻力加速度剖面的求解问题转化为：在热流率、过载、动压以及准平衡滑翔边界确定的再入走廊内，确定一条满足实际轨迹长度的阻力加速度剖面。考虑到式(22)诸多假设的存在，实际阻力加速度剖面应以式(22)为基础微调。之后即可求解出σ的幅值。由于阻力加速剖面的规划方法已经很成熟，本文不再赘述。值得注意的是，由于Rflight是由机动系数直接确定的真实轨迹长度，与基于大圆弧假设的轨迹规划方法相比规划的阻力加速度剖面更长，也更符合实际，避免了再入过程中因为末端航程偏差过大造成的阻力加速的剖面频繁更新，减轻了制导系统的负担。

4.3 轨迹曲率子问题

在确定轨迹长度之后，还需要确定倾侧角的符号，求解轨迹的曲率问题。本文采用一种动态的航向偏差走廊求解轨迹曲率问题。与确定的反转次数相比，该方法具有较大的灵活性，在轨迹长度确定后，动态调节走廊边界即可较好满足末端约束，动态航向偏差走廊如图2所示[16]。

图2 动态航向偏差走廊

图2中，Δψ0与Δψf提前给定，ΔψB为可调节边界。在航向偏差的许可范围之内，轨迹曲率问题转化为搜索合适的ΔψB，使末端约束得以满足。采用牛顿割线法进行搜索：

(24)

式中：i为迭代次数。当确定ΔψB边界之后，即可确定倾侧角的符号：

(25)

式中：Δψup与Δψdown分别为航向偏差的上下界。

在解决轨迹曲率子问题之后，即完成了再入标称轨迹的设计。相比于EAGLE等传统方法，本文提出的方法从宏观角度出发，省去了轨迹长度的迭代求解过程，能在量化机动性的同时，更快速获得制导所需的合理标称轨迹。算法的流程如图3所示。

图3 基于机动系数的轨迹生成算法逻辑

5 仿真校验

本文采用CAV-H这种典型的高超声速滑翔飞行器作为研究对象，对轨迹生成算法进行仿真校验。CAV-H的详细参数由文献[17]给出。假设初始航向对准目标，采用表1数据进行初始化。

表1 初始和终端要求

再入过程中的最大热流率取4 MW/m2，最大过载取3g，最大动压取70 kPa，末端能量管理区域半径取50 km。CAV-H的固定攻角剖面参考文献[18]，并结合实际的气动参数表，确定为：

(26)

表1任务对应的纵程Rl=6671.7 km。许多优化方法都可用于机动系数的求解，本文采用文献[19]提出的优化方法确定Pm的上下界。可得：Pm∈[1.0109,1.1461]，对应的地面轨迹如图4所示。

图4 给定攻角剖面对应的机动系数边界

值得注意的是，攻角剖面对Pm的取值范围有较大影响，当攻角作为控制量参与优化过程时，有Pm∈[1.0109,1.2824]。相比攻角剖面给定的情形，最大机动航程增加了约900 km。但攻角作为控制量在高超声速飞行中存在诸如热防护、稳定性、舵效等问题，因此在任务设计初期需依据机动性与攻角可控边界折中选择攻角的给定方式。由于规划的轨迹长度代表了飞行器的飞行能力，规划的轨迹长度越长，表明飞行能力越强，因此在确定机动系数取值范围之后，应尽量远离下边界取值，以防止在实际飞行过程中由于随机干扰、倾侧角翻转次数过多等影响造成能量的额外耗散。同时，选择较大的机动系数能够保证有足够的应急飞行能力，为故障后的应急返回或目标重新选择等任务场景预留较多的飞行能力。在给定攻角剖面的前提下，取Pm=1.08作为当前任务的机动系数，对应的需求轨迹长度Rflight=Pm·Rl。需要指出的是，由于为末端能量管理段预留了少量的剩余航程，更为精确的真实轨迹长度应该减去该部分。

实际阻力加速剖面根据参考阻力加速度剖面D(E)计算出(L/D)cosσ及|(L/D)cosσ|并积分动力学方程获得。EAGLE方法在侧向简化动力学中采用开环制导指令，在本文中不适用。为了更精确地跟踪参考阻力加速度剖面，在实际飞行过程中考虑偏差的影响，增加反馈环节以提高跟踪效果。忽略对过渡过程的积分项，采用pd控制进行参考轨迹跟踪，可得参考轨迹跟踪控制器如下：

(27)

(28)

式中：取ω0=0.09，ξ=0.1。

满足真实轨迹长度的参考阻力加速度剖面以及在实际制导中跟踪的阻力加速度剖面如图5所示。

相应的速度-高度、航迹倾角与航向角以及航向偏差如图7～图10所示。其中终端高度、速度落入指定区间；航迹倾角基本保持较小值；终端航向偏差较小。

跟踪参考阻力加速剖面的真实倾侧角与攻角指令如图6所示。

图5 再入走廊及参考、实际阻力加速度剖面

图6 攻角及倾侧角

图7 航迹倾角与航向角

不同机动系数对应的三维再入轨迹与地面投影如图11所示。

可知，跟踪通过指定机动系数Pm方式规划的再入轨迹能够到达指定的末端区域。不同机动系数对应的终端误差如表2所示。

图8 高度随归一化能量变化曲线

图9 速度随归一化能量变化曲线

图10 航向偏差随归一化能量变化曲线

图11 再入轨迹与地面投影

表2 终端偏差

可见由机动系数出发的标准轨迹再入制导方法能够满足再入段与能量管理段的交班要求。

6 结 论

本文研究了一种基于机动系数的标准轨迹制导方法，提出了改进机动系数的概念，并采用该系数从宏观角度出发研究了标准轨迹的快速生成方法。不同于EAGLE等方法，轨迹长度作为与机动性相关的参数直接指定，有目的地规划满足需求的阻力加速度剖面。本文提出的轨迹生成方法由于不含大圆弧假设，因此不需要通过逐次迭代的方式获得所需的轨迹长度。机动系数的区间只与任务剖面有关，对于确定的任务可以快速给出具有不同机动程度的飞行轨迹。机动系数直接定义为轨迹长度与纵程的比例关系，直观描述了再入飞行过程中横向机动的程度，给出了一种再入过程中机动性的量化方法。机动系数的取值空间大小，对应了攻击指定目标可选机动性的大小，也反映了备选弹道的多样性程度，后续可根据这一特点开展多枚再入导弹协同攻击、协同突防等问题的研究。用于轨迹跟踪的控制器形式简单，在实际中还可采用LQR等方法进行轨迹跟踪形成新的制导框架。