一种考虑状态灵敏度的火星进入轨迹优化方法

彭玉明 尤伟 马彦涵

（上海卫星工程研究所，上海 200240）

1 引言

大气进入过程的运动方程形式复杂，进入轨迹对控制变量高度敏感且进入过程的非线性约束较强，因此进入轨迹优化问题一直是研究的热点。目前火星探测器大都沿用“海盗号”的气动外形，依靠滚转角的变化实现轨迹控制，控制能力较弱，因此在设计标称轨迹时还需要考虑控制系统的跟踪能力。传统的轨迹优化方法首先根据约束条件确定一条轨迹，然后再设计控制系统，若控制系统无法跟踪标称轨迹，则重新规划，如此反复直到满足要求。在对轨迹进行优化设计时，完全不考虑控制系统的跟踪能力，设计过程完全依靠经验费时费力，且不一定能取得满意的效果。因此本文提出一种考虑控制饱和的进入轨迹优化算法，把末端状态的灵敏度加入到目标函数中，并与末端高度最高、燃料最省通过罚函数法结合到一起形成最终的优化目标[1-3]。

2 Simpson 法

直接配点法的核心思想是通过离散化把最优控制问题转化成非线性规划问题，直接配点法是将整个时间过程划分成N 段，每一段的两端称为节点，然后用多项式逼近每一段的状态和控制变量[4-5]。根据多项式阶次的不同，直接配点法又可以分为低阶的梯形法、Simpson 法和高阶的四阶、五阶方法，本文采用三阶Simpson 法。

在时间[t0,tf]内将连续时间分成N 段（tf为未端时刻），每一个子区间为[ti,ti+1]，记Ti为相邻2个节点的时间间隔，s=（t -ti）Ti（s为归一化变量；t为时间），在每个子区间上状态变量x 可以用三次Hermite多项式表示

式中 c0、c1、c2、c3为多项式系数。

其边界条件为

式中 x1和x2分别表示相邻两个节点的状态。

把边界条件代入式(1)中得

解上述线性方程组得

在子区间的中点处，即s=0.5 时，将式(3)代入式(1)得

式中 xmi为中点处的状态量；xi为第i个节点的状态量；xi+1为第i+1个节点的状态量；fi，fi1+分别表示函数 f（x,u,t）在第i个子区间两端点处的函数值，（其中，u为控制量）即

式中iu为第i个节点的控制量。

式中 ui1+为第i+1个节点的控制量。

当由公式(5)近似得到的中点处的导数与动力学方程在中点处的值足够接近时，就认为多项式很好的逼近动力学微分方程，如图1所示。

图1 Simpson 法Fig.1 Simpson method

原来的动力学方程由N×m个Hermite-Simpson 积分形式的约束代替

式中 umi为中点处的控制变量，由两端线性插值得到；tmi为中点对应的时刻。

定义决策变量Z为

式中 x为各节点状态变量构成的状态向量；mu为中点处控制量构成的控制向量。

定义Δi为各节点处状态量逼近误差

这样动力学微分方程约束转变成N×m个等式约束，问题转化为有约束的非线性规划问题，搜索决策变量最优值，使得Δi趋于零。

3 考虑状态灵敏度的进入段最优轨迹

3.1 优化目标

由于火星大气特别稀薄，降落伞减速效果有限，为了使探测器能够安全着陆，希望开伞高度越高越好，所以目标函数中应包含末端高度；另一方面，探测器上携带的燃料有限，过多的消耗燃料是不可取的，因此燃料最省也是优化的目标之一。然而，如果简单的以末端高度最高或燃料最省为性能指标进行优化，得到的控制变量的形式很可能是最大–最小–最大的形式，虽然控制量留取一定余量，但是在存在较大不确定性扰动的情况下，会导致控制饱和，更加不利的是，像“海盗号”这一类小升阻比的探测器，轨迹控制能力较弱，加之进入段恶劣的气动加热环境以及大量不确定性因素，探测器很可能会因为控制饱和而无法跟踪标称轨迹，因此在选取目标函数时还要考虑控制系统能否跟踪上标称轨迹[6-7]。末端高度、燃料消耗以及跟踪控制性能等几个指标是相互对立的，因此采用罚函数法把3个性能指标统一到一起，本文选取目标函数J为

式中 cS,cu为加权系数；h（tf）为探测器末端高度；Js为末端状态对进入轨迹的灵敏度函数；Ju为进入过程中关于燃料消耗的目标函数

式中 σ为滚转角。

3.2 约束条件

进入轨迹除了动力学方程的约束外，还考虑了初始条件约束、终端条件约束、过程约束和控制约束。初始状态约束包括初始高度0r 、速度0v 、飞行路径角0γ和航程 s0，即

探测器最后要到达指定开伞区域，并且满足开伞条件，由于末端高度是目标函数之一，因此这里不做约束，末端条件约束为

式中fs为末端状态灵敏度；fv为末端速度。末端时刻 tf自由无约束。

进入过程中，为了结构和设备的安全，需要考虑过载约束[8]。鉴于本文的研究对象为弹道升力式飞行器，其轴向和法向均有可能产生较大过载，因此本文过载约束na取总过载约束

式中 L为升力加速度；D为阻力加速度；na为过载；nmax为最大过载。

对动压的限制应从其峰值和末端时刻来考虑，进入过程中动压峰值q 应不超过给定的最大值qmax。

式中 ρ为大气密度；V为飞行速度。

为了安全打开降落伞，末端时刻动压应不超过给定的最大值

式中 c为常数；b 是与探测器半径相关的常数[9]。

由于只考虑纵向运动，所以假设σ ∈ [ 0,π ]，为了进行侧向制导，控制量需要留取一定余量，控制约束为

式中 σmin为最小滚转角；σmax为最大滚转角。同时由于探测器上姿态执行机构的限制，滚转角角速度和角加速度还必须满足以下条件

式中 σ˙min为最小滚转角速率；为最大滚转角速率；为最小滚转角加速度；为最大滚转角加速度。

3.3 状态灵敏度计算

控制系统的抗扰动性能体现在末端状态对进入轨迹的灵敏度上，灵敏度越小，说明末端状态对进入轨迹的变化越不敏感，也就是说抗扰动性能越强，反之亦然。因此，若是把所有时间点上状态扰动对末端状态的影响用灵敏度函数表示出来，然后加入到目标函数中，这样在设计最优轨迹时就考虑了控制系统的抗扰动性能。

根据灵敏度的定义，任意轨迹X（t t0,x0）相对初始时刻状态 x0的灵敏度矩阵为

为简化表达把S（t t0,x0）表示成 S（t,t0），它满足以下微分方程

则反馈控制律变为

灵敏度矩阵微分方程变为

根据灵敏度矩阵的定义，容易证明上述灵敏度矩阵具有如下性质

式中 I7×7为7 维单位矩阵。

在进入过程中，主要关注状态扰动对末端状态的影响，因此定义ξ（x（tf）,tf）函数

根据定义ξ（x（tf）,tf）相对t 时刻的状态扰动的灵敏度矩阵可以表示为

所以ξ（x（tf）,tf）相对t 时刻的状态扰动的灵敏度矩阵等于S(tf,t)。很显然只有矩阵中的每一个元素同时最小才能使灵敏度最小[10]。通过积分得到末端状态对所有时刻状态扰动的灵敏度，因此灵敏度函数JS为

式中Si,j（tf,t ）表示灵敏度矩阵中第(i,j)个元素；ci为加权系数。

4 仿真分析

4.1 仿真参数

探测器进入火星大气层的初始约束条件为：[r0,v0,γ0,s0]=[3 522 000,6 000,-0.2,0]，末端约束条件v（tf）≤350m/s ；s（tf）=935km ；最大过载值nmax=8gn；最大动压值q=10 0 00Pa，末端动压 qf≤750Pa；最大加热率=70kW/m2；控制变量约束为，11.5°≤ σ≤168.5°，-20(°) s≤ σ˙≤20(°) s，-5 (°) s2≤≤5(°) s2。

不考虑扰动对末端飞行路径角的影响，灵敏度函数的加权系数为 c1=1,c2=1,c3=0,c4=0.01；反馈增益为k=-[0.01 0.005 50 0.001]。

4.2 仿真结果与分析

图2、3 是加权系数 cs和cu取不同值时所得到的最优控制曲线。

图2 进入段最优控制轨迹Fig.2 Entry optimal control trajectory

图3 cu=10 时进入段最优控制轨迹Fig.3 Entry optimal control trajectory while cu=10

从图中可以看出，当cs=0，即目标函数中不包括灵敏度，最优控制轨迹大部分时间处于控制约束边缘，当存在扰动时，控制系统很可能会因为控制饱和而无法跟踪标称轨迹。当uc=10，即考虑燃料消耗以后，控制轨迹变得平滑，但是仍有大部分时间处于饱和状态。随着 cs的不断增大，最优控制轨迹离约束边界越来越远，抗扰动性能逐渐增强，但代价是末端高度越来越低，因此实际应用时需要根据具体情况采取不同方案。表1 给出了cs和cu取不同值时的末端高度。最优轨迹确定以后采用蒙特卡洛仿真方法进一步验证考虑控制饱和的火星进入轨迹优化方法的有效性。具体仿真参数见表2。

表1 不同参数情况下的末端高度Tab.1 Terminal altitude for different parameters

表2 进入段蒙特卡洛仿真参数Tab.2 Entry Monte Carlo simulation parameters

图4 是采用线性反馈控制律跟踪最优轨迹得到的实际控制轨迹，如前文所料，cs越大，控制量饱和的时间越少，当cs=0时，即不考虑状态灵敏度，在200s 以后控制量一直处于饱和状态，而考虑了状态灵敏度以后，在150s 以后基本不会出现控制量饱和。图5 是 cs取不同值时末端状态误差分布。从图中可以看到仅考虑末端高度时，末端误差高度、速度、航程误差分别约为2000m，100m，10km，随着 cs的逐渐增大，目标函数中状态灵敏度的权重越来越大，末端状态误差逐渐减小，当cs=0.005时末端高度、速度和航程误差分别减小到100m，5m/s，1km，进入着陆精度显著提高，这说明考虑控制饱和的进入轨迹优化算法是有效的。

图4 反馈控制Fig.4 Feedback control

图5 末端状态误差分布Fig.5 Terminal state error distribution

5 结束语

本文采用直接配点，对火星大气进入轨迹进行优化设计，并考虑了状态灵敏度对进入轨迹的影响，仿真算例表明，在目标函数中加入状态灵敏度以后，进入轨迹不再靠近边界，控制饱和的现象明显减少，进而减小了末端状态估计误差。

References)

[1]Leavitt J,Mease K.Feasible Trajectory Generation for Atmospheric Entry Guidance[J].Journal of Guidance,Control,and Dynamics,2007,30(2)：473-481.

[2]高滨.火星探测器着陆技术[J].航天返回与遥感,2009,30(1)：1-9,20.GAO Bin.Mars Exploration Entry,Descent and Landing Technologies[J].Spacecraft Recovery & Remote Sensing,2009,30(1)：1-9,20.(in Chinese)

[3]于莹潇, 田佳林.美国新型载人火星探测技术方案[J].航天返回与遥感,2009,30(3)：1-7,15.YU Yingxiao,TIAN Jialin.USA New Mars Exploration Technology Concept[J].Spacecraft Recovery & Remote Sensing,2009,30(3)：1-7,15.(in Chinese)

[4]Shen H J,Seywald H,Powell R W.Desensitizing the Minimum-fuel Power Descent for Mars Pinpoint Landing[J].Journal of Guidance Control and Dynamics,2010,33(1)：108-115.

[5]Seywald H,Kumar R R.Desensitized Optimal Trajectories[R].AMA report No 03-16,2003.

[6]包进进, 荣伟.火星探测器进入阶段稳定性分析[J].航天返回与遥感,2011,32(4)：6-13.BAO Jinjin,RONG Wei.Mars Explorer Stability Analysis in the Entry Phase[J].Spacecraft Recovery & Remote Sensing,2011,32(4)：6-13.(in Chinese)

[7]Seywald H.Desensitized Optimal Trajectories with Control constraints[R].AMA report No.03-28,2003.

[8]Shen H J,Seywald H,Powell R W.Desensitizing the Pin-point Landing Trajectory on Mars[C].AAS/AIAA Astrodynamics Specialist Conference,Honolulu,HI,2008.

[9]Tang S,Conway B A.Optimization of Low-thrust Interplanetary Trajectories Using Collocation and Nonlinear Programming[J].Journal of Guidance Control and Dynamics,1995,18(3)：599-604.

[10]陈晓, 张伟, 彭玉明.基于器间测量的火星进入过程实时高精度导航[J].航天返回与遥感,2012,33(6)：17-23.CHEN Xiao,ZHANG Wei,PENG Yuming.Mars Entry Real-time Navigation Based on Orbiter Tracking Data[J].Spacecraft Recovery & Remote Sensing,2012,33(6)：17-23.(in Chinese)