等再入航程返回轨道模糊混合优化设计

陈伟跃,王国军,王阳,陈蒙,张治国

中国空间技术研究院,北京 100094

1 引言

长期运行的航天器的轨道高度是变化的,如国际空间站的运行轨道高度在340～483km的范围内。对于长期运行的航天器变高度返回的工程问题,通过返回轨道设计使得再入段过载、驻点热流密度等弹道特性参数类似,对返回再入的工程实施具有重要意义。等航程再入弹道具有再入角一致,再入段航程近似不变、再入段过载和驻点热流密度等特性类似的优点,有利于返回再入控制,因此可以利用等航程再入弹道来确保制动高度变化时返回确定着陆点的再入弹道特性一致。

在再入控制过程中不论是标称轨道制导法还是预测校正制导算法,都需要用到标称返回轨道的先验信息[1-2],标称返回轨道的优化对于返回再入飞行任务设计和制导及控制系统系统设计具有重要意义[3-5]。目前标称返回轨道的设计方法多采用工程试算调试法[2]、有界试位迭代法[6]、微分修正法[7],此类方法依赖于设计经验确定初始值,需要进行反复调试才能确定初值。Gauss伪谱法[8]、局部配点法[9]、凸优化[10]、序列凸规划[11]等离散节点类优化算法则是标称返回轨道优化使用较多的方法。采用Gauss伪谱法、局部配点法进行标称返回轨道优化,为了达到较高的计算精度,需要增加离散点的数量,计算量较大。凸优化、序列凸规划对求解的问题有凸性的要求,需要将弹道优化问题中原始非凸问题转化为凸问题,转化过程较为复杂。上述离散节点类方法节点数量多,处理过程较复杂,计算量大。

对于标称返回轨道优化问题,序列二次规划算法(SQP)[12-13]对初值敏感且容易陷入局部最优解,遗传算法[14]不依赖于设计经验确定初始值,但为了获得收敛精度较高的结果,种群数量大,迭代次数多,计算量大。为了克服梯度类算法对初值的敏感性,并提高遗传算法的计算精度,遗传与序列二次规划混合算法(GA-SQP)被应用于弹道优化设计[15-16]。文献[15]在遗传算法初始种群生成时采用混沌搜索,在执行选择操作前选取符合单目标的优秀个体作为初值。文献[16]在遗传算法适应度计算时采用模拟退火罚函数处理约束。现有的遗传与序列二次规划混合算法(GA-SQP)[15-16]以遗传算法给出的准全局最优解作为SQP的输入,对遗传算法的求解精度要求较高,种群大小较大,迭代次数较多。

如何减小种群大小,减少迭代次数,提高搜索速度,降低对遗传算法搜索全局最优解的精度要求,同时保证遗传算法得出的结果对于SQP而言是一个“较好”的初值,是改进遗传与序列二次规划混合算法并将其应用于工程实践的关键问题。针对此问题,本文基于模糊评价[17]提出了改进的稳态遗传算法[18]与SQP混合的求解方法,并在等航程返回轨道设计中进行了应用。此方法不需要提供初始值,有效地提高了计算效率,适用于近地及深空返回再入轨迹优化,具备良好的工程应用前景。

2 等再入航程返回轨道

从不同高度返回时,若再入段航程相等,则称之为等再入航程。当制动点的初始条件变化时,可通过恰当的选择制动点和制动速度增量,使得再入段的航程偏差控制在一定范围内,如±2km,在工程实际中可以认为再入段航程相等。即使制动点的初始条件产生较大变化,再入段航程仍相等,再入环境类似,有利于再入控制。

令航程为S,速度为v,当地速度倾斜角为θ,地球平均半径为Rp,再入过程中高度为h,开伞点的高度为hf,航天器质量为m,升力为L。采用一阶近似解进行分析,可以忽略重力和离心力对再入弹道的影响[19],则有:

(1)

对于近地航天器再入弹道,h≪Rp,则有

(2)

式中:CD为阻力系数;s为参考面积;D为阻力;ρ为大气密度,ρ=ρ0e-βh,β=1.3889×10-5,ρ0为高度为0处的大气密度。

忽略重力和离心力,根据二阶近似解有[19]

(3)

式中:θe为再入角;ρe为再入点的大气密度。将式(3)代入式(2)可得

(4)

其中

对于近地轨道的半弹道式和弹道式返回,|θe|<5°,因此有cosθe≈1。再入点的大气稀薄,可以认为ρe≈0,因此b≈-1,则积分(4)式可得

S=(θf-θe)/β+[cot(θf/2)-cot(θe/2)]/β

(5)

由(3)式可知,当开伞点的高度hf、升阻比、阻力系数、参考面积、质量和大气模型确定的情况下,开伞点的速度倾角θf由再入角θe确定。再结合(5)式可知,再入段航程S由再入角θe确定。等再入航程返回轨道的基本特征是再入角一致,再入航程一致,再入段的弹道特性保持近似不变。若设计的返回轨道对某项特性具有最优性,在初始条件变化时,采用等再入航程返回轨道该项特性仍可保持最优性,有利于再入控制。

一般出现最大加速度时,Ma较大,可认为升力系数CL和阻力系数CD为常数。忽略重力和离心力对再入弹道的影响[19],最大过载Nmax可表述为

(6)

采用经验公式[20]计算驻点热流密度:

(7)

式中:Rn为驻点曲率半径;vc为航天器近地轨道环绕速度;c为经验系数。

在升阻比、阻力系数、参考面积、质量、大气模型确定的情况下,再入角决定最大热流密度出现时的速度倾角、高度以及速度与再入速度之比。等航程再入弹道最大热流密度出现的高度一致。

3 等再入航程返回轨道设计方法

3.1 等再入航程返回轨道设计问题的数学描述

将航天器星下点轨迹经过指定着陆区域的圈次称为返回圈。航天器返回圈升交点的轨道六根数采用σ表示,航天器的轨道运动表示如下:

(8)

航天器的轨道运动方程可根据需要考虑非球形引力摄动、大气阻力摄动等轨道摄动项[21-22]。

以返回圈升交点时刻为时间起点,设制动开机时刻相对于返回圈升交点时刻的相对时间为Tb。Tb>0表示在返回圈升交点之后制动,Tb<0表示在返回圈升交点之前制动。

在制动开机点建立返回坐标系o0x0y0z0。原点是航天器质心在制动时刻的地心矢量与标准地球椭球体表面的交点o0,o0y0轴沿o0点与制动时刻航天器质心连线的方向,o0x0在返回制动时刻航天器运行轨道平面内且与o0y0相互垂直,并指向航天器的运动方向,o0z0与o0x0、o0y0轴组成右手直角坐标系;在再入点(高度为100km)建立再入坐标系oRxRyRzR,原点在再入时刻的地心矢量与标准地球椭球体表面的交点oR,oRyR轴沿oR点与再入时刻航天器质心连线的方向,oRxR在返回制动时刻航天器运行轨道平面内且与oRyR相互垂直,并指向航天器的再入运动方向,oRzR与oRxR、oRyR轴组成右手直角坐标系。再入运动方向由再入点地心矢量与返回瞄准点地心矢量构成的大圆在再入点的地心方位角描述;半速度坐标系为o1xhyhzh,该坐标系的原点o1为航天器质心,o1xh轴沿速度方向,与速度坐标系o1xv方向重合,o1yh在返回坐标系x0o0y0平面内垂直于o1xh轴,o1xhyhzh构成右手直角坐标系。体坐标系为o1xByBzB,该坐标系的原点o1为航天器质心,o1xB为航天器的纵轴,指向头部,o1yB在航天器主对称面内,该平面在返回制动时刻与返回坐标系x0o0y0平面重合,o1yB垂直于o1xB,o1zB垂直于主对称面,沿运动方向看o1zB指向右方。坐标系关系如图1所示。

图1 坐标系关系

图1中右边子图的角度定义为:θ=arctg(vy/vx),σ=arcsin(-vz/v)。vx,vy,vz为返回坐标系下航天器相对于地球的飞行速度分量,v为航天器相对于地球的飞行速度大小,θ为速度倾角,σ为航迹偏航角。

航天器返回再入质心动力学方程描述如下[1]:

(9)

(10)

(11)

(12)

(13)

(14)

(15)

(16)

(17)

(18)

等再入航程返回轨道设计的数学描述为,求解制动开机时刻Tb、制动开机时长Tp、再入段倾侧角序列γi(i=1,2,3,4),优化目标函数为开伞高度处的航程偏差[1]ΔRL与横程偏差[1]ΔRZ的平方和最小:

(19)

并使得再入角θe、再入过载N、再入段航程S和过载超过过载界限的累计时间tg满足如下约束[23]:

Smin≤S≤Smax,θe=θc,N≤Nmax,tg≤tmax

(20)

式中:Smin为允许的再入段航程最小值;Smax为允许的再入段航程最大值;θc为满足再入航天器气动热设计要求的目标再入角;tmax为过载超过过载界限G0的累计时间的最大允许值。

3.2 等再入航程返回轨道设计问题的求解方法

本文给出一种基于模糊评价的改进稳态遗传算法与序列二次规划算法的混合求解方法。首先采用稳态遗传算法搜索等再入航程问题的粗略解,采用模糊评价方法确定搜索结果的好坏。将改进的稳态遗传算法搜索的粗略解作为序列二次规划算法的初始值,采用序列二次规划算法获得最终解。此方法不需要设计者给出初值,可自行搜索到较好的粗略解,同时避免陷入局部最优解,并可以获得较高的求解精度。

3.3 等再入航程返回轨道设计问题的求解步骤

(1)采用改进的稳态遗传算法求解初值

采用稳态(Steady State)遗传算法[18]求解等再入航程返回轨道设计问题时,根据优化变量Tb、制动开机时长Tp、再入段倾侧角序列γi(i=1,2,3,4)的数值范围采用二进制编码随机产生初始种群。此后克隆上一代的原始种群中一定比例或一定数量的个体以产生临时种群,临时种群按照一定的概率交叉、突变后与原始种群合并,合并后按照评价分数依次移除评分最低的个体,种群中个体数量恢复至原始种群的数量时,完成新一代种群的产生。以每一代及之前qc代的最好分数与上一代及其之前qc代的最好分数的比值作为收敛度,当收敛度达到指定的值时,遗传算法搜索停止。

使用稳态(steady state)遗传算法求解等再入航程返回轨道设计问题的粗略解时有3个关键的步骤:①随机产生初始种群,②根据最差评分移除劣质个体,③根据评分计算收敛度,依据收敛度停止迭代。由于初始种群是随机产生的,在第一代种群中产生的劣质基因需要在步骤②中根据评分移除,稳态遗传算法新一代种群与上一代种群存在一定比例的重叠,完全剔除劣质基因的影响需要一定的遗传代数,第一代种群的劣质基因会减缓收敛速度。步骤②和③依赖于遗传算法适应度函数的形式,对于多目标优化问题,鲜见采用模糊思想建立评分函数。以下针对等再入航程返回轨道设计问题,给出初始种群的基因检测方法、基于模糊隶属度的评分函数建立方法和搜索过程中劣质基因的评分方法。

为适应较大的运行轨道高度变化范围,等再入航程返回轨道设计的优化变量数值存在较大的变化范围。随机生成的初值中,制动开机时刻初值取值不合适可能导致再入点在瞄准点之后;制动开机时长初值取值不合适可能导致制动完成后滑行轨道的近地点高度较高,无法形成再入弹道。初始基因检测方法用于对随机产生的初始种群的健康状态进行检查,存在上述两种情况中的任何一种均标记为不健康个体,对于不健康的个体需要重新随机生成。基因检查方法如下:对于随机产生的初始种群个体X0l(l=1,2,…,l,l为种群大小),计算式(8)～(18),检查式(21)(22)。其中Hp为制动后滑行轨道的近地点高度,H0为制动后滑行轨道的近地点高度的允许上限,Ae为再入点地心方位角。若如下2式之一成立

Hp≥H0

(21)

Ae<0

(22)

则标记为不健康个体,剔除不健康个体并重新随机生成。

对于式(19)～(20)描述的优化问题,基于模糊隶属度[17]建立遗传算法的评分函数。

对于式(19)采用正态分布描述开伞点偏差的好坏程度:

μ(d)=e-(d/a)2

(23)

对于再入段航程约束采用正态分布描述约束满足的好坏程度:

(24)

对于再入角等式约束采用“Λ”型模糊隶属度函数[17]描述再入角约束满足的好坏程度:

(25)

式中:θ0根据需要取值,如θ0=0.1°;θc为再入角目标值,认为再入角偏离目标值±0.1°,再入角的优化结果不满足,隶属度为0。

对于最大过载约束N≤Nmax和过载超限时长约束tg≤tmax,采用戒上型隶属度函数[17]描述。对于最大过载约束满足程度的好坏描述为:

μ(N)=

(26)

式中:N0根据需要取值,如N0=0.1g,即认为最大过载超过允许值的数值达到0.1g时,最大过载约束不满足,隶属度为0。

对于过载超限时间约束满足的好坏程度如下:

(27)

式中:t0根据需要取值,如t0=10s,即认为最大过载超过过载界限G0累计时间达到10s,过载超限时长约束不满足,隶属度为0。

遗传算法的评分函数建立如下:

V=100[λ1μ(d)+λ2μ(S)+λ3μ(θe)+

λ4μ(N)+λ5μ(tg)]

(28)

式中:λj(j=1,2,3,4,5)为模糊综合评分的权重[17],须满足如下约束:

(29)

式(28)将遗传算法解的隶属度转化成了百分制评分函数,利用式(28)采用稳态遗传算法求解等再入航程返回轨道设计问题的粗略解。

在初始种群生成过程中进行了基因检测,对于不健康个体重新随机生成。由于在搜索过程中每一代种群个体的基因存在变异,仍然可能产生劣质基因。以下给出搜索过程中劣质基因的评分方法。

制动开机点基因变异可能导致再入点在瞄准点之后,令d、S、θe、N、tg的取值分别为:

d=1000km,S=-Smax,θe=2θc,

N=2Nmax,tg=2tmax

(30)

将式(30)代入式(28),带有劣质基因的个体评分结果为0,将在搜索过程中作为最差的个体被移除。

制动开机时长基因变异可能导致制动完成后滑行轨道的近地点高度较高无法形成再入弹道,令d、S、θe、N、tg的取值分别为:

d=1000km,S=0,θe=0,

N=2Nmax,tg=2tmax

(31)

将式(31)代入式(28),带有劣质基因的个体评分结果为0,将在搜索过程中作为最差的个体被移除。

(2)采用序列二次规划算法求解

将遗传算法计算得到的粗略解作为SQP的初始值,采用文献[24]给出的序列二次规划算法计算等再入航程返回轨道设计问题的最终解。对于(20)式中再入航程约束、最大过载约束和过载超限时间约束,需要转化为g(x)≥0的形式,g1(S)与g2(S)为航程约束,g3(θe)为再入角约束,g4(N)为最大过载约束,g5(tg)为过载超限时间约束,约束函数的形式如下:

g1(S)=S-Smin,g2(S)=-(S-Smax),

g3(θe)=θe-θc,g4(N)=-(N-Nmax),

g5(tg)=-(tg-tmax)

(32)

对于制动点不合适导致再入点在瞄准点之后的特殊情况,优化目标函数(19)和约束函数(32)中的S、θe、N、tg按照式(30)取值,d=2000km。

对于制动时长过短导致无法形成再入弹道的特殊情况,按照约束函数连续可微分的原则,d、S、θe、N、tg取值如下:

d=1000km,S=0,θe=0,N=0,tg=0

(3)等再入航程返回轨道设计问题的求解流程

等航程再入弹道设计问题的求解流程包括两个主要步骤,首先采用改进的稳态遗传算法获得粗略解,然后将其作为序列二次规划求解的初始值,采用序列二次规划获得等再入航程弹道设计问题的精确解。改进的稳态GA-SQP的具体计算步骤如下:

1)稳态遗传算法参数设置。设置种群大小、突变概率、交叉概率、收敛度界限、计算收敛度使用的代数。

2)随机生成初始种群。

3)初始种群个体的基因检查,若个体基因不健康则重新随机生成个体。

4)采用轮盘赌选择算法按照一定比例或指定个体数量从上一代种群中选择临时种群的父母基因。

5)交叉、变异。按照概率交叉、变异,生成临时种群。

6)将临时种群与上一代种群合并。

7)依次剔除评分最差的个体至种群数量恢复至原始值。

8)检查收敛度是否满足,若不满足继续生成新一代种群,否则输出等再入航程弹道设计问题的初略解。

9)序列二次规划算法参数设置。以稳态遗传算法得出的粗略解作为优化变量的初始值,定义优化目标函数、约束函数。

10)将带约束的优化问题转化为拉格朗日扩展的优化问题。

11)在当前解Xk处将约束线性化获得拉格朗日扩展优化问题的二次逼近。

12)求解二次逼近子问题。

13)以扩展的拉格朗日函数作为优化目标函数线性搜索下一步解Xk+1。

14)判断是否满足收敛精度要求,如果满足精度要求则输出等再入航程弹道设计问题的精确解,否则采用BGFS方法更新二次逼近子问题的二阶偏导数矩阵[24],转至步骤11)。

GA-SQP算法计算流程如图2所示。

图2 GA-SQP算法计算流程

4 等再入航程返回轨道设计实例

以CEV从国际空间站ISS轨道以半弹道式返回至爱德华空军基地(W 117°50′06.77",N 34°55′53.02")为例,说明等再入航程返回轨道的设计。国际空间站ISS的初始轨道根数如表1所示,CEV与返回相关的参数如表2所示[25]。

表1 ISS返回圈升交点J2000坐标系下的初始轨道根数

表2 CEV与返回相关的参数

为对本文提出的改进的稳态GA-SQP与未改进的稳态遗传-序列二次规划算法(UGA-SQP)进行比较,选取UGA-SQP的适应度函数为[26]:

采用UGA-SQP优化时,根据所得的优化解采用式(28)得出评分值,在改进的稳态GA-SQP与UGA-SQP评分值偏差小于0.5分的情况下进行对比。为使得UGA-SQP给出的解的评分值与改进的稳态GA-SQP给出的解的评分值偏差小于0.5分,需要调整UGA-SQP的种群大小和收敛代数,其他参数两种方法取值一致。两种方法的SQP部分均考虑对再入点在瞄准点之后和制动时长过短两种特殊情况的处理。

采用VC++编程语言实现本文提出的改进的稳态GA-SQP和UGA-SQP,运行环境为Win7系统,8核CPU,主频3.40GHz。半长轴为6715.8km时,采用改进的稳态GA-SQP和UGA-SQP分别计算5次的结果如表3和表4所示,表3和表4中p为种群大小,q为遗传算法的遗传代数、SQP算法的迭代次数。采用改进的稳态GA-SQP优化时,种群大小为20,收敛代数为40。在改进的稳态GA-SQP与UGA-SQP评分值偏差小于0.5分的前提下,采用UGA-SQP优化时,种群大小为80,收敛代数为80,采用UGA-SQP优化时的种群大下与文献[15]一致。

表3 改进的稳态GA-SQP等航程再入弹道的设计结果

表4 未改进的GA-SQP等航程再入弹道的设计结果

采用改进的稳态GA-SQP算法计算半长轴变化时的结果如表5所示,图3与图4为表5中编号为1的工况对应的再入弹道设计结果,图4中的经度范围为0°～360°。

图3 横向偏差及倾侧角曲线

图4 CEV返回过程的星下点轨迹

从表3可以看出,在参数设置不变的情况,由于遗传算法随机搜索的特性,每次计算遗传算法(GA)得出的收敛程度是不一致的,模糊综合评分结果V的百分数值可以表示遗传算法搜索结果的好坏程度。对于改进的稳态遗传算法的搜索结果,SQP均可找到满足各项约束的设计结果。表3中5次计算的制动开机时刻与平均值的最大偏差为146ms,制动开机时长与平均值的最大偏差为5ms,开伞点最大偏差1.141km,再入角偏差小于0.001°,再入航程约束的上限偏差为1.234km,满足最大过载约束,超过3.0g的时间最大为9.2s。采用改进的稳态GA-SQP计算最小耗时为169s,最大耗时为261s。从表4可以看出,表4中5次计算的制动开机时刻与平均值的最大偏差为720ms,制动开机时长与平均值的最大偏差为3ms,除第二次计算开伞点偏差略大外其他几次计算开伞点偏差均满足要求,再入角偏差小于0.001°,再入航程约束的下限偏差为-5.652km,过载约束存在未满足的情况,最大过载偏差为0.017g,超过3.0g的时间最大为15.4s。采用UGA-SQP计算最小耗时1161s,最大耗时为2027s,采用UGA-SQP计算遗传代数会明显增加,导致计算耗时显著增长。上述对比仿真表明,在不降低计算精度的前提下改进的稳态GA-SQP种群数量减小,迭代次数减少34.7%,计算耗时明显缩短。改进的稳态GA-SQP与工程试算调试法[2]、有界试位迭代法[6]相比,不需要设计者给出初值,可自行搜索到较好的解,计算精度和速度均达到工程适用程度。

从表3可以看出倾侧角γ1～γ4的大小整体呈现依次减小的规律,不同次数倾侧角的计算结果存在一定程度的变化,为再入控制系统的设计留有一定的自由度。从表5可以看出,在固定再入角的条件下,当ISS轨道高度增大时,制动点的位置沿着返回圈的星下点轨迹后退,制动点高度增加,再入弹道的再入点位置不变,制动时长和速度增量随轨道高度的增高而增大。再入段航程与平均值的最大偏差为2km,实现了制动高度变化时再入段航程一致的设计目标,且再入段最大过载基本一致。

在ISS轨道高度发生变化的情况下,再入弹道的特性保持不变,对再入控制与防热是十分有利的。

5 结论

本文推导了再入航程与再入角的关系,分析了等航程再入弹道的优点,提出了基于模糊综合评分的稳态遗传与序列二次规划混合求解算法。针对多目标等航程再入弹道设计问题,设计了稳态遗传算法初始种群的基因检测方法,给出了基于模糊隶属度的评分函数建立方法、搜索过程中的劣质基因评分方法,针对序列二次规划求解过程中再入点在瞄准点之后和制动时长过短的特殊情况,给出了约束的处理方法。考虑了再入过载超限时间约束、倾侧角翻转的最大角速度限制等实际工程约束。

以CEV从ISS轨道返回为例,说明了等航程再入弹道设计方法和设计结果的特性。采用改进的稳态遗传算法的初始搜索结果,SQP算法均可找到满足各项约束的设计结果,计算精度达到工程适用程度,为控制系统的设计留有一定的自由度,实现了制动高度变化时再入弹道特性一致的设计目标,计算速度达到工程适用程度。

针对轨道高度发生变化的情况,等航程再入弹道能保证再入弹道的基本特性不变,为返回再入控制、防热创造良好的条件。本文的结果可为中国空间站阶段的返回轨道及返回方案设计提供参考。