-
潜在威胁的小行星与地球轨道在接近或交叉的情况下,存在着撞击地球的风险。历史上地球曾多次遭到小行星撞击,部分撞击对地球造成了严重破坏,6 500万年前的希克苏鲁伯撞击事件曾导致全球60% ~ 80%的物种灭绝[1],因此围绕潜在威胁的小行星开展行星防御对人类长远生存发展具有重要的意义。评估核爆、动能撞击等防御手段需要精确的小行星星历与小行星质量、密度、自转等信息[2]。受时空与天体运动的限制,现有的地基、天基观测手段所获取的小行星信息具有一定的局限性,对于庞大的小行星群体,通过实施小行星探测任务,人类仅对少数小行星的物理化学性质有过详细的探查[3]。
Vetrisano等[4]的研究表明,近距离抵近飞越是获取小行星特性信息的一种有效手段。飞越任务通常利用光学、光谱与红外远程成像系统对小行星进行探测,飞越探访可获取小行星表面的图像,测量小行星光谱,反演其物质组成、大小、体积、质量、密度和结构等来了解小行星的物理化学性质。针对有潜在撞地威胁的小行星,获取其精确的星历与质量、形状数据可准确计算撞击、偏转效应,对于制定有效的小行星防御策略具有重要的意义[5]。联合国空间任务规划咨询组(Space Mission Planning Advisory Group,SMPAG)也在2021年呼吁加强对小行星飞越探测任务的支持,以获取小行星详细的目标特性信息[6]。
“隼鸟号”(Hayabusa)、“黎明号”(Dawn)[7-8]探测器对特定小行星进行了探测。小行星飞越探测可在单次任务中飞越数量可观的小行星,开展多目标特性普查,通过合理规划访问顺序与优化轨道可显著地降低其任务成本。小行星探测任务也由单一目标向多目标发展,美国国家航空航天局(National Aeronautics Space and Administration,NASA)于2021年10月发射了“露西号”(Lucy)探测器,计划在12 a的任务周期内飞越探访木星特洛伊族的8颗小行星,将创下单个航天器访问小行星数目最多的记录[9]。日本宇宙航空研究开发机构(Japan Aerospace eXploration Agency,JAXA)开发的“命运﹢”(DESTINY+)探测器计划飞越多颗近地小天体并最终于2028年飞越双子座流星雨母体Phaethon 3200(TB 1983)[10]。欧洲航天局(European Space Agency,ESA)计划于2029年发射彗星拦截器(Comet Interceptor),探测器将停泊在日地L2点周期轨道以期择机飞越新发现的长周期彗星[11]。
既往小行星探测轨道设计多聚焦于交会探访任务,飞越探测仅作为整个任务的衍生。Gao等[12]研究了以地球特洛伊族小行星2010TK7为主要探测目标的交会轨道设计,在抵达2010TK7前探测器飞越其它2颗小行星;夏炎等[13]研究了主带小行星的可达性,提出了探测器最终伴飞小行星2004VZ60并在途中飞越4颗主带小行星的轨道设计;Qiao等[14]研究了对近地小行星进行采样返回并在途中对不同光谱类型小行星进行飞越探测的轨道设计。McNutt等[15]提出利用连续推力探测器对近地小行星进行近距离飞越探测的侦查兵(Near Earth Asteroid Scout,NEA Scout)计划。Greco等[16]提出利用立方星双星系统对选定的近地天体飞越进行探访的任务,认为在近地天体过黄道面时对其进行飞越探访可有效节约改变轨道倾角、偏心率所导致的燃料消耗。Giuseppe等[17]提出了针对Apollo型小行星的黄道面内多目标飞越序列规划方法,采取迭代扩大搜索区域的飞越序列拓展策略,一旦搜索区域中出现可行候选飞越目标就终止搜索,可导致解集的过度削减并遗漏最优解,需研究全局性更强、范围更广的多目标小行星飞越探测任务轨道优化设计方法。
多目标探测任务轨道优化的主要难点在于探测序列搜索。目前常用的路径搜索[18-19]有:分支–限界算法、智能优化算法与混合优化方法。分支–限界算法通常以树结构表示解序列,通过不同策略或约束条件对树结构剪枝,可在较低计算成本下得出较好的优化结果。贪婪算法可视为分支–限界算法的一种特殊情况,在树结构的每层均选择局部最优值的策略易使其陷入局部最优。近年来,波束搜索算法作为分支–限界算法的一种,其易于部署、计算复杂度可由波束宽度控制的特点,在国际轨道优化竞赛当中展示了良好的性能[20-21]。以蚁群算法为代表的智能优化算法需要设计者对信息素矩阵、启发式信息矩阵及算法参数设计有较深的认识,该方法在多目标交会任务当中应用广泛。混合优化方法,通常将多种路径搜索算法、人工筛选方法与结合专家知识的设计方法混合使用,部分还涉及到求解最优序列后进一步的局部优化[22],要求具有较强的专家经验。
不同于多目标交会任务,飞越探测任务在飞越时刻探测器与目标的速度并不相等,这使得两条不同飞越序列的相同子序列如
$ C\to A\to B\to D $ 与$ E\to A\to B\to D $ 当中,由$ A\to B $ 的转移速度增量不一致。因此,在交会类任务当中被广泛使用的多种轨道转移速度增量近似计算方法,Edelbaum等[23]估算方法无法在飞越类任务当中使用。转移速度增量的不恒定也会导致以启发式信息矩阵、信息素矩阵为群体智能算法进化依据的蚁群算法在多目标飞越问题上表现不佳。在交会任务中常用的分割轨迹[24]方法,将完整的飞越任务分为多个子飞越任务,会遇到无法确定子轨迹速度初始值的问题。本文本文的潜在威胁小行星数据来自喷气推进实验室(Jet Propulsion Laboratory,JPL)小行星数据库,候选目标为全体2 073个潜在威胁小行星。在Greco等[16]工作的基础上,将黄道面内飞越基本策略进一步应用于全体潜在威胁小行星,在分析潜在威胁小行星的轨道特点的情况下,设计了在黄道面内多目标飞越潜在威胁小行星的轨道优化模型;针对任务时间长、候选目标多导致的序贯飞越解空间巨大的问题,将黄道面内飞越策略与波束选择树搜索算法相结合,对多目标飞越任务的轨道进行了快速求解。研究了不同策略的性能表现,分析了该类任务发射时间窗口的普遍存在性;对比了波束选择树搜索算法与迭代扩大搜索区域算法[17]在相同任务场景下的性能,并进一步对比分析了蚁群算法在黄道面内飞越探测模型下的性能表现及原因;结合运载火箭发射能力、卫星燃料质量比等工程约束,对2024—2028年不同发射时间进行优化计算,得到探测器飞越数颗小行星的序列集合,就有代表性的序列进行结果分析与展示,可为未来飞越探测任务提供设计参考。
-
近地天体是轨道近日点在1.3 AU以内的太阳系小行星。包含小行星和彗星两类,近地小行星占绝大多数。与地球轨道的最近距离在0.05 AU以内且绝对星等小于22等(等效直径大于140 m)的近地小行星定义为潜在威胁小行星,其轨道可近距离接近地球,并且大到足以在发生撞击时造成严重的区域破坏。潜在威胁小行星倾角主要分布在0° ~ 30°,轨道半长轴决定其轨道周期分布在0.5~4 a,如图1所示。
图 1 潜在威胁小行星轨道半长轴、倾角分布
Figure 1. Distribution of semi-major axis and inclination angle of potentially hazardous asteroids' orbit
由于探测器改变轨道倾角的机动需要消耗的燃料远高于轨道面内机动,且无论何种轨道构型的潜在威胁小行星其运动轨迹总会经过黄道面,可得到一种确定性的多目标小行星飞越基本策略为探测器仅在固定轨道平面内进行轨道面内机动,调整自身轨道与穿过该面的小行星进行序贯飞越的基本策略。考虑到探测器需从地月系逃逸与潜在威胁小行星过黄道面的位置分布情况,黄道面是该固定轨道平面的理想选择,确定潜在威胁小行星在特定时间段穿越黄道面时刻与位置是该算法的基础。
-
计算潜在威胁小行星飞越黄道面的时刻和位置时,只考虑太阳中心引力体的引力作用,小行星轨道状态采用开普勒轨道根数表示
$ \left\{a,e,i,\varOmega ,\omega ,\nu \right\} $ ,分别为半长轴、偏心率、轨道倾角、升交点赤经、近日点幅角与真近点角。潜在威胁小行星过黄道面位置、时刻求解利用Giuseppe等[17]在飞越任务中使用的计算方法,借助近焦点坐标系向日心–黄道坐标系转换可以快速求取。$$ \left\{\begin{array}{l}x=r\left(\nu \right)[\left(\mathrm{cos}\varOmega \mathrm{cos}\omega -\mathrm{sin}\varOmega \mathrm{sin}\omega \mathrm{cos}i\right)\mathrm{cos}\nu +\\ \quad \left(-\mathrm{cos}\varOmega \mathrm{sin}\omega -\mathrm{sin}\varOmega \mathrm{cos}\omega \mathrm{cos}i\right)\mathrm{sin}\nu ]\\ y=r\left(\nu \right)[\left(\mathrm{sin}\varOmega \mathrm{cos}\omega +\mathrm{cos}\varOmega \mathrm{sin}\omega \mathrm{cos}i\right)\mathrm{cos}\nu +\\ \quad \left(-\mathrm{sin}\varOmega \mathrm{sin}\omega +\mathrm{cos}\varOmega \mathrm{cos}\omega \mathrm{cos}i\right)\mathrm{sin}\nu ]\\ z=r\left(\nu \right)\left[\mathrm{sin}\omega \mathrm{sin}i\mathrm{cos}\nu +\mathrm{cos}\omega \mathrm{sin}i\mathrm{sin}\nu \right]\end{array}\right. $$ (1) 潜在威胁小行星过黄道面位置计算见式(1),其中
$ r\left(\nu \right)=a(1-{e}^{2})/(1+e\mathrm{cos}\nu ) $ ,令$ z=0 $ ,有$ {\mathrm{tan}\nu }_{\mathrm{tr}}=-\mathrm{tan}\omega $ ,用$ {\nu }_{\mathrm{t}\mathrm{r}}^{1} $ ,$ {\nu }_{\mathrm{t}\mathrm{r}}^{2} $ 代表小行星两次穿过黄道面对应真近点角($ {\nu }_{\mathrm{t}\mathrm{r}}^{1} < {\nu }_{\mathrm{t}\mathrm{r}}^{2} $ )。将求得的
$ {\nu }_{\mathrm{t}\mathrm{r}}^{1} $ ,$ {\nu }_{\mathrm{t}\mathrm{r}}^{2} $ 带回式(1)可得到小行星从2个方向穿越黄道面位置$r_\mathrm{tr}^1$ ,$r_\mathrm{tr}^2$ 、速度$ {\nu }_{\mathrm{t}\mathrm{r}}^{1} $ ,$ {\nu }_{\mathrm{t}\mathrm{r}}^{2} $ ,由真近点角可求取潜在威胁小行星飞越黄道面的时刻${t}_{\mathrm{t}\mathrm{r}}^{\mathrm{I}}$ ,${t}_{\mathrm{t}\mathrm{r}}^{\mathrm{II}}$ 。全体潜在威胁小行星穿越黄道面位置如图2所示,其主要分布在距日心3 AU的距离内。图 2 潜在威胁小行星穿越黄道面位置分布
Figure 2. Location distribution of potentially hazardous asteroids crossing the ecliptic plane
为简化卫星转移弧段的计算,本文设定卫星在飞越上一颗小行星后立即变轨,至此在潜在威胁小速度增量可以利用兰伯特定理[25]解出探测器在出发时刻、位置
$Sat \langle t,r,v\rangle$ 到目标小行星过面时刻,位置$As{t}_{\mathrm next}\langle t,r,v\rangle$ 的转移速度增量$ \Delta v $ 来表示。但需要注意的是,由于重复计算兰伯算法消耗大量程序时间,在计算时进行合理的矩阵化与并行化可以显著提高程序运算效率。飞越黄道面位置时刻算法最终输出满足探测器最大机动速度增量约束的候选小行星即$$ \Delta v=\left|{v}_{\mathrm{i}\mathrm{n}}-{v}_{\text{af}}\right| < \Delta {v}_{\mathrm{max}} $$ (2) 其中:
$ {v}_{\mathrm{i}\mathrm{n}} $ 为机动前速度;$ {v}_{\text{af}} $ 为机动后速度。探测器在单次飞越完成后即开始对全体潜在威胁小行星进行单圈兰伯特转移计算,搜索下一次黄道面内飞越的候选小行星,除此之外不再对单次转移时间做其他约束。因为小行星过黄道面位置、时刻存在周期关系避免了重复计算,所以确定的黄道面内多目标小行星飞越的基本策略比传统依靠时间窗口确定兰伯特转移始末时刻的大规模计算大幅降低了复杂度。
-
确定黄道面内飞越小行星的基本原则与兰伯特转移弧段表示的单次飞越轨迹后,转化为如何在给定任务的约束条件下,飞越尽可能多的潜在威胁小行星。小行星飞越问题的可行解是一个典型的序贯结构,可利用决策树生成。树节点表示探测器从出发小行星经过兰伯特转移到飞越目标小行星,每个节点包括的信息:出发小行星编号
${A}_{\text {Ns}}$ 、目标小行星编号${A}_{\text {Nt}}$ ,出发时刻$ {t}_{\text{Ns}} $ 、飞越时刻$ {t}_{\text{Nf}} $ ,本弧段速度增量$ {v}_{\text{Narc}} $ ,剩余任务时间$ {t}_{\mathrm{Re}} $ ,剩余机动速度增量$ {v}_{\text{Re}} $ 。对于树结构的解探索,有基于宽度或深度优先搜索的全局遍历方法与基于波束生成的选择性方法。全局遍历方法在构造序列解受到多重约束的特定条件下,存在着在可接受的时间、空间复杂度下完成遍历的情况。本文将转移搜索区域调整为整个黄道面,搜索区域内的候选目标较多,且候选解随着序列长度的增加呈指数增长,难以在可接受的时间内完成对可行解的完全遍历,因此采用基于宽度优先的波束树搜索算法[21]作为解优化的方法。
波束选择树搜索算法的子树拓展策略对于算法整体表现有着决定性的作用。树节点值函数
$ r $ 作为子树拓展策略的直接参考,一般可由探测器飞越下一个子节点所需的速度增量$ {v}_{\text{Narc}} $ 或子节点飞越时刻$ {t}_{\text{Nf}} $ 定义。考虑到直接定义下的$ {v}_{\text{Narc}} $ 与$ {t}_{\text{Nf}} $ 在各弧段互相独立缺乏全局性,对以上两种树节点值函数做无量纲化[26]的定义:剩余速度增量比$ {r}_{v}={v}_{\mathrm{R}\mathrm{e}}/{v}_{\mathrm{a}\mathrm{l}\mathrm{l}} $ 与剩余任务时间比$ {r}_{t}={t}_{\mathrm{Re}}/{t}_{\mathrm{a}\mathrm{l}\mathrm{l}} $ ,其中${v}_{\mathit{\text {Re}}}$ 为剩余速度增量,$ {v}_{\mathrm{a}\mathrm{l}\mathrm{l}} $ 为总速度增量,$ {t}_{\mathrm{Re}} $ 为任务剩余时间,$ {t}_{\mathrm{a}\mathrm{l}\mathrm{l}} $ 为任务总时间。在只考虑速度增量或任务时间的拓展策略基础上,融合两种评价指标的拓展策略及其变体在波束树搜索处理空间序贯决策任务当中也有出色的效果[21]。波束选择树搜索算法的子树拓展策略值函数可有表1的不同定义方式,其中
$ {r}_{1} $ 、$ {r}_{2} $ 分别为仅考虑燃料的燃料最省策略与仅考虑剩余任务时间的时间最省策略,$ {r}_{3} $ 、$ {r}_{4} $ 、$ {r}_{5} $ 为将剩余速度增量与剩余任务时间综合考虑,分别求取平均值、加权值与最小值的拓展策略。表 1 节点拓展策略
Table 1. Node expansion strategy
策略名 策略定义 燃料最省
策略$ {r}_{1}=\text{max}\text{(}{r}_{v1},{r}_{v2},\cdots \cdots \text{,}{r}_{vN}\text{)} $ 时间最省
策略$ {r}_{2}=\text{max}\text{(}{r}_{t1},{r}_{t2},\cdots \cdots {,r}_{tN}\text{)} $ 燃料、时间平均策略 $ {r}_{3}=\text{max}{(r}_{\text{ave}1},{r}_{\text{ave}2},\cdots \cdots {,r}_{\text{ave}N}),{r}_{\text{ave}}=\dfrac{1}{2}\left({r}_{v}+{r}_{t}\right) $ 燃料、时间加权策略 ${r}_{4}=\text{max}\text{(}{r}_{ {\rm{w} }1},{r}_{ {\rm{w} }2},\cdots \cdots ,{r}_{{\rm{w}}N}),{r}_{{\rm{w}}}=\dfrac{ {r}_{v}{r}_{t} }{ {r}_{v}+{r}_{t} }$ 燃料、时间最大最小策略 ${r}_{5}=\mathrm{max}\{ {r}_{\rm m1},{r}_{\rm m2},\dots \dots ,{r}_{{\rm{m}}N}\},{r}_{\rm m}=\mathrm{min}({r}_{v},{r}_{t})$ 不同的子树拓展规则本质是对剩余资源的选择分配,策略的优劣需要结合任务场景分析,优先选择在当前任务节点中更为稀缺的资源,避免出现一种资源耗尽而另一种资源剩余较多,从而影响最终的结果。
波束选择树搜索通过限制生成树拓展子树时的宽度来达到剪枝的效果,本文采用的波束宽度
${bw}$ 为10~10 000,在构造每一层的解序列后,根据节点的$ r $ 值排列所有可行解,若满足约束的单次转移兰伯特弧段数量小于波束宽度上限则拓展全部子弧段,若单次转移兰伯特弧段数量大于${bw}$ 则依据子树拓展规则对特定的子弧段进行宽度优先拓展,其余节点被剪枝删去。波束选择树搜索算法如表2所示。表 2 黄道面内多目标小行星飞越算法
Table 2. Multi-target asteroid flyby algorithm in ecliptic plane
黄道面内多目标小行星飞越算法 1. 计算潜在威胁小行星星历时刻对应过面位置$r_{\rm tr}^1 $、$r_{\rm tr}^2$与过面时刻$t_{\rm tr}^I $,$t_{\rm tr}^{II} $; 2. 根据探测器状态Sat<t, r, v>,利用兰伯特算法计算探测器飞越各小行星的速度增量Δν; 3. 筛选满足Δν<Δνmax 的目标集合N<Astnext1,……,AstnextN>; 4. 利用波束选择树搜索决策下一层子树的飞越小行星; 5. 更新探测器状态Sat<t, r, v>;更新飞越序列Sol; 6. 重复计算1~5直到N<Astnext1,……,AstnextN> = φ或违反任务约束。 1)生成根节点,根节点为地球E,设定发射时间与任务约束等任务指标;
2)拓展生成子节点,从根节点出发计算向目标小行星飞越弧段,如果飞越目标小行星满足速度增量约束则生成子节点,所有可行的子节点构成了搜索树的新一层解集;
3)排列与剪枝,对新一层子节点解集按照
$ r $ 值进行排列,保留${bw}$ 个子节点,对该层其余子节点进行剪枝删除;4)重复步骤2与3,直到违反剩余机动速度增量约束、任务总时间约束或无可行目标小行星。
波束选择树搜索算法流程如图3所示,最终根据总速度增量约束与总任务时间约束来终止子树的生成,待完整搜索树生成后,根据目标函数筛选出最优子树分支,本例目标函数设定为
$$ q=\mathrm{m}\mathrm{a}\mathrm{x} $$ (3) $$ q=\mathrm{m}\mathrm{a}\mathrm{x}\Bigg(N+\frac{{V}_{\mathrm{s}\mathrm{u}\mathrm{m}}}{{V}_{\mathrm{a}\mathrm{l}\mathrm{l}}}\Bigg)或q=\mathrm{m}\mathrm{a}\mathrm{x}\Bigg(N+\frac{{T}_{\mathrm{l}\mathrm{a}\mathrm{s}\mathrm{t}\_\mathrm{f}\mathrm{l}\mathrm{y}\mathrm{b}\mathrm{y}}}{{T}_{\mathrm{t}\mathrm{a}\mathrm{s}\mathrm{k}\_\mathrm{a}\mathrm{l}\mathrm{l}}}\Bigg) $$ (4) 即选择最长的子树分支作为最优解,
$ N $ 为探测器飞行单条轨迹飞越的小行星总数,但在单次空间任务生成的决策树中可能得到数条飞越小行星数目相同的轨迹,因此可用总速度增量或任务总时长占比作为除飞越目标个数外的评价指标供多个解飞越数量一致时筛选使用。 -
假设卫星总速度增量3.5 km/s,单次机动最大脉冲速度增量不大于1 km/s,任务周期10 a,两次飞越间隔时间不做约束。任务目标为在一次任务中飞越探测小行星数目最多,多个解同数目情况下按照式(4)选择剩余燃料多者,飞越最后一颗小行星时间作为任务终止时间。
-
假定发射时间选择在2025年1月1日,波束宽度设置为1 000。采用5种不同策略最终形成的解飞越个数、目标函数值、任务终止时间与剩余燃料的对比结果如表3所示。
表 3 不同树节点拓展规则结果
Table 3. Results of different tree node expansion rules
策略
名称飞越
个数目标函数
平均值最优序列任务
终止时间最优序列剩余速度
增量/(km∙s−1)$ {r}_{1} $ 19 19.0177 28th Dec 2034 0.1782 $ {r}_{2} $ 15 15.0187 22nd Aug 2031 0.1825 $ {r}_{3} $ 19 19.0105 14th Apr 2034 0.0864 $ {r}_{4} $ 19 19.0135 14th Apr 2034 0.0864 $ {r}_{5} $ 19 19.0122 28th Dec 2034 0.1782 结果显示,在该问题场景下以飞越目标个数为评价指标,相比任务时间,燃料无疑是更加“稀缺”的资源,策略2时间优先原则下的结果明显差于考虑燃料的其他原则,策略1、5的最优飞越结果一致也证明了这种观点;同时在波束宽度取值较大时,出现了不同策略最优序列相同的情况,说明该波束宽度对整个决策树的探索较为完全,决策树波束宽度
${bw}$ 取值合理。同时以波束宽度为10、50、200、500、5 000、10 000等为例在策略3下对比说明波束宽度的选择对结果的影响,发射时间设置为2025年1月1日,结果见表4。当波束宽度选择较小时,对于决策树的解探索不完全,得到的最优飞越个数相较于大波束宽度有一定差距;波束宽度增大,在飞越数量一致情况下最优解序列剩余速度增量变大,最终解序列则有一定的收敛趋势。
表 4 不同波束宽度结果
Table 4. Results of different beam widths
bw 飞越
个数目标函数
平均值最优序列任务
终止时间最优序列剩余速度
增量/(km∙s−1)10 17 17.038 6 23rd Feb 2034 0.286 3 50 18 18.027 5 06th Nov 2034 0.323 3 200 19 19.006 5 25th Aug 2034 0.040 7 500 19 19.010 8 14th Apr 2034 0.086 4 1 000 19 19.010 8 14th Apr 2034 0.086 4 2 000 19 19.013 5 28th Dec 2034 0.178 2 5 000 19 19.012 5 28th Dec 2034 0.178 2 10 000 19 19.012 5 28th Dec 2034 0.178 2 由于不同时刻潜在威胁小行星位置分布的差异,发射时间对于黄道面内多目标飞越小行星任务决策树的解生成与优化也有着重要影响,本文自2024年1月1日—2029年1月1日,以每6个月为一个发射窗口,选择策略
$ {r}_{4} $ ,波束宽度设置为1 000,对比验证发射时间对于解影响,结果如表5所示。在不同发射窗口下的最优序列飞越小行星个数在18 ~ 21颗,该任务对发射窗口要求不敏感,不同发射时间的该类任务均可以飞越数量可观的小行星。表 5 不同发射时间结果
Table 5. Results of different launch times
发射时间 飞越
个数目标
函数q
平均值最优序列任务
终止时间最优序列剩余
速度增量/
(km∙s−1)01stJan2024 19 19.023 3 25th Nov 2033 0.101 5 30thJun 2024 19 19.040 0 02nd May 2034 0.316 4 01stJan 2025 19 19.012 2 28th Dec 2034 0.178 2 30thJun 2025 18 18.025 6 24th May 2035 0.246 1 01stJan 2026 19 19.009 1 11th Aug 2035 0.031 9 30thJun 2026 18 18.013 7 15th May 2036 0.105 3 01stJan 2027 21 21.006 1 11th Dec 2036 0.021 3 30thJun 2027 20 20.016 1 25th Jun 2037 0.197 7 01stJan 2028 20 20.015 2 19th Dec 2037 0.053 3 30thJun 2028 20 20.024 2 14th May 2038 0.084 7 在表5的基础上,对发射窗口按月份进一步密集化搜索,发现于2027年9月1日发射的任务序列
$ L $ 最终也可飞越21颗小行星,但剩余燃料多,为0.194 5 km/s。该条轨道飞越各个小行星时刻与各次机动速度增量见表6,其中小行星直径、小行星轨道距地球轨道最近距离由ESA小行星数据库获得[27],探测器轨迹与飞越小行星位置如图4所示。表 6 多小行星飞越序列表
Table 6. Multiple asteroids’ flyby sequence
小行星编号 飞越时刻 速度增量/
(km∙s−1)小行星
直径/m轨道最小
交会距离
MOID/AU2001 SG286
2006 SU49
2015 UR51
2019 TN
2015 TV144
2017 HT2
2003 RB5
2012 UU136
2004 RQ10
Zephyr
2014 VL6
2006 HV5
2020 FR3
2011 LT17
2001 TC45
2017 NM6
2017 SL17
2004 AS1
2006 CY10
2012 KU12
2010 DA18th May 2028
24th Jan 2029
08th Nov 2029
09th Apr 2030
08th Jan 2031
01st Aug 2031
19th Sept 2031
05th May 2032
05th Sept 2032
20th Sept 2032
15th Nov 2032
14th May 2033
08th Oct 2033
08th Jul 2034
29th Nov 2034
12th Aug 2035
28th Nov 2035
31st Mar 2036
09th Jun 2036
28th Jan 2037
27th May 20370.136 5
0.224 4
0.075 8
0.194 6
0.102 4
0.056 6
0.188 5
0.207 6
0.302 3
0.380 5
0.096 5
0.121 0
0.146 2
0.052 6
0.020 3
0.081 0
0.173 5
0.138 2
0.340 4
0.143 2
0.123 6230
400
300
140
400
800
220
200
240
2060
210
300
260
170
500
600
150
290
300
170
3000.004 68
0.000 02
0.026 64
0.009 37
0.017 68
0.009 95
0.024 16
0.034 70
0.034 01
0.021 70
0.018 18
0.015 00
0.011 57
0.002 03
0.025 21
0.007 92
0.049 27
0.022 72
0.009 78
0.040 64
0.005 54本文以序列
$ L $ 为例考虑“长征三号乙”(CZ-3B)运载火箭发射探测器,探究单次任务中飞越小行星总数与剩余质量的关系,结果如图5所示。发射能量为零时“长征三号乙”运载火箭可提供约3 284 kg发射质量[28],结合工程需要,本方法可快速估计干重与飞越总数的关系,若在此发射条件下,干重为1 000 kg的探测器可飞越21颗小行星,干重为1 500 kg则可飞越16颗小行星。 -
为比较本文方法与Giuseppe等[17]迭代扩大搜索区域算法性能,在原方法搜索得到最优序列的同等任务条件下利用本文方法仿真,结果对比如表7所示。
表 7 本文方法与迭代扩大搜索区域算法对比
Table 7. Comparison between beam tree search and iterative extended region search algorithm
方法名称 飞越
个数最优序列总速度
增量/(km∙s−1)最优序列任务
终止时间本文方法 25 5.12 12th Dec 2028 迭代扩大搜索区域算法 21 5.30 14th May 2029[2] ——任务时间:2024年1月1日—2029年1月1日;
——单次机动最大脉冲速度增量:
$ \Delta {v}_{\mathrm{max}}=0.5\; \mathrm{k}\mathrm{m}/\mathrm{s} $ ;——单次任务总速度增量:
${v}_{\mathrm{t}\mathrm{o}\mathrm{t}}\text{ = 5}.5 \; \mathrm{k}\mathrm{m}/\mathrm{s}$ ;利用本文方法搜索,设置波束宽度
$ bw $ 为1 000,选择策略$ {r}_{4} $ 得到的针对Apollo型小行星的最优多目标飞越轨道共可飞越25颗小行星,多于其方法的21颗小行星,所需总速度增量为5.12 km/s也少于其所需的5.30 km/s,本文方法在个人电脑单次仿真用时为1 277.81 s[3]。仿真结果表明本文方法在保证快速性的同时,相比原方法搜索性能更强,更具全局性。其次本文参考西安飞行控制中心在GTOC9(Global Trajectory optimization Competition 9th)中的方法[29],针对飞越场景改造的蚁群算法对黄道面内多目标小行星飞越模型进行优化。本文所用方法与燃料最优的贪婪算法、蚁群算法优化结果的对比如表8所示。
表 8 本文方法与贪婪算法、蚁群算法对比
Table 8. Comparison table of beam tree search with greedy algorithm and ant colony algorithm
方法名称 最优飞越
个数最优序列总速度
增量/(km∙s−1)最优序列任务
终止时间贪婪–燃料最优 14 3.418 7 06th Jun 2035 蚁群算法 21 3.511 2 27th May 2037 本文方法 21 3.305 5 27th May 2037 以2027年9月1日发射的序列
$ L $ 为对照,多数情况下蚁群算法生成的序列飞越总数在17颗左右,需要重复调参并多次运行才得到飞越21颗小行星的解序列。而波束树搜索作为一种确定性算法,在确定节点拓展策略与波束宽度后结果稳定收敛。因此,对于多目标小行星飞越问题,本算法收敛稳定性优于蚁群算法。此外,蚁群算法收敛至飞越21颗小行星的序列所需总速度增量也在3.5~3.9 km/s,远大于本文方法。蚁群算法在多目标小行星飞越问题上表现不佳,原因可能为飞越探访转移速度增量不恒定导致启发式信息矩阵不断变化、目标群体数量较大导致信息素矩阵过于稀疏更新迭代慢及其swap、2-opt等局部搜索算子作用有限。
-
本文研究了在单次任务中对多颗潜在威胁小行星进行多目标飞越探测的轨道优化设计。针对该任务背景解空间大、任务轨迹计算耗时长的问题,利用黄道面内飞越潜在威胁小行星的基本策略简化模型,通过波束选择树搜索算法建立了一种快速有效求解多目标小行星飞越任务轨道的优化方法。仿真结果显示,波束宽度与树节点选择策略共同影响算法优化结果,在时间、空间复杂度允许的前提下,加大波束宽度减少剪枝范围有助于算法寻优;在本文任务场景下,燃料相比任务时间为更稀缺资源,仅考虑任务时间的树节点选择策略结果差于综合考虑燃料、飞越时间的树节点选择策略。不同发射窗口最优序列飞越潜在威胁小行星的个数有差异,仿真结果显示2024—2028年发射的任务至少可飞越探测18颗潜在威胁小行星,进一步细分时间窗口优化可得到数条最优飞越21颗小行星的序列。与迭代扩大区域搜索算法、蚁群算法的对比证明了本方法对多目标小行星飞越问题良好的求解性能。黄道面内多目标小行星飞越轨道设计方法,可为未来实施对潜在威胁小行星物理化学性质飞越普查任务设计提供参考。
未来将从黄道面内飞越小行星模型入手,将探测器飞越小行星后即刻变轨改为其可在原轨道运行一段时间后变轨;不再限定探测器仅在黄道面内进行飞越,允许探测器变倾角机动对在单次任务当中飞越探测更多潜在威胁小行星进行下一步的研究。
Trajectory Optimization Design for Multiple-Target Asteroid Flyby Mission in Ecliptic Plane
-
摘要: 近距离飞越小行星可获取小行星的表面图像、测量小行星的光谱、反演小行星物理化学性质。对潜在威胁小行星实施多目标飞越探测任务进行轨道优化设计,分析了潜在威胁小行星穿越黄道面时刻与位置分布,以小行星穿越黄道面时刻作为探测器飞越小行星的时刻,通过波束选择树搜索算法优化求解序贯飞越序列,建立了一种快速有效求解小行星序贯飞越任务轨道的优化方法模型。仿真结果表明,2024—2028年发射的小行星探测器可在黄道面内飞越探测至少18颗潜在威胁小行星,2027年9月发射窗口在10 a任务周期内可飞越21颗潜在威胁小行星。Abstract: Closely flying by asteroids can help to capture asteroid surface images, measure asteroid spectra, and obtain physical and chemical properties of asteroids. In particular, flying by multiple asteroids with potential hazards to the earth in one mission will significantly improve the understanding of the characteristics of potentially hazardous asteroids, and it is also of great significance to asteroid defense missions. In this paper, the trajectory of the multiple asteroid flyby mission of potentially hazardous asteroids was optimized. Firstly, the time and position distribution of asteroids passing through the ecliptic plane were analyzed, and the basic strategy of asteroids flyby in the ecliptic was determined. The time of asteroids crossing the ecliptic was taken as the time of asteroids’ flyby. Secondly, the sequential flyby sequence was optimized via beam selection tree search algorithm, and an optimization model for fast and effective solution of asteroid sequential flyby mission trajectory was established. Simulation results show that missions launched from 2024 to 2028 can fly by at least 18 potentially hazardous asteroids, especially the launch window in September 2027, which can fly by 21 potentially hazardous asteroids within a ten-year mission duration.Highlights
● The trajectory of multi-target asteroid flyby mission was solved fast and effectively by using the basic strategy of ecliptic in-plane flyover and beam search algorithm. ● Simulation results show that the mission is not sensitive to the launch window and the mission launch window from 2024 to 2028 can fly by at least 18 potentially hazardous asteroids. ● Combined with the constraints of launch vehicle and fuel to mass ratio of satellite,this method can quickly determine the number of sequential flyby asteroids in a single mission. -
表 1 节点拓展策略
Table 1 Node expansion strategy
策略名 策略定义 燃料最省
策略$ {r}_{1}=\text{max}\text{(}{r}_{v1},{r}_{v2},\cdots \cdots \text{,}{r}_{vN}\text{)} $ 时间最省
策略$ {r}_{2}=\text{max}\text{(}{r}_{t1},{r}_{t2},\cdots \cdots {,r}_{tN}\text{)} $ 燃料、时间平均策略 $ {r}_{3}=\text{max}{(r}_{\text{ave}1},{r}_{\text{ave}2},\cdots \cdots {,r}_{\text{ave}N}),{r}_{\text{ave}}=\dfrac{1}{2}\left({r}_{v}+{r}_{t}\right) $ 燃料、时间加权策略 ${r}_{4}=\text{max}\text{(}{r}_{ {\rm{w} }1},{r}_{ {\rm{w} }2},\cdots \cdots ,{r}_{{\rm{w}}N}),{r}_{{\rm{w}}}=\dfrac{ {r}_{v}{r}_{t} }{ {r}_{v}+{r}_{t} }$ 燃料、时间最大最小策略 ${r}_{5}=\mathrm{max}\{ {r}_{\rm m1},{r}_{\rm m2},\dots \dots ,{r}_{{\rm{m}}N}\},{r}_{\rm m}=\mathrm{min}({r}_{v},{r}_{t})$ 表 2 黄道面内多目标小行星飞越算法
Table 2 Multi-target asteroid flyby algorithm in ecliptic plane
黄道面内多目标小行星飞越算法 1. 计算潜在威胁小行星星历时刻对应过面位置$r_{\rm tr}^1 $、$r_{\rm tr}^2$与过面时刻$t_{\rm tr}^I $,$t_{\rm tr}^{II} $; 2. 根据探测器状态Sat<t, r, v>,利用兰伯特算法计算探测器飞越各小行星的速度增量Δν; 3. 筛选满足Δν<Δνmax 的目标集合N<Astnext1,……,AstnextN>; 4. 利用波束选择树搜索决策下一层子树的飞越小行星; 5. 更新探测器状态Sat<t, r, v>;更新飞越序列Sol; 6. 重复计算1~5直到N<Astnext1,……,AstnextN> = φ或违反任务约束。 表 3 不同树节点拓展规则结果
Table 3 Results of different tree node expansion rules
策略
名称飞越
个数目标函数
平均值最优序列任务
终止时间最优序列剩余速度
增量/(km∙s−1)$ {r}_{1} $ 19 19.0177 28th Dec 2034 0.1782 $ {r}_{2} $ 15 15.0187 22nd Aug 2031 0.1825 $ {r}_{3} $ 19 19.0105 14th Apr 2034 0.0864 $ {r}_{4} $ 19 19.0135 14th Apr 2034 0.0864 $ {r}_{5} $ 19 19.0122 28th Dec 2034 0.1782 表 4 不同波束宽度结果
Table 4 Results of different beam widths
bw 飞越
个数目标函数
平均值最优序列任务
终止时间最优序列剩余速度
增量/(km∙s−1)10 17 17.038 6 23rd Feb 2034 0.286 3 50 18 18.027 5 06th Nov 2034 0.323 3 200 19 19.006 5 25th Aug 2034 0.040 7 500 19 19.010 8 14th Apr 2034 0.086 4 1 000 19 19.010 8 14th Apr 2034 0.086 4 2 000 19 19.013 5 28th Dec 2034 0.178 2 5 000 19 19.012 5 28th Dec 2034 0.178 2 10 000 19 19.012 5 28th Dec 2034 0.178 2 表 5 不同发射时间结果
Table 5 Results of different launch times
发射时间 飞越
个数目标
函数q
平均值最优序列任务
终止时间最优序列剩余
速度增量/
(km∙s−1)01stJan2024 19 19.023 3 25th Nov 2033 0.101 5 30thJun 2024 19 19.040 0 02nd May 2034 0.316 4 01stJan 2025 19 19.012 2 28th Dec 2034 0.178 2 30thJun 2025 18 18.025 6 24th May 2035 0.246 1 01stJan 2026 19 19.009 1 11th Aug 2035 0.031 9 30thJun 2026 18 18.013 7 15th May 2036 0.105 3 01stJan 2027 21 21.006 1 11th Dec 2036 0.021 3 30thJun 2027 20 20.016 1 25th Jun 2037 0.197 7 01stJan 2028 20 20.015 2 19th Dec 2037 0.053 3 30thJun 2028 20 20.024 2 14th May 2038 0.084 7 表 6 多小行星飞越序列表
Table 6 Multiple asteroids’ flyby sequence
小行星编号 飞越时刻 速度增量/
(km∙s−1)小行星
直径/m轨道最小
交会距离
MOID/AU2001 SG286
2006 SU49
2015 UR51
2019 TN
2015 TV144
2017 HT2
2003 RB5
2012 UU136
2004 RQ10
Zephyr
2014 VL6
2006 HV5
2020 FR3
2011 LT17
2001 TC45
2017 NM6
2017 SL17
2004 AS1
2006 CY10
2012 KU12
2010 DA18th May 2028
24th Jan 2029
08th Nov 2029
09th Apr 2030
08th Jan 2031
01st Aug 2031
19th Sept 2031
05th May 2032
05th Sept 2032
20th Sept 2032
15th Nov 2032
14th May 2033
08th Oct 2033
08th Jul 2034
29th Nov 2034
12th Aug 2035
28th Nov 2035
31st Mar 2036
09th Jun 2036
28th Jan 2037
27th May 20370.136 5
0.224 4
0.075 8
0.194 6
0.102 4
0.056 6
0.188 5
0.207 6
0.302 3
0.380 5
0.096 5
0.121 0
0.146 2
0.052 6
0.020 3
0.081 0
0.173 5
0.138 2
0.340 4
0.143 2
0.123 6230
400
300
140
400
800
220
200
240
2060
210
300
260
170
500
600
150
290
300
170
3000.004 68
0.000 02
0.026 64
0.009 37
0.017 68
0.009 95
0.024 16
0.034 70
0.034 01
0.021 70
0.018 18
0.015 00
0.011 57
0.002 03
0.025 21
0.007 92
0.049 27
0.022 72
0.009 78
0.040 64
0.005 54表 7 本文方法与迭代扩大搜索区域算法对比
Table 7 Comparison between beam tree search and iterative extended region search algorithm
方法名称 飞越
个数最优序列总速度
增量/(km∙s−1)最优序列任务
终止时间本文方法 25 5.12 12th Dec 2028 迭代扩大搜索区域算法 21 5.30 14th May 2029[2] 表 8 本文方法与贪婪算法、蚁群算法对比
Table 8 Comparison table of beam tree search with greedy algorithm and ant colony algorithm
方法名称 最优飞越
个数最优序列总速度
增量/(km∙s−1)最优序列任务
终止时间贪婪–燃料最优 14 3.418 7 06th Jun 2035 蚁群算法 21 3.511 2 27th May 2037 本文方法 21 3.305 5 27th May 2037 -
[1] SCHULTE P,ALEGRET L,ARENILLAS I,et al. The chicxulub asteroid impact and mass extinction at the cretaceous-paleogene boundary[J]. Science,2010,327(5970):1214-1218. doi: 10.1126/science.1177265 [2] 龚自正,李明,陈川,等. 小行星监测预警、安全防御和资源利用的前沿科学问题及关键技术[J]. 科学通报,2020,65(5):346-372. doi: 10.1360/TB-2019-0425GONG Z Z,LI M,CHENG C,et al. The frontier science and key technologies of asteroid monitoring and early warning,security defense and resource utilization[J]. Chinese Science Bulletin,2020,65(5):346-372. doi: 10.1360/TB-2019-0425 [3] 廖慧兮,王彤,贾晓宇. 小行星探测进展及技术特点分析[J]. 国际太空,2017(7):2-9. doi: 10.3969/j.issn.1009-2366.2017.07.001 [4] VETRISANO M,VASILE M. Autonomous navigation of a spacecraft formation in the proximity of an asteroid[J]. Advance in Space Research,2016,57(8):1783-1804. doi: 10.1016/j.asr.2015.07.024 [5] WALKER L,CARLO D M,GRECO C,et al. A mission concept for the low-cost large-scale exploration and characterisation of near Earth objects[J]. Advances in Space Research,2020,67(11):3880-3908. [6] SMPAG. SMPAG statement of support for small-class, high-velocity flyby missions to small bodies for planetary defence[EB/OL]. (2021-03-25). https://www.cosmos.esa.int/web/smpag/meeting-16-mar-2021-. [7] WATANABE S I,TSUDA Y,YOSHIKAWA M,et al. Hayabusa2 Mission Overview[J]. Space Science Reviews,2017,208(1-4):3-16. doi: 10.1007/s11214-017-0377-1 [8] RUSSELL T C,CAPACCIONI F,CORADINI A,et al. Dawn mission to Vesta and Ceres[J]. Earth,Moon,and Planets,2007,101(1-2):65-91. [9] MARCHI S,OLKIN C B. Lucy in the sky with Trojan asteroids[J]. Nature Astronomy,2021,5(11):1178-1178. doi: 10.1038/s41550-021-01534-6 [10] SARLI B V,HORIKAWA M,YAM C H,et al. DESTINY+ trajectory design to(3200)Phaethon[J]. The Journal of the Astronautical Sciences,2018,65(1):82-110. doi: 10.1007/s40295-017-0117-5 [11] GATER W. Comet mission given green light by European Space Agency[J]. Physics World,2019,32(8):13. doi: 10.1088/2058-7058/32/8/22 [12] GAO Y,LU X,PENG Y,et al. Trajectory optimization of multiple asteroids exploration with asteroid 2010TK7 as main target[J]. Advances in Space Research,2019,63(1):432-442. doi: 10.1016/j.asr.2018.08.047 [13] 夏炎,罗永杰,赵海斌,等. 主带小行星深空探测可接近性与多目标探测轨道的实现[J]. 天文学报,2010,51(2):163-172.XIA Y,LUO Y J,ZHAO H B,et al. Accessibility for main belt asteroid exploration and trajectory design for multiple asteroids[J]. Acta Astronomica Sinica,2010,51(2):163-172. [14] QIAO D,CUI P Y,CUI H H. Proposal for a multiple-asteroid-flyby mission with sample return[J]. Advances in Space Research,2012,50(3):327-333. doi: 10.1016/j.asr.2012.04.014 [15] MCNUTT L, JOHNSON L, CLARDY D, et al. Near-earth asteroid(NEA)scout[C]//AIAA Space 2014 Conference and Exposition. San Diego, CA: AIAA, 2014. [16] GRECO C, DI CARLO M, WALKER L, et al. Analysis of NEOs reachability with nano-satallites and low-thrust propulsion[C]//4S Symposium 2018-Small Satellites Systems and Services. [S. l. ]: AIAA, 2018. [17] CATALDI G,MARCUCCIO S. A 2-D Trajectory design algorithm for multiple asteroid flyby missions[J]. Aerotecnica Missili & Spazio,2020,99(4):287-295. [18] LI S,HUANG X Y,YANG B. Review of optimization methodologies in global and China trajectory optimization competitions[J]. Progress in Aerospace Sciences,2018,102:60-75. doi: 10.1016/j.paerosci.2018.07.004 [19] SHIRAZI A,CEBERIO J,LOZANO J A. Spacecraft trajectory optimization:a review of models,objectives,approaches and solutions[J]. Progress in Aerospace Sciences,2018,102:76-98. doi: 10.1016/j.paerosci.2018.07.007 [20] JIANG F,CHEN Y,LIU Y,et al. GTOC5:results from the Tsinghua University[J]. Acta Futura,2014,8(1):37-44. [21] IZZO D, HENNES D, SIMÕES L F, et al. Designing complex interplanetary trajectories for the global trajectory optimization competitions[M]. Switzerland: Springer, Cham, 2016. [22] LUO Y, ZHU Y, ZHU H, et al, GTOC9: results from the national university of defense technology[J]. Acta Futura, 2018(11): 37-47. [23] CASALINO L,PASTRONE D,SIMEONI F,et al. GTOC5:results from the Politecnico di Torino and Università di Roma Sapienza[J]. Acta Futura,2014,8(1):29-36. [24] PETROPOULOS A E,BONFIGLIO E P,GREBOW D J,et al. GTOC5:results from the Jet Propulsion Laboratory[J]. Acta Futura,2014,8(1):21-27. [25] HD· 柯蒂斯 (美). 轨道力学[M]. 周建华, 徐波, 冯全胜, 译. 北京: 科学出版社, 2009. [26] LI H Y,BAOYIN H X. Sequence optimization for multiple asteroids rendezvous via cluster analysis and probability-based beam search[J]. Science China Technological Sciences,2021,64(1):122-130. doi: 10.1007/s11431-020-1560-9 [27] ESA. Near-Earth objects coordination centre [EB/OL]. [2021-11-20]. https: //neo.ssa.esa.int/search-for-asteroids. [28] LIU J,ZHENG J,LI M. Dry mass optimization for the impulsive transfer trajectory of a near-Earth asteroid sample return mission[J]. Astrophysics and Space Science,2019,364(12):1-14. [29] SHEN H X,ZHANG T J,HUANG A Y,et al. GTOC9:results from the Xi’an Satellite Control Center(team XSCC)[J]. Acta Futura,2018,11:49-55. -