下载:
-
我国探月工程按照“绕”“落”“回”逐步实施，“嫦娥五号”（Chang'E-5，CE-5）采样返回为我国探月工程三步走画上圆满句号。在采样返回任务中，登月舱完成采样任务后，需要从月面起飞，进入环绕器所在轨道面与之完成交会对接后，一起返回地球。采样返回任务除了涉及月面下降及软着陆技术外，还涉及到月面上升制导技术。月面动力上升段可以看成软着陆动力下降段的逆过程。20世纪“阿波罗”（Apollo）计划创造了人类首次登陆月球并返回的历史，文献[1]中描述了Apollo登月舱上升和下降的过程。
随着我国探月计划的进一步推进，对月面上升制导技术的需求越来越突出。为了进一步提高任务的成功率，需要研究月面应急情况下的上升制导问题，即在出现紧急情况时，上升器能够快速从月面上升，直接进入环绕器停泊轨道，并与之完成交会对接。在应急返回过程中，存在着终端速度位置状态不确定和异面上升等难题。
目前，针对月面下降和软着陆技术研究的成果很多[2-4]，而研究月面上升制导技术的文献相对较少。“嫦娥三号”在动力下降段中，成功应用了自适应动力显式制导技术，并通过终端状态的预测，对制导目标进行修正，以满足接近段的初始条件[5]。对于月面上升制导技术的研究，则主要集中在迭代制导、多项式制导、间接法以及直接法优化求解上升轨迹。巩庆海等[6]研究了迭代制导在月面上升中的应用，对比分析了月面上升段与运载火箭上升段的异同。迭代制导所需的轨道参数需要通过推算或射前注入，迭代初值也需要在起飞前确定，对于应急起飞情况适应性较差。李鑫等[7]针对载人登月任务中登月舱上升入轨的制导问题，分别推导了显式制导律和燃料最优制导律，并进行了对比分析。该显式制导律能够综合考虑节省燃料、抗干扰和入轨精度，只能够适应小角度异面发射问题。李桃取等[8]针对月球探测器大角度异面上升入轨问题，基于E制导和轨道机动原理，提出了一种共面间接上升策略，可以适应大角度异面任务，但存在等待耗时问题。马克茂等[9]以燃料消耗为最优指标，利用极小值原理，将问题转化为时间自由的两点边值问题。采用了一种基于初值的预估方法和向前扫描法对该问题进行求解。但是在建模过程中，只考虑了二维平面。邱丰等[10]针对航天器月面上升在线轨迹规划问题，提出了一种求解最优轨迹的联立框架。通过内点法求解离散后的非线性问题，采用收纳深度控制策略从平衡解的精度和计算效率的角度来改进优化算法的实时性，但仍然存在求解计算量大的问题。
传统的迭代制导、多项式制导都需要提供终端目标位置和速度信息，而对于应急上升任务需要地面及时测算目标参数，过程较为复杂。而间接法则可以将目标轨道根数直接放入终端约束中，解的精度高且满足一阶最优性必要条件，对于一些动力学模型简单的问题，间接法具备一定的优势。国内外学者在间接法求解上升轨迹方面做了大量工作，其中，Lu 等[11-12]基于间接法研究了大气层内的上升制导，采用有限差分和密度同轮技术来求解两点边值问题，国内崔乃刚、黄盘兴等[13-15]也对该方法进行了研究。目前该方法只推导了目标轨道为圆轨道的五约束条件，而没有考虑椭圆轨道的五约束条件；此外，文中采用割线法来调整飞行时间，在工程实现中存在着迭代不收敛的情况。李超兵等[16]针对传统迭代制导小角度修正假设的不足，基于间接法研究了一种多终端约束的最优制导方法，并通过引入权重因子来提高制导方程数值求解精度，该方法仍然需要给出终端的速度和位置，不适用应急上升情况。
本文针对月面应急上升问题，基于间接法研究了一种自适应显式制导算法，通过工程化的近似处理和双层迭代求解策略，使得该制导算法能够满足实时制导需求。该方法目标参数装订简单，能够同时满足应急状况下的共面和异面上升任务需求。
-
月面最优上升制导问题可以描述为根据上升器当前的飞行状态，计算最优制导指令，使得上升器在能够精准进入预定轨道的同时，也能满足燃料消耗最少的性能指标，对于恒定推力的上升器，燃料消耗最少代表着飞行时间最短。对于具有终端速度约束的最优上升制导问题，终端飞行时间自由、以飞行时间最短为性能指标的优化求解过程，与终端时间固定以末端能量最大为性能指标、再通过调节终端飞行时间使得终端速度与目标速度相等的求解过程具有等价性，而且后者在问题求解的收敛性上更优于前者[17-18]。因此，本文首先构建给定终端时间、入轨点能量最大的两点边值问题。
上升器在月球表面处于真空状态，没有大气的影响，主要受到月球引力和发动机的恒定推力，在月面发惯系中，上升器的质心运动模型可以描述为
$$\left\{ {\begin{array}{*{20}{l}} {{\dot{{r}}} = {{V}}} \\ {{\dot{{V}}} = {{g}} + \dfrac{{{T_{{\rm{vac}}}}{{{I}}_{\rm{b}}}}}{m}} \\ {m = {m_0} - \dot mt} \end{array}} \right.$$ (1) 其中：
${{r}} \in {R^3}$ 为发惯系下月心矢径；${{V}} \in {R^3}$ 为发惯系速度；$m$ 为上升器质量；${{g}} \in {R^3}$ 为月面引力加速度；${{{I}}_{\rm{b}}} \in {R^3}$ 为发动机推力方向矢量；${T_{{\rm{vac}}}}$ 为真空发动机推力；${m_0}$ 、$\dot m$ 为初始质量和秒耗量。对于月面上升制导，一般入轨高度在15 km左右，月球引力加速度变化较小，因此对引力加速度做出如下简化
$${{g}}({{r}}) = - \frac{\mu }{{r_0^3}}{{r}}$$ (2) 其中：
$\mu $ 为月球引力常量；${r_0}$ 为当前月心距离大小。随着距离目标越近，引力加速度简化带来的误差越小。另外，为了增加后续数值计算的稳定性，对质心运动模型进行无量纲化处理。取当前月心距离
${r_0}$ 为位置无量纲化系数，当前的引力加速度大小${g_0} = \mu /r_0^2$ 为加速度无量纲化系数，则速度无量纲系数为${v_0} = \sqrt {{g_0}{r_0}} $ ，时间无量纲系数为${t_0} = \sqrt {{r_0}/{g_0}} $ 。无量纲化后的动力学模型为$$\left\{ {\begin{array}{*{20}{l}} {{\dot {\bar{{ r}}}} = {\bar{{V}}}} \\ {{\dot{\bar{{ V}}}} = - {\bar{{r}}} + {{\bar T}_{{\rm{vac}}}}{{{I}}_{\rm{b}}}} \end{array}} \right.$$ (3) 其中：
${\bar T_{{\rm{vac}}}} = {T_{{\rm{vac}}}}/m{g_0}$ ；${\bar{{r}}} \in R{}^3$ 表示无量纲位置矢量；${\bar{{V}}} \in {R^3}$ 表示无量纲速度矢量。入轨能量最大性能指标为
$$\max {\rm{ }}J = \phi ({{{r}}_{\rm{f}}},{{{V}}_{\rm{f}}}) = \frac{1}{2}{{\bar{{V}}}_{\rm{f}}}^{\rm{T}}{{\bar{{V}}}_{\rm{f}}} - \frac{1}{{\left\| {{{{\bar{{r}}}}_{\rm{f}}}} \right\|}}$$ (4) 其中，
${{\bar{{V}}}_{\rm{f}}} \in {R^3}$ 为末端速度矢量；${{\bar{{r}}}_{\rm{f}}} \in {R^3}$ 为末端位置矢量。对于半长轴、偏心率、轨道倾角、入轨点真近点角和升交点赤经约束的五约束问题，可以先考虑终端约束为
$$\left\{ {\begin{array}{*{20}{l}} {\dfrac{1}{2}{{{\bar{{r}}}}_{\rm{f}}}^{\rm{T}}{{{\bar{{r}}}}_{\rm{f}}} - \dfrac{1}{2}\bar {{r}}_{\rm{f}}^{*2} = 0} \\ {{{{\bar{{r}}}}_{\rm{f}}}^{\rm{T}}{{{\bar{{V}}}}_{\rm{f}}} - {{\bar {{r}}}_{\rm{f}}}{{\bar {{V}}}_{\rm{f}}}\sin \gamma _{\rm{f}}^* = 0} \\ {{{{\bar{{r}}}}_{\rm{f}}}{}^{\rm{T}}{{{I}}_{\rm{h}}} = 0} \\ {{{{\bar{{V}}}}_{\rm{f}}}^{\rm{T}}{{{I}}_{\rm{h}}} = 0} \end{array}} \right.$$ (5) 其中：
${{{I}}_{\rm{h}}} \in {R^3}$ 为目标轨道角动量方向矢量。式（5）中第3、第4个等式可确定轨道倾角和升交点赤经，第1、第2个等式配合调整飞行时间使得终端速度满足要求，可确定轨道的半长轴、偏心率和真近点角。根据最优控制理论，构建如下哈密尔顿函数
$$H = \lambda _{\rm{r}}^{\rm{T}}\bar {{V}} + {{\lambda}} _{\rm{V}}^{\rm{T}}\left( { - {\bar{{r}}} + {{\bar T}_{{\rm{vac}}}}{{{I}}_{\rm{b}}}} \right) + \nu \left( {{{I}}_{\rm{b}}^{\rm{T}}{{{I}}_{\rm{b}}} - 1} \right)$$ (6) 其中：
${{{\lambda}} _{\rm{r}}} \in {R^3}$ 、${{{\lambda}} _{\rm{V}}} \in {R^3}$ 为协态变量；$\nu $ 为乘子变量。采用极大值原理有$$H\left( {{{{\lambda}} _{\rm{r}}},{{{\lambda}} _{\rm{V}}},{{\bar {{V}}}^*},{{\bar {{r}}}^*},{{I}}_{\rm{b}}^*} \right) = \mathop {\max }\limits_{{{{I}}_{\rm{b}}}} H\left( {{{{\lambda}} _{\rm{r}}},{{{\lambda}} _{\rm{V}}},{{\bar {{V}}}^*},{{\bar {{r}}}^*},{{{I}}_{\rm{b}}}} \right)$$ (7) 根据式（6）和式（7）可知，当哈密尔顿函数取得最大值时，最优体轴方向
${{I}}_{\rm{b}}^*$ 与${{{\lambda}} _{\rm{V}}}$ 方向重合，因此可得$${{I}}_{\rm{b}}^* = \frac{{{{{\lambda}} _{\rm{V}}}}}{{\left\| {{{{\lambda}} _{\rm{V}}}} \right\|}}$$ (8) 通常
${{{\lambda}} _{\rm{V}}}$ 也被称为主矢量。协态变量方程为
$$\left\{ {\begin{array}{*{20}{c}} {{{\dot {{\lambda}} }_{\rm{r}}} = - \dfrac{{\partial H}}{{\partial \bar {{r}}}}} \\ {{{\dot {{\lambda}} }_{\rm{V}}} = - \dfrac{{\partial H}}{{\partial \bar {{V}}}}} \end{array}} \right.$$ (9) 求解可得
$$\left\{ {\begin{array}{*{20}{l}} {{{\dot {{\lambda}} }_{\rm{r}}} = {{\dot {{\lambda}} }_{\rm{V}}}} \\ {{{\dot {{\lambda}} }_{\rm{V}}} = - {{{{\bf{\lambda}} }}_{\rm{r}}}} \end{array}} \right.$$ (10) 边值条件为
$$\left\{ {\begin{array}{*{20}{l}} {{{{\lambda }}_{\rm{rf}}} = \dfrac{{\partial \phi ({{{\bar{{r}}}}_{\rm{f}}},{{{\bar{{V}}}}_{\rm{f}}})}}{{\partial {{{\bar{{r}}}}_{\rm{f}}}}} + {{\left( {\dfrac{{\partial {{\psi }}}}{{\partial {{{\bar{{r}}}}_{\rm{f}}}}}} \right)}^{\rm{T}}}{{\xi }}} \\ {{{{\lambda }}_{\rm{Vf}}} = \dfrac{{\partial \phi ({{{\bar{{r}}}}_{\rm{f}}},{{{\bar{{V}}}}_{\rm{f}}})}}{{\partial {{{\bar{{V}}}}_{\rm{f}}}}} + {{\left( {\dfrac{{\partial {{\psi }}}}{{\partial {{{\bar{{V}}}}_{\rm{f}}}}}} \right)}^{\rm{T}}}{{\xi }}} \end{array}} \right.$$ (11) 其中：
${{{\lambda }}_{\rm{rf}}} \in {R^3}$ 、${{{\lambda }}_{\rm{Vf}}} \in {R^3}$ 为协态变量末端值；${{\psi }}$ 表示终端约束；${{\xi }}$ 为乘子变量。边值条件的求解具体可参考文献[15]中的推导过程，消去乘子变量可得
$$\left\{ {\begin{array}{*{20}{l}} {{{{\bar{{V}}}}_{\rm{f}}}^{\rm{T}}{{{\lambda }}_{{\rm{Vf}}}} - {{{\bar{{V}}}}_{\rm{f}}}^{\rm{T}}{{{\bar{{V}}}}_{\rm{f}}} = 0} \\ {\left( {{{{\bar{{V}}}}_{\rm{f}}}^{\rm{T}}{{{\lambda }}_{{\rm{rf}}}}} \right){{\bar {{r}}}_{\rm{f}}}^2 - \left( {{{{\bar{{r}}}}_{\rm{f}}}^{\rm{T}}{{{\lambda }}_{{\rm{Vf}}}}} \right){{\bar {{V}}}_{\rm{f}}}^2 - \left( {{{{\bar{{r}}}}_{\rm{f}}}^{\rm{T}}{{{\lambda }}_{{\rm{rf}}}} - {{\bar {{V}}}_{\rm{f}}}^2} \right){{{\bar{{r}}}}_{\rm{f}}}^{\rm{T}}{{{\bar{{V}}}}_{\rm{f}}} = 0} \end{array}} \right.$$ (12) 结合式（5）和式（12）可以得到五约束问题下的6个等式约束
$$\left\{ {\begin{array}{*{20}{l}} {\dfrac{1}{2}{{{\bar{{r}}}}_{\rm{f}}}^{\rm{T}}{{{\bar{{r}}}}_{\rm{f}}} - \dfrac{1}{2}\bar {{r}}_{\rm{f}}^{*2} = 0} \\ {{{{\bar{{r}}}}_{\rm{f}}}^{\rm{T}}{{{\bar{{V}}}}_{\rm{f}}} - {{\bar {{r}}}_{\rm{f}}}{{\bar {{V}}}_{\rm{f}}}\sin \gamma _{\rm{f}}^* = 0} \\ {{{{\bar{{r}}}}_{\rm{f}}}{}^{\rm{T}}{{{I}}_{\rm{h}}} = 0} \\ {{{{\bar{{V}}}}_{\rm{f}}}^{\rm{T}}{{{I}}_{\rm{h}}} = 0} \\ {{{{\bar{{V}}}}_{\rm{f}}}^{\rm{T}}{{{\lambda }}_{{\rm{Vf}}}} - {{{\bar{{V}}}}_{\rm{f}}}^{\rm{T}}{{{\bar{{V}}}}_{\rm{f}}} = 0} \\ {\left( {{{{\bar{{V}}}}_{\rm{f}}}^{\rm{T}}{{{\lambda }}_{{\rm{rf}}}}} \right){{\bar r}_{\rm{f}}}^2 - \left( {{{{\bar{{r}}}}_{\rm{f}}}^{\rm{T}}{{{\lambda }}_{{\rm{Vf}}}}} \right){{\bar V}_{\rm{f}}}^2 - \left( {{{{\bar{{r}}}}_{\rm{f}}}^{\rm{T}}{{{\lambda }}_{{\rm{rf}}}} - {{\bar V}_{\rm{f}}}^2} \right){{{\bar{{r}}}}_{\rm{f}}}^{\rm{T}}{{{\bar{{V}}}}_{\rm{f}}} = 0} \end{array}} \right.$$ (13) 至此，给定终端时间条件下，入轨点能量最大的两点边值问题构建完毕。
-
为了提高计算效率，满足可实时迭代求解的要求，需要对求解两点边值问题进行一些工程化近似处理。要求解终端参数，需要求解微分方程式（3），给定无量纲飞行时间
${\tau _{\rm{f}}}$ ，可得上升器位置和速度矢量为$$ \begin{aligned} & {\left[ {\begin{array}{*{20}{l}} {{\bar{{r}}}({\tau _{\rm{f}}})} \\ {{\bar{{V}}}({\tau _{\rm{f}}})} \end{array}} \right] = \left[ {\begin{array}{*{20}{l}} {\cos ({\tau _{\rm{f}}})}&{\sin ({\tau _{\rm{f}}})} \\ { - \sin ({\tau _{\rm{f}}})}&{\cos ({\tau _{\rm{f}}})} \end{array}} \right]\left[ {\begin{array}{*{20}{l}} {{{{\bar{{r}}}}_0}} \\ {{{{\bar{{V}}}}_0}} \end{array}} \right] +} \\ & {\left[ {\begin{array}{*{20}{l}} {\sin ({\tau _{\rm{f}}})}&{ - \cos ({\tau _{\rm{f}}})} \\ {\cos ({\tau _{\rm{f}}})}&{\sin ({\tau _{\rm{f}}})} \end{array}} \right]\left[ {\begin{array}{*{20}{l}} {{{{I}}_c}({\tau _{\rm{f}}})} \\ {{{{I}}_s}({\tau _{\rm{f}}})} \end{array}} \right]} \end{aligned}$$ (14) $$\left\{ {\begin{array}{*{20}{c}} {{{{I}}_{\rm{c}}}({\tau _{\rm{f}}}) = \int_0^{{\tau _{\rm{f}}}} {{{{I}}_{\rm{b}}}{{\bar T}_{{\rm{vac}}}}\cos (\sigma ){\rm{d}}\sigma } } \\ {{{{I}}_{\rm{s}}}({\tau _{\rm{f}}}) = \int_0^{{\tau _{\rm{f}}}} {{{{I}}_{\rm{b}}}{{\bar T}_{{\rm{vac}}}}\sin (\sigma ){\rm{d}}\sigma } } \end{array}} \right.$$ (15) 其中：
${{\bar{{r}}}_0} \in {R^3}$ 、${{\bar{{V}}}_0} \in {R^3}$ 分别初始位置和速度矢量。对于式（15）中的定积分，可以采用高斯–勒让德求积公式进行快速求解，对于n个节点的高斯–勒让德求积公式，具有2n+ 1次的代数精度[19]。因此可得$$\left\{ {\begin{array}{*{20}{l}} {{\sigma _{\rm{k}}} = \dfrac{{{\tau _{\rm{f}}}}}{2} + \dfrac{{{\tau _{\rm{f}}}}}{2}{x_{\rm{k}}}} \\ {{{{I}}_{\rm{c}}}({\tau _{\rm{f}}}) = \displaystyle\sum\limits_{{\rm{k}} = 0}^5 {{A_{\rm{k}}}{{{I}}_{\rm{b}}}({\sigma _{\rm{k}}})} {{\bar T}_{{\rm{vac}}}}({\sigma _{\rm{k}}})\cos ({\sigma _{\rm{k}}})} \\ {{{{I}}_{\rm{s}}}({\tau _{\rm{f}}}) = \displaystyle\sum\limits_{{\rm{k}} = 0}^5 {{A_{\rm{k}}}{{{I}}_{\rm{b}}}({\sigma _{\rm{k}}})} {{\bar T}_{{\rm{vac}}}}({\sigma _{\rm{k}}})\sin ({\sigma _{\rm{k}}})} \end{array}} \right.$$ (16) 其中：
${{{I}}_{\rm{b}}}({\sigma _{\rm{k}}})$ 按式（8）计算，${x_{\rm{k}}}$ 和${A_{\rm{k}}}$ 分别为5节点求积公式的节点和系数，取表1中的值。表 1高斯–勒让德求积公式节点和系数
Table 1.Gauss-Legendre quadrature formula nodes and coefficients
序号 节点x 系数A 0 –0.932 469 5 0.171 324 5 1 –0.661 209 4 0.360 761 6 2 –0.238 619 2 0.467 913 9 3 +0.238 619 2 0.467 913 9 4 +0.661 209 4 0.360 761 6 5 +0.932 469 5 0.1713 24 5 如果按照
${\bar T_{{\rm{vac}}}} = {T_{{\rm{vac}}}}/m{g_0}$ 计算${\bar T_{{\rm{vac}}}}({\sigma _{\rm{k}}})$ ，则需要显知实际推力和质量，当实际推力和质量存在偏差时，会导致计算不准确。为了提高制导算法对推力和质量偏差的自适应性，需要对${\bar T_{{\rm{vac}}}}({\sigma _{\rm{k}}})$ 的计算做出如下改进$${\bar T_{{\rm{vac}}}}({\sigma _{\rm{k}}}) = \frac{{{V_{{\rm{isp}}}}f_0^{{\rm{ins}}}}}{{({V_{{\rm{isp}}}} - f_0^{{\rm{ins}}}{t_0}{\sigma _{\rm{k}}}){g_0}}}$$ (17) 其中：
${V_{{\rm{isp}}}}$ 为发动机喷气速度大小，$f_0^{{\rm{ins}}}$ 为当前导航视加速度大小；${t_0}$ 、${g_0}$ 为无量纲系数。求解协态变量方程可得
$$\left\{ {\begin{array}{*{20}{l}} {{{{\lambda }}_{\rm{r}}}(\tau ) = {{{\lambda }}_{{\rm{r0}}}}\cos (\tau ) + {{{\lambda }}_{{\rm{V0}}}}\sin (\tau )} \\ {{{{\lambda }}_{\rm{V}}}(\tau ) = - {{{\lambda }}_{{\rm{r0}}}}\sin (\tau ) + {{{\lambda }}_{{\rm{V0}}}}\cos (\tau )} \end{array}} \right.$$ (18) 其中：
${{{\lambda }}_{{\rm{r0}}}} \in {R^3}$ 、${{{\lambda }}_{{\rm{V0}}}} \in {R^3}$ 为协态变量初值，$\tau $ 为无量纲化时间。按照式（8）、式（14）~（18），可以快速计算上升器的终端状态。 -
本文采用收敛性较好的双层迭代求解策略，来计算制导指令。首先求解入轨点能量最大的两点边值问题，再通过一定策略调整飞行时间，使得入轨点速度等于目标速度，即可获得最优的制导程序角。
给定飞行时间
${\tau _{\rm{f}}}$ ，不同的协态变量初值，对应着不同的终端状态，只有满足等式约束的协态变量初值才是需要的。令$ {Y}=\left[ {{{\lambda }}_{\rm{r0}}}^{\rm{T}};{{{\lambda }}_{\rm{V0}}}^{\rm{T}} \right]$ ，$ {{E}} \in {R^6}$ 表示式（13）描述的6个等式约束，则协态变量初值迭代求解如下$$ {{{Y}}_{{\rm{k}} + 1}} = {{{Y}}_{\rm{k}}} - \beta {\left( {\frac{{\partial {{E}}\left( {{{{Y}}_{\rm{k}}}} \right)}}{{\partial {{{Y}}_{\rm{k}}}}}} \right)^{ - 1}}{{E}}\left( {{{{Y}}_{\rm{k}}}} \right) $$ (19) 其中：
$ 0 < \beta \leqslant 1$ ；当$ \left\| {{{E}}\left( {{{{Y}}_{\rm{k}}}} \right)} \right\| \leqslant {\varepsilon _{\rm{\psi }}}$ ，停止迭代；$ {\varepsilon _{\rm{\psi }}}$ 为相应的收敛精度；雅克比矩阵$\partial {{E}}({{{Y}}_{\rm{k}}}){\rm{/}}\partial {{{Y}}_{\rm{k}}}$ 可以通过数值差分计算得到。最后需要调整
${\tau _{\rm{f}}}$ ，使得入轨点的速度大小等于目标值。在文献[15,17~18]中，都采用了如下所示的割线法来迭代求解飞行时间$${\tau _{{\rm{f,k + 1}}}} = {\tau _{{\rm{f,k}}}} - \frac{{{\tau _{{\rm{f,k}}}} - {\tau _{{\rm{f,k - 1}}}}}}{{{{\bar V}_{{\rm{f,k}}}} - {{\bar V}_{{\rm{f,k - 1}}}}}}(\bar V_{\rm{f}}^* - {\bar V_{{\rm{f,k}}}})$$ (20) 在实际仿真中发现，远离终端目标值时，较大范围的
${\tau _{\rm{f}}}$ 能够使入轨点能量最大的两点边值问题存在解，随着上升器越来越接近终端目标，两点边值问题存在解的${\tau _{\rm{f}}}$ 范围越来越小，割线法在迭代求解过程中使得${\tau _{\rm{f}}}$ 不受控，容易超出使得问题有解的范围，导致问题求解失败。为了提高求解的收敛性和可靠性，本文设计了下所示的迭代策略来调整${\tau _{\rm{f}}}$ $$\left\{ \begin{array}{*{20}{l}} {\Delta {{\bar V}_{{\rm{f}},{\rm{k}}}} = \bar V_{\rm{f}}^* - {{\bar V}_{{\rm{f}},{\rm{k}}}}}\\ {\Delta {\tau _{{\rm{f}},{\rm{k}}}} = \left( {\dfrac{{{V_{{\rm{isp}}}}}}{{f_0^{{\rm{ins}}}}} - {t_0}{\tau _{{\rm{f}},{\rm{k}}}}} \right)\left( {1 - \exp \left( { - \dfrac{{\Delta {{\bar V}_{{\rm{f}},{\rm{k}}}}{v_0}}}{{{V_{{\rm{isp}}}}}}} \right)} \right)}\\ {{\tau _{{\rm{f}},{\rm{k}} + 1}} = {\tau _{{\rm{f}},{\rm{k}}}} + \rho \Delta {\tau _{{\rm{f}},{\rm{k}}}},0 < \rho \leqslant 1} \end{array}\right. $$ (21) 其中：
${t_0}$ 、${g_0}$ 、${v_0}$ 为无量纲系数；$\rho $ 为修正因子，用来控制迭代的收敛性；当$\left| {\Delta {{\bar V}_{{\rm{f,k}}}}} \right| < {\varepsilon _{\rm{v}}}$ 时，停止迭代，${\varepsilon _{\rm{v}}}$ 为速度收敛精度。一般对于共面上升问题，$\rho $ 可取1，对于异面上升制导，$\rho $ 需要随着异面角度的增大，适当减小，保证迭代的收敛性。整个制导指令双层迭代求解流程如图1所示。其中，第一次启动算法时，对于${\tau _{\rm{f}}}$ 的初值选取需要选的偏大，对于一般设计好的任务，直接装订额定工作时间，而在后续迭代中，可直接使用上一次迭代出的结果。
图 1制导指令双层迭代求解策略
Figure 1.Two-layer iterative solution strategy for guidance instruction
-
对于一般的两点边值问题，协态变量初值猜测是困难的，但是在本文研究的问题中，协态变量主矢量
${{{\lambda }}_{\rm{V}}}$ 代表了推力的最佳方向，间接赋予了协态变量${{{\lambda }}_{{\rm{V0}}}}$ 一定的“物理含义”。通过仿真分析发现，在协态变量初值中，影响两点边值问题收敛性较大的是主矢量${{{\lambda }}_{{\rm{V0}}}}$ ，通过对实际飞行任务的分析可知，推力方向决定了速度的方向，进而决定了上升器运动的方向。因此，在迭代过程中，主矢量${{{\lambda }}_{{\rm{V0}}}}$ 初始猜测值可按照下式计算获得，而协态变量${{{\lambda }}_{{\rm{r0}}}}$ 的初始猜测值直接取${{\bf{\lambda }}_{{\rm{r0}}}} = {\bf{0}}$ 。$${{{\lambda }}_{{\rm{V0}}}} = {{{I}}_{\rm{h}}} \times {\rm{unit(}}{{{r}}_0}{\rm{)}}$$ (22) 其中：
${\rm{unit(}}{{{r}}_0}{\rm{)}}$ 表示当前月心位置方向矢量；${{{I}}_{\rm{h}}} \times {\rm{unit(}}{{{r}}_0}{\rm{)}}$ 表示沿着目标轨道面的运动方向。 -
为了验证制导算法的性能，本文分别对月面共面上升问题和异面上升问题进行了仿真，并考核了制导算法对秒耗量和比冲偏差的适应性。其中，仿真参数如表2所示，暂不考虑发动机后效的影响。
表 2仿真参数
Table 2.Simulation parameters
参数 数值 初始质量/ kg 5 070 推力/ N 24 500 秒耗量/ （kg·s-1） 8.167 比冲/（m·s-1） 3 000 发射点方位角/（°） 90 发射点纬度/（°） 18 发射点经度/（°） 56.784 轨道半长轴/m 1 795 582 轨道倾角/° 20 轨道偏心率 0.023 669 2 入轨高度/m 15 000 对共面任务仿真，设目标轨道面升交点赤经为300°，外层迭代修正因子
$\rho = 1.0$ 。图2为额定工况下共面上升制导仿真结果。先垂直上升，10 s后加入本文的自适应制导律。从图2中可以看出，俯仰、偏航程序角满足线性正切规律。由表3可知，半长轴、偏心率、轨道倾角和升交点赤经具有相当高的精度，自适应制导律对秒耗量偏差和比冲偏差均具有良好的适应性。求解共面上升制导指令过程中，内层求解协态变量初值平均迭代收敛次数为3次，外层调整${\tau _{\rm{f}}}$ 平均迭代收敛次数为2次。3.4 GHz PC机Matlab仿真环境下，共面制导任务单次计算耗时约20 ms。
图 2月面共面上升制导仿真结果
Figure 2.Simulation results of lunar coplanar ascending guidance
表 3共面上升入轨精度
Table 3.Orbit entry accuracy in coplanar ascent
仿真工况 半长轴△a/km 偏心率△e 轨道倾角△i/（°） 升交点赤经△Ω/（°） 飞行时间/s 额定 0.005 21 6.685 8×10-6 –1.725 1×10-5 0.800×10-4 280.301 秒耗量（+10%） 0.00848 4.042 3×10-6 –2.226 4×10-5 1.072 4×10-4 255.145 秒耗量（–10%） 0.003 00 6.086 2×10-6 –1.134 4×10-5 0.518 2×10-4 311.430 比冲（+5%） –0.003 46 3.082 4×10-6 –1.116 7×10-5 0.529 3×10-4 270.577 比冲（–5%） 0.006 22 1.178 2×10-5 –1.877 8×10-5 0.875 7×10-4 290.805 对异面任务仿真，设目标轨道面升交点赤经为286.823°，外层迭代修正因子
$\rho = 0.75$ 。图3为异面上升情况下的仿真结果。从图3（d）可以看出，异面上升时，存在较大的偏航角。图3（f）为异面上升的时，垂直轨道面侧向距离的变化曲线。入轨高度为15 km，垂直轨道面侧向机动距离高达到72 km。仿真结果表明，本文设计的自适应制导算法能够适应一定程度的异面上升任务。表3为异面上升入轨精度统计，对于异面上升情况，半长轴和偏心率仍然具有较高的精度，但是轨道倾角精度、升交点赤经精度不如共面上升的高，其中主要原因是在上升过程中，同时存在垂直轨道面的侧向机动，会对入轨精度产生影响。在求解异面上升制导指令过程中，内层求解协态变量初值平均迭代收敛次数为4次，外层调整${\tau _{\rm{f}}}$ 平均迭代收敛次数为3次。3.4 GHz PC机Matlab仿真环境下，异面制导任务单次计算耗时约21.8 ms。
图 3月面异面上升制导仿真结果
Figure 3.Simulation results of lunar non-coplanar ascent guidance
表 4异面上升入轨精度
Table 4.Orbit entry accuracy in non-coplanar ascent
仿真工况 半长轴△a/ km 偏心率△e 轨道倾角△i/（°） 升交点赤经△Ω/（°） 飞行时间/s 额定 0.004 55 6.694×10-6 0.002 42 –0.007 154 304.586 秒耗量（+10%） –0.002 79 1.223×10-5 0.003 20 –0.009 676 281.021 秒耗量（–10%） 0.004 55 5.713×10-6 0.002 53 –0.006 156 318.409 比冲（+5%） 0.004 20 5.737×10-6 0.000 88 –0.002 625 295.794 比冲（–5%） –0.007 07 1.496×10-5 0.003 75 –0.011 031 314.131 本文设计的制导算法能够适应0~9°异面角度上升任务。其中，直接异面上升入轨和先共面上升再轨道转移入轨的燃料消耗如图4所示，可以看出异面角在3°以内时，两者燃耗相差不大，随着异面角度进一步增大，直接异面上升入轨相比轨道转移入轨燃耗显著增大，此时直接异面上升入轨不再适用。文献[20]中针对载人登月短期访问任务，考虑应急返回需求，给出了月面调整角度计算方法和任务期间返回上升轨道与目标轨道平面夹角最小化的条件。按照该方法设计目标交会轨道的轨道倾角，对于非极地地区的6 d以内的短期月面探测任务，可使得异面角度最大不超过4.86°；对于低纬度地区探测任务，可通过设计目标轨道倾角使得异面角度最大不超过3°。本文设计的制导方法在任务周期内能够全程适用。
图 4异面入轨和轨道转移入轨燃料消耗仿真结果
Figure 4.Simulation result of fuel consumption for non-coplanar ascent and orbital transfer
-
针对月面应急上升任务，本文研究了一种自适应显式制导律。该自适应制导律目标参数装订简单，能够同时适应月面共面和异面上升任务。从对偏差的适应性来看，在存在大秒耗量偏差和比冲偏差的情况下，该自适应制导律仍然具有较高的入轨精度。从收敛性上来看，本文设计的双层求解策略，改进的外层时间调整策略和内层协态变量初值的选取策略，能够有效保证制导律迭代求解收敛。
Research on Adaptive Guidance Technology for Lunar Emergency Ascent
-
摘要:针对月面应急上升问题，研究了一种显式自适应制导方法。根据极大值原理，推导建立了5个约束条件下最大能量入轨的两点边值问题。设计了一种双层迭代求解策略，内层通过牛顿迭代求解最大能量入轨的两点边值问题，外层通过调整时间使得速度满足目标速度。在外层迭代中，设计了一种时间迭代调整策略；在内层求解两点边值问题中，根据任务中推力方向变化规律，设计一种协态变量主矢量初值的选取策略。仿真结果表明，设计的制导律能可靠收敛，目标参数装订简单，可适应共面上升任务和异面上升任务。在存在秒耗量和比冲偏差的情况下，仍具有较高的制导精度。Abstract:In this paper, an adaptive guidance method is studied for the emergency ascent of the lunar surface. Firstly, according to the principle of maximum value, the two-point boundary value problem of maximum energy entering orbit under five constraints is derived and established. Secondly, a two-layer iterative solution strategy is designed. The inner layer uses Newton iteration to solve the two-point boundary value problem of maximum energy into orbit, and the outer loop adjusts the time to make the speed meet the target speed. In the outer loop iteration, a time iterative adjustment strategy is designed. In the inner loop, according to the change rule of the thrust direction in the task, a strategy for selecting the initial value of the principal vector of the covariant variable is designed to solve the two-point boundary value problem. The simulation results show that the guidance law designed in this paper can converge reliably, and the target parameter binding is simple. It can adapt to coplanar ascending tasks and different-plane ascending tasks. In the presence of second consumption and specific impulse deviation, it still has high guidance accuracy.Highlights
● Five constraints on semi-major axis，eccentricity，orbit inclination，true anomalous angle of entry point，and right ascension of ascending node are established. The five constraints are more general and can satisfy both circular and elliptical orbits. ● The secant method adjusting the flight time strategy is improved by engineering，and the improved adjustment strategy can make the flight time converge to the target value quickly and monotonously，which improves the convergence of the guidance algorithm. ● The guidance method designed in this paper can be adapted to coplanar and non-planar lunar ascent missions at the same time. The target parameter binding is simple，the solving efficiency is high，and the real-time guidance needs can be met. ● The adaptive guidance law still has high orbit accuracy under the condition of large second consumption deviation and specific impulse deviation. -
图 4异面入轨和轨道转移入轨燃料消耗仿真结果
Fig. 4Simulation result of fuel consumption for non-coplanar ascent and orbital transfer
表 1高斯–勒让德求积公式节点和系数
Table 1Gauss-Legendre quadrature formula nodes and coefficients
序号 节点x 系数A 0 –0.932 469 5 0.171 324 5 1 –0.661 209 4 0.360 761 6 2 –0.238 619 2 0.467 913 9 3 +0.238 619 2 0.467 913 9 4 +0.661 209 4 0.360 761 6 5 +0.932 469 5 0.1713 24 5 
下载:
导出CSV
表 2仿真参数
Table 2Simulation parameters
参数 数值 初始质量/ kg 5 070 推力/ N 24 500 秒耗量/ （kg·s-1） 8.167 比冲/（m·s-1） 3 000 发射点方位角/（°） 90 发射点纬度/（°） 18 发射点经度/（°） 56.784 轨道半长轴/m 1 795 582 轨道倾角/° 20 轨道偏心率 0.023 669 2 入轨高度/m 15 000 下载:
导出CSV
表 3共面上升入轨精度
Table 3Orbit entry accuracy in coplanar ascent
仿真工况 半长轴△a/km 偏心率△e 轨道倾角△i/（°） 升交点赤经△Ω/（°） 飞行时间/s 额定 0.005 21 6.685 8×10-6 –1.725 1×10-5 0.800×10-4 280.301 秒耗量（+10%） 0.00848 4.042 3×10-6 –2.226 4×10-5 1.072 4×10-4 255.145 秒耗量（–10%） 0.003 00 6.086 2×10-6 –1.134 4×10-5 0.518 2×10-4 311.430 比冲（+5%） –0.003 46 3.082 4×10-6 –1.116 7×10-5 0.529 3×10-4 270.577 比冲（–5%） 0.006 22 1.178 2×10-5 –1.877 8×10-5 0.875 7×10-4 290.805 下载:
导出CSV
表 4异面上升入轨精度
Table 4Orbit entry accuracy in non-coplanar ascent
仿真工况 半长轴△a/ km 偏心率△e 轨道倾角△i/（°） 升交点赤经△Ω/（°） 飞行时间/s 额定 0.004 55 6.694×10-6 0.002 42 –0.007 154 304.586 秒耗量（+10%） –0.002 79 1.223×10-5 0.003 20 –0.009 676 281.021 秒耗量（–10%） 0.004 55 5.713×10-6 0.002 53 –0.006 156 318.409 比冲（+5%） 0.004 20 5.737×10-6 0.000 88 –0.002 625 295.794 比冲（–5%） –0.007 07 1.496×10-5 0.003 75 –0.011 031 314.131 下载:
导出CSV
-
[1] BENNETTT F V. Apollo experience report mission planning for lunar module descent and ascent: NASA Technical Note D-6846[R]. [S. l.]: NASA, 1972. [2] 薛志飞. 月球飞船动力下降制导技术研究[D]. 南京: 南京航空航天大学, 2018.XUE Z F. Guidance for the powered descent of lunar spacecraft[D]. Nanjing: Nanjing University of Aeronautics and Astronautics, 2018. [3] 李茂登. 月球软着陆自主导航、制导与控制问题研究[D]. 哈尔滨: 哈尔滨工业大学, 2007.LI M D. Study the guidance, autonomous navigation and control of lunar soft landing[D]. Harbin: Harbin Institute of Technology, 2007. [4] 薛志飞,韩艳铧. 一种月球飞船动力下降预测制导方法的研究[J]. 航天控制,2017,5(6):23-30.XUE Z F, HAN Y H. A predictive guidance method for the powered descent of lunar spacecraft[J]. Aerospace Control,2017,5(6):23-30. [5] 黄翔宇,张宏华,王大轶,等. “嫦娥三号”探测器软着陆自主导航与制导技术[J]. 深空探测学报(中英文),2014,1(1):52-59.HUANG X Y, ZHANG H H, WANG D Y, et al. Autonomous navigation and guidance for Chang’e-3 soft landing[J]. Journal of Deep Space Exploration,2014,1(1):52-59. [6] 巩庆海,宋征宇,吕新广. 迭代制导在月面上升段的应用研究[J]. 载人航天,2015,21(3):33-38.GONG Q H, SONG Z Y, LV X G. Study on application of iterative guidance in lunar ascent[J]. Manned Spaceflight,2015,21(3):33-38. [7] 李鑫,刘莹莹,周军. 载人登月舱上升入轨段的制导律设计[J]. 系统工程与电子技术,2011,33(11):2480-2484.doi:10.3969/j.issn.1001-506X.2011.11.26LI X,LIU Y Y,ZHOU J. Design of guidance law for lunar ascent phase of manned lunar module[J]. Systems Engineering and Electronics,2011,33(11):2480-2484.doi:10.3969/j.issn.1001-506X.2011.11.26 [8] 李桃取,董长虹,王峰波. 大角度异面动力上升制导策略研究[J]. 航天控制,2015,33(2):32-38.doi:10.3969/j.issn.1006-3242.2015.02.006LI T Q, DONG C H, WANG F B. Research on guidance strategy of lunar non-coplanar powered ascent with large angle[J]. Aerospace Control,2015,33(2):32-38.doi:10.3969/j.issn.1006-3242.2015.02.006 [9] 马克茂,陈海朋. 登月舱上升段最优轨迹设计[J]. 中国空间科学技术,2013,33(2):54-60.doi:10.3780/j.issn.1000-758X.2014.02.007MA K M, CHEN H P. Optimal trajectory design for the lunar module in ascent stage[J]. Chinese Space Science and Technology,2013,33(2):54-60.doi:10.3780/j.issn.1000-758X.2014.02.007 [10] 邱丰,宋征宇. 联立法求解月面上升段最优轨迹的快速收敛控制技术[J]. 载人航天,2015,21(1):6-18.doi:10.3969/j.issn.1674-5825.2015.01.002QIU F,SONG Z Y. Fast convergence control on simultaneous approach based trajectory design during lunar ascent[J]. Manned Spaceflight,2015,21(1):6-18.doi:10.3969/j.issn.1674-5825.2015.01.002 [11] LU P,PAN B F. Highly constrained optimal launch ascent guidance[J]. Journal of Guidance Control,and Dynamics,2010,33(2):404-414.doi:10.2514/1.45632 [12] LU P,GRIFFIN B J,DUKEMAN G A,et al. Rapid optimal multi burn ascent planning and guidance[J]. Journal of Guidance,Control,and Dynamics,2008,31(6):1656-1664. [13] 崔乃刚,黄盘兴,韦常柱,等. 基于混合优化的运载器大气层内闭环制导方法[J]. 中国惯性技术学报,2015,23(3):328-333.CUI N G,HUANG P X, WEI C Z,et al. Closed-loop endo-atmospheric guidance of launch vehicle based on hybrid optimization approach[J]. Journal of Chinese Inertial Technology,2015,23(3):328-333. [14] HUANG P,WEI C,GU Y,et al. A symplectic optimisation method for rapid endo-atmospheric ascent trajectory planning[J]. International Journal of Modelling,Identification and Control,2015,24(3):196-205. [15] 黄盘兴. 运载器大气层内上升段轨迹快速优化方法研究[D]. 哈尔滨: 哈尔滨工业大学, 2015.HUANG P X. Study on endo-atmospheric ascent trajectory rapid optimization for launch vehicles[D]. Harbin: Harbin Institute of Technology, 2015. [16] 李超兵,王晋麟,李海. 一种基于多终端约束的最优制导方法[J]. 中国空间科学技术,2016,36(5):9-17.LI C B, WANG J L,LI H. An optimal guidance method based on multiple terminal constraints[J]. Chinese Space Science and Technology,2016,36(5):9-17. [17] 傅瑜. 升力式天地往返飞行器自主制导方法研究[D]. 哈尔滨: 哈尔滨工业大学, 2012.FU Y. Autonomous guidance method for lift transportation vehicle[D]. Harbin: Harbin Institute of Technology, 2012. [18] LU P,SUN H,TSAI B. Closed-loop endoatmospheric ascent guidance[J]. Journal of Guidance,Control,and Dynamics,2003,26(2):283-294. [19] 李庆杨, 王能超, 易大义. 数值分析（第5版）[M]. 北京: 清华大学出版社, 2008: 228-229. [20] 陈海朋,余薛浩,黄飞. 载人登月应急返回轨道倾角优化设计[J]. 中国空间科学技术,2017,37(4):69-74.CHEN H P,YU X H, HUANG F. Optimization design of inclination orbit of human lunar landing emergency return trajectory[J]. Chinese Space Science and Technology,2017,37(4):69-74. [21] 柳嘉润,巩庆海,翟雯婧. 智能自主系统及其航天控制应用[J]. 飞控与探测,2018,1(1):059-062.LIU J R, GONG Q H,ZHAI W J. Intelligent autonomous system and application in aerospace[J]. Flight Control & Detection,2018,1(1):059-062. [22] 刘付成. 人工智能在航天控制中的应用[J]. 飞控与探测,2018,1(1):016-025.LIU F C. Application of artificial intelligence in spacecraft[J]. Flight Control & Detection,2018,1(1):016-025. -
点击查看大图
图(4)/
表 (4)
计量
- 文章访问数:368
- HTML全文浏览量:142
- PDF下载量:39
- 被引次数:0
