三体轨道动力学研究进展1)

李翔宇 乔栋 程潏

(北京理工大学宇航学院,北京 100081)

(飞行器动力学与控制教育部重点实验室,北京 100081)

(深空探测导航与控制工信部重点实验室,北京 100081)

引言

物理学和经典力学中,三体问题是指:已知三质点的初始位置和速度,根据牛顿运动定律和牛顿万有引力定律求解其后续运动的问题[1].三体问题动力学模型是天体力学中继传统二体模型、受摄二体模型后更为精确的动力学模型,也是当前航天动力学中广泛采用的动力学模型之一.研究三体问题具有丰富的理论与工程意义.理论上,三体问题研究能推动天体力学的发展,为探索更为复杂的N体问题奠定基础,进而帮助人类更好地理解自然界天体运动,丰富天文学理论;工程上,三体问题研究有助于建立精确的动力学模型,为空间探测任务提供有力的理论支撑.

三体问题可分为一般三体问题和限制性三体问题两类.一般三体问题的研究可以追溯到17 世纪80 年代.1687 年,Newton 在其《自然哲学的数学原理》中首次研究了3 个大质量物体在相互扰动的引力作用下的运动问题,并将其应用于日地月系统.18世纪40 年代,d’Alembert 等曾尝试在更普遍的意义上研究此问题,并首次使用名词“三体问题”.Bruns和Poincar´e分别于1887 年和1889 年指出,由代数表达式与积分式表示的三体问题不存在通解.1913 年,Sundman[2]证明了三体问题存在一类幂级数解,但其收敛速度极其缓慢.此后,一般三体问题的研究重点由寻求通解转向特解,并找到了三族周期性特解:Lagrange-Euler 族,即3 个天体等间距地在圆轨道上运动;Broucke-H´enon 族,即两个天体在一个区域内运动,第3 个天体在它们外围做环绕运动;“8”字形族,即3 个天体运动在一条“8”字形的轨道上互相追赶.2013 年,Suvakov 等[3]借助计算机模拟发现了新的13 族周期性特解,将三体问题的周期解扩充至16 族.研究表明除了Lagrange-Euler 族轨道适用于任意质量比的三体系统外,其他族周期轨道要求三天体的质量相等.

以上研究针对一般的三体问题,它与现实航天动力学的情况存在较大的差异.为了将三体模型应用于航天动力学,科学家对模型进行了简化与改进,将三体问题简化成限制性三体问题RTBP,即认为系统中某一天体的质量相对于另两天体可忽略.根据两个质量天体的相互运动,可以将限制性三体问题分为圆型限制性三体问题CRTBP 和椭圆限制性三体问题ERTBP[4].此外,为了研究类似日地月系统这类存在卫星的行星系统,学者们又建立了四体问题模型,包括双圆模型[5]和拟双圆模型[6]等.

基于简化的三体问题模型,科学家对运动形式进行了深入研究,并发现了平衡点这一特殊解集.但受限于数学认知,简化模型下的全局运动仍然无法完全求解.随着航天技术的发展,20 世纪60 年代开始,科学家将研究重点转向平衡点附近的运动,并讨论了平衡点在航天领域的潜在应用.平衡点附近周期轨道的近似解析解[7]和数值计算方法[8]被相继发展并应用于地月系统的轨道设计.随后,平衡点附近不同类型的周期运动及其转移轨道成为科学家研究的热点,基于周期轨道的应用设想也如雨后春笋般出现.

进入21 世纪,随着人们对非线性动力学理论研究的深入,不变流形理论[9]被引入三体系统研究中,并用于三体系统内的低能量转移轨道设计.随后,三体系统的研究也由平衡点附近区域的运动进一步拓展至三体系统中质量较小天体附近区域,提出了弱稳定边界理论[10]和稳定集理论[11],并应用于三体系统下的弹道捕获轨道设计.弱稳定边界理论与不变流形间的联系也得到了学者的关注[12],有望建立起三体系统中统一的轨道描述方法.

三体系统动力学的理论研究也推动了航天工程应用的发展,尤其在月球探测和深空探测领域得到了广泛的应用.三体系统平衡点可对主天体形成全面的覆盖,决定了其在深空导航和通信中的巨大应用价值.就地月系统而言,L2点附近Halo 轨道上运行的飞行器能够提供对地球和月球背面连续不间断的通信,因而是月球背面探测中通信卫星的理想飞行轨道;L4点和L5点的卫星能形成非常大的基线,是天基VLBI 基站的理想部署位置.1978 年,ISEE-3成为首颗进入日地平衡点轨道的人造卫星[13].此后1991 年,Belbruno 等[14]创造性地将弱稳定转移方式应用于日本失控卫星Hiten 的任务救援中,在太阳引力摄动下,经过两次月球借力飞行,成功实现月球捕获.2018 年,我国的嫦娥4 号中继星也成功进入地月L2点Halo 轨道中[15],成为首个在地月周期轨道长期稳定运行的通讯卫星,为月球背面探测提供通信中继服务.截止目前共有15 颗探测器到达地月或日地平衡点开展探测,NASA,ESA,JAXA 等均有相关任务计划,平衡点探测任务如表1 所示.

表1 平衡点探测任务概述Table 1 Summary of libration point exploration missions

随着各国小行星探测和月球探测计划的稳步推进,三体系统动力学的应用将更加深入与广泛.本文以当前和未来的深空探测任务为背景,对航天动力学中三体问题相关研究进展进行全面的回顾与总结.

1 三体系统的动力学模型

三体系统的动力学建模是研究三体系统内运动的基础,长期以来,科学家对三体系统动力学模型开展了深入研究,根据运动的特性进行了不同程度的简化与假设,提出了圆型限制性三体模型、椭圆限制性三体问题模型、Hill 模型、双圆模型和拟双圆模型等动力学模型.下面将对各动力学模型的特性与应用范围进行综述.

1.1 圆型限制性三体模型

圆型限制性三体问题CRTBP 由于其形式简单且能较好的描述平衡点附近的动力学特性而成为广泛采用的三体系统动力学模型.圆型限制性三体问题中两个主天体(P1和P2) 绕它们的公共质心做圆周运动,第三体的质量可忽略,即不考虑第三体对两个主天体运动的影响.通过归一化,质心旋转系下的圆型限制性三体动力学可以表示为

1.1.1 Jacobi 常数与零速度曲面

1836 年,Jacobi 发现圆型限制性三体模型在旋转坐标系下的能量积分为常数,即Jacobi 常数C,也是目前为止CRTBP 中唯一存在的积分常数

令方程(3) 中的速度v=0,可得曲面方程2Ω(x,y,z)=C,即为零速度曲面.零速度曲面的结构随Jacobi 常数C的变化而变化,如图1 所示.通过计算平衡点对应的Jacobi 常数,可以确定探测器在系统内的可能运动区域.当C>C1时,探测器只能在两个主天体各自附近的空间内运动;C1>C>C2时,探测器可以经L1点实现主天体附近空间的转移;C2>C>C3时,探测器可以经L2点进入外部空间,因而通常称日地L2点为通往深空的大门;进一步增大探测器能量C3>C时,探测器可以从L3点进入外部空间,因而Jacobi 常数和零速度曲面是表述三体系统内轨道运动的重要指标.

图1 三体系统中的零速度曲面Fig.1 Zero-velocity curve in the three-body system

1.1.2 Hill 三体模型

为了研究日地月系统中月球的运动,Hill 对CRTBP 进行了进一步简化,将旋转坐标的原点由系统质心平移至第二个主天体中心,同时假设R/R1≪1,从而对方程(1)右侧的引力项进行线性化,从而得到方程

若忽略探测器在Z方向的运动,方程(4)即为Hill 模型,此后Brown[16]做了改进,形成了Hill-Brown 理论,用于描述月球的运动.Hill 模型同样可以近似描述质量比较小系统中探测器在第二主天体附近的运动.

1.1.3 圆型限制性三体问题的局限性

圆型限制性三体模型是研究三体问题的基础,但其存在两个很重要的假设:(1)主天体的运行轨道为圆;(2)系统中不存在其他天体对探测器的运动产生摄动作用.圆型限制性可以较好地应用于太阳-地球(日地) 三体系统的运动分析,但在研究地球-月球(地月)三体系统和太阳-火星(日火)、太阳-水星(日水)三体系统时遇到了问题.

对于地月系统内的运动,太阳引力的影响不可忽略.太阳引力对地月系统内探测器的影响分为直接和间接影响两方面:一方面,Farquhar[17]对太阳引力的影响作了全面评估并得出结论:太阳引力对探测器的直接影响不可忽略;另一方面,太阳引力还会影响地月系统的运动[18],从而间接影响探测器的运动.太阳引力对地月系统运动的影响体现在:(1)月球绕地球的运行轨道不再是理想的圆或椭圆.月球轨道的实际偏心率在1/15～1/23 之间变化(周期为31.8 d)[19];(2)月球绕地球运行的轨道平面不再是空间中的固定平面,而是存在着章动(平均周期为173 d)和进动(平均周期18.6 年)[20].因此,对地月系统,圆型限制性三体模型的理论计算与实际结果存在较大误差.

此外,水星、火星的轨道偏心率分别约为0.205 6,0.093 4,日-火和日-水系统内的运动若采用圆型限制性三体模型,主天体的运动与实际椭圆轨道存在较大差异,因而相关结果也无法直接应用于轨道设计.

为了进一步提高动力学模型在处理上述问题中的精度,科学家分别考虑主天体的椭圆运动和太阳引力摄动作用,建立了椭圆型限制性三体模型与限制性四体模型.

1.2 椭圆限制性三体问题

椭圆型限制性三体模型ERTBP 中考虑了两主天体轨道运动的偏心率,如图2 所示[21],动力学方程被改写为

图2 椭圆型限制性三体模型[21]Fig.2 Elliptic restricted three-body problem[21]

式中,ω 为旋转坐标系相对惯性系的瞬时角速度,˙ω为瞬时角加速度.半长轴,e为该椭圆运动的偏心率,ν 为真近点角.同时两主天体P1和P2在旋转坐标系x轴上的位置也不再固定不变.

ERTBP 为时变系统,系统中不存在Jacobi 积分,进而不存在CRTBP 中的零速度面和平衡点,仅存在任一瞬时方程(5)左端项为0 对应的瞬时平衡点.针对ERTBP 的研究主要集中于周期轨道及其稳定性.研究发现,在平面内有一族稳定轨道和六族不稳定轨道[4].ERTBP 模型被应用于分析探测器在水星附近的运动及稳定性[22-23].2018 年10 月,ESA 和JAXA合作研制的“Bepi-Colombo”水星探测器于圭亚那太空中心发射升空,该任务基于ERTBP 模型开展设计转移与捕获轨道设计[24].

1.3 四体模型

为了考虑太阳的引力影响,从20 世纪50 年代开始,国内外的学者对受摄三体系统进行了大量的研究[5,20,25].针对日-地-月系统提出了双圆四体动力学模型,即认为地-月系统和日-地-月质心系统为两个独立的圆型限制性三体系统,在此基础上研究探测器的运动.Nicholson[18]忽略了月球轨道相对于黄道面的轨道倾角,讨论了太阳扰动下的地月共线平衡点附近轨道.Lewallen 和Tapley[20]考虑了月球轨道倾角并分析了地月三角平衡点的变化,但没有考虑月球轨道的进动和偏心率.G´omez 等[26]同样基于双圆四体模型,分析了三角平衡点的运动稳定性.Castelli[27]基于惯性系对双圆模型和三体动力学模型下运动的误差进行了对比分析,并给出了各模型的最佳适用范围.由于双圆四体动力学模型中两个三体模型相互独立,实际上并不符合牛顿第二定律,模型非自洽.为了满足自洽性,Andr´ea 等[28]提出了拟双圆四体动力学模型,如图3 所示[29].

图3 拟双圆四体动力学模型[29]Fig.3 Quasi-bicircular four-body problem[29]

对于日-地-月系统,拟双圆四体模型的哈密顿函数可表示为

式中,θ=ωst,ωs为太阳的归一化运动频率,方程(6)中系数的具体值可参见文献[29].

拟双圆四体模型保证日-地-月的运动满足三体问题的解.基于该模型,Andr´ea[29]讨论了平衡点的不变流形与低能量转移.此后Guzman[30]通过多步打靶法实现拟双圆四体动力学模型下地月平衡点附近的周期轨道/拟周期轨道的设计,结果与真实星历情况接近,验证了拟双圆模型在轨道设计中的可行性.

此外,在研究火星-火卫系统以及双小行星系统内的运动时,也有学者将小天体的形状摄动力加入三体动力学模型中[31-33],用来分析不规则小天体系统内的运动,扩展了三体动力学模型的应用范围.

从以上研究可以看出,对三体系统的动力学建模经历了由繁到简再到繁的过程,当前的动力学模型已能够满足多数太阳-行星,行星-卫星系统内的运动分析,目前的研究集中在将除天体质心引力外的其他摄动力融入三体系统动力学,如天体非球型摄动、太阳光压、洛伦兹力等,从而更精确地描述真实空间环境下三体系统内的运动特征.

2 三体系统平衡点附近的局部运动

三体系统内的运动始终是三体问题研究的重点.根据运动分布的区域,三体系统内的运动可分为在三体系统平衡点附近的局部运动、主天体附近的局部运动以及三体系统内的全局运动,其中主天体附近的局部运动与受摄二体轨道较为相似,但在三体动力学下将表现出不同的运动形式.而平衡点及其附近的周期运动由于其独特的动力学特性得到了广泛的关注,并在一系列的深空探测任务中得到了应用.本节将着重介绍平衡点附近的局部运动,首先给出三体系统的平衡点,然后回顾圆型限制性模型下平衡点附近周期运动的求解及轨道保持方法,最后给出拟周期轨道与特殊周期轨道的求解.

2.1 三体系统平衡点

若令方程(1)中¨x=¨y=¨z=0,˙x=˙y=˙z=0,对应的解即为三体系统的平衡点,位于平衡点的质点将与主天体保持相对不变的位置和速度.圆型限制性三体系统中共存在5 个平衡点,其中包括3 个共线平衡点和两个三角平衡点.Euler 于1767 年发现了的3 个共线平衡点L1,L2和L3,其中L1点位于两个主天体之间,L2点位于P2外侧,而L3点位于P1一侧,L1和L2点距离P2较近,而L3距离主天体较远.此后,Lagrange 于1772 年求解了三角平衡点L4和L5,它们与主天体构成等边三角形.随着人们在太阳-木星(日木)系统L4和L5点附近发现大量的小行星群落(特洛伊族小天体),三体动力学模型平衡点的存在得以证实.图4 给出了三体系统平衡点分布示意图.

将方程(1)在平衡点附近线性化,通过求解对应的特征方程可以判断平衡点的稳定性.共线平衡点为不稳定平衡点,施加扰动后将无法继续保持在平衡点位置;三角平衡点为条件稳定平衡点,当两主天体质量满足以下条件时,三角平衡点是稳定的

图4 三体系统平衡点分布示意图Fig.4 Distributions of five libration points in the three-body system

探测器位于稳定的三角平衡点,受到扰动后仍维持在平衡点附近,不会发生逃逸.目前太阳系常见的三体系统,如日地系统、地月系统、日木系统等均满足上述质量约束条件,因此它们的三角平衡点均是稳定的,但部分质量接近的双小行星系统具有不稳定的平衡点,此外研究表明小行星的不规则形状摄动也会影响三角平衡点的稳定性[32].地月系统三角平衡点的卫星能形成非常大的基线,是天基VLBI基站的理想部署位置.

2.2 平衡点附近的周期与拟周期运动

共线平衡点(L1,L2和L3) 的拓扑类型为“中心×中心×鞍点”,鞍点型的分量使其成为不稳定的平衡点[34],而“中心×中心”的分量则形成为四维中心流形,分别对应着平衡点附近的两族周期轨道:平面Lyapunov 轨道和垂直Lyapunov 轨道.从平面Lyapunov 轨道族中还分岔出另一族周期轨道:Halo 轨道.由KAM 定理可知[35],四维中心流形内存在两个频率不共振的周期运动耦合而成的无数拟周期轨道,也称为Lissajous 轨道.需要指出的是,Halo 轨道对应着两个频率相等的特殊情况.当频率的比值可以表示成两个整数比时,拟周期轨道将退化为共振轨道.圆型限制性三体系统中平衡点附近的周期轨道求解经历了解析近似解和精确数值解两个阶段,在周期解的基础上,学者又通过延拓和频闪映射得到了拟周期轨道以及特殊的周期运动族.

2.2.1 周期轨道的近似解析解

解析法的基础是摄动理论中求解常微分方程近似周期解的经典方法——LP 方法(linstedt-poincare),构建平衡点运动的扰动方程并将扰动项进行幂级数展开,消除长期项后得到正余弦形式的扰动运动周期解.利用LP 法,20 世纪70 年代,Farquhar 等[36]推导了日地系和地月系点附近的大幅值Halo 轨道以及小幅值的Lissajous 轨道的线性解析表达式.Richardson[7]进一步推导了Halo 轨道的三阶展开式,成为求解Halo 轨道近似解析解的经典解法.此后,尽管也有学者求解了Halo 轨道和Lissajous 轨道的高阶解析表达式[37],但计算形式繁琐,且三阶解析解已可满足数值计算的精度要求.总体来看,解析法的形式较复杂,且由于受到展开式级数的限制,解析法的精度有限,无法直接用于周期轨道设计.随着轨道距离平衡点的距离增大,解析解的拟合误差也显著增大.但其可以为数值求解提供较好的初值猜测,为更大范围相空间内周期轨道的求解奠定了基础.

2.2.2 数值法求解精确周期轨道

数值法求解周期轨道的本质为微分修正法,其求解过程可描述为:将解析法求得的周期轨道作为初值猜测,迭代修正初始点的状态直至满足周期性条件,也即轨道在经过一个周期的积分后返回初始点.Howell[38]为周期轨道的数值法求解做出了突出的贡献,首次提出了微分修正法求解Halo 轨道的思路,取三阶解析解作为初值,通过修正初值来获得满足周期性条件的轨道,该方法也称为单点打靶法,如图5(a) 所示.此外,利用CRTBP 系统的对称性,选择初始点从对称轴(或平面) 出发,可降低问题的维度并缩减积分时间.然而,上述方法并不适用于Lissajous 轨道,这是由于Lissajous 轨道拟周期的属性使其永远无法返回同一初始点.针对这一问题,Howell等[8]提出了用于求解Lissajous 轨道的两级微分修正法(图5(b)):首先,利用解析法得到一系列的离散点;然后,利用单点打靶法得到一条位置连续的轨道,此为第一级修正;最后,再利用微分修正法抵消中间拼接点的速度误差,从而得到一条位置与速度都连续的轨道.Folta 等[39]进一步将多点打靶法扩展到了高精度星历模型下,成功设计了ARTEMIS 探测器的轨道,这也是人类历史上第一颗地月平衡点探测器.

图5 微分修正法的示意图Fig.5 Illustration of the differential correction method

基于微分修正的数值方法只适用于单个周期轨道的求解.若要获得整个轨道族,需要求解不同的解析解初值进行修正,过程繁琐且对于振幅较大的周期轨道失效.而利用数值延拓技术可以较快地获得轨道族.延拓法从已知的轨道出发,选择合适的延拓参数来猜测出新的周期轨道的初值,然后利用微分修正法得到一条新的周期轨道.延拓参数一般为表征周期轨道的特征参数,比如雅可比积分或轨道周期.针对单参数延拓法容易遇到奇点的情况,Simo[40]提出了稳定且高效的变步长弧长延拓法,沿着特征曲线方向进行延拓,猜测初值的精度较高且延拓速度快.利用变步长弧长延拓法,图6 给出了几类典型的平衡点轨道族,包括平面与垂直Lyapunov 轨道族、Halo 轨道族、以及Lissajous 轨道[41].

2.2.3 平衡点周期轨道的保持

图6 几类典型的平衡点轨道示意图Fig.6 Typical orbits about libration points

由于平衡点轨道具有弱稳定的特性,受到复杂的自然扰动力的影响,航天器将会偏离标称轨道,因此必须定期执行轨道保持机动将其控制在平衡点轨道附近的一定区域内.1961 年Colombo[42]首次评估了稳定地月L1点附近平衡点轨道所需的能量,受限于当时对平衡点轨道的认识,控制目标仅仅是将航天器维持在平衡点附近.1966 年Farquhar[17]评估了地月L2点的中继卫星轨道保持每年所需的能量约为10 m/s.随着有关于平衡点附近的周期、拟周期轨道理论的完善,轨道保持理论也得到了充分的发展.

平衡点轨道保持策略的研究可以分为3 个阶段.1974 年Breakwell 等[43]基于最优控制理论,详细论述了地月L2点附近周期轨道跟踪标称轨道的方案,该方法可看作是靶点法的雏形.Howell 和Pernicka[44-45]基于线性控制思想,对靶点法进行了改进,给出了靶点模式的脉冲控制策略,可作为研究的第一阶段.在这一阶段的研究通常假设导航数据足够精确、修正脉冲(或者小推力)可以频繁施加.在该假设下,对于一个日地系平衡点探测器(周期约为6 m),每一周或两周施加一次微小的脉冲修正,甚至可以得出轨道保持能量接近于零的结论.Dunham 等[46]则指出:研究无工程应用价值的轨道保持方案是无意义的.

随着第一次平衡点飞行任务ISEE-3 的执行,来自工程执行方面的约束和导航测量的约束被考虑到轨道保持策略的设计中[47].自此,平衡点轨道保持策略的研究进入第2 个阶段.脉冲控制方案由于其在工程上易实现、可靠性高、计算效率高而被广泛采用,其中文献[46]中记载的X轴速度控制法被用于日地系统平衡点飞行器SOHO 以及ACE 的轨道保持中.该算法根据三体条件下周期轨道穿越X-Z平面(日地连线与地球公转轨道角动量方向确定的平面) 时刻,X方向速度为零的特征,对探测器每次穿越X-Z平面的X轴方向的速度进行限制.文献[48-49]中的方法都延续了ISEE-3 的控制思想,对日地系统平衡点轨道使用多次的脉冲来建立一种“松”控制的轨道保持策略,将探测器保持在标称轨道附近.

随着平衡点探测研究的不断深入以及对于地月系统平衡点探测成为热点,平衡点轨道保持的研究进入第3 个阶段,即星历模型条件下的轨道保持研究阶段.第2 阶段的研究中,研究的对象通常是日地系统平衡点附近轨道上的探测器,或者是忽略了太阳引力的地月平衡点探测器.对应的目标轨道由圆型或者椭圆型限制性三体模型生成,对于日地系统(或者忽略第四体引力的地月系统),模型较为合理,但对于考虑太阳引力影响的地月系统误差较大.相应的轨道保持方案将不能胜任地月系统的平衡点任务.例如,在星历模型中,X轴速度控制法在穿越X轴时X方向速度为零的假设将不成立.

第3 阶段的研究延续了第2 阶段中对于工程实用性以及导航测量的考虑,但仅要求探测器的实际轨迹在标称轨迹的附近运行,并不要求其严格按照标称轨迹飞行.Pavlak 等[50]在靶点法的基础上提出了基准轨迹靶点法,轨道保持仅要求探测器能朝向下一目标靶点的位置靠近即可.Folta 等[51]针对ARTEMIS任务,从工程实际的角度对Floquet 法、最优轨迹延续法和全局最优搜索法进行了验证,并指出最优轨迹延续法是工程上可实施控制的最优方法,该方法不严格依赖于标称轨道,但需要进行轨迹优化,计算量较大.

2.2.4 拟周期轨道的参数化求解方法

由于拟周期轨道不满足周期性条件,拟周期轨道族的求解比周期轨道族复杂.拟周期轨道上的点经过一个周期的积分后不会返回同一点,而是返回与初始点相差固定相位角的另外一点,使得微分修正法中目标函数构造成为难题.因此,除了文献[8]中两级微分修正算法,Lissajous 拟周期轨道的求解一直没有实质性的突破.直至21 世纪初,G´omez 等[52]首次提出了Lissajous 轨道的参数化求解方法,利用关于两个不共振频率的Fourier 级数来参数化表示Lissajous轨道,从而可以利用数值延拓法,通过延拓Fourier 级数表达式的系数来求解Lissajous 轨道族.

参数化方法的思路可将Lissajous 轨道表示为两个频率(ξ,η)的函数ψ(ξ,η)

式中,ω1和ω2为两个周期性运动的频率,φt为CRTBP 微分方程的映射流,φt[ψ(ξ,η)]则为ψ(ξ,η)在映射φt下经过时间t的象.

通过引入频闪映射,也即时间固定的映射φT,可以将Lissajous 轨道的求解问题简化为其频闪映射——不变曲线的求解问题.通过固定两个频率中的任意一个(如令η=0),得到不变曲线φ(ξ)=ψ(ξ,0),其满足不变方程

式中,T=2π/ω2为由频率ω2确定的返回时间,ρ=2πω1/ω2为旋转数.

将不变曲线φ(ξ)表示成Fourier 级数的形式

式中,Fourier 系数C0∈R6,Ck∈R6,Sk∈R6为待求解变量.

经过以上处理,Lissajous 轨道的求解问题就转化为其频闪映射下不变曲线的Fourier 系数的求解问题.通过对ξ 在[0,2π] 之间进行离散,从而得到由不变方程(9)构成的方程组,然后利用类牛顿法来迭代求解[53].

基于频闪映射的参数化方法可以用于数值求解Lissajous 轨道并进行延拓,且具有较强的普适性,可以很容易扩展到高维(频率≥ 3) 的拟周期轨道求解,也可以应用于其他系统中周期轨道的求解.但是该方法每一次牛顿迭代都需要求解6(2Nf+1) ×[6(2Nf+1)-1] 个方程,运算量相对较大.Kolemen等[54]基于庞加莱映射进行降维处理,然后利用参数化方法来求解不变曲线,但是受限于拟周期轨道的几何形状,需要谨慎选取庞加莱截面,以免出现庞加莱截面上的曲线不闭合的情况.

2.2.5 椭圆型限制性三体问题中的周期与拟周期轨道

以上的研究大都是基于CRTBP 模型展开的,然而在实际的三体系统中,小主天体的轨道通常为绕大主天体的椭圆轨道,轨道偏心率将对三体系统内航天器的运动产生不可忽略的影响.如何在更精确的动力学模型下,比如ERTBP 或者星历模型,高效并快速地求解周期轨道或者拟周期轨道,成为国内外的热点研究问题.

由于小主天体轨道偏心率的周期性扰动,CRTBP模型下的周期轨道将无法一直满足其周期性条件,只存在周期为小主天体轨道的周期整数倍的周期轨道,也即共振轨道,因此也不再存在连续的周期轨道族.Campagnola 等[55]将CRTBP 模型下Halo 轨道在ERTBP 下对应的周期轨道称为椭圆Halo 轨道(elliptic Halo orbit),通过以Halo 轨道为初值并基于偏心率的数值延拓来求解.Peng 等[56]将ERTBP 下的周期轨道命名为多圈Halo 轨道(multi-revolution Halo orbit),采用了弧长延拓法来求解多圈Halo 轨道,有效解决了基于偏心率单参数延拓出现奇点的问题.Ferrari等[57]则研究了ERTBP 模型下与质量较小的主天体轨道周期呈整数倍的周期轨道,通过以CRTBP 下的共线平衡点L1,L2和L3附近的Halo 轨道、Lyapunov轨道和垂直Lyapunov 轨道为初值,沿着偏心率进行延拓求解.图7 给出了ERTBP 模型下的多圈周期轨道及其对应于CRTBP 模型下的周期轨道[57].

图7 ERTBP 下的多圈周期轨道Fig.7 Multi-revolutions periodic orbits in the ERTBP

2.2.6 特殊类型的局部周期运动

除了以上介绍的几类常见的平衡点周期轨道,如水平Lyapunov 轨道、垂直Lyapunov 轨道、Halo 轨道和Lissajous 轨道,通过数值延拓和微分修正,学者们还发现了几类独特的、具有特殊应用前景的三体系统内局部周期运动,包括近垂直Halo 轨道(nearrectilinear Halo orbit,NRHO)、蝴蝶形轨道(butterfly orbit)[58]和大幅值逆行轨道(distant retrograde orbit,DRO)[59].其中,NRHO 作为NASA 未来月球门户的参考任务轨道,最近几年得到了一系列深入的研究,包括以其为出发或到达轨道的转移轨道设计、从NRHO 逃逸或者在轨航天器的抛弃轨道设计等等.DRO 也在转移轨道设计中得到了关注[60],由于其可以延展到地球附近,并与Lyapunov 轨道相切,可以与后者拼接实现地月空间的转移.Butterfly 轨道的形状较为特殊,同时途径L1点、L2点和质量较小的主天体,且轨道在空间中有较大的弯折.几类特殊的平衡点周期轨道如图8 所示.

图8 几类特殊的平衡点周期轨道示意图Fig.8 Several types of special periodic orbits around libration points

从平衡点轨道研究的发展历程可以看出,其经历了从解析法到数值法、从单个轨道到整个轨道族、从简单的CRTBP 到更复杂的ERTBP 或四体模型、从周期轨道到多个频率耦合的拟周期轨道的过程.随着周期轨道种类与数量的增加,对三体系统相空间的认知也在加深,可以从平衡点的小范围邻域扩展到了有界运动的边界.从已开展的平衡点任务中可以看出,实际任务实施中受到测控精度、空间复杂摄动干扰和模型不确定性的影响,理想的周期运动并不存在或无法实现,拟周期运动更具有实际应用价值.在未来的研究中,高精度的动力学模型下的拟周期轨道设计与特性分析将成为研究的难点.如何提出更高效的算法求解,并探索拟周期轨道在航天任务中的应用潜能,已成为亟待解决的问题.

3 三体系统中的全局周期运动

除平衡点附近轨道外,三体系统下还存在着全局周期轨道,由于不存在类似平衡点的解析初值,只能采用纯数值的方法进行求解.由于CRTBP 问题状态空间的维度为六维,采用遍历搜索方法所需的计算量时间和存储空间巨大.因此,学者们多利用CRTBP系统的对称性来研究周期轨道,通过选择对称轴或者对称面,可将需要遍历的状态空间减半,且积分时间也可减半,从而极大的提高了计算效率.Russell 在三体问题全局轨道的遍历求解中做出了重要的贡献,在文献[61]中对木星-木卫二CRTBP 模型下木卫二附近的对称周期轨道进行了遍历求解,利用微分修正法和网格搜索法得到了大量的周期轨道,并进行了稳定性分析,部分轨道的示意图如图9 所示[61].

图9 木星-木卫二系统下的全局周期轨道Fig.9 Global periodic orbits in the Jupiter-Europa system

在三体系统全局周期轨道中,存在几类特殊的轨道,包括共振轨道、循环轨道、自由返回轨道等,他们具有较好的稳定性或途径主天体,因此在航天任务设计中得到广泛的应用.

3.1 共振轨道

共振轨道(resonance orbit)指轨道周期与主天体的轨道周期呈整数比关系的轨道.如果定义B 与A形成p:q的共振轨道,那么B 绕中心天体飞行p圈所需的时间与A 绕中心天体飞行q圈相同.靠近主天体的共振轨道通常具有较好的稳定性,受扰后不会发散或与天体发生碰撞,是探测任务理想的任务轨道.Escribano 和Vaquero[62-63]从动力学结构方面,对地月系统共振轨道进行了深入分析,基于庞加莱截面求解了平面以及空间的共振轨道,并将共振轨道与不变流形相结合,提出了低能量转移轨道的设计方法.地月三体系统下的平面以及空间共振轨道如图10 所示[63].此后,Antoniadou 等[64]研究了地外行星系统中共振比为4:3,3:2,5:2,3:1 及4:1 的平面及空间共振轨道并讨论了轨道稳定性,发现稳定共振轨道对应的偏心率较大.Anderson 等[65]进一步利用网格搜索方法对木卫系统中共振轨道进行了研究,利用连续法在稳定共振轨道的计算得到不稳定的共振轨道.

图10 地月三体系统下的共振轨道示意图Fig.10 Resonance orbits in the Earth-Moon three-body system

3.2 循环轨道

循环轨道(cycler trajectory)为连接两颗天体的自然周期轨道,可以看作是一种特殊的共振轨道.通过循环轨道可以实现天体之间的连续往返,然而无法在其中任意一颗天体上逗留,典型的地月循环轨道如图11 所示[66].Aldrin[67]于1985 年首次提出了地月循环轨道的概念,为载人登月、开采月球资源等任务提供了一种经济可行的方案,该方案通常包括两类飞行器,第一类转移飞行器利用循环轨道实现地月往返,另一类运输飞行器则负责在转移飞行器与月球表面之间输送人员和物资.由于转移飞行器无需被地球或者月球捕获,并能多次重复使用,从而能显著的节省入轨燃料.随着月球探测任务的增加和“地月经济带”等概念的提出,地月循环轨道有望得到实际的应用.但目前针对该类轨道的研究较少且模型仅基于圆型限制性三体问题,考虑月球轨道偏心率和其他摄动因素,研究地月循环轨道的稳定性,转移周期和与近地/近月飞行器的交会窗口将更有现实意义.

图11 地月循环轨道Fig.11 Earth-Moon cycler trajectory

3.3 自由返回轨道

自由返回轨道(free return orbit)则是另外一类特殊的周期轨道,主要应用于载人任务的备用轨道.自由返回轨道能保证在近月点变轨时点火失败的情况下,即使不进行月球捕获机动,也可经自由返回轨道返回地球附近.图12 为地月自由返回轨道的示意图.阿波罗13 号探测器利用了自由返回轨道,在发生动力故障后安全返回地球[68].

图12 地月自由返回轨道Fig.12 Earth-Moon free return orbit

早期的关于地月转移轨道的研究大多是在双二体模型下展开的.黄文德等[69]给出了自由返回轨道的设计方法,并对影响返回轨道的主要参数,比如飞行时间、近月距、轨道倾角等进行了分析.彭祺擘等[70]则将优化算法引入自由返回轨道的设计过程中,大大加快了轨道的搜索速度.然而,以上文献都没有考虑载人登月中固定近月点高度的工程约束.王丹阳和邓辉[71]针对该问题,提出了固定近月点高度的搜索策略,通过调节其他约束参数来设计更符合工程实际需求的自由返回轨道.候锡云等[72]则在限制性三体模型的框架下,研究了白道面内自由转移轨道的性质,并分析了近月段与近地段轨道倾角的变化规律.

由于人类对三体系统认知仍不全面,特别是引力常数对周期运动存在性和稳定性的影响尚不清楚,因此,除以上几类周期轨道外,三体系统中可能仍有大量具有良好性能的全局周期运动未被发现.相比于平衡点附近的局部运动,对三体系统内全局周期运动的认识仍是一个漫长的过程,一方面可以结合机器学习和人工智能[73],开发更高效和高精度的轨道搜索算法,另一方面,采用先进的数学工具(如辛几何、李群[74]等)重新认识三体系统,探寻其中周期运动的共同特征和特殊描述方式,可能为三体系统全局周期运动的研究带来新的思路.同时,如何将工程约束、复杂的动力学环境与特殊轨道的应用潜能结合起来,提出快速高效的参考轨道设计方法,也是未来的研究中需要思考与解决的问题.

4 三体系统中的转移轨道设计

针对文中讨论三体系统中的周期运动,实现不同周期运动间的转移也是三体系统轨道动力学研究的重点.根据不同的初始与终端状态,转移轨道包括主天体至平衡点的转移、同一系统不同平衡点轨道间的转移、不同系统平衡点轨道之间的转移等.其中主天体至平衡点的转移又包括较大主天体至平衡点的转移(如地球至地月平衡点转移)和较小主天体至平衡点的转移(如地球至日地平衡点转移).从设计思路及其所依赖的动力学模型来讲,转移轨道又可以分为直接转移、基于不变流形理论的低能量转移、考虑主天体引力影响的弱稳定边界弹道转移等.下面将分别针对这3 种设计思路进行总结.

4.1 直接转移方式

直接转移方式与传统二体模型下的双脉冲转移方式类似,以地月系统下地球至地月L1/L2点转移为例,探测器的第一次为月球转移脉冲,从地球停泊轨道离轨进入地球转移段;第2 次脉冲为月球附近轨道(比如Halo 轨道)的入轨脉冲.直接双脉冲转移如图13 所示[75],目标轨道分别为地月L1或者L2点附近的Halo 轨道.Parker 等[75]对地月系统下的直接转移方式进行详尽的研究与分析,并指出直接转移所消耗的时间通常较短,约5 d～2 m 之间.针对L1和L2的目标Halo 轨道,月球转移脉冲对应的能量分别为-2.6～-2.0 km2/s2和-2.6～-1.0 km2/s2;而第2个Halo 入轨脉冲的大小通常为500～600 m/s 之间,最大与最小值则分别约为950 m/s 与430 m/s.

图13 地月L1/L2 Halo 轨道直接转移示意图Fig.13 Direct transfer trajectory to the Earth-Moon L1/L2 Halo orbit

4.2 基于不变流形的低能量转移轨道设计

随着对平衡点轨道认知的深入,学者们开始将平衡点轨道的不变流形引入到低能量转移轨道的设计中.20 世纪90 年代,G´omez 等[76]首次将动力系统理论应用到CRTBP 的动力学分析中,发现Halo 轨道的稳定流形能够延展到月球,并利用这类自然的轨道实现了从地球到地月L1点Halo 轨道的转移轨道设计.研究发现基于不变流形的低能量转移轨道所消耗的能量要远小于直接转移方式.图14 给出了典型的地月低能量转移的示意图.低能量转移轨道沿着Halo 轨道的流形方向切入Halo 轨道,因此无需施第2 次入轨脉冲.

图14 地月L1/L2 Halo 轨道低能量转移示意图Fig.14 Low-energy transfer trajectory to Earth-Moon L1/L2 Halo orbits

在此之后,很多学者开始广泛研究不变流形理论在低能量轨道转移中的应用.Barden 等[77]展开了一系列的研究,包括利用稳定流形实现地球至Halo轨道的转移,利用不稳定流形实现Halo 轨道返回至地球的转移,不同Halo 轨道之间的转移等.Born 和Parker[78-79]则全面研究了地月低能量转移的问题,考虑了不同的近地出发轨道,不同的环月目标轨道,以及不同的转移方式,并从发射能量C3、转移时间、转移总脉冲大小等方面对以上各种情况进行了比较与分析,并将不同路线下低能量地月转移的各种性能指标总结如图15 所示.

图15 不同路线的低能量地月转移概览图Fig.15 Flowchart of different low-energy lunar transfer itineraries

Parker 等[75]在研究中发现,当Halo 轨道的稳定流形经由月球附近时,转移所需消耗的脉冲更少,但是转移时间会增长.Gordon[80]对这一现象进行了深入的研究,将不变流形与月球近旁飞越相结合,设计了地球至地月L2点Halo 轨道的双脉冲转移轨道.然而以上研究中只考虑了沿着Halo 轨道的自然流形进入Halo 轨道的情况,转移时间通常较长.Li 等[81]则提出了基于月球飞越的三脉冲转移方式,通过增加Halo 轨道入轨脉冲这一额外的机动,可极大缩减在流形上滑行的时间.Zeng 等[82]对文献[81]提出的方法进行了深入的研究,详细分析了不同的借力高度、飞行角以及Halo 轨道入轨点等对转移总脉冲的影响.图16 为考虑月球旁近飞越的低能量转移方式示意图.Cheng 等[83]研究了利用流形实现平衡点轨道向其他轨道(环月轨道) 的转移问题,考虑了月球自转轴与白道面之间的夹角,提出了适用于任意倾角、任意高度环月轨道的双脉冲通用转移方案.

图16 基于月球旁近的地月L2 点Halo 轨道低能量转移示意图Fig.16 Low-energy transfer trajectory to Earth-Moon L2 Halo via lunar flyby

另外,学者们还利用不变流形研究了平衡点轨道间的转移,包括同一平衡点之间的转移(同宿连接)以及不同平衡点之间的转移(异宿连接).G´omez等[84]首先考虑能量相等的两个Halo 轨道,发现利用异宿连接可以实现零燃耗转移;此后,G´omez 等[85]针对能量不相等的情况,又提出了基于双脉冲转移方法,也即分别沿着出发Halo 轨道的切向和目标Halo轨道稳定流形的方向施加一次机动.Davis 等[86]进一步研究了不同能量的平衡点轨道之间的转移问题,通过拼接出发Halo 轨道的不稳定流形与目标Halo轨道的稳定流形,设计了双脉冲低能量转移轨道,并利用主矢量法优化了施加脉冲的位置和大小[87].

不变流形理论还可以用于两个或多个三体系统之间的转移[88].Howell 等[89]研究了出发轨道为月球轨道或者地月L1/L2Halo 轨道,目标轨道为日-地/月系统L2Halo 轨道的转移轨道设计问题,并考虑了不变流形以及穿越轨道两种不同的转移思路.Kakoi等[90]利用地-月系统L1/L2Halo 轨道的不变流形,研究了向日-火系统L1/L2Halo 轨道转移的低能量设计问题.Peng 等[91]则研究了主天体轨道偏心率对低能量转移轨道的影响,以日地ERTBP 和地月ERTBP 拼接双椭圆三体问题模型为动力学模型,研究了从日地平衡点轨道至地月ERTBP 模型下的多圈Halo 轨道的转移轨道.

除了Halo 轨道以外,近年来受到美国重返月球计划的驱动,NASA 提出了以NRHO 为参考轨道的月球门户任务.以NASA 与普渡大学的Howell 团队为主的学者们,对NRHO 在低能量转移轨道中的应用进行了一系列全面的探索,包括NRHO 分别作为出发轨道、中转轨道或目标轨道的转移轨道设计问题[92-94],以及航天器从NRHO 逃逸或者废弃处理等任务.未来对低能量转移的研究将会向着实际工程任务发展,结合探测器的可达性[95],如何处理复杂工程约束,探索更多类型的平衡点轨道在空间探测任务中的应用前景,必将成为未来研究的趋势.

4.3 三体/四体系统的弱稳定边界理论与低能量捕获轨道设计

在三体系统中,除了利用不变流形的低能量转移方式外,还存在一种利用弱稳定边界理论的低能量转移.Belbruno[96]在研究地月低能转移轨道时发现,若转移过程中探测器相对地球的远地点距离远大于月球轨道半径时,探测器可以利用太阳扰动力实现转移,到达月球附近无需施加捕获机动即可被月球临时捕获,形成月球环绕轨道.该现象反映在月球附近存在一个区域,区域内的运动经多天体的引力作用可实现自然的逃逸或捕获.Belbruno[97]将该边界定义为弱稳定边界(WSB).此后学者们从理论角度对三体系统中较小天体附近的运动稳定性进行了详细的研究,并发展成为弱稳定边界理论.利用弱稳定边界理论可以设计低能量的弹道捕获轨道,相比传统的脉冲捕获形式大大节省了燃料消耗.

4.3.1 弱稳定边界理论的发展

弱稳定边界最初基于限制性三体问题的模型确定,表示质量较小天体附近二体运动稳定与逃逸的边界.随后,Belbruno[97]给出了更准确的数学描述,如图17 所示.

图17 弱稳定边界的定义Fig.17 Definition of the weak stability boundary

弱稳定边界将三体系统中平衡点附近的研究拓展到主星附近,并利用二体运动的参数对运动进行表述,进而建立了主天体捕获与逃逸运动的联系.同时,尽管弱稳定边界在三体问题模型下定义,但其运动稳定性分析并不受限于三体动力学.采用四体模型甚至考虑其他摄动力的动力学模型同样可以得到相应的弱稳定边界[98],进而使得弱稳定边界的研究可由三体系统拓展到四体和多体系统中.基于限制性三体与四体系统模型,地-月、日-火系统的弱稳定边界结构和运动特性被详细研究[99-100],其中日-火系统的弱稳定边界如图18 所示.

图18 日-火系统的弱稳定边界结构Fig.18 Geometry structure of the weak stability boundary of the Sun-Mars system

此外,Garcia 等[11]发现质量较小主天体附近的稳定区域并非为以天体为中心的闭合区域,随着探测器在射线上初始半径的增大,对应的运动稳定区域和不稳定区域交替出现.因而Garcia 和G´omez将运动稳定的区域定义为稳定集.稳定集是若干个闭合的区域,其内部可能不包含主天体,稳定集可看作弱稳定边界理论的拓展与精细化描述.此后,Topputo 等[101]改进了稳定集的求解算法并应用于日-木系统附近的稳定区域计算.利用圆型限制性模型和椭圆型限制性模型,地-月系统、日-木系统、木星-木卫系统、日-水系统、日-火系统的稳定集结构相继被研究[102-104].逆向稳定集的概念也被提出[105],用于描述探测器沿时间逆向积分对应的运动稳定性,对稳定集结构的研究,为低能量的捕获轨道设计奠定了重要的基础.

弱稳定边界和稳定集与三体系统内不变流形的联系是弱稳定边界理论研究的另一热点.文献[99]指出,利用Jacobi 常数C(L1)的等势面可以非常粗略的估计弱稳定边界的范围.Fantino 等[105]定义了固定Jacobi 常数时的弱稳定边界,发现相同Jacobi 常数下L1和L2附近周期轨道的稳定流形上速度和位置垂直的点和弱稳定边界基本一致.通过几何证明和数值仿真,验证在一定Jacobi 常数范围内,流形上径向速度为零的点就是弱稳定边界上的点[12],为统一弱稳定边界理论和不变流形理论做出了重要的贡献.

4.3.2 弱稳定边界理论的应用

弱稳定边界理论的重要应用是设计低能量的捕获轨道.日本的Hiten 卫星是第一个成功应用利用弱稳定边界理论实现月球弹道捕获的任务[106].由于发动机故障,Hiten 卫星采用传统的霍曼转移方式无法实现月球捕获,Belbruno 利用月球借力增大探测器的远地点,并利用太阳摄动力使探测器进入月球的弱稳定边界,实现临时捕获.计算发现相比霍曼转移可以节省15%的速度增量.

此后针对地-月系统的弹道捕获被深入研究[107-108].Sweetser[109]总结了利用不同方法从地球停泊轨道转移至环月圆极轨所需的速度脉冲,并发现弹道捕获所需的速度增量最小,且接近CRTBP 模型下理论计算的最小值.弹道捕获轨道设计的关键是构造恰当的日-地-月构型使得太阳引力可以有效提高轨道的近地点.Koon 等[110]创造性的采用流形拼接的方法构造了地-月系统弹道捕获轨道,解释了弹道捕获转移轨道的原理,使得月球弹道捕获轨道可以通过双三体模型研究来进行初步设计,提高了转移轨道的设计效率.

此外,针对其他天体系统的弱稳定边界转移与弹道捕获轨道也进行了广泛的研究.Hyeraci 等[103]提出捕获集的概念.Circi 等[111]分析了行星公转轨道的偏心率对捕获集的影响,发现偏心率越大越有利于生成捕获轨道.在日-火系统的弹道捕获轨道的研究中,学者们发现由于地-火转移对应的火星交会速度较大,较难得到弹道捕获的机会.一种方式是采用小推力轨道与相应的火星弱稳定边界状态匹配[112];另一种策略是将探测器在距离火星较远的星际空间捕获进入弹道轨道,然后自然的进入火星弱稳定边界,与霍曼转移结果相对比,Topputo 等[113]发现该方法在近火点高度较高的捕获任务中有优势且增加了火星转移的灵活性.Li 等[21,114]利用逆向稳定集的思想,将弱稳定边界轨道应用于火星高轨道的捕获中,作为星际转移轨道和任务轨道间的过渡轨道,尽管捕获过程中需要施加两次机动,但相比直接捕获可以明显降低速度增量.

除Hiten 卫星外,NASA 的GRAIL 月球探测器也采用了弹道捕获方式进入月球轨道[115-116],如图19所示.Bepi-Colombo 水星探测器同样计划借助水星的弱稳定边界逐渐降低探测器的能量[24],实现低能量的弹道捕获.弱稳定边界因其较小的燃料消耗将在未来深空探测中有更广泛的应用前景,但其主要劣势在于转移时间较长.深入了解弱稳定边界内的运动机理,在不显著增大燃料消耗的情况下通过施加人为扰动,改变弹道捕获轨道的运动轨迹,缩短转移时间将是一个值得深入研究的动力学问题.

图19 GRAIL 地月转移轨道Fig.19 Translunar trajectory for GRAIL

5 三体系统轨道动力学的应用

随着对三体系统轨道动力学研究的深入和对运动特性了解,学者们也将三体系统轨道动力学与航天动力学中其他技术和概念相结合,将三体系统的运动特性应用于相关问题的研究中.较有代表性的应用包括三体系统下的编队飞行与导航定位星座设计.

5.1 三体系统中的编队飞行

在空间观测任务中,为了提高观测分辨率,需要增大光学设备的基线长度,但单个探测器尺寸受到运载器尺寸和运载能力的限制存在极限,无法实现几十米甚至上百米的基线长度,因此,编队飞行的概念应运而生.通过多个探测器组成飞行阵列,探测器按一定的构型保持相对位置和速度,从而形成一个虚拟的探测器整体,通过协同控制共同实现观测任务,可以大大增加探测器的观测能力.编队飞行早期主要针对近地轨道,随着三体系统动力学研究的深入和平衡点附近轨道独特的观测条件,针对三体系统平衡点附近的编队动力学研究也成为近年来的研究热点之一.

日地平衡点距离地球较远,受到地球的电磁干扰较少,因此适合开展对深空天体或太阳的高分辨率观测.美国国家航空航天局(NASA)、欧空局(ESA)先后启动了类地行星探测器(terrestrial planet finder,TPF)计划、Darwin 计划、MAXIM(X 射线干涉测量)任务等空间观测计划.中国也提出了觅音计划,计划在日地平衡点开展干涉测量,而三体系统中的编队飞行动力学研究是任务开展的基础.

5.1.1 平衡点编队轨道动力学

与近地轨道的编队飞行不同,平衡点轨道编队飞行所面临的动力学环境更加复杂.轨道编队在平衡点附近的运动演化是研究的关键问题.基于圆型限制性三体模型,H´eritier 等[117]研究了探测器编队的自然演化过程和稳定性,发现线性化无控自然编队的相对距离将逐渐发散.H´eritier 和Howell 定义了编队误差较小的稳定域.探测器在稳定域内可以形成松散编队,相对位置将在一定范围内波动但不会发散.研究发现,Halo 轨道的稳定域为二次曲面.日地L2点大尺度Lissajous 轨道下的低自然漂移区域也被研究[118],可实现6.0×105km 振幅的Lissajous轨道100 d 内的误差<5.0×104km.以上研究给星间相对位置精度要求较低的编队提供了构型设计依据.此外,利用平衡点邻域自然存在的中心流形[119-120],可以将编队设计在一定幅度的周期或拟周期轨道上,各探测器间的距离在一定范围内波动,可以形成松散的卫星编队.

针对平衡点高精度编队的动力学问题,Out[121]研究了日地L2Halo 轨道上100 m 量级的双星编队,分别考虑了圆型和椭圆型限制性三体问题、椭圆限制性四体问题、考虑太阳光压的四体问题下的编队演化情况.研究发现,太阳光压是影响最大的扰动因素.太阳光压受航天器的面质比和表面材质属性和探测器姿态影响,如探测器相对于太阳光压的法向旋转1°,将引起1µm/s2的相对加速度变化,是引力相对加速度的1000 倍.在高精度的星历模型下的TPF 任务仿真表明[122],高精度平衡点编队对速度增量需求较小,但对推进器的控制精度要求高,且轨道维持频繁,同时也对相对测量精度提出了较高的要求.

5.1.2 编队轨控策略

由于环境摄动力和测控与执行误差的存在,想要实现高精度的平衡点编队,编队控制必不可少,即形成非自然编队.Howell 和Marchand 等[123-125]对自然编队和非自然编队做了比较,并对非自然编队的控制策略给出概述.目前针对平衡点编队的控制策略主要分为两类,一类是编队的每颗星分别进行标称轨道的跟踪控制,另一类是编队主星进行周期轨道保持,子星只考虑与主星的相对距离控制.仿真分析表明采用相对控制所需的燃料消耗较省.基于线性化模型的LQR 方法被用于编队保持[126],但是由于线性模型为实际模型在平衡点附近的近似,因此此类基于线性模型设计的控制算法难以到达较高的精度要求.非线性自适应控制器[127]和自抗扰控制[128](active disturbance rejection control,ADRC)等非线性编队飞行控制方法被应用于编队控制中,在理论模型下可以得到较好的位置保持精度,误差优于0.5 mm.此外,Infeld 等[129]还将卫星编队问题处理为多智能体、非线性、有约束的最优控制问题,同时解决了轨道设计和编队控制.但以上控制方法的研究采用的动力学较简单,且未考虑测量与控制误差,实际的编队控制情况仍有待验证.

从以上研究可以看出目前针对三体平衡点编队飞行的动力学与控制研究相比单航天器的研究还较为基础,三体系统中的相对运动规律仍未探明,影响相对运动稳定性的关键参数仍不清晰,同时多星编队的姿轨高精度测量与协同控制也有待解决.三体系统下的多星编队具有潜在的应用价值和空间探测需求,相关任务也仍在研究论证中,解决上述问题将为任务的精确实施奠定基础.

5.2 三体系统中的导航与定位问题

5.2.1 三体系统导航方法研究

截至目前,所有在轨的平衡点探测器均采用地面测控站测定轨的导航方式来实现.采用两次轨控之间所有的数据或部分数据进行交叉定轨,定轨时间段较长,定轨的精度在千米量级,外推的精度在10 km 量级.典型的有:(1)ISEE-3 使用地基的S 波段雷达进行测控,测量的弧长不固定,一般由多个5 min左右的短弧段构成;(2)ACE 和SOHO 都主要使用26 m 的美国深空网(deep space network,DSN)天线数据;(3) MAP,WIND,Genesis 以及2010 年成功到达地月平衡点的ARTEMIS 都使用34 m 的DSN 天线数据为主,辅助以34 m 和70 m 天线数据或者26 m 的天线数据;(4)HERSCHEL 和PLANCK 等欧空局的飞行器则是使用欧空局deep space antenna-2,DSA-2 天线数据.随着平衡点任务的增多,地面测控站将难以兼顾所有的测定轨任务.此外,地基测控导航本身也存在着时延以及精度不够高的缺陷[130].对于需要频繁机动的平衡点飞行器,研究能够提高导航能力和在轨生存能力的自主导航方式是重要的任务.

尽管国内外的学者已经根据不同测量方式提出了几种适用于行星际飞行的自主导航方法,包括基于太阳和地球方位矢量、基于相对测量信息、天体方位测量与多普勒雷达组合等.但以上自主定轨在三体系统下的应用还较少.2005 年Hill 等[131]提出一种称为LiAISON (linked autonomous interplanetary satellite orbit navigation) 的导航概念,仅利用星间伪距测量信息即可进行卫星的绝对导航[132],其可行性就建立在平衡点轨道的非对称性和局部唯一性之上.通过利用交联测距(crosslink range),Hill 等[133]还提出了平衡点自主导航方法,实时在轨测量精度可达到10 m 量级,大大提高了三体系统下轨道的定轨精度.但当前针对三体系统的自主导航研究主要针对平衡点轨道卫星,且采用的动力学模型较为简单,当考虑高精度的星历模型和摄动因素后,探测器的运动将变得复杂,星间相对测量数据的解算将相应变得复杂,且可能需要引入其他观测量,如何在降低系统复杂度的同时利用三体动力学特性提高定轨精度仍需进一步探讨.

5.2.2 三体系统导航星座设计

除了三体系统内自主导航方法的研究,学者们还研究了利用三体系统特性的导航星座设计,特别是服务于月球探测的中继与导航星座设计.由于月球公转与自转同步,月球背面无法与地面直接通信,在20 世纪60 年代末Apollo 计划期间,Farquhar 等[134]即提出利用地月平衡点L2点Halo 轨道实现地球和月球背面通信的方案.此后William 等[135]提出了通过载人进入地月L2点轨道,对月背表面探测器进行遥操作的月球背面探测方案,并分析了不同平衡点轨道的月球覆盖率,发现小振幅的北向Halo 轨道对月球南极附近着陆点的覆盖性较好.中国的嫦娥4 号中继星是首个在地月L2点Halo 轨道实现月背中继通信的通讯卫星.结合工程约束,学者们在日-地-月系统下开展了任务轨道选择[136-137],月背通讯覆盖率[138]等方面的研究,发现月背通讯覆盖与轨道类型和轨道振幅密切相关,振幅较小的周期轨道对月面整体的覆盖性较好而较大振幅的任务轨道对月球极地的覆盖较好.

此外,为实现未来载人登月和月球基地的建设,需连续无间断和有冗余的地月通讯,甚至实现类似GPS 功能实时定位能力.相比月球低轨导航星座,放置在平衡点轨道的导航星座在相同的覆盖性下所需的卫星数量较少,且导航卫星受月影影响小,可以长期与地球保持通信.轨道上的光照时长也比月球低轨要长.Zhang 等[139-141]提出了地月L1/L2点的双星星座方案,并分析了导航性能.而根据伪距方程,至少需要采用4 颗卫星实现对月面的导航定位.Daniele等[142]提出在L1和L2点附近的轨道周期约为1 年的Lissajous 轨道上各放置4 颗卫星的星座方案,且XY平面和Z方向上的相角平均分布.Ren 等[143]在此基础上讨论了该星座方案的定位精度因子.研究发现,由于Lissajous 轨道星座靠近平衡点,从月面观测的相角接近,定位精度因子较大,对应的定位精度较差.Ren 等还提出了其他平衡点轨道组成的最小星座构型,包括Halo 构型和水平和垂直Lyapunov 构型.研究同时发现,在真实星历下,受太阳扰动的作用Halo轨道和Lissajous 轨道变成了不变环面,理想模型下的全覆盖星座仅能覆盖83.25%的区域,为保证可见性,需要在任两个轨道上各添加一颗卫星,星座的卫星个数为14 颗,但仍小于Walker 星座构型(18 颗).

从以上分析可以看出,当前针对地月三体系统的导航星座研究已较为成熟,多种星座构形被提出并分析了其观测性能.同时嫦娥4 号中继星任务的成功实施也为中继星座的轨道设计与轨道维持提供了宝贵的经验.随着月球表面和附近探测活动的增多,未来地月三体系统下的导航星座必将成为现实.针对导航星座的研究也应从理论化的方案设计转向考虑更多工程约束的任务分析,如星座系统的发射规模分析,星座维护成本评估,以及冗余备份设计等,同时真实星历模型下多星的精确部署与构形初始化是保证导航星座性能的关键,也是三体系统轨道动力学应用研究中值得深入探讨的问题.

6 总结与展望

随着人类航天任务由近地走向深空,三体系统轨道动力学在任务设计中将扮演越来越重要的角色.例如,美国计划在地月L2点近垂直Halo 轨道上建造“深空之门”空间站,用于长期低成本的地月往返任务.我国也计划在日地L2点开展空间探测活动.同时,随着商业航天的迅速发展,利用低成本小尺寸的立方星开展太空活动成为未来发展的趋势,基于微小卫星进行月球和火星探测也已被验证.低成本深空探测发展趋势给三体系统的轨道动力学提出了新的命题和挑战.结合未来的月球和深空探测活动,未来三体系统的轨道动力学方面的研究可在以下3 个方向展开.

6.1 地月高效可重复利用的转移与全月面到达轨道设计

未来的月球探测需满足全月面覆盖、高效且可重复利用的要求.利用三体系统的不变流形可以从平衡点轨道实现弹道着陆,但自然流形所能覆盖的月球表面区域有限.地月三体系统存在多种转移方式,但存在转移时间和燃料消耗间的矛盾.如何设计着陆轨道,实现尽可能少的燃料消耗下全球到达？如何实现货物运输与人员往返？如何针对不同的探测任务需求,快速对地月往返轨道,特别是地球与平衡点周期轨道间的往返轨道,进行设计与重规划？这些问题将是未来建立月球基地开展长期月球探测的关键,是未来三体系统轨道动力学值得系统研究的问题.

6.2 三体系统下的小推力轨道优化问题

小推力相比化学推进具有比冲高,效率高,工作时间长的优点,适合开展长期的深空探测任务.对于微小卫星,小推力推进器更具有质量和体积上的优势,将是未来重要的推进方式.然而小推力转移轨道优化始终是轨道设计的难点,经过科学家多年的研究与实践,基于惯性系下的二体与二体摄动模型的小推力转移优化设计取得较大的研究进展.但在三体动力学模型下的小推力转移优化研究有限,优化算法收敛性差、计算效率低.小推力与三体系统相结合的低能量转移轨道仍处在理论研究阶段,未得到很好的解决,也是未来值得深入探讨的问题之一.

6.3 三体系统中L4 和L5 点附近的运动特性与开发利用

相比三体系统共线平衡点,三角平衡点的稳定性更好,所需的轨道保持消耗更省,具有较大的开发应用潜力.但三角平衡点附近的动力学特性研究相比共线平衡点较少,对其附近的运动行为认知有限.三角平衡点稳定的特性使其无法通过简单的流形拼接实现低能量转移,也给转移轨道的设计提出了挑战.在日木、日火和日地的三角平衡点均发现了小天体族群.对三角平衡点的研究可以解释小行星族群的迁徙与演化,也为开展小天体探测提供理论支持.

7 结论

三体系统轨道动力学是近代轨道动力学领域研究的重点.国内外对三体系统的轨道动力学与控制问题已经进行了大量的理论研究和工程应用,构建了成熟的理论体系,提出了一系列服务于工程应用的任务轨道设计与控制方法.但随着月球与深空探测活动的日益频繁、工程约束与任务需求的增加,对三体系统轨道动力学研究提出了新的要求,分析三体系统内运动状态和运动特性还面临着诸多基础性的理论难题.因此,未来三体系统轨道动力学研究有必要结合工程需求,对三体系统轨道动力学机理等问题展开更为深入和系统的研究,为未来深空探测任务奠定理论基础.