低成本多小天体并行交会技术

王鹏，武小宇，张立华

（1.航天东方红卫星有限公司，北京 100094；2.北京理工大学 宇航学院，北京 100081；3.飞行器动力学与控制教育部重点实验室，北京 100081）

引 言

通过对小天体探测可以了解太阳系演化、生命起源等基本规律，并有可能在未来实现小天体资源利用。随着航天技术的发展，小天体探测逐渐成为人类深空探测任务的热点目标之一[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轴构成右手坐标系。质量比 µ为行星质量与M∗之 比。则在旋转坐标系oxyz中，探测器的动力学方程可以描述为

其中：Ω 为探测器的伪势能函数，具有如下表达式

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

圆型限制性三体模型中存在5个动平衡点，如图1所示。其中 L1、L2和L3位于太阳与行星的连线上，称为共线平衡点。L4和L5分别与太阳和行星构成等边三角形，称为三角平衡点。平衡点及其附近周期和拟周期轨道存在与之相连的稳定和不稳定流形。不变流形的渐近特性为设计低能量转移轨道提供了可能，探测器可以沿着稳定或不稳定流形逐渐地趋近或者脱离（拟）周期轨道。

1.2 三体问题下的流形结构

为获得halo轨道的稳定和不稳定流形，首先在halo轨道上选取一参考点X0=X(t0)（本文选取halo轨道正向穿越xy平面的点为参考点），计算halo轨道对应的单值矩阵M(τ=1,t0)，τ∈[0,1]为halo轨道周期归一化后的时间参数，则τ对应的时间为

求取单值矩阵M(τ=1,t0)的特征值及对应的特征向量，得到稳定特征向量ηS(X0)和不稳定特征向量ηU(X0)。在[0,1]区 间内对τ进行等分离散化，则各离散点对应的稳定和不稳定特征向量可以用下式计算得到

其中：M(τ,t0)为 τ时刻相对参考点的状态转移矩阵。

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

其中：ε 为微小的摄动，一般选取ε 使位置摄动处于20～50 km范围内。

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

1.3 小天体探测交会模型

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

图2 太阳–地球系统halo轨道不变流形结构Fig.2 The invariant manifolds of halo orbits in the Sun-Earth system

图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×10−6～3.0×10−1km/s。扰动流形的速度初值可以写成

其中：∆V为 施加的速度摄动，χS和 χU分别为稳定和不稳定特征向量的速度分量。

图4给出了太阳—地球系统 L1点 附近z方向幅值为1.6×105km北向halo轨道的扰动流形结构，速度摄动大小为 ∆V=100m/s。 由图中可以看出，L1点附近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−µ，计算L1点附近halo轨道右侧不稳定流形第一次到达该截面的时间，如图5所示，图中横轴为一个halo轨道周期内halo轨道上不同流形分支对应的出发时间，纵轴为流形从halo轨道到庞加莱截面的飞行时间。

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

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

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

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

其中：t0为探测器从初始halo轨道出发的时间；τD表征了出发时探测器在halo轨道上的位置；∆VD为出发不稳定扰动流形的速度摄动大小；∆vD为在出发行星近拱点处施加的轨道机动脉冲；tS为由出发行星近拱点到目标行星近拱点的飞行时间；∆vA为在目标行星近拱点处施加的轨道机动脉冲； ∆VA为目标稳定扰动流形的速度摄动大小； τA表征了探测器进入目标halo轨道时在halo轨道上的位置。

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

其中：tD为出发halo轨道到行星近拱点的飞行时间；tD为目标行星近拱点到halo轨道的飞行时间；β>0为罚因子。

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

其中：ρ为探测器的轨道状态；XD和XA分别为出发和目标halo轨道上的状态；tf=t0+tD+tS+tA。

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

其中：ξD和ξA分别为探测器相对于出发和目标行星的位置矢量；tPD=t0+tD；tPA=tPD+tS。

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

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

1）在日心二体模型中，对于给定的出发行星逃逸时间段和日心转移时间段，利用二体Lambert算法对日心转移轨道进行遍历搜索，求得各遍历点对应的探测器逃逸出发行星时的时间tD、逃逸双曲线超速、日心转移段飞行时间tS以及到达目标行星时的双曲线超速。

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

3）在以行星为中心的二体模型中，分别完成不变流形近拱点轨道状态与步骤1）各遍历点对应和的匹配，计算得到在行星近拱点需要施加的速度脉冲∆vD和∆vA，确定对应最小的日心转移轨道参数tD、、tS和。

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

5）在以行星为中心的二体模型中，分别完成扰动流形近拱点轨道状态与步骤3）计算的和的匹配，并根据性能指标（9）确定最佳的机动脉冲∆vD和∆vA，进一步确定halo出发和入轨段轨道参数t0、τD、∆VD、∆VA和τ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

选取日地 L1点z方 向幅值为1 .6×105km的北向halo轨道作为出发轨道，微纳飞行器逃逸地球影响球的时间区间为2018年1月1日至2018年12月31日，日心转移段最大飞行时间为500天。速度摄动 ∆V搜索范围为10−6～300m/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 km/s以下的转移机会，但飞行时间较长，约为250天。利用该等高线图，可以容易得到总速度增量最小解，约为 2 .72km/s。与之对应的地球逃逸时间为2018年6月7日，日心转移段飞行时间263天，总速度增量为0.40 km/s。根据该转移机会，采用轨道匹配算法对扰动流形参数 τ和 ∆V进行搜索，分别绘制了在地球和小行星附近需要施加的速度脉冲和探测器在扰动流形上飞行时间的能量等高线图，如图6～7所示。

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

由图7和图8可以看出，对于日地halo轨道出发段，在近拱点需要施加的速度增量存在多个局部极小和极大值，而由出发halo轨道到地球近拱点的飞行时间则随着摄动速度∆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为扰动流形解，tD和tA分别为探测器在出发扰动流形和捕获扰动流形上的飞行时间。可以看出，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轨道的方式，提出通过在线选取探测目标，逐个释放微纳飞行器，与目标小行星进行交会探测的方式，有效解耦多小天体探测任务、降低任务成本并且达到探测目标的多样化，实现一次发射，多目标探测的“并行”探测方式；②不变流形和行星借力是实现低能量深空探测任务的重要天体力学机理，本文基于不变流形的演化规律将两种机理结合起来，构建了大幅度降低燃料消耗的转移轨道方案；③针对传统不变流形转移时间长的问题，引入扰动流形思想，通过施加较大的速度摄动减小周期轨道的稳定系数，缩短了探测器在流形上的飞行时间，并增加了解空间的多样性；④提出了一种多目标小行星转移轨迹设计搜索方法，可对星历模型下结合扰动流形与行星借力机制的转移轨道进行大范围全局搜索，为探测任务设计提供可行的快速低能量转移方案初始解集。相比传统方案，本文获得的转移轨道在燃料消耗和飞行时间上都具有较大的改进，可以为未来的小天体商业勘查提供有益参考。