文章快速检索    
  深空探测学报  2019, Vol. 6 Issue (1): 73-81  DOI: 10.15982/j.issn.2095-7777.2019.01.011
0

引用本文 

王鹏, 武小宇, 张立华. 低成本多小天体并行交会技术[J]. 深空探测学报, 2019, 6(1): 73-81. DOI: 10.15982/j.issn.2095-7777.2019.01.011.
WANG P, WU X Y, ZHANG L H. Low-Cost and Multi-Objective Asteroid Transfer Trajectory Design[J]. Journal of Deep Space Exploration, 2019, 6(1): 73-81. DOI: 10.15982/j.issn.2095-7777.2019.01.011.

基金项目

国家自然科学基金资助项目(11773004)

作者简介

王鹏(1974– ),男,硕士,研究员,主要研究方向:深空探测任务设计、深空小型探测器设计等。电话:13911574704 E-mail:13911574704@139.com;
武小宇(1991– ),男,博士,主要研究方向:航天器轨道动力学与控制、轨道设计与优化等。电话:18810684138 E-mail:987002661@qq.com;
张立华(1970– ),男,博士,研究员,主要研究方向:小型深空探测器总体设计、引力波空间探测任务分析与设计等。电话:13910314318 E-mail:zlh70717@sina.com

文章历史

收稿日期:2017-12-28
修回日期:2018-04-25
低成本多小天体并行交会技术
王鹏1,2,3, 武小宇2,3, 张立华1    
1. 航天东方红卫星有限公司,北京 100094;
2. 北京理工大学 宇航学院,北京 100081;
3. 飞行器动力学与控制教育部重点实验室,北京 100081
摘要: 小天体资源开发方兴未艾,为降低开发风险和成本,需要发射无人探测器交会观测多个待选目标小天体。传统多目标探测方案存在无法多次交会、成本高,周期长等不足。提出的低成本多小天体并行交会技术,能够对多小天体探测任务进行解耦,实现基于微纳飞行器的低成本多目标并行交会勘查,从而降低探测成本,缩短探测周期。该技术结合行星借力与不变流形机制构建了低能量星际转移方案。然后引入扰动流形思想,使微纳飞行器能够实现与目标小天体的快速交会。进一步,提出了一种多目标小天体探测全局搜索方法,该方法基于在日地halo轨道上停泊的微纳飞行器集群,逐次确定小天体探测目标,并利用上述方法完成了多微纳飞行器与多小天体的交会。数值仿真结果表明,该方案能够大幅度降低转移过程的燃料消耗,并缩短转移时间。
关键词: 小天体探测    多目标交会    不变流行    引力甩摆    
Low-Cost and Multi-Objective Asteroid Transfer Trajectory Design
WANG Peng1,2,3, WU Xiaoyu2,3, ZHANG Lihua1     
1. DFH Satellite Co.,Ltd.,Beijing 100094,China;
2. School of Aerospace Engineering,Beijing Institute of Technology,Beijing 100081,China;
3. Key Laboratory of Dynamics and Control of Flight Vehicle,Ministry of Education,Beijing100081,China
Abstract: To reduce the total cost and risk of asteroid exploration, it is necessary to select multi-objective to avoid the disadvantages such as single target, high-cost and long-period. In this paper, a multi-target transfer trajectory design method is proposed based on the cubesats with the cooperation of gravity flyby and invariant manifolds. Introduced by the perturbation invariant manifolds, the cubesat can rapidly transfer to the target and carry out a global launch chance searching for asgteroids based on the Sun-Earth Halo parking orbit and determine the next target. The results indicated that the proposed method is able to significantly reduce the fuel consumption of trajectory transfer as well as transfer time.
Key words: asteroid exploration    multi-objective transfer    invariant manifold    gravity flyby    

Hight lights:

  1. ● A multi-objective transfer trajectory strategy is proposed which significantly reduces the mission cost.
  2. ● The low-energy transfer trajectory is constructed based on the mechanism of invariant manifolds and gravity flyby to construct.
  3. ● The perturbation invariant manifolds is introduced in this paper to reduce the transfer time by adding the velocity increment to the periodic orbit and their manifolds.
引 言

通过对小天体探测可以了解太阳系演化、生命起源等基本规律,并有可能在未来实现小天体资源利用。随着航天技术的发展,小天体探测逐渐成为人类深空探测任务的热点目标之一[1-2]。1997年, 第一个专门探测小行星的太空计划NEAR[3]成功交会了433 Eros 小行星[3-4],并在途中飞越了253 Mathilde小行星。“深空1号”探测器于1999年造访了9969 Braille 小行星[5]。2005年,日本“隼鸟号”完成了“丝川”小行星的采样任务并返回地球,使人类第一次从小行星上获取了样本[6]。美国“黎明号”探测器于2007年飞往 Vesta和 Ceres 小行星并进行探测[7]。美国国家航空航天局(NASA)宣布了近地小行星采样返回任务,计划于2022年发射探测器,2028年到达101955 1999RQ36 小行星采集土壤并于2032年送回地球[8]

迄今为止,通过天文观测和探测任务,人类已经基本确定了太阳系小天体的化学组成和内部结构,小天体中蕴藏有丰富资源已成为共识[7]。2015年美国通过的《2015外空资源探索和利用法》标志着小天体商业开发热潮的开始。小天体抵近资源勘查是满足小天体资源商业开发风险控制的可行途径,但传统抵近探测模式又很难满足商业开发对成本的要求。以欧洲最为成功的小天体抵近任务ROSETTA为例,花费了7亿欧元,10年时间,利用甩摆轨道技术[9]实现了一颗彗星的附着和2颗小天体的飞越[10]。其采用的一颗探测器串行探测多颗小天体的方式,使得多个小天体探测任务紧密耦合,任务复杂度高、风险大,并不完全适合商业勘查的要求。

随着微纳飞行器相关技术的不断进步,50 kg级小天体交会探测任务正成为现实。微纳飞行器能够在相对短的周期(12~18个月)实现从任务概念设计到任务实施,其尺寸、重量和功耗都相对较小,成本低廉。如果多个微纳小天体探测器能够以一箭多星的方式发射至日地拉格朗日点进行停泊与轨道保持,而后对多个小天体进行一对一交会探测,则能够有效缩短探测任务的时间并降低经济成本、提高探测效能。利用一颗CZ3A或CZ4C火箭可以实现6~10颗50 kg级的微纳飞行器小天体交会勘查任务,经济成本有可能降到30万人民币/千克,任务全寿命周期有可能缩短到5年以内,是实现“一探多”的一种可行途径。

以飞行时间为约束,降低微纳飞行器的星际转移成本是多小天体并行探测的关键技术之一。本文在圆型限制性三体模型下构建了结合行星借力与不变流形机制的低能量转移方案。引入了扰动流形思想,实现与目标小天体的快速交会。进一步,提出了一种基于多小天体探测的全局搜索方法,该方法基于在日地轨道上停泊的微纳飞行器集群,在线确定探测目标,逐次离轨进行探测,从而达到发射一次、探测多个目标的“并行”探测方式。本文数值仿真采用不同类型的近地小天体(2009BD,2008EA9,99942 Apophis)作为示例进行数值仿真分析,结果表明,该方案能够大幅度降低转移任务的燃料消耗,并缩短转移时间。

1 动力学原理 1.1 动力学模型

在近地小天体探测过程中,微纳飞行器主要受到太阳、出发行星两个天体的引力作用,为便于分析,本文采用圆型限制性三体模型对探测器的轨道动力学进行建模,忽略探测器对太阳和行星运动的影响,并且假设太阳和行星绕两者的共同质心作圆周运动,具体如图 1所示。

图 1 圆型限制性三体问题 Fig. 1 The sketch of circular restricted three-body problem

定义旋转坐标系 $oxyz$ :原点为太阳与地球的共同质心, $x$ 轴为原点指向地球方向, $y$ 轴在太阳—行星运动平面内指向地球运动方向, $z$ 轴构成右手坐标系。质量比 $\mu $ 为行星质量与 ${M^ * }$ 之比。则在旋转坐标系 $oxyz$ 中,探测器的动力学方程可以描述为

$ \left\{ \begin{split} \ddot x - 2\dot y = \dfrac{{\partial \varOmega }}{{\partial x}}\\ \ddot y + 2\dot x = \dfrac{{\partial \varOmega }}{{\partial y}}\\ \!\!\!\! \ddot z = \dfrac{{\partial \varOmega }}{{\partial z}} \end{split} \right. $

其中: $\varOmega $ 为探测器的伪势能函数,具有如下表达式

$\varOmega = \displaystyle\frac{1}{2}\left( {{x^2} + {y^2} + {z^2}} \right) + \displaystyle\frac{{\left( {1 - \mu } \right)}}{{{r_1}}} + \frac{\mu }{{{r_2}}}$ (2)

其中: ${r_1} = \sqrt {{{\left( {x + \mu } \right)}^2} + {y^2} + {z^2}} $ 为探测器与太阳的距离, ${r_2} = \sqrt {{{\left( {x - 1 + \mu } \right)}^2} + {y^2} + {z^2}} $ 为探测器与行星的距离。

在圆型限制性三体模型中,探测器的运动只存在一个独立积分,即雅可比积分

$C = 2\varOmega - \left( {{{\dot x}^2} + {{\dot y}^2} + {{\dot z}^2}} \right)$ (3)

圆型限制性三体模型中存在5个动平衡点,如图 1所示。其中 ${{\rm L}_1}$ ${{\rm L}_2}$ ${{\rm L}_3}$ 位于太阳与行星的连线上,称为共线平衡点。 ${{\rm L}_4}$ ${{\rm L}_5}$ 分别与太阳和行星构成等边三角形,称为三角平衡点。平衡点及其附近周期和拟周期轨道存在与之相连的稳定和不稳定流形。不变流形的渐近特性为设计低能量转移轨道提供了可能,探测器可以沿着稳定或不稳定流形逐渐地趋近或者脱离(拟)周期轨道。

1.2 三体问题下的流形结构

为获得halo轨道的稳定和不稳定流形,首先在halo轨道上选取一参考点 ${{X}_0}{ = X}\left( {{t_0}} \right)$ (本文选取halo轨道正向穿越 $xy$ 平面的点为参考点),计算halo轨道对应的单值矩阵 ${M}\left( {\tau = 1,{t_0}} \right)$ $\tau \in \left[ {0,1} \right]$ 为halo轨道周期归一化后的时间参数,则 $\tau $ 对应的时间为

${t_\tau } = {t_0} + \tau {T_{_{\rm H}}}$ (4)

求取单值矩阵 ${M}\left( {\tau = 1,{t_0}} \right)$ 的特征值及对应的特征向量,得到稳定特征向量 ${{\eta }^{{\rm {\tiny{S}}}}}\left( {{{X}_0}} \right)$ 和不稳定特征向量 ${{\eta }^{\rm U}}\left( {{{X}_0}} \right)$ 。在 $\left[ {0,1} \right]$ 区间内对 $\tau $ 进行等分离散化,则各离散点对应的稳定和不稳定特征向量可以用下式计算得到

$ \left\{ \begin{array}{l} {{ \eta} ^{\rm{S}}}\left( {X\left( \tau \right)} \right) = { M}\left( {\tau ,{t_0}} \right){{ \eta} ^{\rm{S}}}\left( {{X_0}} \right)\\ {{ \eta} ^{\rm{U}}}\left( {X\left( \tau \right)} \right) = { M}\left( {\tau ,{t_0}} \right){{ \eta} ^{\rm{U}}}\left( {{X_0}} \right) \end{array} \right. $ (5)

其中: ${M}\left( {\tau ,{t_0}} \right)$ $\tau $ 时刻相对参考点的状态转移矩阵。

进一步,采用摄动法可以求得各离散点对应的稳定和不稳定不变流形初值如下

$\left\{ \begin{array}{l} {X^{{\rm{S}} \pm }}\left( {X\left( \tau \right)} \right) = X\left( \tau \right) \pm \varepsilon \dfrac{{{\eta ^{\rm{S}}}\left( {X\left( \tau \right)} \right)}}{{\left\| {{\eta ^{\rm{S}}}\left( {X\left( \tau \right)} \right)} \right\|}}\\ {X^{{\rm{U}} \pm }}\left( {X\left( \tau \right)} \right) = X\left( \tau \right) \pm \varepsilon \dfrac{{{\eta ^{\rm{U}}}\left( {X\left( \tau \right)} \right)}}{{\left\| {{\eta ^{\rm{U}}}\left( {X\left( \tau \right)} \right)} \right\|}} \end{array} \right.$ (6)

其中: $\varepsilon $ 为微小的摄动,一般选取 $\varepsilon $ 使位置摄动处于20~50 km范围内。

以式(6)为初值,分别通过逆向和正向递推轨道可以得到halo轨道的稳定和不稳定流形。图 2给出了太阳—地球系统中 ${{\rm L}_1}$ ${{\rm L}_2}$ 点附近的halo轨道及与其相连的不变流形。由图中可以看出,halo轨道周围存在4个不同的不变流形分支。在局部区域内,halo轨道的左端和右端各存在两个分支,一支为稳定流形,另外一支为不稳定流形。

图 2 太阳–地球系统halo轨道不变流形结构 Fig. 2 The invariant manifolds of halo orbits in the Sun-Earth system
1.3 小天体探测交会模型

微纳飞行器从地球停泊轨道出发转移至目标小天体的过程(如图 3所示)可以描述为:地球停泊微纳飞行器集群中的一飞行器施加机动 $\Delta {{v}_{\rm D}}$ 到达日地/地月halo轨道,进而继续施加速度机动 $\Delta {{v}_{\rm M}}$ 到达停泊halo轨道(只针对地月halo轨道出发情况)。当探测器进入目标影响范围时,在目标小天体附近施加机动 $\Delta {{v}_{\rm A}}$ 使探测器进入目标附近,实现目标观测轨道的入轨。从微纳飞行器的飞行历程来看,整条转移轨道可以分为halo轨道出发、出发行星逃逸、星际转移、目标行星捕获和观测轨道入轨5个阶段。

图 3 微纳飞行器从地球停泊轨道出发转移至目标小天体的过程 Fig. 3 The process of a cubesat transfering from the Earth parking orbit to the target asteroid
2 多小天体并行探测关键技术 2.1 扰动流形技术

设计小天体行星际探测任务转移轨道时,燃料消耗和转移时间是通常考虑的两个标准。合理利用halo轨道的不变流形可以降低转移过程中的燃料消耗,同时不变流形引入也必将增加探测器的转移时间。例如探测月球的低能量转移轨道,其本质是利用太阳的长时间引力摄动增大轨道的近地点直至月球轨道。影响探测器在不变流形上飞行时间的是halo轨道的稳定系数,若减小稳定系数则可使探测器能够较快地脱离halo轨道。计算传统不变流形时本质上就是通过施加微小摄动改变其稳定系数。如果适当地增大摄动,进一步减小其稳定系数,减少流形在halo轨道附近的运动时间,则既可以利用不变流形轨道的低能量特性,又能够使探测器快速脱离或抵达halo轨道。为此,本文引入扰动流形进行探测器转移轨道的设计。本文定义的扰动流形具有以下3个特点:①仅对速度矢量施加摄动;②速度摄动方向仍沿单值矩阵稳定和不稳定特征矢量方向;③速度摄动范围为 $1.0 \times {10^{ - 6}}$ $ 3.0 \times {10^{ - 1}}\,{\rm{km/s}}$ 。扰动流形的速度初值可以写成

$\begin{gathered} {V^{{\rm S} \pm }} = {V_0} \pm \Delta {V_\chi }^{\rm S} \\ {V^{{\rm U} \pm }} = {V_0} \pm \Delta {V_\chi }^{\rm U} \\ \end{gathered} $ (7)

其中: $\Delta V$ 为施加的速度摄动, ${{\chi }^{\rm S}}$ ${{\chi }^{\rm U}}$ 分别为稳定和不稳定特征向量的速度分量。

图 4给出了太阳—地球系统 ${{\rm L}_1}$ 点附近 $z$ 方向幅值为 $1.6 \times {10^5}\;{\rm{km}}$ 北向halo轨道的扰动流形结构,速度摄动大小为 $\Delta V = 100\;{\rm{m/s}}$ 。由图中可以看出, ${{\rm L}_1}$ 点附近halo轨道的扰动流形结构仍然由4个分支构成,并且扰动流形与传统不变流形的演化趋势基本一致。

图 4 日地系统L1点Halo轨道扰动流形结构 Fig. 4 The perturbation invariant manifolds of the Halo orbits near the L1 point in the Sun-Earth system

与传统不变流形不同的是,扰动流形是在较大的速度摄动下产生的,其可以更快地脱离halo轨道区域,从而节省探测器的飞行时间。为考察扰动流形对飞行时间的影响,设定如下庞加莱截面 $x = 1 - \mu $ ,计算 ${{\rm L}_1}$ 点附近halo轨道右侧不稳定流形第一次到达该截面的时间,如图 5所示,图中横轴为一个halo轨道周期内halo轨道上不同流形分支对应的出发时间,纵轴为流形从halo轨道到庞加莱截面的飞行时间。

图 5 扰动流形对飞行时间的影响 Fig. 5 The influence of the perturbation invariant manifolds on the fly time

图 5可以看到,对于 ${{\rm L}_1}$ 点附近 $z$ 方向幅值为 $1.6 \times {10^5}\;{\rm{km}}$ 北向halo轨道,其右端的传统不稳定流形到达庞加莱截面的飞行时间为200天左右,而扰动流形到达庞加莱截面的时间则大大减少,例如 $\Delta V = 50\;{\rm{m/s}}$ 扰动流形最小飞行时间为73.79天, $\Delta V = 300\;{\rm{m/s}}$ 扰动流形最小飞行时间仅为37.98天。由此可见,扰动流形能够以较小的燃料消耗代价大大缩短探测器在流形上的飞行时间。

2.2 停泊halo轨道期间小天体发射机会搜索技术

对给定的出发和目标halo轨道,利用行星借力和扰动流形的低能量转移轨道关键参数包括

$Z = \left[ {{t_0},{\tau _{\rm D}},\Delta {V_{\rm D}},\Delta {{v}_{\rm D}}^{\rm{T}},{t_s}, \Delta {v_{\rm A}}^{\rm{T}}, \Delta {V_{\rm A}},{\tau _{\rm A}}} \right]$ (8)

其中: ${t_0}$ 为探测器从初始halo轨道出发的时间; ${\tau _{\rm D}}$ 表征了出发时探测器在halo轨道上的位置; $\Delta {V_{\rm D}}$ 为出发不稳定扰动流形的速度摄动大小; $\Delta {{v}_{\rm D}}$ 为在出发行星近拱点处施加的轨道机动脉冲; ${t_{\rm S}}$ 为由出发行星近拱点到目标行星近拱点的飞行时间; $\Delta {{v}_{\rm A}}$ 为在目标行星近拱点处施加的轨道机动脉冲; $\Delta {V_{\rm A}}$ 为目标稳定扰动流形的速度摄动大小; ${\tau _{\rm A}}$ 表征了探测器进入目标halo轨道时在halo轨道上的位置。

性能指标可综合考虑燃料消耗和转移时间,选取为

$J\left( Z \right) = \Delta {V_{\rm D}} + \left\| {\Delta {v_{\rm D}}} \right\| + \left\| {\Delta {v_{\rm A}}} \right\| + \Delta {V_{\rm A}} + \beta \left( {{t_{\rm D}} + {t_{\rm S}} + {t_{\rm A}}} \right)$ (9)

其中: ${t_{ D}}$ 为出发halo轨道到行星近拱点的飞行时间; ${t_{ D}}$ 为目标行星近拱点到halo轨道的飞行时间; $\beta > 0$ 为罚因子。

转移轨道必须满足出发和目标halo轨道的轨线约束。

${{\varPhi} _1}\left( Z \right) = \left[ {\begin{array}{*{20}{c}} {\rho \left( {{t_0}} \right)}&{{{X}_{\rm D}}\left( {{\tau _{\rm D}},{t_0}} \right)} \\ {\rho \left( {{t_f}} \right)}&{{{X}_{\rm A}}\left( {{\tau _{\rm A}},{t_f}} \right)} \end{array}} \right] = 0$ (10)

其中: ${\rho }$ 为探测器的轨道状态; ${{X}_{\rm D}}$ ${{X}_{\rm A}}$ 分别为出发和目标halo轨道上的状态; ${t_f} = {t_0} + {t_{\rm D}} + {t_{\rm S}} + {t_{\rm A}}$

另外,由于需要在行星近拱点施加轨道机动,转移轨道还应该满足如下内点约束

${{\varPhi} _2}\left( Z \right) = \left[ {\begin{array}{*{20}{c}} {{\xi _{\rm D}}^{\rm T}{{\left( {{t_{P{\rm D}}}} \right)}}}&{{{\dot \xi }_{\rm D}}\left( {{t_{P{\rm D}}}} \right)} \\ {{\xi _{\rm A}}^{\rm T}{{\left( {{t_{P{\rm A}}}} \right)}}}&{{{\dot \xi }_{\rm A}}\left( {{t_{P{\rm A}}}} \right)} \end{array}} \right] = 0$ (11)
${{\varPhi} _3}\left( Z \right) = \left[ {\begin{array}{*{20}{c}} {{{\dot \xi }_{\rm D}}^{\;\;\, \rm T}{{\left( {{t_{P{\rm D}}}} \right)}}{{\dot \xi }_{\rm D}}\left( {{t_{P{\rm D}}}} \right)}&{{\xi _{\rm D}}^{ \rm T}{{\left( {{t_{P{\rm D}}}} \right)}}{{\ddot \xi }_{\rm D}}\left( {{t_{P{\rm D}}}} \right)} \\ {{{\dot \xi }_{\rm A}}^{\;\;\, \rm T}{{\left( {{t_{P{\rm A}}}} \right)}}{{\dot \xi }_{\rm A}}\left( {{t_{P{\rm A}}}} \right)}&{{\xi _{\rm A}}^{ \rm T}{{\left( {{t_{P{\rm A}}}} \right)}}{{\ddot \xi }_{\rm A}}\left( {{t_{P{\rm A}}}} \right)} \end{array}} \right] \geqslant 0$ (12)

其中: ${{\xi }_D}$ ${{\xi }_A}$ 分别为探测器相对于出发和目标行星的位置矢量; ${t_{PD}} = {t_0} + {t_{ D}}$ ${t_{PA}} = {t_{PD}} + {t_S}$

上述模型构成了一个复杂的多点边值问题,理论上可采用非线性规划等技术直接求解。然而,由于此类型转移轨道途经多个引力场,动力学非线性强,轨道约束对设计参数异常敏感,导致数值优化算法很难收敛,提供合理的设计初值是解决该问题的关键。针对该问题,本文提出一种基于模型拼接的初始转移轨道构建方法。由前面的讨论可知,整条转移轨道由5个轨道段构成,5个分段的轨道特性各不相同。综合探测器考虑飞行时间和受力特点,可分别在不同动力学模型下对各分段轨道进行讨论。对于halo轨道出发和入轨段,在圆型限制性三体模型中讨论。对于出发行星逃逸和目标行星捕获段,在以行星为中心的二体模型中讨论,且忽略其飞行时间。对于日心转移段,则只考虑太阳的引力影响。需要指出的是,以上讨论均考虑星历约束下的行星相位。

根据以上假设,本文提出的初始转移轨道搜索方法可以描述为

1)在日心二体模型中,对于给定的出发行星逃逸时间段和日心转移时间段,利用二体Lambert算法对日心转移轨道进行遍历搜索,求得各遍历点对应的探测器逃逸出发行星时的时间 ${t_{\rm D}}$ 、逃逸双曲线超速 ${v}_{\rm D}^\infty $ 、日心转移段飞行时间 ${t_{\rm S}}$ 以及到达目标行星时的双曲线超速 ${v}_{\rm A}^\infty $

2)在圆型限制性三体模型中,分别计算出发halo轨道不稳定流形(传统不变流形)和目标halo轨道稳定流形对应的第一次行星近拱点庞加莱映射,并将其转化到以行星为原点的惯性坐标系中。

3)在以行星为中心的二体模型中,分别完成不变流形近拱点轨道状态与步骤1)各遍历点对应 ${v}_{\rm D}^\infty $ ${v}_{\rm A}^\infty $ 的匹配,计算得到在行星近拱点需要施加的速度脉冲 $\Delta {{v}_{\rm D}}$ $\Delta {{v}_{\rm A}}$ ,确定对应最小 $\Delta v = \left\| {\Delta {{v}_{\rm D}}} \right\| + \left\| {\Delta {{v}_{\rm A}}} \right\|$ 的日心转移轨道参数 ${t_{\rm D}}$ ${v}_{\rm D}^\infty $ ${t_{\rm S}}$ ${v}_{\rm A}^\infty $

4)在圆型限制性三体模型下,分别计算出发和目标halo轨道对应的不同速度摄动 $\Delta V$ 和不同出发位置 $\tau $ 扰动流形的第一次近拱点庞加莱映射,并将其转化到以行星为原点的惯性坐标系中。

5)在以行星为中心的二体模型中,分别完成扰动流形近拱点轨道状态与步骤3)计算的 ${v}_{\rm D}^\infty $ ${v}_{\rm A}^\infty $ 的匹配,并根据性能指标(9)确定最佳的机动脉冲 $\Delta {{v}_{\rm D}}$ $\Delta {{v}_{\rm A}}$ ,进一步确定halo出发和入轨段轨道参数 ${t_0}$ ${\tau _{\rm D}}$ $\Delta {V_{\rm D}}$ $\Delta {V_{\rm A}}$ ${\tau _{\rm A}}$

该搜索方法通过对完整转移轨道进行分段处理,分别在不同的动力学模型下对关键参数进行搜索,有效降低了转移轨道对关键参数的敏感度。该搜索算法本质上分为两个层次,第一层是确定日心转移轨道参数,第二层是确定其它4个轨道段的参数。通过第一层的搜索,出发行星和目标行星附近的轨道实现了解耦,可独立进行关键参数搜索,降低了搜索难度。

3 数值仿真分析

本文分别选取了具有较高探测价值的近地小行星2009BD(Aollo型,S类)、2008EA9(Amor,C类)以及99942 Apophis(Aten型,S类),上述目标同时考虑到探测的可行性与科学价值的多样性,可以有效地覆盖当今小天体探测所关注的热点问题,以2009BD为例,小行星基本轨道参数如表1

表 1 目标小行星轨道根数 Table 1 The orbital elements of the target asteroid

选取日地 ${{\rm L}_1}$ $z$ 方向幅值为 $1.6 \times {10^5}{\rm{km}}$ 的北向halo轨道作为出发轨道,微纳飞行器逃逸地球影响球的时间区间为2018年1月1日至2018年12月31日,日心转移段最大飞行时间为500天。速度摄动 $\Delta V$ 搜索范围为 ${10^{ - 6}}$ $ 300\;{\rm{m/s}}$ 。采用本文所提搜索方法,可以获得采用传统不变流形时在地球和小行星交会施加机动脉冲之和的等高线图(图6)。

图 6 地球–目标小行星(2009BD为例)halo轨道间转移近拱点速度增量等高线图 Fig. 6 The contour map of the periapsis velocity increment to the target 2009 BD asteroid based on the halo orbit

图 6可以看出,等高线图将解空间分成两个区域,并且上半区的总速度增量更小,可以寻找到总速度增量在 $3\;{\rm{km/s}}$ 以下的转移机会,但飞行时间较长,约为250天。利用该等高线图,可以容易得到总速度增量最小解,约为 $2.72\;{\rm{km/s}}$ 。与之对应的地球逃逸时间为2018年6月7日,日心转移段飞行时间263天,总速度增量为0.40 km/s。根据该转移机会,采用轨道匹配算法对扰动流形参数 $\tau $ $\Delta V$ 进行搜索,分别绘制了在地球和小行星附近需要施加的速度脉冲和探测器在扰动流形上飞行时间的能量等高线图,如图 6~7所示。

图 7 地球近拱点速度增量等高线图 Fig. 7 The contour map of the periapsis velocity increment near the Earth

图 7图 8可以看出,对于日地halo轨道出发段,在近拱点需要施加的速度增量存在多个局部极小和极大值,而由出发halo轨道到地球近拱点的飞行时间则随着摄动速度 $\Delta V$ 的增大而增大。选择合适的转移轨道需要综合考虑机动速度增量和飞行时间,对于图 7中的3个局部极小值,解1的飞行时间过长,解2具有最小的速度增量,飞行时间约为140天,解3虽然速度增量有所增加,但飞行时间则大大减少,约为74天。因此,可以选取解2和解3作为可行的出发halo轨道离轨机会。通过以上分析,可以获得可行的出发halo轨道离轨机会,如表 2所示。

图 8 出发halo轨道到目标小行星的飞行时间 Fig. 8 The flight time from the halo orbit to the target asteroid
表 2 出发halo轨道离轨机会 Table 2 The departure chances from halo orbits

表 2中IM为传统不变流形解,PM为扰动流形解, ${t_D}$ ${t_A}$ 分别为探测器在出发扰动流形和捕获扰动流形上的飞行时间。可以看出,EPM1不仅比EIM大大减少了微纳飞行器在流形上的飞行时间,并且在地球近拱点需要施加的速度增量更小,这主要是由于速度摄动改变了halo轨道流形的近拱点分布,增加了解的多样性。同样,MPM1和MPM2相比EIM也同时节省了飞行时间和速度增量,这充分证明了扰动流形在构建快速低能量转移轨道方面的优势。由表 2可以看到,由于引入了行星借力机制,使得转移轨道的燃料消耗大幅降低。另一方面,速度扰动有效地减少了探测器在流形上的飞行时间,获得的转移轨道兼具快速低能的特性。

通过上述获取的离轨机会参数,微纳飞行器转移至目标小行星(以2009BD为例)轨迹如图 9~10所示。

图 9 从日地halo轨道转移至目标小行星轨迹 Fig. 9 The trajectory from the halo orbit in the Sun-Earth system to the target asteroid
图 10 从地月halo轨道转移至目标小行星轨迹前段 Fig. 10 The trajectory from the Earth-Moon system to the target asteroid

通过对上述3颗目标小行星进行转移轨迹设计,可以发现总速度增量明显小于前人轨迹设计方案[11],如表 3所示。

表 3 转移至目标小行星总速度增量 Table 3 The total velocity increment during the transfer process to the target asteroids
4 结 论

本文首次提出了基于微纳飞行器的多小天体并行探测轨迹设计方案,所做的研究工作具有如下几个特点:①针对传统多小天体串行探测任务设计复杂、周期长,成本高的特点,本文利用搭载微纳飞行器的母舰集群停泊在日地/地月halo轨道的方式,提出通过在线选取探测目标,逐个释放微纳飞行器,与目标小行星进行交会探测的方式,有效解耦多小天体探测任务、降低任务成本并且达到探测目标的多样化,实现一次发射,多目标探测的“并行”探测方式;②不变流形和行星借力是实现低能量深空探测任务的重要天体力学机理,本文基于不变流形的演化规律将两种机理结合起来,构建了大幅度降低燃料消耗的转移轨道方案;③针对传统不变流形转移时间长的问题,引入扰动流形思想,通过施加较大的速度摄动减小周期轨道的稳定系数,缩短了探测器在流形上的飞行时间,并增加了解空间的多样性;④提出了一种多目标小行星转移轨迹设计搜索方法,可对星历模型下结合扰动流形与行星借力机制的转移轨道进行大范围全局搜索,为探测任务设计提供可行的快速低能量转移方案初始解集。相比传统方案,本文获得的转移轨道在燃料消耗和飞行时间上都具有较大的改进,可以为未来的小天体商业勘查提供有益参考。

参考文献
[1]
李俊峰,曾祥远. 不规则小行星引力场内的飞行动力学[J]. 力学进展,2017,47:429-451.
LI J F,ZENG X Y. Flight dynamics in the gravitational fields of irregular asteroids[J]. Advances In Mechanics,2017,47:429-451.   http://d.old.wanfangdata.com.cn/Periodical/lxjz201701012 (0)
[2]
李俊峰,曾祥远,张韵. 小行星的奇特动力学1)[J]. 力学与实践,2016,38(6):603-611.
LI J F,ZENG X Y,ZHANG Y. Unique dynamics of asteroids 1)[J]. Mechanics in Engineering,2016,38(6):603-611.   http://d.old.wanfangdata.com.cn/Periodical/lxysj201606001 (0)
[3]
CANALIAS E, GOMEZ G, MARCOTE M, et al. Assessment of mission design includingutilization of libration points and weak stability boudaries [R]. ESOC Contract18142/04/NL/MV, Final Report, 2004. (0)
[4]
FARQUHAR R,MUHONEN D,RICHARDSON D. Mission design for a halo orbiter of the Earth[J]. Journal of Spacecraft and Rockets,1977,14:170-177. DOI:10.2514/3.57176 (0)
[5]
张振江,崔祜涛,崔平远. 基于三阶解析解的小行星平衡点附近 halo 轨道确定方法研究[J]. 宇航学报,2011,32(2):277-283.
ZHANG Z J,CUI H T,CUI P Y. Research on third-order analytical solution determination of halo orbits around equilibrium point of asteroid[J]. Journal of Astronautics,2011,32(2):277-283. DOI:10.3873/j.issn.1000-1328.2011.02.007 (0)
[6]
FOLTA D, RICHON K. Libration orbit mission design at L2: A MAP and NGST perspective[C]//AIAA/AAS Astrodynamics Conference. Boston: AIAA/AAS, 1998. (0)
[7]
于洋,宝音贺西. 小天体附近的轨道动力学研究综述[J]. 深空探测学报,2014,1(2):93-104.
YU Y,BAOYIN H X. Review of orbital dynamics in the vicinity of solar system small celestial bodies scientific vision for future missions[J]. Journal of Deep Exploration,2014,1(2):93-104.   http://d.old.wanfangdata.com.cn/Periodical/sktcxb201402002 (0)
[8]
吴伟仁,崔平远,乔栋,等. 嫦娥二号日地拉格朗日 L2 点探测轨道设计与实施[J]. 中国科学:科学通报,2012,57(21):1987-1991.
WU W R,CUI P Y,QIAO D,et al. Design and performance of exploring trajectory to Sun-Earth L2 point for Chang’E-2 mission[J]. Chin Sci Bull,2012,57(21):1987-1991.   http://www.wanfangdata.com.cn/details/detail.do?_type=perio&id=kxtb201221005 (0)
[9]
曾祥远,李俊峰,刘向东. 不规则小天体引力场内的广义甩摆轨道[J]. 深空探测学报,2016,3(1):29-33.
ZENG X Y,LI J F,LIU X D. Generalized flyby trajectories over irregular-shaped small bodies[J]. Journal of Deep Space Exploration,2016,3(1):29-33.   http://d.old.wanfangdata.com.cn/Periodical/sktcxb201601004 (0)
[10]
姜宇,宝音贺西. 强不规则天体引力场中的动力学研究进展[J]. 深空探测学报,2014,1(4):250-261.
JIANG Y,BAOYIN H X. Research trend of dynamics in the gravitational field of irregular celestial body[J]. Journal of Deep Space Exploration,2014,1(4):250-261.   http://d.old.wanfangdata.com.cn/Periodical/sktcxb201404002 (0)
[11]
HOWELL K C,PERNICKA H J. Numerical determination of Lissajous orbits in the restrictedthree-body problem[J]. Celestial Mechanics,1988,41:107-124. (0)