一种基于解析方法的滑翔再入轨迹快速规划方法

曾 亮,薛 辉,刘 欢,涂国勇

(酒泉卫星发射中心 技术部,甘肃 酒泉 735000)

为了提升飞行器的任务适应性、可靠性及性价比,新一代再入飞行器正逐渐向自主化、智能化方向发展。对各类航天器飞行事故的调查表明,如果飞行器具备更强的自主性和智能性,很多飞行事故是可以避免的。在NASA的空间发射行动计划(Space Launch Initiative,SLI)中,再入飞行器的自主性和自适应性已经作为衡量飞行器的重要指标。要使飞行器实现自主化、智能化,则飞行器应具备在线进行轨迹规划的能力。

对于在线轨迹规划,数值方法被广泛采用。LIN等为解决飞行器临时侦查任务的轨迹设计问题,将轨迹规划问题转化为一个两级优化问题,实现了轨迹的快速规划。SHEN等引入一套侧向控制算法,基于平衡滑翔理论,发展了一种三维轨迹在线规划算法。尽管这些方法在一定程度上提升了轨迹规划速度,采用数值方法也能够实现轨迹的高精度计算,但其庞大的计算量会大幅加重飞行器自载计算机的计算负担,其要实现工程应用还有待于计算机技术的快速发展。

相对于数值方法,解析方法的计算量要更小。在再入运动解析理论的研究上,文献[7]推导出了飞行器再入的解析解,拓展了解析解的适用范围,并给出了航程和时间的解析解。相较于直接再入飞行器,升力式再入飞行器的气动系数不为零,且随着飞行器自身的状态变化。气动力升力的存在使得飞行器的再入轨迹变得可控可设计,但也加大了解析理论的研究难度。MITITELU推导了升阻比为常值情况下升力再入2-D轨迹的解析解。文献[10]推导了适用于高超声速滑翔问题的解析解。相对于数值方法,这些解析理论的精度较低且前提假设较多,难以满足再入轨迹在线规划的需求。

本文研究了一种基于解析方法的再入轨迹在线快速规划方法大幅减小轨迹预测的计算量。为提升预测轨迹的精度,分阶段推导了速度和当地速度倾角相对于高度的解析解,并给出了阶段划分的方法和依据。为简化轨迹规划过程,基于推导得到的解析解,将再入相关约束转化为轨迹参数约束范围。在轨迹参数的约束范围内,根据轨迹参数与待飞航程之间一一对应的关系,通过迭代校正轨迹参数获取满足任务需求的再入轨迹,实现再入轨迹的在线规划。文中通过与数值仿真对比,分析解析方法的精度并验证轨迹规划方法的适用性。

1 问题描述

1.1 运动方程

再入飞行器的动力学方程的矢量表达式为

(1)

式中:为从地球中心到飞行器质心的矢量,为气动力矢量,为重力加速度矢量,为科氏力加速度矢量,为离心惯性力加速度矢量。由于离心惯性力加速度对飞行器的影响量级仅为10,因此忽略该项产生的误差也相对较小。对应于矢量形式的动力学方程,基于圆球旋转地球假设的三维运动可以表示为

(2)

(3)

(4)

(5)

(6)

(7)

式中:和分别为经度和纬度,为相对速度,为当地速度倾角,为速度方位角,为从地球中心到飞行器质心的矢量长度,为倾侧角,为地球自转角速度,和分别为升力系数和阻力系数,为飞行器的参考面积,为飞行器的质量,为大气密度,重力加速度=981 m/s。根据圆球地球假设,地球矢径长度与高度的关系为

=+

(8)

式中:为地球平均半径。为便于分析推导,大气密度采用指数大气模型:

=e-

(9)

式中:为大气系数,可认为是常数,为海平面的标准大气密度。

1.2 再入过程约束条件

(10)

(11)

(12)

若上述3种过程约束描述在速度-高度坐标系下,则这3种过程约束的表达式为

(13)

(14)

(15)

2 再入阶段划分及运动特性分析

2.1 飞行器滑翔再入阶段划分方法

(16)

在飞行器再入过程中,这2种状况都有可能出现。第一种状况一般出现在以高于第一宇宙速度再入的初期阶段,而第二种状况一般出现在以低于第一宇宙速度再入的初期阶段和再入末期阶段。根据式(16)中2种状况的定义,结合当地速度倾角微分式(6),在忽略科氏力条件下,可得到速度和高度的关系式为

(17)

式(17)中的不等式与当地速度倾角的余弦相关,而在整个再入阶段中,当地速度倾角的幅值都比较小,因此当地速度倾角的余弦可近似为1。由此,式(17)中的不等式变为

(18)

结合过程约束式(13)～(15)及不能满足平衡滑翔条件的状态区间式(18),可在速度-高度坐标系下得到可实现平衡滑翔飞行区域,如图1所示。图中,Case 1限制区为第一种状况中因能力过剩而不能实现平衡滑翔飞行的区域;Case 2限制区为第一种状况中因能力不足而不能实现平衡滑翔飞行的区域,约束限制区为过程约束(如热流约束、过载约束、动压约束等)而形成的限制区域。

图1 某飞行器再入可平衡滑翔飞行能力区域图

此外,飞行器能否平衡滑翔飞行,还与再入飞行器的实际当地速度倾角相关。所以,可以将再入飞行器能否实现平衡滑翔飞行和再入飞行器的实际当地速度倾角是否满足滑翔条件,作为飞行器实现滑翔飞行的判断条件。对于再入速度小于第一宇宙的升力式再入飞行器,其再入过程一般可划分为2个阶段:初始下降段和滑翔飞行段。其中,初始下降段为飞行器初始进入大气层的阶段,该阶段大气密度稀薄,飞行器在各力的作用下不能实现平衡滑翔飞行;滑翔飞行段是在飞行器高度降低到一定程度,当地速度倾角接近于零后,气动升力足以使飞行器实现平衡滑翔飞行的阶段。如果飞行器开始能够保持平衡滑翔飞行状态,飞行器状态则由初始下降段过渡到滑翔飞行段。判断飞行器是否进入滑翔飞行段的条件,如下式所示:

(19)

式中:,f0为偏滑翔飞行段起始点处的当地速度倾角。

结合图1中的某飞行器再入可平衡滑翔飞行能力区域图,图2中给出了该飞行器某条再入轨迹的阶段划分。

图2 某飞行器器升力式再入阶段划分

对于探月飞船这类再入速度大于第一宇宙速度依靠升力再入的飞行器,其整个再入过程可以划分为3个阶段:初始下降段、滑翔飞行段和末段,如图3所示。其中,在初始下降段,由于探月飞船再入的速度大(接近第二宇宙速度)、能力强,其前期不能保持平衡滑翔飞行主要受第一种状况限制;滑翔飞行段与再入速度小于第一宇宙的再入飞行器相似,末段探月飞船由于飞行能力不足,不再能保持平衡滑翔飞行。对于探月飞船再入,飞船开始由初始下降段过渡到滑翔飞行段的状态应满足:

(20)

滑翔飞行段过渡到末段的状态应满足:

(21)

图3 探月飞船滑翔再入阶段划分

2.2 飞行器滑翔再入轨迹特性

结合式(6),可求得以高度为自变量的当地速度倾角的二阶微分方程在忽略科氏力后为

(22)

式中:

(23)

设函数:

(24)

则式(22)转变为

(25)

由式(25)可知,若(,)始终小于零,且当地速度倾角的初值小于零,则当地速度倾角将始终小于零,即若,f0<0,则在整个滑翔飞行段,再入飞行器的高度将单调持续下降。

在滑翔飞行段,气动升力占据主导作用,结合式(24),式中气动升力项的表达式为主导项,通过设计倾侧角幅值函数,可将(,)构造为

(,)=-e-

(26)

式中指数函数恒为正,则轨迹参数的符号即决定了函数(,)的正负。若轨迹参数取正值,则飞行器的高度在滑行飞行段将单调下降。此外,结合式(24)和式(26)可知,构造式(26)所示函数对应的倾侧角幅值函数为

(27)

从式(27)可以看到,轨迹参数和倾侧角幅值具有一一对应的函数关系,根据式(27)所示可确定倾侧角任意取值下对应的轨迹参数。结合式(25)～式(26),则可得到:

(28)

引入轨迹参数,即可用简单的参数描述一条高度单调的滑翔再入轨迹,改变轨迹参数的大小即可改变飞行器高度升降的快慢,进而直接改变再入轨迹航程的长短。在滑翔飞行段,相较于倾侧角,能更好描述滑翔飞行类轨迹的特性。

3 再入运动解析理论推导

为快速获取飞行器滑翔再入的纵向轨迹,可以求解飞行器再入轨迹速度和高度的解析表达式。根据前面所述的再入阶段划分,一般飞行器滑翔再入都包含初始下降段和滑翔飞行段。因此,本节以这两段为例进行推导。

3.1 前提和假设

再入运动方程非线性强,涉及系统复杂,直接求解运动方程得到解析解难度较大。为得到初始下降段和滑翔飞行段的解析解,作如下假设:

①假设升力系数和阻力系数在一定范围内不变。虽然飞行器实际再入过程中气动系数随着飞行器状态变化,但变化速率相对较小,在一定范围内认为其不变是合适的;

②假设攻角剖面由再入过程约束和再入任务的需求确定,倾侧角是质心运动唯一的控制参数;

③在整个再入过程中,当地速度倾角数值较小,能够满足sin≈和cos≈cos的假设。

3.2 初始下降段解析解推导

根据式(2)、式(6)和式(23),有:

(29)

式中:

(30)

由于在初始下降段对、和不敏感,因此,可被看做一个初始常值,表达式为

(31)

式中:,,,,为再入点处初始状态量。由此,求解式(29)得:

(32)

微分方程的积分常数为

(33)

根据式(29)的解,有:

(34)

根据式(2)和式(5),有:

(35)

在初始下降段,高度变化对阻力加速度的影响要远大于速度变化的影响。因此,在该阶段可以忽略速度对阻力加速度的影响。将式(34)代入式(35)中,可以得到:

(36)

尽管式(36)是一个一阶微分方程,但仍难以求解得到它的解析解。但是,微分方程(36)的解却可以通过数值方法快速求解得到。式(36)的数值方法求解表达式为

(37)

积分常数,的表达式为

(38)

(39)

因此,速度的数值积分表达形式的为

(40)

由于当地速度倾角在初始下降段的变化范围较小,式(35)中的当地速度倾角可近似看作常值,由此式(35)变为

(41)

式(41)的解析解为

(42)

3.3 滑翔飞行段解析解推导

在滑翔飞行段引入轨迹参数,求解式(28),可得当地速度倾角的表达式为

(43)

积分常数的表达式为

(44)

式中:分别是滑翔飞行段起始点处的高度。

在滑翔飞行段,阻力加速度要远大于重力项在速度方向上的分量。因此,重力项的分量可以忽略。将当地速度倾角的表达式代入式(35)中,式(35)可简化为

(45)

式(45)的解为

(46)

积分常数的表达式为

(47)

式中:为滑翔飞行段起始点处的速度。速度的表达式为

(48)

3.4 解析解有效性分析

引入数值解和罗赫二阶解,来验证初始下降段和滑翔飞行段解析解的有效性和准确性。在仿真对比分析中,飞行器的质量为800 kg,参考面积为0.5 m。和分别为0.966 2和0.707 4。飞行器的再入速度为7 500 m/s,再入高度为120 km,再入角为-1.0°,再入方位角为50°,再入点的经纬度都为30°。初始再入段的当地速度倾角为45°,滑翔飞行段的轨迹参数=0.1。

图4和图5为初始下降段当地速度倾角和速度的计算对比分析图。从图4和图5可以看到,初始下降段的解析解得到的结果与数值解间的误差较小,但要明显优于罗赫二阶解。这是由于在罗赫二阶解的推导过程中忽略了科氏力,而在初始下降段科氏力在纵向分量的量级和纵向上其他力合力的量级相当,因此不能被忽略。

图4 初始下降段中当地速度倾角解的对比图

图5 初始下降段中速度解的对比图

图6和图7中分别给出了滑翔飞行段速度和当地速度倾角的计算对比分析图。从图中可以看到,滑翔飞行段解析解的与数值解的结果比较接近,整体优于罗赫二阶解。但相较于初始下降段,罗赫解的准确性有了较大的提升,这是由于滑翔飞行段气动升力成为主导力,科氏力的影响基本可以忽略。

图6 滑翔飞行段当地速度倾角各解比较

图7 滑翔飞行段速度各解比较

从图4～图7的对比图中可以看到,初始下降段和滑翔飞行段解析解与数值解的差别较小。与经典罗赫二阶相比较,分阶段推导的初始下降段和滑翔飞行段解析解的精度明显更高,尤其是在初始下降段。此外,初始下降段和滑翔飞行段解析解的表达式更为简单,更便于求解。

4 再入轨迹规划方法

基于上节推导得到的飞行器滑翔式再入解析解,进行再入轨迹快速规划。对于一般再入任务,满足再入约束条件和再入终端指标是规划再入轨迹的基本要求。其中,再入约束条件有倾侧角幅值约束、热流约束、过载约束和动压约束;再入终端纵向指标有再入终端的高度、速度和待飞航程。

4.1 再入约束条件处理

在滑翔飞行段,轨迹构型通过轨迹参数设计。根据式(27)中轨迹参数与倾侧角间的关系及倾侧角本身的取值空间,首先可以得到控制变量约束条件下的轨迹参数的取值范围为

(49)

结合式(26)定义轨迹参数时确定其为正值,因此,在控制变量约束条件下轨迹参数的取值范围为

(50)

图8 控制量约束下的轨迹参数取值范围

由式(10)～式(12)可知,在再入过程中,其他状态不变的情况下,飞行器的热流密度、过载和动压都与大气密度的大小呈正相关关系。初始下降段整个阶段都处于大气密度稀薄区域,3种过程约束的峰值一般都出现在滑翔飞行段。因此,主要在滑翔飞行段考虑再入过程约束处理。

上一节中推导的滑翔飞行段速度相对于高度的解析解及式(13)～式(15)所示的过程约束表达式,可以得到由轨迹参数和高度表示的过程约束表达式:

(51)

(52)

(53)

根据式(51)～式(53)中所示的不等式,可进一步转化为轨迹特征参数的取值区间:

(54)

(55)

(56)

过程约束条件下轨迹参数的取值区间即为由式(54)～式(56)组成的不等式方程组的解。而上述3个不等式方程都有如下结构:

(57)

式(57)所示不等式方程的解为

(58)

由式(58)可知,取值越大,的取值区间越小;取值越小,的取值区间越大。因此,在过程约束条件下轨迹参数的取值区间为

(59)

(60)

结合式(50)中轨迹参数在控制量约束下的取值范围和式(59)中轨迹参数在过程约束下的取值范围,轨迹参数的实际取值范围为

0<≤

(61)

(62)

4.2 再入终端指标满足

再入终端纵向指标有再入终端的高度、速度和待飞航程,其中终端的高度作为再入段的截止条件;终端的速度可通过攻角剖面的设计来满足。在终端的待飞航程上,初始下降段气动力小,对飞行器的再入轨迹影响有限,可将初始下降段倾侧角设置为零,通过调整滑翔飞行段的轨迹参数来满足待飞航程。

(63)

根据再入运动微分方程,在不考虑侧向运动对纵向运动影响的条件下,飞行器的待飞航程对高度的微分为

(64)

在滑翔飞行段,当地速度倾角为小量,式(64)可近似为

(65)

(66)

(67)

(68)

通过引入有界试位法,迭代校正轨迹参数。有界试位法的基本原理为提前试位给出2个轨迹参数和,使得:

()·()<0

(69)

并通过下式求得新的制导参数:

(70)

并通过轨迹预测进一步求得(),并根据()的符号判断用(,())替代(,())还是(,()),最终使得()能够满足:

|()|≤

(71)

式中:为终端待飞航程指标精度要求。

4.3 轨迹规划方法有效性分析

在某一再入任务下,基于上述轨迹规划方法规划出满足再入任务要求的轨迹并获得该轨迹对应的轨迹参数。基于式(27),利用数值方法得到该指令剖面下的再入轨迹,期间不调整攻角剖面。以数值方法结果作为实际再入轨迹,通过对比分析2条再入轨迹,验证该基于解析方法的滑翔再入轨迹规划方法的有效性。

在仿真验证中,飞行器的质量为800 kg,参考面积为0.5 m,飞行器的气动系数表与X-34一致。飞行器的再入速度为7 500 m/s,再入高度为120 km,再入角为-1.0°,再入方位角为50°,再入起始点的经纬度都为30°,再入任务航程为7 263.3 km,再入终端点高度为20 km。飞行器再入的攻角剖面为

(72)

基于解析方法的滑翔再入轨迹规划方法规划出的再入滑翔轨迹的当地速度倾角、速度、航程和倾侧角,并与实际再入轨迹进行比较,结果如图9～图12所示。从图9中可以看到,利用解析解规划的速度倾角轨迹与数值方法计算得到实际轨迹总体吻合,其结果与图4和图6一致。图10中重合的航程高度对比图进一步验证了由解析解推导得到的航程计算式(65)和式(66)的正确性。图11给出了过程约束边界。从图中可以看到,在对再入约束条件进行处理后,规划得到的再入轨迹都能够满足过程约束要求。图11给出了将解析解代入式(26)得到的倾侧角与数值方法的对比图。从图中可以看到,由于解析方法规划得到的再入轨迹只考虑飞行器再入的纵向运动,因此得到的倾侧角曲线只有幅值没有符号;数值方法既考虑了飞行器纵向运动也考虑了飞行器侧向运动,其倾侧角曲线有幅值和符号,但2种方法计算得到的倾侧角幅值相当。从图9～图12可以看到,在给定的再入任务条件下,基于该轨迹规划方法不仅可以快速规划出满足任务航程的再入滑翔轨迹,并且基于解析方法规划出的再入轨迹与数值方法规划出的轨迹之间误差较小,具有较高的精度。

图9 当地速度倾角比较

图10 航程-高度曲线比较

图11 速度-高度曲线比较

图12 倾侧角曲线比较

5 结束语

本文推导了一种基于解析方法的滑翔再入轨迹规划方法。该轨迹规划方法能基于推导的解析方法快速预测轨迹,并将过程约束、航程、倾侧角等转化为解析表达式,实现了满足约束条件的再入轨迹的快速规划。该轨迹规划方法能快速规划出较高精度的再入轨迹,使得高精度的再入轨迹规划由离线走向在线,推动再入飞行器向自主智能方向发展。