一种月球借力约束的地月平动点转移轨道设计方法

摘要

本发明涉及一种月球借力约束的地月平动点转移轨道设计方法，属于航天器轨道设计与优化领域。选择目标Halo轨道入轨点位置与期望的月球借力约束大小，利用方程(1)对入轨点的状态量进行积分，至航天器到达预期的月球借力位置。通过相关的优化算法对入轨点状态量进行修正调整，使得借力位置的状态满足方程(3)的约束条件，进而确定航天器进入Halo轨道所需要的机动速度增量大小。通过合理选择月球借力约束条件及Halo轨道入轨点等参数，能够有效地设计出满足任务要求的低耗能转移轨道。

说明书

技术领域

本发明涉及一种月球借力约束的地月平动点转移轨道设计方法，属于航天器 轨道设计与优化领域。

背景技术

近年来，地月系统L2平动点(EML2)附近Halo轨道成为航天工程及应用 领域的研究重点。通过EML2点Halo轨道，探测器能够观测月球背面并为月球 着陆任务提供导航策略。同时，围绕EML2点能够设计出到达月球、火星与其 它行星低耗能的转移轨道，可以为探测深空环境提供最佳的落脚点(徐明，徐 世杰.地-月系平动点及Halo轨道的应用研究[J].宇航学报.2006,27(4)： 695-699)。因此，为了降低平动点转移轨道设计过程中的燃料消耗，有必要对 探测器采用的飞行方式进行研究。

对于EML2点转移轨道设计问题，轨道转移方式通常采用霍曼直接转移，间 接转移，弱稳定转移与月球借力转移。例如Zazzera采用遗传算法与序列二次规 划的混合优化方法，对不考虑月球借力的多条间接转移轨道进行分析(ZazzeraF B,TopputoF,MassariM.Assessmentofmissiondesignincludingutilizationof librationpointsandweakstabilityboundaries[R].ESTECContractNo. 18147/01/NL/MV,2004)；Parker基于弱稳定边界与不变流形转移策略，通过连 接日地系统与地月系统的不变流形，设计了从低地球轨道到达EML2点Halo轨 道的长时间飞行转移轨道(ParkerJS.Low-EnergyBallisticLunarTransfers[D]. Colorado:UniversityofColorado,2007)；Gordon利用稳定流形理论与月球借力 技术，构造了两脉冲转移轨道，但降低燃耗的能力有限；Li则在Gordon的基础 上，设计了低燃耗的三脉冲转移轨道方案，但并未给出Halo轨道入轨点(HOI) 的选择策略与分析所考虑的月球借力约束及其影响，在工程应用中存在着一定 的局限性(GordonDP.Transferstoearth-moonL2haloorbitsusinglunarproximity andinvariantmanifolds[D].Indiana:PurdueUniversity,2008；LiMT,ZhengJH. Impulsivelunarhalotransfersusingthestablemanifoldsandlunarflybys[J].Acta Astronautica,2010,66:1481-1492)。

因此，对以EML2点Halo轨道为目标轨道的深空探测器而言，在低燃耗条 件下设计转移轨道的技术成为一个关键性问题，如何选择月球借力约束条件， 并综合考虑Halo轨道入轨点的影响成为突出问题。

发明内容

本发明的目的是为了降低平动点转移轨道设计过程中的燃耗问题。本发明公 开的一种月球借力约束的地月平动点转移轨道设计方法，通过合理选择月球借 力约束条件及Halo轨道入轨点等参数，能够有效地设计出满足任务要求的低耗 能转移轨道。

本发明的目的是通过下述技术方案实现的：

一种月球借力约束的地月平动点转移轨道设计方法，具体步骤如下：

步骤一、航天器施加第一次脉冲，由目标Halo轨道反推至月球借力位置；

选择目标Halo轨道入轨点位置与期望的月球借力约束大小，利用方程(1) 对入轨点的状态量进行积分，至航天器到达预期的月球借力位置。通过相关的 优化算法对入轨点状态量进行修正调整，使得借力位置的状态满足方程(3)的 约束条件，进而确定航天器进入Halo轨道所需要的机动速度增量大小。

步骤二、针对不同的任务轨道，航天器在借力点施加不同的脉冲次数，使 得航天器由月球借力位置进入地-月转移轨道段，最终反推至地球停泊轨道。

当所述任务轨道包含绕月飞行的地-月转移轨道时，(a)选择绕月飞行轨迹 末端满足的约束大小，利用方程(1)对步骤一中确定的月球借力点状态量进行 积分，通过相关的优化算法对状态量进行修正调整，使得借力位置的状态满足 方程(3)的约束条件，进而确定航天器绕月飞行所需的速度增量大小。(b)选 择期望的地球停泊轨道约束大小，利用方程(1)对(a)中确定的绕月飞行轨 迹末端的状态进行积分，至航天器到达预期的地球停泊轨道位置，通过相关的 优化算法对绕月飞行轨迹末端的状态量进行修正调整，进而确定航天器进入地- 月转移轨道段所需的速度增量大小、地球逃逸点速度增量大小与总飞行时间。

当所述任务轨道为地-月转移轨道时，选择期望的地球停泊轨道约束大小， 利用方程(1)对步骤一中确定的月球借力点状态量进行积分，至航天器到达预 期的地球停泊轨道位置。通过相关的优化算法对月球借力点状态量进行修正调 整，进而确定航天器在月球借力点与地球逃逸点的速度增量大小与总飞行时间。

以上两种任务轨道，在地球停泊轨道上都满足高度约束与航迹角约束，并 且逃逸点机动速度增量方向与航天器速度方向共线。

步骤三、整个平动点转移轨道在设计过程中采用的是逆向积分的策略，因 此，步骤一和步骤二中施加的脉冲与实际工程中航天器转移轨道中机动顺序相 反。

步骤一所述方程(1)为在地球-月球质心旋转系下建立的航天器动力学方程。

其中坐标系的原点为系统的质心，x轴由地球指向月球，z轴与地-月系统的 角速度方向一致，y轴与x轴，z轴垂直，构成右手坐标系。

在质心旋转系下，航天器的动力学模型表示为

$\left(\begin{array}{c}\stackrel{··}{x}-2\stackrel{·}{y}=x-\frac{(1-\mu )(x+\mu )}{r_e^3}-\frac{\mu (x-1+\mu )}{r_m^3}\\ \stackrel{··}{y}+2\stackrel{·}{x}=y-\frac{(1-\mu )y}{r_e^3}-\frac{\mu y}{r_m^3}\\ \stackrel{··}{z}=-\frac{(1-\mu )z}{r_e^3}-\frac{\mu z}{r_m^3}\end{array}\right)---\left(1\right)$

其中表征航天器在质心旋转系中的位置、速度与加速度状 态量，μ表示系统的质量系数，为航天器与地球之间的距离， 为航天器与月球之间的距离。

所述月球借力位置的约束条件为：在转移轨道设计过程中，当月球借力位 置考虑不同的约束条件时，将对整个轨道转移过程的速度增量与飞行时间等参 数产生较大的影响。为了达到降低燃耗的目的，考虑月球借力位置的约束条件 分别为：近月高度h_m，相对月球飞行航迹角γ_m与借力方位角δ_m，可表示为：

$\left(\begin{array}{c}{h}_m=\sqrt{{\left(x_{fm}-1+\mu \right)}^2+y_{fm}^2+2_{fm}^2}-R_m\\ {\gamma }_m=\arcsin \left(\frac{\left(x_{fm}-1+\mu \right){\stackrel{·}{x}}_{fm}+y_{fm}{\stackrel{·}{y}}_{fm}+z_{fm}{\stackrel{·}{z}}_{fm}}{\left|\right|P_{fm}\left|\right|\left|\right|V_{fm}\left|\right|}\right)\\ {\delta }_m=\arccos \left(\frac{\left(x_{fm}-1+\mu \right){\stackrel{·}{y}}_{fm}-y_{fm}{\stackrel{·}{x}}_{fm}}{\left|\right|Q_{fm}\left|\right|}\right)\end{array}\right)---\left(2\right)$

其中，r_fm＝[x_fm,y_fm,z_fm]^T与为旋转系中借力点位置与速度矢 量，P_fm＝[x_fm-1+μ,y_fm,z_fm]^T为航天器相对于月球质心的位置矢量，R_m为月球半径， 矢量$Q_{fm}=\left(\begin{array}{c}{y}_{fm}{\stackrel{·}{z}}_{fm}-z_{fm}{\stackrel{·}{y}}_{fm}\\ z_{fm}{\stackrel{·}{x}}_{fm}-(x_{fm}-1+\mu ){\stackrel{·}{z}}_{fm}\\ (x_{fm}-1+\mu ){\stackrel{·}{y}}_{fm}-y_{fm}{\stackrel{·}{x}}_{fm}\end{array}\right).$

同时，为了充分利用目标Halo轨道的特性，采用逆向积分的策略设计地月 平动点转移轨道。在保证高度h_m与航迹角γ_m约束的条件下，改变借力方位角δ_m的 取值大小，能够得到不同约束条件时燃料消耗变化情况，进而确定设计方案的 各项参数。

所述采用逆向积分的策略，即为将步骤一中的公式(1)逆向积分；

步骤一、步骤二中的优化算法为微分修正算法，其作为一种局部优化算法， 具有快速收敛特性，基于此算法对由公式(1)积分确定的轨迹进行调整，可以 保证轨迹末端的状态量满足给定的约束条件。首先，选取自由变量C与约束矢量 F(C)，满足

$\left(\begin{array}{c}C={[V_x,V_y,V_z,T]}^T\\ F\left(C\right)={\left(\begin{array}{c}{h}_m-h_m^{*}\\ {\mathrm{sin\gamma }}_m-{\mathrm{sin\gamma }}_m^{*}\\ {\mathrm{cos\delta }}_m-{\mathrm{cos\delta }}_m^{*}\end{array}\right)}_{moon}\end{array}\right)---\left(3\right)$

其中，[V_x,V_y,V_z]^T为Halo轨道入轨点三轴方向速度矢量，T为Halo轨道至月 球借力位置的飞行时间；(h_m,γ_m,δ_m)由公式(2)确定，为给定的期望约 束值。

然后，推导约束矢量关于控制变量的偏导数关系

$\delta F=\frac{\partial F}{\partial V_x}{\mathrm{\delta V}}_x+\frac{\partial F}{\partial V_y}{\mathrm{\delta V}}_y+\frac{\partial F}{\partial V_z}{\mathrm{\delta V}}_z+\frac{\partial F}{\partial T}\delta T---\left(4\right)$

其中，符号δ表示变分关系，为约束矢量关于速度变量V_x的偏导数。最 后利用最小二乘法对自由变量C进行微分修正，迭代求解能够保证最终的转移轨 道满足给定的约束条件。

有益效果

1、本发明公开的一种月球借力约束的地月平动点转移轨道设计方法，由于 同时考虑了多种月球借力约束条件，相比于不考虑月球借力或少约束情况的速 度增量小，进而降低了任务燃料消耗。

2、本发明公开的一种月球借力约束的地月平动点转移轨道设计方法，根据 不同的飞行时间与燃料消耗等任务要求，选择不同的目标Halo轨道、不同的Halo 轨道入轨位置、航天器在月球借力点满足的约束条件，可以完成相应的转移轨 道设计任务。

3、本发明公开的一种月球借力约束的地月平动点转移轨道设计方法，在月 球借力点多施加一次脉冲机动，能够形成绕月飞行轨道，从而近距离的对月球 背面进行观测，为未来月球背面搭建无线电望远镜时进行物质运输投放任务提 供参考。

附图说明

图1为月球借力约束下的三脉冲与绕月飞行四脉冲轨道设计示意图；

图2为一种月球借力约束的地月平动点转移轨道方法流程图；

图3为公式(1)中地球-月球质心旋转坐标系的示意图；

图4为步骤一中月球借力点约束条件示意图；

图5为本发明实施例1不同借力约束条件的三脉冲转移轨道；

图6为本发明实施例1不同参数下转移轨道所需速度增量对比分析图；

图7为本发明实施例2不同借力约束条件的绕月飞行四脉冲转移轨道。

具体实施方式

为了更好地说明本发明的目的和优点，下面将结合附图对本发明作进一步 的详细说明。

实施例1

如图1、图2所示，针对三脉冲设计方案，一种月球借力约束的地月平动点 转移轨道设计方法，包括如下步骤：

步骤一、航天器施加第一次脉冲，由目标Halo轨道反推至月球借力位置, 月球借力点约束如图4所示；

选择法向幅值Az＝10000km的Halo轨道为目标轨道，期望的月球借力约束 条件满足：近月高度飞行航迹角借力方位角选取 Halo轨道入轨点状态为[1.12507,-0.08292,-0.00461,-0.03415,0.06426,-0.04819]^T，利用 公式(1)对Halo轨道入轨点状态量进行积分，在满足航迹角约束的位置中， 距离月球最近点的高度为h_m＝1.3978×10⁴km，借力方位角δ_m＝16.434°，不满足期望 的借力高度与方位角约束，采用相关的优化算法对入轨点状态量进行修正调整， 使得借力位置的状态量满足公式(3)的约束矢量，进而确定航天器进入Halo 轨道所需的机动速度增量大小为39.62m/s。

步骤二、航天器在月球借力位置施加第二次脉冲，使航天器由月球借力位 置进入地-月转移轨道段，最终反推至地球停泊轨道。

选取期望的地球停泊轨道约束条件满足：轨道高度飞行航迹角 利用公式(1)对月球借力位置状态进行积分，通过相关的优化算法进 一步对月球借力位置状态量修正调整，直至满足给定的停泊轨道约束条件，则 计算出航天器在月球借力点的速度增量大小为190.24m/s，总飞行时间为13.63 天，从地球停泊轨道逃逸的速度增量为3120.33m/s。

当月球借力位置考虑少约束条件，即仅包含高度约束与航迹角约束时，设 计整个平动点转移轨道所需的燃耗为：地球逃逸点为3133.41m/s，月球借力点 为685.01m/s，Halo轨道入轨点为9.67m/s，而总飞行时间为13.25天。比较以上 两种情况可知，在飞行时间相近条件下，考虑方位角约束后，总速度增量显著 降低，两者相差477.9m/s。

针对以上两种情况设计了三脉冲转移轨道，如图5所示，其中虚线为仅考 虑借力高度与航迹角约束的转移轨道，实线为添加借力方位角后的转移轨道。

为了确定Halo轨道入轨点的选择依据，选取不同目标Halo轨道、不同Halo 轨道入轨位置、航天器在月球借力点满足的约束条件进行设计与分析，不同参 数下的轨道所需速度增量的对比分析图如图6所示。由图6可知：对于相同编 号入轨点，Halo轨道入轨点速度增量随着法向幅值的增大而增加，借力点机动 速度增量不断下降，总速度增量逐步增大，范围满足ΔV_total＝(3320,3340)m/s。同时， 能够满足约束的入轨点数目随着借力方位角的不断减小而减少，综合分析数据， 能够较好地满足月球借力高度、航迹角与借力方位角约束的Halo轨道入轨点选 取区间为(90,160)与(240,315)。以上参数能够为平动点轨道的设计与应用、参数 的选择等方面提供依据。

步骤三、整个平动点转移轨道在设计过程中采用的是逆向积分策略，因此， 步骤一和步骤二中施加的脉冲顺序与实际工程中航天器转移轨道中机动顺序相 反。即航天器首先施加一次脉冲从地球停泊轨道逃逸，进入地-月转移轨道段， 到达满足任务要求的月球借力位置，施加第二次脉冲，进入月球-Halo转移轨道， 最终在Halo轨道入轨点施加第三次脉冲，实现目标轨道捕获。

实施例2

如图1、图2所示，针对绕月四脉冲设计方案，一种月球借力约束的地月平 动点转移轨道设计方法，包括如下步骤：

步骤一、航天器施加第一次脉冲，由目标Halo轨道反推至月球借力位置， 月球借力点约束如图4所示；

选择法向幅值Az＝5000km的Halo轨道为目标轨道，期望的月球借力约束 条件满足：近月高度飞行航迹角借力方位角选取 Halo轨道入轨点状态为[1.12135,-0.05816,-0.00740,-0.01368,0.12341,-0.00535]^T，利用 公式(1)对Halo轨道入轨点状态量进行积分，在满足航迹角约束的位置中， 距离月球最近点的高度为h_m＝2.18664×10⁴km，借力方位角δ_m＝7.54890°，不满足期 望的借力不满足期望的借力高度与方位角约束，采用相关的优化算法对入轨点 状态量进行修正调整，使得借力位置的状态量满足公式(3)的约束矢量，进而 确定航天器进入Halo轨道所需的机动速度增量大小为15.61m/s。

步骤二、航天器在月球借力位置施加第二次脉冲，形成绕月飞行轨迹，在 施加第三次脉冲，使航天器由月球借力位置进入地-月转移轨道段，最终反推至 地球停泊轨道。

(a)选择绕月飞行轨迹末端满足的约束条件为：近月高度飞行 航迹角利用公式(1)对步骤一中确定的月球借力点状态量 [0.98420,-0.00307,-0.00032,1.42238,-1.69059,-0.03850]^T进行积分，并结合相关的优化 算法对状态量进行修正调整，直至满足高度与航迹角约束条件，进而确定航天 器绕月飞行所需的速度增量大小为51.38m/s，飞行时间1.33天。(b)选择期望 的地球停泊轨道约束条件：停泊轨道高度航迹角利用公式(1) 对(a)中确定的绕月飞行轨迹末端的状态量进行积分，状态量为 [0.98485,-0.00467,-0.00038,1.67400,-1.07063,-0.04063]^T。通过相关优化算法对状态量 进行修正调整，直至满足给定的地球停泊轨道约束条件。最终确定航天器进入 地-月转移轨道段所需的速度增量大小为250.92m/s，地球逃逸点为3118.35m/s 与总飞行时间为14.72天。

当月球借力位置仅考虑高度约束与航迹角约束时，总速度增量大小为 3490.896m/s，总飞行时间为14.52天，，比较以上两种情况，在飞行时间相近条 件下，考虑方位角约束后，总速度增量降低，两者相差54.636m/s。

针对以上两种情况设计了绕月飞行四脉冲转移轨道，如图7所示，其中虚 线为仅考虑借力高度与航迹角约束的转移轨道，实线为添加借力方位角后的转 移轨道。

实施例1与实施例2所述公式(1)为在地球-月球质心旋转系下建立航天器 的动力学方程。

如图3所示，其中坐标系的原点为系统的质心，x轴由地球指向月球，z轴 与地-月系统的角速度方向一致，y轴与x轴，z轴垂直构成右手坐标系。

在质心旋转系下，航天器的动力学模型表示为

$\left(\begin{array}{c}\stackrel{··}{x}-2\stackrel{·}{y}=x-\frac{(1-\mu )(x+\mu )}{r_e^3}-\frac{\mu (x-1+\mu )}{r_m^3}\\ \stackrel{··}{y}+2\stackrel{·}{x}=y-\frac{(1-\mu )y}{r_e^3}-\frac{\mu y}{r_m^3}\\ \stackrel{··}{z}=-\frac{(1-\mu )z}{r_e^3}-\frac{\mu z}{r_m^3}\end{array}\right)---\left(5\right)$

其中表征航天器在质心旋转系中的位置、速度与加速度状 态量，μ表示系统的质量系数，地-月系统中，取值为航天器与地球之间的距离，为航天器与月球间的距离。

所述月球借力位置的约束条件为：在转移轨道设计过程中，当月球借力位 置考虑不同的约束条件时，将对整个轨道转移过程的速度增量与飞行时间等参 数产生较大的影响。为了达到降低燃耗的目的，考虑月球借力位置的约束条件 分别为：近月高度h_m，相对月球飞行航迹角γ_m与借力方位角δ_m，可表示为：

$\left(\begin{array}{c}{h}_m=\sqrt{{\left(x_{fm}-1+\mu \right)}^2+y_{fm}^2+2_{fm}^2}-R_m\\ {\gamma }_m=\arcsin \left(\frac{\left(x_{fm}-1+\mu \right){\stackrel{·}{x}}_{fm}+y_{fm}{\stackrel{·}{y}}_{fm}+z_{fm}{\stackrel{·}{z}}_{fm}}{\left|\right|P_{fm}\left|\right|\left|\right|V_{fm}\left|\right|}\right)\\ {\delta }_m=\arccos \left(\frac{\left(x_{fm}-1+\mu \right){\stackrel{·}{y}}_{fm}-y_{fm}{\stackrel{·}{x}}_{fm}}{\left|\right|Q_{fm}\left|\right|}\right)\end{array}\right)---\left(6\right)$

其中，r_fm＝[x_fm,y_fm,z_fm]^T与为旋转系中借力点位置与速度矢 量，P_fm＝[x_fm-1+μ,y_fm,z_fm]^T为航天器相对于月球质心的位置矢量，R_m为月球半径， 矢量$Q_{fm}=\left(\begin{array}{c}{y}_{fm}{\stackrel{·}{z}}_{fm}-z_{fm}{\stackrel{·}{y}}_{fm}\\ z_{fm}{\stackrel{·}{x}}_{fm}-(x_{fm}-1+\mu ){\stackrel{·}{z}}_{fm}\\ (x_{fm}-1+\mu ){\stackrel{·}{y}}_{fm}-y_{fm}{\stackrel{·}{x}}_{fm}\end{array}\right).$

同时，为了充分利用目标Halo轨道的特性，采用逆向积分的策略设计地月 平动点转移轨道。在保证高度h_m与航迹角γ_m约束的条件下，改变借力方位角δ_m的 取值大小，能够得到不同约束条件时燃料消耗变化情况，进而确定设计方案的 各项参数。

所述采用逆向积分的策略，即为将步骤一中的公式(1)逆向积分；

步骤一、步骤二中的优化算法为微分修正算法，其作为一种局部优化算法， 具有快速收敛特性，基于此算法对由公式(1)积分确定的轨迹进行调整，可以 保证轨迹末端的状态量满足给定的约束条件。首先，选取自由变量C与约束矢量 F(C)，满足

$\left(\begin{array}{c}C={[V_x,V_y,V_z,T]}^T\\ F\left(C\right)={\left(\begin{array}{c}{h}_m-h_m^{*}\\ {\mathrm{sin\gamma }}_m-{\mathrm{sin\gamma }}_m^{*}\\ {\mathrm{cos\delta }}_m-{\mathrm{cos\delta }}_m^{*}\end{array}\right)}_{moon}\end{array}\right)---\left(7\right)$

其中，[V_x,V_y,V_z]^T为Halo轨道入轨点三轴方向速度矢量，T为Halo轨道至月 球借力位置的飞行时间；(h_m,γ_m,δ_m)由公式(2)确定，为给定的期望约 束值。

然后，推导约束矢量关于控制变量的偏导数关系

$\delta F=\frac{\partial F}{\partial V_x}{\mathrm{\delta V}}_x+\frac{\partial F}{\partial V_y}{\mathrm{\delta V}}_y+\frac{\partial F}{\partial V_z}{\mathrm{\delta V}}_z+\frac{\partial F}{\partial T}\delta T---\left(8\right)$

其中，符号δ表示变分关系，为约束矢量关于速度变量V_x的偏导数。最 后利用最小二乘法对自由变量C进行微分修正，迭代求解能够保证最终的转移轨 道满足给定的约束条件。

尽管已经给出并描述了本发明的实施例，本领域的普通技术人员都能够理 解，因此，在不脱离本发明设计原理与方法的情况下，可以对给定的实施例进 行修改，参数替换。