雙小行星探測軌道動力學研究進展1)

王 悅 伏 韜 張瑞康

(北京航空航天大學宇航學院,北京 102206)

引言

雙小行星系統由在萬有引力作用下彼此環繞運行的兩顆小行星組成,在近地小行星、主帶小行星和柯伊伯帶小天體中普遍存在.據估算,大約16%的近地小行星和直徑小于10 km 的主帶小行星是雙星系統[1-2].自1993 年“伽利略”探測器發現首個雙小行星系統243 Ida 以來,大量雙小行星系統被相繼發現,引起和推動了行星科學和航天動力學等領域對雙小行星系統的關注和研究.

最早引發關注的是雙小行星系統的形成和演化理論.目前的主流理論認為,較小的雙小行星系統是由“碎石堆疊”的單個小行星在熱力學的YORP(Yarkovsky-O'Keefe-Radzievskii-Paddack effect)效應驅動下不斷加速旋轉,超過臨界轉速后分裂而形成[3-4],較大的離心力使得小行星表面高緯度地區的物質向赤道積聚,形成了主星頂部錐臺形、赤道隆起的典型形狀[5-6].分裂后雙星系統的動力學演化取決于次星和主星的質量比[7-8].當雙星質量比較接近時(質量比≥0.2),雙星系統在潮汐力的長期作用下逐步演化為雙同步系統,即兩顆小行星互相潮汐鎖定;當雙星質量懸殊時(質量比 ＜ 0.2),系統是非束縛的,可能在混沌的動力學演化中逃逸并形成小行星對(asteroid pair),亦或次星在引力矩作用下加速旋轉分裂,并最終演化為穩定的三星或雙星系統[9].而擁有較大主星(尺寸＞20 km)的雙小行星系統,通常次星尺寸較小且主星轉速低于旋轉分裂的臨界速度,這表明該類雙星系統可能形成于大型的碰撞過程[10-11].此外,部分雙小行星系統還可能通過主星引力捕獲另外一顆小行星,或在接近大天體時被潮汐力撕裂而形成[12-13].

在航天動力學方面,近二十多年來,美國和日本成功實施了多次小行星探測任務,不斷掀起小行星探測的熱潮.小行星探測對揭示太陽系起源、開展行星防御和推動小行星資源開采利用都具有極為重要的意義.與單小行星相比,雙小行星系統具有更高的研究和探測價值:其一,對雙小行星系統形成和演化過程的研究可為行星系統的形成和演化提供寶貴線索;其次,旋轉分裂形成的雙星系統為研究小行星內部物質組成和構造提供了可能;最后,雙星相對運動的微小改變更容易被觀測,是小行星防御技術驗證的絕佳對象.因而,雙小行星系統正在成為小行星探測與防御備受關注的新目標.但同時,雙小行星系統附近的復雜力學環境為探測器軌道動力學和任務設計帶來了全新的挑戰,促進了航天動力學在多個方向的進展.

對雙小行星系統不規則引力場和雙星間相對運動的準確建模是研究探測器軌道動力學問題的基礎.針對不規則天體的引力場描述,目前主要形成了三種通用方法:質點群法、球諧函數法和多面體方法,以及一些基于特定簡單幾何體的近似方法.它們具有各自的優缺點,例如質點群法簡單、靈活性高,但不具備邊界條件,隨精度要求的提高計算量急劇增加;球諧函數法具有解析形式、計算效率高,但靠近天體表面收斂速度慢甚至發散;多面體法精確、完備,但計算量大,不利于解析研究[14];基于特定幾何體的方法只適用于具有特定形狀特征的天體.因此,通常需要根據具體問題和應用場景選擇合適的方法或將不同方法混合使用.

在雙星不規則引力的相互作用下,它們的相對運動構成了復雜的姿態軌道耦合的二體問題,被稱為完全二體問題(the full two-body problem,F2 BP)[15].這種兩個任意剛體構成的二體問題是天體力學中的一個基本問題,得到了廣泛的研究[16-19].學者們發現采用相對位置和姿態坐標來描述二體的相對運動,可以降低系統自由度,并有效地簡化問題,得到了一般性的系統運動方程[16],而相對運動的具體演變由方程中的相互引力勢模型所決定.

雙星的相互引力勢本質上是將分屬兩星的兩質量元間的萬有引力勢分別對兩個剛體求一次體積分,即是一個雙重體積分.除了一些特殊情況,兩個任意質量分布剛體之間的相互引力勢沒有封閉的表達式,通常采用無窮級數展開的方法.以不規則天體引力場的建模方法為基礎,發展出了多種相互引力勢的建模方法.根據處理雙重體積分的數學思路,可以將其分為三種類型:基于簡單幾何體的模型、基于剛體質量分布參數(如慣性張量/積分、引力場球諧系數等)的模型,以及基于多面體形狀和剛體質量分布參數的混合模型.針對一些特殊的簡化情況,學者們在球-橢球、橢球-橢球和球-任意形狀的雙星模型下,得到雙星相對運動的一些解析結果,如系統穩定性和雙星相對運動的平衡構型及穩定性[8,15,20-22];在更精確的基于質量分布參數模型下,利用數值方法研究了兩顆多面體小行星的相對運動[23].研究所揭示的雙星相對運動的平衡構型與實際天文觀測數據一致,為探測器軌道動力學研究和任務設計奠定了可靠的基礎,同時也檢驗和豐富了行星系統的演化理論.

在不規則引力場建模和雙星相對運動問題得到解決之后,就可以在此基礎上研究探測器的軌道動力學.雙星和探測器構成的三體系統滿足限制性假設,因而可將探測器的運動視為一個限制性三體問題(the restricted three-body problem,R3 BP).與球形主天體的經典限制性三體問題相比,由于雙星對探測器的非球形引力和雙星之間的姿態軌道耦合動力學,該動力學系統呈現出新的豐富的動力學特性,被稱為限制性完全三體問題(the restricted full threebody problem,RF3 BP).2005 年,Scheeres 和Bellerose[24]在附加太陽引力的限制性希爾完全三體問題模型下,給出了系統的運動方程,奠定了這類問題研究的基礎.接著,學者們在同步橢球-橢球(或球-橢球)的簡化模型下,對系統的平動點及其穩定性,周期軌道及其穩定性和基于不變流形的轉移軌道等一系列問題開展了深入的研究[25-29].同時,為了更加真實地模擬雙星附近的動力學環境,學者們采用了更加精確的多面體-橢球或雙多面體模型,或一并考慮太陽光壓力和太陽引力攝動,利用數值方法對上述問題進行了研究[30-33].系統存在的多種類型的周期軌道能夠滿足雙星近距離探測的不同任務需求,但雙星附近力學環境的不確定性等因素將導致實際的飛行軌道偏離名義軌道,甚至撞擊小行星表面.因此,雙星系統附近高效的軌道維持是另一項重要的研究內容,目前主要針對軌道維持策略和實現手段展開了研究[34-36].

近年來,學者們認識到太陽光壓力在雙星系統附近的軌道動力學中也扮演著重要的角色,在精確的任務軌道設計階段需要加以考慮.研究發現在某些情況下,選擇不同的太陽光壓力模型可能會導致顯著不同的結果[37].另一方面,小行星的弱引力場為太陽帆航天器提供了非常適用的場景.與傳統的將太陽光壓力視為攝動不同,太陽帆可以利用太陽光壓的微小推力產生傳統航天器難以實現的平衡點及周期軌道(如懸浮軌道),并可通過調整太陽帆姿態實現軌道的維持和機動.2004 年,McInnes[38]廣泛討論了太陽帆在限制性三體問題中的應用.2018 年,Heiligers 和Scheeres[39]將相關理論推廣到了雙小行星系統中,驗證了懸浮軌道的存在.2019 年Guzzetti 等[36]將太陽帆推進應用于平動點軌道的維持.

上述探測器軌道動力學的研究是在限制性三體問題框架內,重點研究了平動點、周期軌道及其衍生問題,這是對球形主天體的經典限制性三體問題的自然推廣.然而事實上,在雙小行星附近的軌道動力學中,還存在對更為經典的受攝二體問題的推廣,即環繞系統單顆星的軌道動力學.以環繞雙小行星系統主星的軌道為例,在靠近主星的區域,主星的中心引力遠大于次星的引力和其他攝動力,本質上屬于受攝二體問題的范疇,其開普勒軌道要素變化的時間尺度遠大于探測器軌道周期.這時如果仍然采用限制性三體問題模型,其積分效率太低,更糟糕的是掩蓋了受攝二體問題所蘊含的清晰的軌道演化圖景,難以揭示系統固有的演化機制和特性.受攝二體問題是天體力學最為經典的問題,最初被用來研究行星和行星衛星的軌道運動.近幾十年來,二體問題的軌道攝動理論在航天動力學,尤其是人造地球衛星軌道動力學領域得到了極大的發展,形成了豐富的研究成果.2020 年,Wang 和Fu[40]進一步將二體問題的軌道攝動理論拓展到雙小行星系統中,建立了環繞主星的半解析軌道動力學模型,揭示了環繞主星軌道運動的Lidov-Kozai 機制,分析了軌道長期演化的穩定性[41].

進入21 世紀的第三個十年,對雙小行星系統的研究和探測將由理論轉為實踐.目前已有多個雙小行星探測任務被提出和批準:(1)美國和歐洲聯合提出了AIDA (asteroid impact and deflection assessment)任務,包括美國DART (double asteroid redirection test)和歐洲Hera 兩個探測器,計劃以DART 探測器(已于2021 年11 月發射)撞擊雙小行星系統65803 Didymos 的次星,以Hera 探測器(擬2024 年發射)對雙星系統的兩個成員進行詳細的勘察,并測量DART 探測器的撞擊對次星軌道的影響,評估動能撞擊在小行星防御中的可行性[42];(2) 2021 年10 月發射的美國Lucy 探測器將首次造訪特洛伊小天體,它將于2033 年飛掠最后一個探測目標—位于日木系統L5 點的特洛伊雙小行星系統617 Patroclus;(3)美國科羅拉多大學Daniel J.Scheeres 教授擔任首席科學家的Janus 探測計劃(擬2022 年發射)將對兩個近地雙小行星系統1996 FG3 和1991 VH 進行飛掠探測[43].

與前述的軌道動力學理論研究不同,面向探測任務的軌道設計還需要進一步考慮實際的任務需求和來自探測器平臺的約束條件,目前的研究主要包括兩個方面,一是基于不變流形的低能轉移和著陸/起飛軌道研究,二是針對特定探測目標和任務需求的力學環境分析和任務軌道設計,例如在DART 任務中為其攜帶的立方星設計任務軌道,實現對濺射粒子和撞擊坑的長時間、近距離觀測.

本文的主體結構安排及各部分間的關系如圖1所示.第1 節首先論述雙星引力勢建模與相對運動研究,以不規則天體的引力場模型作為起點,其是雙星系統相互引力勢建模的基礎,而相互引力勢又是研究雙星相對運動的基礎.通過雙星相對運動的研究可以確定雙星系統自身的運動狀態,是探測器軌道動力學研究的前提和基礎.探測器在雙星系統附近的軌道動力學按照理論框架、求解思路和關注點的不同,主要可分為限制性三體問題和受攝二體問題兩類,將在第2 節和第3 節分別對其研究現狀進行論述.隨后,第4 節將介紹基于探測器軌道動力學研究的理論成果,進一步考慮實際任務需求和約束的軌道動力學分析與設計.

圖1 本文結構安排及各部分之間的關系Fig.1 Organization of this paper and connections between different parts

1 雙星引力勢建模與相對運動研究

1.1 不規則天體的引力場模型

單顆小行星的引力場模型是雙星引力勢建模的基礎.不規則天體引力場的通用建模方法目前主要有三種:質點群法、球諧函數法和多面體方法,此外還有一些基于特定簡單幾何體的引力場近似方法.

1.1.1 質點群法

質點群法是描述不規則天體引力場的一種直觀方法.該方法將小行星表面所覆蓋的實體離散成一系列體積元,計算引力時用質點代替體積元,最后對所有體積元求和,以近似描述小行星的引力場.

假設將小天體劃分為n個體積元,ri和Mi分別表示第i個體積元的位置和質量,則空間位置為r的單位質量的引力勢表達式為

可以看到,質點群法的優勢在于算法簡單,易于編程實現;容易通過增加體積元的數量來提高引力場建模精度;針對“碎石堆疊”的質量分布不均勻的小天體,可通過調整各體積元的質量權重和離散規則,實現對真實質量分布的模擬,如圖2 所示.

圖2 質點群法對橢球的不同離散方式[14]Fig.2 The different way of filling the ellipsoid with point masses[14]

質點群法在小行星附近軌道動力學研究中得到了成功的應用.Ren 和Shan[44]利用質點群法研究了小行星的著陸軌跡優化.Yu 等[45]將雙星系統(65803)Didymos 中的次星用質點群近似描述,研究了次星受到撞擊后濺射物的軌道演化.

質點群法的缺陷也很明顯.收斂非常慢,隨著對精度要求的提高,所需體積元數量的急劇增加會造成計算量的迅速攀升;其次,離散的體積元數量過多,會造成計算的累計誤差較為顯著;最后,質點群法無法描述小行星的表面邊界[46],也就無法提供直接的碰撞檢測條件.

1.1.2 球諧函數法

球諧函數法是描述行星引力場最成熟和經典的方法,在人造地球衛星軌道動力學中具有重要的應用.該方法將引力勢表示成由正交球諧函數組成的無窮級數形式.引力勢球諧函數展開的一般表達式為

式中,r,φ,λ 為質點在球坐標系下的坐標,即距離,經度和緯度;Plm為連帶勒讓德函數;re表示參考球半徑,可取為天體表面距離質心的最遠距離,該參考球又被稱為布里淵球;各階球諧系數Clm和Slm用以描述天體的質量分布.

球諧函數法有著解析的級數表達式,可根據精度需要在任意階截斷,在動力學的解析研究方面具有巨大的優勢.應用球諧函數法的關鍵是獲取天體的引力場球諧系數,其主要技術途徑是通過天體的外形進行估算或利用軌道飛行數據進行反演計算.

在布里淵球內部,球諧級數將發散,意味著對于形狀不規則的小天體,在靠近表面的較大區域內其引力場無法通過球諧函數法描述.

針對球諧函數法收斂域的局限性,Werner[47]首次提出了內球諧函數法.該方法的收斂域為一個球心在天體外部、與天體表面相切的球體,稱為內布里淵球,據此可利用球諧級數描述該切點上方區域的引力勢.Takahashi 等[48-50]將其成功地應用到小天體表面引力場建模和動力學分析中.結合球諧函數法和內球諧函數法描述小天體表面附近的引力場,具有比下文將要介紹的多面體模型更高的計算效率.Wu 等[51]應用該方法研究了火衛一的表面起飛問題,將起飛軌跡經過的區域分段用內球諧函數法和球諧函數法來描述,如圖3 所示.

圖3 火衛一復合球諧函數模型分區示意圖[51]Fig.3 Regionalization of the compound gravity harmonic approach of Phobos[51]

內球諧函數法表述的場點引力勢函數的一般表達式為

該方法使用的坐標系原點為內布里淵球球心,r,φ 和λ為該坐標系下場點的球坐標;表示參考球半徑,一般取為內布里淵球半徑;是與天體質量分布相關的球諧系數.

橢球諧函數法是球諧函數法的另一個重要的推廣和改進,它將引力勢展開為正交橢球諧函數的級數,其收斂域為包含天體的最小橢球,使之更適合描述扁長和不規則形狀小天體的引力場.如圖4,橢球諧函數法可以極大地拓寬球諧函數法的有效范圍,其代價是正交橢球諧函數計算的復雜性.而且,該方法仍然不具有全局收斂性,無法精確計算小天體表面的引力場.另外,與質點群法類似,球諧函數法和橢球諧函數法也無法提供直接的天體表面碰撞檢測條件.

圖4 小行星4769 Castalia 的球諧(左)和橢球諧(右)引力場模型的收斂域[14]Fig.4 The convergence domain of spherical harmonics (left) and ellipsoidal harmonics (right) gravity models of the asteroid 4769[14]

1.1.3 多面體法

針對球諧函數法的局限性,20 世紀90 年代Werner 等[46,52]提出了利用勻質多面體模型計算不規則小天體的引力勢、引力和引力梯度的方法.小行星的多面體模型中,每一個面均為三角形,如圖5所示為小行星433 Eros 的多面體模型.

圖5 小行星433 Eros 多面體模型Fig.5 The polyhedron model of the asteroid 433 Eros

利用所有三角形頂點、邊和面的信息,就可以得到場點引力勢的一般表達式

式中,σ 為多面體的密度;re和rf為場點到多面體棱上和面上任意一點的矢量;Ee和Ef是基于二元并矢的與多面體形狀相關的常數矩陣;Le和 ωf是描述場點與邊、面相對關系的幾何量.對引力勢分別計算一次和二次梯度,可得到引力和引力梯度的解析表達式(詳細的推導和表達式可參考文獻[46]).

多面體法可以精確模擬天體不規則的表面形狀,建模精度取決于對天體觀測數據的精度.此外,多面體法具有封閉的表達式,不存在收斂性和截斷誤差的問題,能夠精確描述小天體表面引力場,因而非常適合小行星著陸和起飛等表面動力學問題的研究.但多面體法的高精度以大量計算為代價,因而在遠離小行星表面的區域常切換為球諧函數法以提高計算效率.多面體法的另一個優勢在于提供了天體表面碰撞檢測的直接指標.

多面體法憑借其對不規則引力場的精確建模能力,被廣泛用于小行星附近的軌道動力學研究中,并被成功應用到人類首個小行星環繞探測器NEAR(near Earth asteroid rendezvous).

1.1.4 基于特定簡單幾何體的近似方法

除了上述三種通用的不規則天體引力場建模方法外,還存在一些基于特定簡單幾何體的引力場近似方法,適用于具有特定幾何形狀的天體.例如勻質橢球體模型、鏈接質點模型和質點連桿模型等.

勻質橢球體的引力勢可表示為基于橢圓積分的封閉表達式[53]

Ms為球體質量; ρ 為橢球密度; α ,β 和 γ 分別表示橢球三軸半徑,存在關系 0 ＜γ ≤β ≤α ; λ 滿足

對引力勢分別計算一次和二次梯度可得到引力和引力梯度.該模型引力場表達式中涉及的無窮積分可通過卡爾松(Carlson)橢圓積分算法進行快速計算.

鏈接質點模型由兩個或兩個以上的質點用無質量的桿元件(質點間相對位置固定) 連接構成.1987 年,Chermnykh[54]首次提出了軸對稱剛體的兩質點旋轉啞鈴模型,Kokoriev 和Kirpichnikov[55-56]利用該模型研究了二體的平面運動.該啞鈴模型附近的限制性軌道運動可以視為圓型限制性三體問題的推廣,因此能夠對諸如平動點和穩定性的動力學問題進行解析的研究[57-60].2015 年,Zeng 等[61]利用質點偶極子模型描述形狀狹長的小行星,建立了該模型與真實小行星系統參數之間的關系.2017 年,Lan 等[62]利用三質點鏈接模型建立了拱形外形天體的引力模型.Qian 等[63]利用該模型研究了小行星附近考慮太陽引力的參數共振軌道.兩質點和三質點的鏈接模型如圖6 所示.

圖6 小行星216 Kleopatra 兩質點鏈接模型(上)與243 Ida 三質點鏈接模型(下)[62]Fig.6 The double-particle-linkage model of 216 Kleopatra (up) and the triple-particle-linkage model of 243 Ida (down)[62]

為了提高簡單模型對形狀狹長小行星的引力場近似精度,2017 年Zeng 等[64]在鏈接質點模型和直桿模型的基礎上進一步提出了質點連桿模型,即質點間用考慮質量的桿連接.對于兩質點的連桿模型,在其質心固連坐標系下,引力勢表達式為

式中,μ1和 μ2分別為質點和連桿質量參數;l為兩質點之間距離;r1和r2分別表示場點與兩質點之間的距離.

除了上述的幾種簡單幾何模型,在一些研究中還建立了直棒模型[65-66]、勻質圓環和圓盤模型[67-68]、勻質立方體模型[69-70]、非質點啞鈴體模型[71-72]等.雖然簡單幾何體模型對引力場刻畫的精度比多面體和質點群模型差,但其表達式具有簡單和解析的形式,易于分析動力學特性與模型參數之間的關系,并得到一般性的結論.

需要注意的是,上述引力場模型都具有各自的優勢和適用場景,在實際研究中,需要根據研究目標、精度要求和計算量約束進行權衡,選擇最為適合的模型,或將不同模型相互結合.

1.2 雙星系統相互引力勢模型

雙星系統的相互引力勢模型是研究雙星相對運動的基礎,其本質是將分屬兩星的兩質量元間的萬有引力勢對兩顆星各求一次體積分,是一個雙重體積分.除特殊情況外,兩個任意質量分布小行星之間的相互引力勢是沒有封閉表達式的,通常采用無窮級數展開的方法[73].根據處理雙重體積分的數學思路,可以將其分為三類:基于簡單幾何體的模型、基于剛體質量分布參數的模型和基于多面體形狀和剛體質量分布參數的混合模型.

1.2.1 基于簡單幾何體的模型

典型的雙小行星系統通常包含一個快速自旋接近球形的主星和一個被主星潮汐鎖定的同步旋轉的接近橢球形的次星[74].此時,通常可將其中一顆小行星用勻質球體代替,其萬有引力勢等價為一個質點的引力勢,那么雙星的相互引力勢就退化為質點在另一顆不規則小行星引力場中的引力勢,可利用1.1 節中的不規則引力場模型進行建模計算.

基于這類模型,學者們對雙小行星系統的平衡構型和穩定性問題進行了深入的研究.尤其當另一顆不規則小行星也用簡單幾何體代替時,其相互引力勢具有簡單的形式,由此可對雙星相對運動的平衡構型及穩定性與系統模型參數(如形狀尺寸參數、相對距離參數、質量參數等)之間的關系進行定量分析.例如,在研究這類問題時,Scheeres 等[21-22]采用球-橢球模型,李旭等[75]采用了旋轉偶極子-球/質點模型.

對該模型另一類重要的拓展,是將其中一顆小行星用質點群法描述,另一顆小行星可采用任意引力場模型,那么雙星相互引力勢的計算轉化為不規則小行星對每個體積元質點的引力勢之和.該模型能夠大幅提高對雙星相互引力的近似精度,尤其是當第二顆小行星采用高精度模型時,但是此時計算也更加復雜,適用于數值研究.例如,Yu 等[45]使用多面體-質點群模型研究了65803 Didymos 的相對運動,并分析了次星受到撞擊后拋射物的軌跡和動力學特性.設mi為每個體積元的質量,則雙星系統的多面體-質點群引力勢為

式中,UPoly為多面體的引力勢,見式(4).雙星的相互引力和引力矩也具有類似的結構[45].

1.2.2 基于剛體質量分布參數的模型

基于剛體質量分布參數的相互引力勢模型是一種經典的描述方法,在研究雙星姿態軌道耦合運動中具有廣泛的應用,其基本思路是將相互引力勢展開為剛體各階質量分布參數的無窮級數形式.剛體質量分布參數有多種具體的形式,如剛體引力場的球諧系數、剛體慣性積分(二階慣性積分即為慣性張量)和STF (symmetric trace free)張量等.雖然這些參數形式不同,但本質都來源于剛體質量元坐標的特定表達式對剛體的體積分,也即是反映了剛體的質量分布.

基于剛體質量分布的二體相互引力勢的一般形式為

式中,r表示兩剛體的質心相對位置矢量; ρ 表示分屬兩個剛體任意質量元之間的相對位置矢量,滿足ρ=r-(h1-h2),Δh=h1-h2,如圖7 所示.上式的計算基于兩質量元相對位置倒數的級數展開

圖7 兩剛體相互引力勢的幾何構型Fig.7 Geometrical representation of mutual gravitational interaction

1971年,Giacaglia 和Jefferys[76]基于剛體引力場的球諧系數推導了兩個剛體間的四階相互引力勢,但僅考慮到二階的球諧系數(等價于二階慣性張量/慣性積分).1979 年,Schutz[77]根據剛體的質量分布推導了兩個剛體的四階相互引力和二階相互引力矩,但同樣也只考慮了剛體的二階慣性積分.1988 年,Paul[78]推導了兩個任意質量分布剛體按慣性積分展開的相互引力勢,具有很強的普適性,但是其表達式包含多重復雜的嵌套求和,導致計算復雜.針對Schutz 模型的不足,2007 年Ashenberg[79]推導了完整的四階相互引力和引力矩,包含了剛體的二階、三階和四階慣性積分.2008 年,Tricarico[80]在Paul 工作的基礎上,進一步推導了兩個任意質量分布剛體的相互引力和引力矩,同時證明了基于剛體引力場球諧系數和基于剛體慣性積分的兩類模型本質上的一致性.2014 年,Compère 和Lema?tre[81]推導了基于STF 張量的二體相互引力勢,并用于雙星姿態軌道耦合動力學的建模和仿真.

另一類基于剛體質量分布參數的模型,由Werner 和Scheeres[73]于2005 年提出,他們將兩個剛體均用多面體描述,而且將每個多面體分解為眾多四面體的組合,通過質量分布參數計算兩個四面體之間的相互引力勢,再求和從而得到最終的相互引力模型.2006 年Fahnestock 和Scheeres[23]在此基礎上推導了兩個多面體之間的相互引力和引力矩,并建立了二體姿態軌道耦合的動力學模型.2013 年,Hirabayashi 和Scheeres[82]進一步建立了各階級數中形狀相關項的遞推關系.2017 年,Hou 等[83]利用類似的遞推關系建立的模型能大幅提高計算效率.之后Hou[84]對該方法進行了改進,定義了新的廣義積分慣量,由此可將積分效率提升至與球諧函數法相當的水平.2018 年,Jiang 等[85]基于上述雙多面體相互引力的建模方法計算了兩個八面體間四至六階的引力勢,并對二體姿態軌道耦合運動進行了仿真研究.2019 年,Yu 等[86]采用基于雙線性四面體元素的有限元方法建立了雙星之間的相互引力和引力矩模型.

Werner 和Scheeres[73]在計算兩個多面體之間引力勢時,首先將每個多面體完全分割成簡單四面體組成的實體,以此可將計算相互引力勢的二重積分轉化為計算A和B兩個多面體包含的所有簡單四面體對 (a∈A,b∈B) 間的引力勢之和,最后得到基于各階慣性積分/張量的引力勢模型[23]

式中,a和b代表兩個多面體的面與多面體質心構成的四面體,密度分別為 ρa和 ρb;Ta和Tb為每個四面體局部坐標到慣性坐標變換的雅可比行列式;Q,Qi,Qij和Qijk為常數慣性張量;w和 Δr為與四面體的位置和姿態相關的矢量.該模型能通過調整四面體體積元的密度參數 ρa和 ρb來適應雙星密度的不均勻性.

基于剛體質量分布參數的模型能夠計算兩個任意質量分布剛體間的相互引力,進而建立二體相對運動的姿態軌道耦合動力學模型,具有廣泛的適用性.可以通過保留低階項,對二體相對運動進行理論定性研究;也可以利用高階的級數展開,得到高精度的相互引力勢模型,進行精確的數值研究.

1.2.3 基于多面體形狀和剛體質量分布參數的混合模型

2017年,Shi 等[87]提出了基于多面體形狀和剛體質量分布參數的混合模型,其基本思路為將其中一顆星建模為勻質多面體,其引力場可采用多面體法進行解析描述,從而將二體相互引力勢展開為第二顆星慣性積分的無窮級數.

在多面體小行星的固連坐標系下,將二體相互引力勢展開為第二顆星質量分布參數的級數,形式為

式中,-V(r) 表示多面體引力勢,由式(4)給出;r表示雙星質心的相對位置矢量;h表示第二顆星質量元相對其質心的位置矢量.式中級數的通項可進一步表示為

式中的積分項即為第二顆星在其固連坐標系下的慣性積分,反映了其質量分布

由此,可計算展開到任意階慣性積分的雙星相互引力勢.

當展開到二階項時,雙星的相互引力勢和引力具有簡潔清晰的形式

式中P為與第二顆星二階慣性積分(即慣性張量)和兩星相對姿態相關的矩陣

第二顆星受到的引力矩為

而多面體小行星受到的引力矩可由如下關系得到

與1.2.2 節中對兩個剛體進行慣性積分展開的引力勢模型相比,該模型能夠對其中一顆星進行不規則引力場的精確建模,避免了同時舍去兩顆星的高階慣性積分.當忽略第二顆星對多面體小行星的引力時,該模型非常適用于描述剛體航天器在小行星不規則引力場中的姿態軌道耦合動力學,其引力模型包含了航天器的質量分布和姿態信息,精度高于經典的將航天器視為質點的多面體引力模型,而引力矩模型則完整包含了多面體小行星引力場的不規則性[88].例如,Lei 等[89]基于該引力矩模型研究了小行星216 Kleopatra 附近航天器的姿態穩定性和姿態周期解.

1.3 雙星相對運動的完全二體問題

雙星相對運動的完全二體問題本質是兩個任意質量分布的小行星在相互引力作用下的姿態軌道耦合運動.在早期的研究中,學者們以分析單小行星旋轉分裂后的雙星系統動力學演化為目標,重點研究了分裂后的雙星在潮汐力等作用下的長期演化過程.本節側重介紹航天動力學領域對雙星相對運動的平衡構型與穩定性的研究成果.在一些簡化的相互引力勢模型下,雙星的相對運動方程具有封閉的解析形式,以此為基礎可對系統的平衡構型與穩定性進行解析的研究.

1.3.1 相對運動方程

1995年,Maciejewski[16]在他關于完全二體問題的經典論文中,采用相對位置和姿態坐標,在其中一個剛體B2的固連坐標系中(如圖6 所示),給出了具有一般性的二體相對運動的牛頓-歐拉方程

式中,r和A分別為二體相對位置矢量和相對姿態矩陣; Ω1和 Ω2為自轉角速 度;H1和H2為自轉角 動量;相互引力勢的一般表達式為一個雙重體積分

μ1和 μ2分別表示兩個剛體受到的引力矩,均在B2的固連坐標系中表示,具有形式

式中 α,β 和 γ 滿足

在該模型下,系統存在能量積分和角動量積分.針對雙星的簡化模型,還能進一步降低系統自由度.例如,在球-橢球二體模型中,可以不考慮球的轉動,那么在橢球固連坐標系下,橢球轉動的歐拉方程可簡化為[26]

如果僅考慮球-橢球的平面運動,利用引力場的空間對稱性,可將雙星的相對運動簡化為一個二自由度的哈密頓系統.2018 年,Hou 和Xin[90]考慮雙星相互引力的二階球諧模型,利用對稱性建立了雙星平面內的姿態軌道耦合運動方程,并得到了基于軌道偏心率和雙星非球形引力一階近似的解析解.

1.3.2 相對運動的平衡構型與穩定性

雙小行星系統需要滿足一定的能量和角動量條件才能保證系統的穩定(不發生逃逸或碰撞).雙星的姿態軌道耦合運動表明雙星平動能量和角動量與轉動能量和角動量之間存在轉換和約束關系.Scheeres[15]基于此研究了二體系統一般性的希爾穩定性準則和雙星避免碰撞的穩定性條件.2009 年,Scheeres[8]在具有空間對稱性的二階相互引力勢模型下,給出了量化的系統平面穩定性條件.杜燕如等[20]則研究了雙橢球模型下,系統的物理參數(如橢球形狀參數、質量比等)對系統平衡態和穩定性的影響.

在雙小行星系統穩定(不發生逃逸或碰撞)的前提下,2004 年,Scheeres[21]研究了球-橢球二體模型下系統的平衡構型及其穩定性,在該對稱引力場中存在6 種平衡構型,其中橢球繞其最大慣量主軸轉動的兩種平衡構型如圖8 所示,它們具有不同的穩定域,其平面運動穩定性是二體質量比、相對距離和橢球形狀的函數.針對長軸平衡構型,對于給定的角動量存在兩個穩定性相反(不同能量)的平衡解,Bellerose 和Scheeres[22]研究了這兩個平衡解附近的周期軌道族,并且發現高能量狀態的平衡解經過能量耗散可轉化為低能量狀態的平衡解.圖9 展示了某一球-橢球二體系統的長軸平衡構型隨二體質量比和相對距離變化的穩定性特征,以及給定角動量下長軸構型的穩定和不穩定平衡解(圖中數值經過了歸一化)[22].

圖8 球與橢球二體系統相對運動平衡構型[21]Fig.8 Relative equilibrium configuration of the sphere-ellipsoid two-body problem[21]

圖9 球與橢球二體系統長軸平衡構型的穩定[22]Fig.9 Stable parameter space of the sphere-ellipsoid long-axis equilibrium configuration[22]

對于雙橢球的平面二體模型,在給定的系統角動量下,不同的系統能量對應了不同的系統穩定狀態—決定系統是否發生碰撞和逃逸,同時也對應了不同的雙星平衡構型的穩定狀態—譜穩定性和能量穩定性[8].2020 年,Wang 和Hou[91-92]從周期軌道的角度,研究了圖8 所示的平衡構型下次星(橢球)的天平動,發現天平動振蕩幅度與其繞主星的軌道偏心率大小成反向關系.

2 基于限制性三體問題的動力學研究

2.1 軌道動力學模型

在研究航天器在雙小行星系統附近的軌道動力學問題時,將主星和次星視為主天體,忽略航天器對兩個主天體的引力作用,此時的軌道動力學問題屬于限制性三體問題.如上節所述,需要首先確定兩顆小行星之間的姿態軌道耦合運動,因此該問題被稱為限制性完全三體問題,該問題是對主天體為質點的經典限制性三體問題的推廣.與經典問題相比,小行星顯著非球形/不規則引力場和系統極小的時空尺度引入了嶄新和復雜的動力學特性,為探測器軌道設計帶來了新挑戰和新機遇.本節主要介紹針對雙小行星系統常采用的模型簡化方法以及所建立的軌道動力學模型.

現有對系統的簡化主要包括兩個方面,引力場模型的簡化和雙星相對運動的簡化.在引力場模型方面,可以根據主星和次星的形狀特點將其建模為球、橢球等簡單幾何體,也可以對主星和次星的球諧引力場模型在合適階次截斷.此外,主星和次星可以根據其特點和研究需要選擇不同的引力場模型進行組合,從而獲得適合的雙星系統的引力場模型.Gabern 等[93]將主星簡化為固連的三個共線質點,將次星簡化為質點/球,這種簡潔形式有利于理論研究和解析推導.對于主星形狀接近球體的雙小行星系統,Bellerose 和Scheeres[25,94]將主星和次星分別簡化成球和橢球.Woo 和Misra[95-96]將兩顆小行星分別簡化成橢球和橢圓臺,相比于球-橢球模型引入了模型的非對稱性.在雙星相對運動方面,對于偏心率較小的雙星系統可忽略偏心率,將雙星軌道運動簡化為圍繞公共質心的圓周運動.此外,由于主星具有近似軸對稱的特點,對主星快速自旋次星同步轉動的單同步雙星系統可忽略主星的自旋,假定系統是雙同步的.

非同步和單同步雙小行星系統附近的力學環境具有時變性,只存在離散分布的、周期與雙星自身運動周期滿足特定整數比的共振周期軌道.目前的研究大部分是針對雙同步的雙小行星系統,如圖10所示.雙同步雙小行星系統附近圓型限制性完全三體問題的軌道動力學模型具有和主天體為質點的經典圓型限制性三體問題相似的形式,主要區別在于主天體的引力勢表達式

圖10 航天器(質點)在雙小行星系統附近的軌道運動Fig.10 A massless particle moving in the gravity field of a binary asteroid system

其中U1,U2分別為主星和次星的引力勢,下標表示對其求對應的偏導數.

2.2 平動點及平動點周期軌道

對限制性三體問題平動點及其附近周期軌道的研究由來已久,具有很高的理論研究價值和廣闊的應用空間.早在1767 年和1772 年,數學家歐拉和拉格朗日在研究主天體為質點的經典圓型限制性三體問題時發現了5 個平動點.1967 年,Szebehely[97]在其專著《Theory of Orbits》中總結了關于圓型限制性三體問題的研究成果,并通過一些數值結果充實了理論研究.1966 到1974 年間,Hénon[98-99]在希爾限制性三體問題模型下對周期和準周期軌道進行了系統的研究,對軌道進行了分類,并研究了各類軌道的穩定性.

在雙同步雙小行星系統附近的限制性完全三體問題中,同樣存在著平動點及其附近的周期軌道,而且小行星的非球形/不規則引力場對平動點及周期軌道存在顯著的影響.關于雙小行星系統附近平動點的研究主要集中在不同引力場模型下平動點的存在性、位置以及穩定性的差別.將雙星均建模為質點或勻質球體時,問題將退化為經典圓型限制性三體問題.2001 年,Sharma 等[100]將主星和次星均假設為剛體,研究了平動點的存在性和穩定性.2006 年,Gabern 等[101]將主星和次星分別建模為球體和橢球體,研究了非球形引力場對平動點位置和穩定性的影響.2005 年～2008 年,Bellerose 和Scheeres[24-27,94,102]對雙小行星系統附近的限制性完全三體問題進行了深入研究,將雙小行星1999 KW4 建模為球-橢球模型,在長軸平衡情況下計算了平動點,并對平動點的穩定性進行了分析.2018 年,Li 等[103]研究了不同尺寸雙橢球模型附近的平動點,并指出平動點的位置受橢球長軸影響較大.這些研究采用的引力場模型均具有一定的對稱性,與小行星形狀不規則的特點并不完全一致.2014 年,Woo 和Misra[95-96]設計了橢圓臺形和花生形等不對稱的引力場模型,研究了不同尺寸雙星系統的平動點位置與穩定性.球諧函數引力場模型也被一些關注平動點的研究所采用,一些新的現象被發現.2015 年,Hou 等[104-105]采用球諧函數模型近似雙橢球系統的引力場,發現系統附近會出現額外的平動點,而直接采用橢圓積分精確計算橢球引力勢時則不會產生這種現象,分析發現,額外出現的平動點是球諧函數模型在小行星表面附近誤差過大而錯誤引入的.多面體模型作為一種高精度的小行星引力場模型也被用于雙小行星系統附近平動點的研究.2018 年,Shi 等[32]采用雙多面體模型對雙小行星1999 KW4 進行引力場建模,計算了雙星系統附近的平動點和零速度曲面,如圖11 所示.可以看到雙星系統附近的平動點和零速度面與經典圓型限制性三體問題相似,存在5 個平動點,其中L1,L2 和L3 點對應于經典問題的共線平動點,L4 和L5 點則對應于三角平動點.由于小行星不規則引力場的影響,L1,L2 和L3 點不再嚴格位于x軸上,而L4 和L5 點的位置也不再與主星和次星構成嚴格的等邊三角形.Bu 等[106]在2017 年研究了雙橢球模型的雙小行星系統在連續小推力作用下的人造平動點及其穩定性,并借助零速度曲面研究了小推力對可達域的影響.此外,對于非同步的雙星系統,由于系統的時變性,不存在嚴格意義上的平動點,只存在瞬時的“平動點”,可被稱為動力學替代(dynamical substitutes).Shi 等[32]計算了雙小行星系統1999 KW4 在主星快速自旋的單同步情況下的瞬時平動點,通過對比發現瞬時平動點在z方向的變化小于在x和y方向的變化,而且主星的自旋并未改變平動點的穩定性.2018 年,Hou 和Xin[90]將雙星建模為兩個剛體,雙星相互引力勢保留至二階,在平面完全二體問題一階解析解的基礎上構建了限制性三體問題模型,并求解了其中的共線平動點.

圖11 雙多面體模型下雙星系統1999 KW4 附近的零速度面與平動點[32]Fig.11 Zero velocity curves and libration points of the binary asteroid system 1999 KW4[32]

經典限制性三體問題的平動點周期軌道已經有了非常深入和廣泛的研究,并應用于實際飛行任務.雙小行星系統附近限制性完全三體問題的平動點周期軌道也吸引了不少學者的關注.

2015年,Chappaz 和Howell[28]利用多段打靶法得到了球-橢球和橢球-橢球雙同步雙星系統附近的平面Lyapunov 軌道,垂直Lyapunov 軌道,南向和北向的Halo 軌道,如圖12 所示.數值計算結果表明,雙同步雙星系統平動點周期軌道的結構和經典圓型限制性三體問題類似.由于球和橢球模型的引力場具有對稱性,平動點周期軌道也具有和經典問題一樣的對稱性,具體體現在平動點周期軌道關于雙星連線的對稱性,南向和北向Halo 軌道的對稱性,以及L4 和L5 平動點軌道的對稱性.而在更高精度的非對稱引力場模型下,平動點周期軌道雖然和對稱模型下的軌道相似,但是軌道的對稱性和各族軌道之間的連接關系會發生顯著變化.2017 年,Dell’Elce等[31]采用雙橢球模型研究了雙小行星系統65803 Didymos 的二體運動,之后通過多面體-橢球模型研究了雙小行星系統附近的周期軌道.2018 年,Shi等[32]將雙小行星系統1999 KW4 建模為雙多面體模型,在雙同步情況下求解了平動點和平動點周期軌道,發現受到不規則引力場的影響,平動點周期軌道族的分岔情況與經典圓型限制性三體問題完全不同.圖13 為經典問題和雙小行星系統平動點周期軌道族連接方式的對比,可以看到在雙小行星附近,Lyapunov 軌道族發生了斷裂,分別與南北向Halo 軌道融合.Capannolo 等[30]在2019 年采用多面體-橢球模型研究了雙小行星系統65803 Didymos 附近的周期軌道族,包括平動點L1 和L4 附近的Lyapunov 軌道,Halo 軌道以及短周期軌道(short period orbit,SPO),并分析了不同主星姿態對軌道的影響.Zhang 等[107]在2020 年借助多面體-橢球模型研究了次星形狀的不確定性對平動點及平動點周期軌道的影響,發現次星形狀的變化將導致平動點位置的變化,穩定性可保持不變,而周期軌道族的分岔情況會發生顯著改變,南北向Halo 軌道不再與平面Lyapunov 軌道融合,而是與垂直Lyapunov 軌道融合.此外,軸向軌道也與垂直Lyapunov 軌道斷裂后的不同部分相融合,而不像經典問題中那樣兩端分別連接平面Lyapunov 軌道和垂直Lyapunov 軌道,如圖14 所示.

圖12 雙橢球模型中的平動點軌道族[28]Fig.12 Libration point orbits (LPO) in the ellipsoid-ellipsoid model[28]

圖13 不同模型下平動點周期軌道族連接方式[32]Fig.13 Connections of libration point orbits in different dynamical model[32]

圖14 不同次星形狀時平動點軌道族連接情況[107]Fig.14 Connections of libration point orbit families with different shapes of the secondary[107]

單同步雙小行星系統主星的快速自旋導致系統是非自治的,動力學模型具有時變的特點.此時雙星系統附近雖然不再存在完整的周期軌道族,但仍存在一些與主星在會合坐標系中自旋周期滿足特定整數比的共振周期軌道.Shi 等[32]在單同步模型下通過微分修正求解得到了滿足特定周期條件的共振周期軌道,如圖15 所示.Shang 等[108]在2017 年借助拉格朗日擬序結構(LCS)研究了非同步球-橢球雙星模型附近的軌道,將發現的15 種軌道分為了四類:短暫軌道、逃逸軌道、撞擊軌道和飛掠軌道.

圖15 單同步雙星系統L1 點和L2 點附近的周期軌道[32]Fig.15 Periodic orbits near L1 and L2 about the single synchronous binary asteroid system[32]

在雙星系統附近,除共線平動點外,學者們對三角平動點L4 和L5 附近的軌道動力學也進行了研究.Liang 等[109]研究了雙小行星系統三角平動點區域的動力學和穩定區域,將雙小行星系統1999KW4簡化為經典圓型限制性三體問題,發現在雅克比常數在特定區間時存在兩個穩定區域,后續在球-橢球模型下發現了可以圍繞L4 點14 天以上的有界軌道.Hou 等[110]在2020 年將雙星系統建模為兩個三軸橢球,研究了受太陽光壓力時三角平動點附近的軌道運動,發現即使對較小面質比的航天器,受迫周期軌道也可以有很大的振幅.

2.3 其余類型的周期軌道

由于平動點軌道往往位于天體的一側,適用于對天體某一區域的持續觀測,但如果需要對雙小行星系統進行全面探測,就需要考慮其他類型的周期軌道.除了受到廣泛關注的平動點周期軌道外,雙星系統附近還存在大量其他類型的周期軌道,比如圍繞主星和次星的周期軌道,與雙星軌道運動具有共振關系的共振周期軌道,以及圍繞雙星系統的大范圍周期軌道等.不同類型的周期軌道各具特點,可以根據探測任務的需要選作名義軌道,對雙小行星探測任務設計具有重要意義.

對于不規則引力場中的周期軌道求解問題,數值搜索是一種有效且全面的手段,Yu 等[111]在2015 年提出了一種單小行星附近的周期軌道搜索算法,并將該算法應用于小行星216 Kleopatra,得到了29 族周期軌道.2019 年,Shi 等[33]對該算法進行了改進,得到了適用于雙小行星系統的搜索算法,搜索得到了雙小行星系統1999 KW4 附近大量的周期軌道,其中就包括圍繞主星和次星的周期軌道,如圖16 所示.此外,Shang 等[29]在2015 年借助雙橢球同步雙星系統中的對稱性,通過搜索算法得到了雙小行星809 Lundia 附近的30 族周期軌道和雙小行星3169 Ostro 附近的28 族周期軌道,并分析了軌道的穩定性以及潛在的應用前景.

圖16 雙小行星系統1999 KW4 中部分圍繞主星和次星的周期軌道[33]Fig.16 Periodic orbits about the primary and the secondary in the binary asteroid system 1999 KW4[33]

在圍繞次星的周期軌道中存在一種逆行軌道,即遠距離逆行軌道(DRO).該類軌道在1968 年由Broucke[112]在地月系統中發現,屬于穩定周期軌道,因其在探測任務中具有廣闊的應用前景而受到了眾多學者的關注.雙小行星系統附近的DRO 也得到了特別關注.2015 年,Chappaz 和Howell[28]將雙星系統建模為球-橢球系統,并對DRO 進行了數值計算以及穩定性分析.Capannolo 等[113]在2017 年研究雙小行星65803 Didymos 附近的探測器編隊飛行問題時考慮了DRO,在2018 年利用多面體-橢球模型研究雙小行星系統65803 Didymos 附近的周期軌道時考慮了圍繞次星的DRO,并分析了不同主星姿態角對DRO 的影響.Jean 等[114]在2019 年研究太陽光壓力對雙星系統附近軌道影響時也將DRO 作為主要研究對象之一.

在限制性三體問題中還存在一種具有較高應用價值的周期軌道,即與主天體軌道運動具有共振關系的共振周期軌道.DRO 在幅度較大時也屬于共振周期軌道的一類.2014 年,Vaquero 和Howell[115]對地月系統中共振軌道的動力學結構進行了深入細致的分析,得到了平面及三維共振軌道,并提出利用共振軌道實現低能轉移.2016 年,Feng 等[116]對相連雙小行星附近動力學環境進行了分析,發現了三類不同共振比的共振軌道.Chappaz 和Howell[28]在2015 年借助球-橢球模型和雙橢球模型研究了雙星系統附近的平面及三維共振軌道,并分析了這些軌道的穩定性.圖17 給出了Chappaz 和Howell[28]得到的平面共振軌道,其中1:1 共振軌道即是前文提到的DRO.2019 年,楊雅迪等[117]將雙小行星1999 KW4 建模為雙橢球模型,計算了共振比1:1,1:2,1:3,1:4 和2:3 的平面和三維共振軌道,分析了共振軌道周期和能量的變化規律.

圖17 雙小行星雙橢球模型下的平面共振軌道[28]Fig.17 Planar resonant orbits in the ellipsoid-ellipsoid model[28]

除共振周期軌道外,還存在著一些圍繞雙星系統的大范圍周期軌道,可用作雙小行星探測的初始停泊軌道或者長期遠距離探測軌道.Shi 等[33]在2019 年借助雙多面體模型搜索雙小行星系統1999 KW4 附近的周期軌道時,得到了許多圍繞整個系統的大范圍周期軌道,如圖18 所示.

圖18 圍繞雙小行星系統1999 KW4 的大范圍周期軌道[33]Fig.18 Period orbits with wide range around the binary asteroid system 1999 KW4[33]

2.4 轉移軌道研究

在雙小行星探測任務中,轉移軌道是必不可少的重要環節,高效安全的轉移軌道可以為雙星系統的全面探測提供有力保障.目前,對探測器在雙小行星系統內的軌道轉移研究較少,主要針對同步雙小行星系統進行了轉移軌道設計,此時系統與經典圓型限制性三體問題類似,可以類比采用經典問題中的方法解決.Ferrari 等[118]使用不同三體問題的拼接研究了探測器的轉移軌道,將雙小行星系統的主星用固連的兩個質點代替,將次星視為質點,將主星附近的軌道運動和遠離主星的軌道運動分別用兩個圓型限制性三體問題模型來近似,根據兩個三體問題不變流形的分布特征,采用流形拼接方式設計了全局的低能轉移軌道和小行星表面的起飛或著陸軌道.Chappaz 和Howell[119]求解了球-橢球模型下平動點軌道的穩定與不穩定流形,研究了同步雙星系統中的同宿和異宿轉移軌道,將2:3 共振軌道以及圍繞主星和次星的軌道作為探測軌道,通過微分修正設計轉移軌道將探測軌道連接起來,從而實現雙星系統的全局探測,如圖19 所示.

圖19 球-橢球模型下的全局探測路徑[119]Fig.19 Exploration trajectories in the sphere-ellipsoid model[119]

2.5 周期軌道維持

由于任務名義軌道的不穩定性以及入軌偏差、未建模攝動力、導航誤差及推力誤差等各類擾動,探測器的飛行軌跡通常會偏離所設計的名義軌道,需要進行軌道維持.除了面臨經典限制性三體問題的共性問題外,雙小行星系統附近的軌道維持需要考慮更多的約束,面臨更大的困難.首先,由于觀測條件和精度的限制,小行星運動狀態的觀測比行星要困難許多,有很強的不確定性,名義軌道無法在更高精度模型中做進一步的優化,與實際飛行軌跡相差更大,對軌道維持能力提出更高要求.其次,小行星的引力場非常微弱,探測器受到的擾動相對更為顯著,更容易偏離名義軌道.最后,雙小行星系統的時間和空間尺度都遠小于地月系統或日地系統,偏離名義軌道后若沒有及時修正,將很快發生逃逸或撞擊.因此,雙小行星系統附近的軌道維持需要更高的控制精度和頻率.2014 年,Woo 和Misra[96]在圓錐臺-橢球模型下計算了平動點L2 和L4 附近的軌道,并通過設計線性反饋控制器實現了軌道維持.2018 年,Li 等[34]通過雙橢球模型計算了雙星的相對運動,從而建立了非同步情況下雙橢球系統附近的軌道動力學模型,將同步模型下的周期軌道作為初值,計算了非同步模型下的有界軌道,并將其作為任務名義軌道,設計了一種自適應的軌道維持控制律.上述研究均采用了連續推力方式.2020 年,Shi 等[35]將主星建模為多面體,次星視為具有二階質量分布參數的剛體,對雙星系統的相對運動進行了仿真,計算了雙星處于激發態和平穩態兩種情況下圍繞主星、次星的軌道和L1 點Halo 軌道的演化,發現在激發態情況下,即使穩定的名義軌道也會迅速發散,如圖20 所示.接著設計了基于目標點法和離散LQR 法的兩種脈沖式軌道維持策略,發現活躍態情況下維持消耗要顯著大于平穩態,改進后的目標點法在脈沖消耗和脈沖間隔上的表現要優于離散LQR 法,更適合雙小行星系統附近的軌道維持.圖21中給出了L1 點Halo 軌道在平穩態與活躍態系統中采用目標點法實現軌道維持時的情況.2019 年,Guzzetti 等[36]研究了太陽帆航天器在雙小行星系統中的人造平動點及擬周期軌道,并設計了一種通過調整太陽帆姿態角實現L4 點區域軌道維持的控制策略,可以實現數周的維持.

圖20 名義軌道在平穩態(上)與活躍態(下)系統中的軌跡[35]Fig.20 Trajectories of nominal orbits in the dynamical systems with relaxed mode (up) and excited mode (down)[35]

圖21 L1 點Halo 軌道在平穩態(上)與活躍態(下)系統中采用目標點法實現軌道維持[35]Fig.21 The trajectories during the stationing-keeping in the dynamical systems with relaxed mode (up) and excited mode (down) by using the target point method with the L1 Halo orbit as nominal orbits[35]

3 基于受攝二體問題的動力學研究

上一節著重介紹了在限制性三體問題框架內,平動點及平動點周期軌道、大范圍周期軌道、轉移軌道和周期軌道維持等方面的研究現狀.這些研究是對球形主天體的經典限制性三體問題的自然推廣.然而事實上,在雙小行星附近的軌道動力學中,還存在對更為經典的受攝二體問題的推廣,即環繞系統單顆星的軌道動力學.本節將針對環繞雙星系統單顆星的環繞型軌道動力學問題,即受攝二體問題的研究進行綜述,以受攝二體問題為主線,首先引入經典的軌道攝動理論和太陽系行星系統中受攝二體問題的相關研究成果,并以此為基礎介紹環繞雙小行星系統單顆星的軌道動力學研究.

3.1 受攝二體問題的動力學模型及基本解法

歷史上對受攝二體問題的研究始于17 世紀末對太陽系行星和行星天然衛星軌道運動的研究.在跨越近兩個世紀的時間里,涌現了諸如牛頓、歐拉、拉格朗日、拉普拉斯等一大批的天才科學家,他們致力于解決月球軌道理論和行星軌道穩定性等難題;到了近現代,1957 年第一顆人造地球衛星的發射,極大地拓展和豐富了軌道攝動理論的發展和應用.長久以來,對攝動問題的求解方法形成了三種基本的類型:特殊攝動法,一般攝動法和半解析法.

3.1.1 特殊攝動法

特殊攝動法(也稱直接法)是指借助數值法求解軌道運動的微分方程,得到基于軌道初始條件的特定解.數值法以Encke 法和Cowell 法為代表.

19世紀前期,Encke 提出了建立密切軌道與參考軌道(無攝動的開普勒軌道)偏差的運動微分方程并進行數值積分的方法,而不直接求解包含中心引力和攝動力的完整運動方程[120].基于該方法建立的攝動方程為

式中,ρ 表示參考軌道的位置矢量; Δr為真實位置與參考位置的偏差矢量;f表示攝動加速度.Encke 法的優勢在于,由于偏差的變化緩慢,數值積分可采用更大的步長,易于獲得高精度的數值解.但這種優勢被現代計算機的強大計算能力所掩蓋,因此該方法在如今的研究中應用較少.

20世紀之前,受計算機能力的限制,數值法沒有受到太多的關注,直到Cowell 對攝動軌道方程(形如式(28),被稱為Cowell 公式)采用直接數值積分,成功定位了八顆木星衛星的軌道,并且對哈雷彗星在1759 年～1910 年三次回訪中的兩次進行了成功的預測[120]

利用上式,可以方便地將任意的攝動力添加到軌道的運動微分方程中(線性疊加).而且,得益于現代計算機的強大性能,基于Cowell 公式進行數值求解仍然是現代攝動軌道計算的常用方法.

特殊攝動法雖然形式和原理上都比較簡單,而且可以得到精確的數值解,但是獲得的解依賴于軌道初始條件和具體的模型參數,難以推廣到更廣泛的問題中.

3.1.2 一般攝動法

一般攝動法是指用解析方法求解攝動軌道的運動微分方程,從而得到具有普遍性的解.解析法是最早被用于軌道攝動研究的方法,可以追溯到17 世紀牛頓在《自然哲學的數學原理》中對月球軌道運動的研究和描述.在此后的兩個多世紀里,克萊羅、漢森、德洛納、希爾、布朗和埃克特等先后研究了月球的軌道運動理論,最終建立了非常精確的月球運動模型[121].

在航天動力學中,一般攝動法大多以軌道要素的攝動方程為出發點,而軌道要素變化微分方程的建立則是以常數變易法(variation of parameters,VOP)為基礎.VOP 的數學方法最早由歐拉提出,其基本思想是用無攝動系統的解(常量)來表示相應攝動系統的解,不過這些未攝動常量需要替換為隨時間變化的量.1782 年拉格朗日將VOP 應用到軌道的攝動問題中,建立了針對經典軌道要素的拉格朗日行星方程[122]

式中,n為平均軌道運動;R為攝動函數(攝動勢能的相反數).拉格朗日行星方程僅適用于攝動力為保守力的情況.

對于攝動力為非保守力的系統,或者即便是保守力,也可采用攝動力分量的形式來建立相應的攝動運動方程,這種形式的方程稱為高斯攝動方程.為方便研究,可將攝動力在不同的坐標系下分解,例如將其沿軌道徑向、橫向和軌道面(副)法向分解(S,T,W型);或沿切向、主法向和副法向分解(U,N,W型),從而能夠得到兩組不同的軌道攝動方程,具體表達式可參考文獻[122]的第三章.

在軌道攝動方程中除了采用經典軌道要素,還可采用米蘭科維奇軌道要素:角動量矢量H、拉普拉斯矢量L(可轉換為偏心率矢量e)和一個表示沿軌道運行時間的標量,分別建立拉格朗日型和高斯型的攝動運動方程[123].矢量形式的方程具有結構緊湊和幾何意義直觀的特點,能夠清晰地反映軌道運動的某些動力學特征.

對于哈密頓系統,還可建立正則形式的軌道攝動運動方程.在軌道力學中通常采用的正則共軛變量是德洛納(Delaunay)變量:L,G,H(角動量),l,g,h(角變量),它們與經典軌道要素之間的關系為[122]

針對上述正則變量的攝動運動方程為[122]

式中 F 為哈密頓函數.

按照上述方式建立的軌道攝動運動方程都是復雜非線性的,難以得到其準確的解析解.但是由于攝動均為小量,可將攝動力展開為小參數的冪級數形式,由此可構造攝動運動方程的小參數冪級數解,這種求解方法稱為攝動法[122].在攝動法的基礎上,日本天文學家古在(Kozai) 最早給出了改進的攝動法——平均根數法[124].該方法首先將攝動變化項拆分為長期項、長周期項和短周期項,構造攝動冪級數解采用的參考軌道不再是初始軌道對應的無攝動橢圓軌道,而是采用包含長期變化的橢圓軌道,即長期進動橢圓,相應的軌道根數稱為平均根數.平均根數法在求解地球扁率、中心天體非球形引力、第三體引力和太陽光壓力等作用下的攝動冪級數解中具有重要的應用.

3.1.3 半解析法

鑒于數值法計算的復雜性和尋找解析解的巨大困難,20 世紀60 年代平均化方法被應用到受攝軌道研究中,基于該方法建立的針對攝動長期項和長周期項的軌道運動模型稱為半解析模型.平均化方法消除了軌道攝動的短周期項,一方面在數值積分時可采用更大的積分步長,提高計算效率;另外,半解析模型在分析軌道長期演化的動力學特征和揭示軌道演化的一般性理論方面具有巨大的優勢,有著廣泛的應用[125-127].

在早期的研究中,通過分離出攝動函數中與快變量(平近點角) 無關的項[128-129],或者通過高斯(環) 平均化方法[130]來實現短周期攝動項的分離.1961 年,Bogoliubov 和Mitropolsky[131]提出了用于非線性動力學系統漸近分析的Krylov-Bogoliubov-Mitropolsky 平均化方法,利用該方法能夠以更加簡單直接的方式消除方程中的短周期項.例如,對攝動函數的平均化可以表述為

式中,M為軌道平近點角; α 代表其他軌道要素和變量.將平均化后的攝動函數代入到拉格朗日行星方程中,就得到了軌道的長期攝動方程;或者對攝動力進行平均化,再利用高斯攝動方程相應能得到高斯型的長期攝動方程.

以矢量要素形式的長期攝動方程為例,其表達式具有簡潔對稱的結構[132]

式中,h和e分別表示歸一化的角動量矢量和偏心率矢量;表示平均化和歸一化的攝動函數.同樣,將平均化后的哈密頓函數代入到正則形式的攝動方程(30)中,也能得到基于Delaunay 正則變量的軌道長期運動方程.

基于上述消除短周期項的軌道長期攝動運動方程,借助一些分析方法和數值方法,就能夠對軌道攝動運動問題展開研究,下一節將具體介紹.

3.2 受攝二體問題的動力學特性研究

在復雜攝動力作用下,航天器環繞中心天體的軌道可能發生長期變化,因此分析軌道變化的規律對于任務設計和長期維持策略都具有重要的意義.一方面要盡量避免攝動力對軌道產生重大的影響,通過分析軌道的穩定性、設計凍結軌道等來保證軌道的長期穩定;另一方面,也能利用軌道攝動的進動特征,設計某些特殊類型的軌道,如太陽同步軌道.在不同的天體系統中,系統參數不同,不同類型攝動力的相對大小也不同,航天器軌道將具有不同的演化特征,這里主要闡述環繞地球軌道和環繞其他行星或行星衛星軌道的相關研究.

3.2.1 近地軌道的長期演化動力學

1957年第一顆人造地球衛星的發射,極大地推動了對軌道攝動問題的研究.1959 年,Brouwer[128]和Kozai[129]幾乎在同一時間研究了低軌道衛星在地球非球形引力攝動下的軌道運動,計算獲得了短周期項、長周期項和長期項的近似解析解,推動了衛星軌道理論的發展.但是他們的研究沒有考慮大氣阻力的作用,而對于低軌衛星大氣阻力有顯著的影響.后來的學者們在他們研究的基礎上不斷發展改進,相關理論被成功應用到工程實踐中.

對于中高軌道的人造地球衛星,不僅受到地球非球形引力的攝動,日月第三體引力的攝動作用也非常顯著.1964 年,Allan 和Cook[132]考慮了地球扁率和日月第三體引力攝動,利用式(33)所示的半解析模型,研究了圓軌道的進動效應,發現各攝動力會致使軌道平面繞攝動軸進動,最終效果為軌道平面繞著某一合成攝動軸進動,如圖22 所示.

圖22 軌道平面法向量在單位球面上的進動軌跡[132]Fig.22 Trajectories of the orbit pole on the unit sphere[132]

在早期的工程實踐中,人們利用軌道演化特性,通過合理配置軌道初始狀態來實現軌道的長期穩定[133].其中,凍結軌道是一類重要的研究對象,它可以最大限度節省軌道維持所需的燃料[134-135].例如,地球低軌衛星在地球引力帶諧項為主要攝動力的環境下,J2臨界傾角(6 3.4°和1 16.6°)附近的軌道能夠保持近地點(偏心率和近地點幅角)的穩定[136-140],但隨著軌道高度的增加,日月第三體引力的顯著作用會導致近地點高度的振蕩和近地點位置的變化[141-142].當軌道半徑位于3 至10 倍地球半徑范圍時,地球扁率和日月第三體引力的攝動都不能忽視,形成了與低軌衛星不同的動力學環境.Ulivieri等[134]研究了這個范圍內圓軌道的凍結軌道方案,Circi 等[135]利用矢量描述的方法研究了類似場景下橢圓軌道的凍結軌道.

3.2.2 環繞其他行星或行星衛星的軌道動力學

自20 世紀70 年代以來,軌道的長期攝動理論還被拓展到環繞月球的軌道和其他行星系統中,用于第三體引力的長期攝動作用[126,143-144]、凍結軌道[135,145-146]、軌道穩定性分析[147-148]和維持軌道長期穩定[149-150]等方面的研究.在不同的系統中,由于具有顯著不同的系統參數,軌道的長期演化呈現出不同的特征.例如,月球和金星等天體球諧引力的扇諧項或高階帶諧項,會導致基于J2項得到的臨界傾角值發生較大改變[151-152].某些系統中,強烈的第三體引力攝動會引起高傾角軌道的極不穩定,致使軌道壽命的大幅縮減[127,147,153].在大多數以中心天體引力J2項和第三體引力為主要攝動力的行星系統中,軌道穩定性差異本質是由這兩種攝動力的相對強弱決定的.Scheeres 等[147]定義了一個參數來表征這兩種攝動力的相對大小,并研究了不同比值下軌道的穩定性.Delsate 等[145]研究了針對該系統參數變化的系統平衡點(凍結軌道)的分岔,如圖23 為平衡點分岔線,將模型參數空間分成了不同的區域,在每個區域內系統的平衡點、相空間結構及穩定性都具有不同的特征,可以發現系統除了 ω =±π/2 的Lidov-Kozai 平衡點外,還存在 ω =0 和 π 的中心和鞍點兩類平衡點.

圖23 系統在參數空間 ( κ,) 中的分岔Fig.23 Bifurcation lines in the parameter space(κ,)

軌道的長期攝動理論在近地空間碎片的軌道分析、碎片減緩和再入預報等方面具有重要的應用.與在役衛星保持長期軌道穩定的目標正好相反,人們希望空間碎片能夠在自然攝動力的作用下盡快實現大氣再入,釋放軌道資源.空間碎片減緩在20 世紀90 年代已有研究[154-155].近年來,Wang 等[156-158]研究了通過設置發射時間調整日月第三體引力攝動的相位和利用日月長期共振等實現空間碎片壽命的減小.

3.3 Lidov-Kozai 動力學機制

Lidov-Kozai 機制是天體力學中的一個重要軌道動力學現象,描述在圓軌道第三天體引力攝動下(也被稱為層級圓型限制性三體問題),當受攝二體問題軌道(即受攝開普勒軌道)傾角大于某一特定值時,軌道偏心率和軌道傾角存在大幅耦合振蕩,且近地點幅角/軌道拱線可在平衡點附近往復振蕩,對軌道全局演化和穩定性產生重要影響.這一機制還被拓展到橢圓軌道第三天體引力攝動下的受攝二體問題,即層級橢圓型限制性三體問題,發現了偏心Lidov-Kozai 機制,將導致軌道在順行和逆行間發生翻轉.另外,在攝動第三天體位于所研究軌道內部時(即所謂外問題),Lidov-Kozai 機制仍然存在,被用于海王星外天體的軌道演化研究.

3.3.1 經典Lidov-Kozai 機制

1961年,Lidov[159]在層級圓型限制性三體問題的構型下研究了環繞地球的軌道動力學,發現經過兩次平均化后,結合系統守恒量(軌道角動量垂直分量和軌道能量),系統自由度可降為一,系統可積.此時,系統呈現出兩個有趣的動力學現象:(1)軌道偏心率和軌道傾角存在長期耦合振蕩;(2)軌道傾角大于某一臨界傾角(3 9.2°)時,近地點幅角可能在±90°附近往復振蕩,不再持續進動.系統的運動可通過兩個積分常量 (c1,c2) 完全確定

其中c2決定了拱線進動的特征.當c2＜0 時,近地點幅角處于振蕩狀態;當c2＞0 時,近地點幅角單調變化,如圖24 所示.

圖24 二次平均化模型下的兩個軌道積分常量[160]Fig.24 The two orbit integral constants under the doubly-averaged model[160]

1962年,Kozai[125]基于相同的層級圓型限制性三體問題構型,在研究木星引力攝動下主帶小行星環繞太陽的軌道運動時得到了類似的結論.因此,這一動力學現象被稱為Lidov-Kozai 機制.在該機制下,高傾角近圓軌道的偏心率將經歷長期大幅振蕩,致使軌道非常不穩定[161].接著Lidov[161]和Kozai[162]進一步考慮了中心天體的扁率攝動,發現了其對Lidov-Kozai 機制的抑制作用.而且,扁率攝動會引起Lidov-Kozai 平衡點的分岔以及新平衡點的出現,從而改變系統的穩定性特征[145,147].Lidov-Kozai 機制被廣泛地應用到行星不規則衛星、彗星、近地小行星和系外行星的形成和演化研究[163-169]以及航天器軌道動力學研究[143,145,147].例如,“熱木星”形成的原因可能是Lidov-Kozai 機制促使其軌道偏心率增長,“近星點”不斷靠近中心恒星,最終在長期的潮汐耗散作用下形成穩定的熱木星[165-166].

3.3.2 偏心Lidov-Kozai 機制

Lidov-Kozai 機制的另一個重要拓展是解除“圓型”限制,在層級橢圓型限制性三體問題的構型下,系統呈現出新的豐富的動力學現象.此時,攝動第三天體位于橢圓軌道上,系統的軸對稱性消失,軌道角動量的垂直分量不再是守恒量,系統不可積.基于二階近似的第三體引力攝動模型不能準確描述這種非對稱性,而必須考慮引力攝動的三階項[170].1969 年,Harrington[171]對三階近似模型進行積分,發現了系統的共振和混沌等動力學現象.2000 年,Ford等[172]考慮了攝動天體軌道的大偏心率(ep=0.9),發現軌道在Lidov-Kozai 長期振蕩的基礎上,還存在周期更長的振蕩,導致軌道偏心率的極大增長.2011 年,Naoz 等[173]首次發現在橢圓軌道第三體引力攝動下,軌道平面可能在順行(i＜90°) 和逆行(i＞90°)之間發生周期性翻轉,并且翻轉時伴隨著偏心率非常接近1,如圖25 所示,這被用來解釋部分逆行“熱木星”的形成.這種動力學行為被稱為偏心Lidov-Kozai 機制.

圖25 偏心Lidov-Kozai 機制下的軌道翻轉[170]Fig.25 Orbit flip of the eccentric Lidov-Kozai mechanism[170]

之后,Litchwick 和Naoz[174]與Katz 等[175]對軌道翻轉發生的條件分別進行了數值和解析研究.Li 等[176]研究了橢圓軌道第三體引力攝動下系統的混沌現象.Naoz 等[170,177]對偏心Lidov-Kozai 機制在系外行星、恒星系統和黑洞系統動力學演化方面的重要作用進行了綜述.

3.3.3 外問題中的Lidov-Kozai 機制

實際上,層級限制性三體問題,即第三天體引力攝動下的受攝二體問題,可分為外問題和內問題,如圖26 所示.因為絕大多數層級三體系統都具有內問題的幾何構型,即攝動天體位于目標天體軌道外側,所以對外問題的研究相對較少.1975 年,Ziglin[178]在層級橢圓型限制性三體問題的構型下,較早地研究了行星環繞雙恒星系統的軌道運動.之后,學者們在二階近似攝動模型下,發現了外問題中不同于傳統內問題Lidov-Kozai 機制的高傾角平衡點[179-180],可以用來解釋HD 98 800 恒星系統極軌道上穩定粒子的存在.2017 年,Naoz 等[181]將偏心Lidov-Kozai機制拓展到外問題中,同樣發現了軌道的翻轉現象.外問題的軌道攝動模型還能用于研究行星外離散盤、奧爾特云和海王星外天體等的軌道演化[182-183].

圖26 層級三體問題中的內問題與外問題Fig.26 The inner and outer problem in the hierarchy three-body system

3.4 環繞雙小行星系統主星的軌道演化

2020年,Wang 和Fu[40]將前述受攝二體問題的理論和方法拓展到雙小行星系統中,研究環繞雙星系統單顆星的軌道長期演化.與上述經典的受攝二體問題不同,雙小行星不僅為中心天體引入了非球形引力,而且為第三天體引力攝動也引入了非球形項.此外,中心天體和攝動天體極近的距離一方面導致第三天體引力攝動的高階項非常顯著,傳統希爾近似不再適用,需保留高階項,另一方面使得中心天體非球形引力攝動的主導區域和第三天體引力攝動的主導區域高度重合,導致了不同攝動力的強烈耦合.在筆者的研究中,系統的幾何構型如圖27 所示,主星為探測器所環繞的中心天體,次星為攝動第三天體,假設次星軌道為圓形,即為層級圓型限制性三體問題.在動力學建模中,首次在第三天體引力攝動函數中包括了主星和次星的非球形引力項.針對第三天體引力攝動高階項顯著的特點,將次星引力攝動勢函數保留至四階,其形式為[40]

圖27 雙小行星的二體構型和環繞主星的軌道[41]Fig.27 The two-body configuration of the binary asteroid system and orbits around the primary[41]

式中,τ0A表示主星的扁率作用項,τ0B和τ2B表示次星的二階球諧引力項.從上式可以發現,主星和次星的非球形引力項分別出現在第三天體引力攝動的三階項和四階項中.

將四階近似的軌道運動方程針對航天器軌道和次星軌道進行兩次平均化,發現由于系統具有軸對稱性,軌道角動量的垂直分量守恒,系統的自由度可降為一,系統可積.基于二次平均化的軌道運動方程具有形式

式中,fe和fω分別表示相應變量的函數,具體的表達式可參見文獻[41];第三式表示軌道角動量的垂直分量守恒,其數值可由初始條件確定.

此外,需要注意的是,在行星系統中較少討論的關于平均化有效條件的問題,在雙小行星系統中則會凸顯:較強的攝動會導致軌道在一個平均化周期內較大程度地偏離開普勒軌道,使得平均化失效.Wang 和Fu[40]利用數值方法研究了保證平均化有效的攝動力大小需滿足的限制條件.在較高精度要求下,可近似認為中心天體扁率和第三體引力攝動大小與中心天體質點引力大小的比值需小于10-2量級.另外,航天器軌道運動與次星軌道運動不發生共振,且航天器軌道進動周期相對次星軌道周期較大.

在基于方程(36)描述的動力學模型下,能夠進一步分析包括軌道穩定性在內的動力學特性,發現在Lidov-Kozai 機制和中心天體扁率攝動的共同作用下,J2臨界傾角附近的軌道具有較強的不穩定性,如圖28 所示[41],但對于這類軌道仍然可以通過將軌道初始近地點幅角置于“長駐留點”(滿足)處,以延長軌道壽命[41,150],如圖29 所示;而高傾角和低傾角軌道近地點高度振蕩的幅度較小,且可通過設置軌道的初始條件來維持軌道的小偏心率狀態.

圖28 J2 臨界傾角附近的 ω -e 相空間軌跡[41]Fig.28 ω -e phase trajectories around the J2 critical inclination[41]

圖29 不穩定軌道的軌道壽命[41]Fig.29 The orbital lifetime of unstable orbits[41]

上述針對環繞主星軌道的研究是在次星圓軌道且次星位于主星赤道面內的條件下進行的,實際的雙小行星系統中次星軌道可能還具有偏心率,且次星軌道平面可能傾斜于主星赤道平面,這兩個因素都將破壞系統的對稱性,使得軌道角動量的垂直分量不再是守恒量,系統不可積,包括混沌和高階共振在內的非線性動力學特征將對軌道運動產生顯著影響.最新的研究表明,在考慮次星軌道偏心率的層級橢圓型限制性三體問題中,探測器軌道平面存在不同于傳統偏心Lidov-Kozai 機制的軌道翻轉,如圖30所示.可以發現,與偏心Lidov-Kozai 機制相比(如圖25),中心天體的J2攝動一方面會導致允許翻轉的軌道傾角范圍縮小,但另一方面導致允許翻轉的軌道偏心率降低,使得軌道可在更低的偏心率下完成軌道平面的翻轉,拓寬了軌道翻轉的場景.該動力學特征不僅可用于解釋雙小行星系統附近粒子特異的軌道演化行為,對于揭示系外行星系統形成和演化機制也具有重要的應用價值.

圖30 J2 攝動的層級橢圓型限制性三體問題下的軌道翻轉Fig.30 Orbit flip in the hierarchical elliptical restricted three-body problem with the J2 perturbation

4 面向任務的軌道動力學分析與設計

前面章節介紹的研究大多是基于簡化力學模型的理論研究.在實際的雙小行星探測任務中,探測器軌道設計不但要考慮更加真實的力學環境,還要滿足任務目標和探測器平臺性能等諸多約束條件.目前,以雙小行星系統為目標的任務中,除了已發射的DART 探測器(AIDA 第一階段),另有處于立項研制階段的Hera 探測器(AIDA 第二階段)和Janus 任務,它們具有各自明確的探測對象、任務目標和探測器屬性.本節主要介紹針對這些任務需求,在實際力學環境和相關的約束下,對探測器軌道動力學的分析和軌道方案的設計.

4.1 任務軌道設計研究

4.1.1 基于周期軌道的任務軌道設計

周期軌道是深空探測活動中具有重要應用價值的軌道,不同類型的周期軌道能夠滿足特定的任務需求,同時能夠以較少的燃料消耗進行軌道維持.在雙小行星系統附近,周期軌道提供了可靠的探測軌道,是探測任務軌道設計的重要參考.但是復雜的攝動力可能增加周期軌道的不穩定性,因而實際的任務軌道需要在接近真實的力學環境中進一步改進計算得到,并進行軌道穩定性的分析.

2017年,Dell’Elce 等[31]針對AIDA 探測任務,找到了兩類具有魯棒性的周期軌道——圍繞主星的內部逆行軌道和圍繞整個雙星系統的晨昏軌道,內部逆行軌道適合對主星的觀測,晨昏軌道提供了對雙星大范圍的觀測,如圖31 所示.Capannolo 等[30]在高精度的引力場模型下,研究了系統附近多種類型的周期軌道,可滿足不同的任務需求,例如多種尺寸的遠距離逆行軌道(DRO)可實現穩定且大范圍的次星環繞,L1和L2平動點軌道能夠實現對次星局部區域的近距離觀測,如圖32 所示.

圖31 繞雙星65803 Didymos 晨昏軌道(從太陽方向看)[31]Fig.31 Terminator orbits around the binary asteroid system 65803 Didymos seen from the Sun direction[31]

圖32 雙小行星系統65803 Didymos 附近的周期軌道族[30]Fig.32 Periodic orbits in the binary asteroid system 65803 Didymos[30]

周期軌道的求解可基于不同精度水平的動力學模型,模型的精度越高,計算得到的周期軌道就越接近真實軌道,軌道維持所需的燃料消耗也往往更少[184].Shi 等[35]以雙小行星系統66391 Moshup(即1999 KW4)為例,比較了各類軌道維持方法的優劣和有效性,計算了各類周期軌道的軌道維持所需的燃料消耗,其結果可為任務軌道設計提供參考.

4.1.2 面向特定任務需求的任務軌道設計

AIDA 任務第一階段是由DART 探測器抵達雙小行星系統65803 Didymos 并對次星實施撞擊.然而若只有單一的探測器,對撞擊后的觀測只能由地面遠程開展.2018 年,Capannolo 等[185]提出由DART 探測器攜帶并釋放一個微納衛星用于撞擊后的觀測,以獲得更多的科學數據.為了滿足對撞擊坑和噴濺物觀測的良好條件,在微納衛星安全性、相機性能等諸多約束條件下,對釋放立方星的時刻、相對速度大小和方向進行了優化設計,驗證了該方案的可行性.最新進展表明,DART 探測器攜帶了一個由意大利航天局(ASI)負責的6 U 立方星LICIACube,其主要任務是實施對撞擊坑和噴濺物的觀測和對次星形狀的測量等[186].為此,Capannolo 等[187]針對LICIACube 的具體參數,結合撞擊坑模型和噴濺物擴散模型,考慮平臺和任務目標的約束條件,將軌道設計問題轉換為以立方星釋放方向為優化變量的優化問題,得到了相應的軌道方案,并且對不確定性因素進行了魯棒性分析.

AIDA 任務第二階段由Hera 探測器于撞擊四年后到達該雙小行星系統,對撞擊坑和軌道偏轉等特征進行觀測和測量.Hera 任務的實施可分為遠距離初步測量階段和近距離詳細探測階段,G i l-Fernandez 等[188]根據任務特定需求,分別分析和設計了滿足任務約束(如脈沖機動頻率、脈沖大小、軌道高度、光照條件和安全裕度等)的遠距離和近距離雙曲線飛掠軌道,如圖33 所示,圖中閉合的軌跡由幾段雙曲線軌道連接構成,可通過不同的組合調整近星點高度和次星的緯度覆蓋范圍.基于類似的方法,Ferrari 等[189]對僅覆蓋雙星系統光照面的兩點和多點往復飛掠軌跡進行了詳細的研究,分析了近星點高度、所需的推進劑消耗和軌道機動頻率,以此為依據可選擇滿足任務需求的軌道.

圖33 Hera 探測器的雙曲線飛掠軌跡設計Fig.33 The hyperbolic fly-by trajectories design of Hera mission

4.2 轉移軌道與彈道著陸軌道研究

在雙小行星探測任務中,由于雙星系統附近強攝動的復雜力學環境,主探測器靠近雙星可能存在較大的風險.為了實現近距離和表面探測,有科學家提出采用成本更低、體積更小的微納衛星方案,例如AIDA 任務中DART 探測器攜帶的6 U 立方星LICIACube;Hera 探測器計劃攜帶的兩個6 U 立方星,其中一個為著陸器.因為微納衛星的軌道機動能力有限,所以基于不變流形的低能轉移和著陸軌道成為必要的選擇.

2013年,Tardivel 和Scheeres[190]在雙星視為質點的經典圓型限制性三體問題模型下,研究了利用平動點附近動力學特征的無控著陸軌道,并將理論應用到MarcoPolo-R 任務(已取消)的備份目標雙小行星系統1996 FG3 中,設計了可靠性較高的L2平動點到次星表面的彈道著陸軌道[191].2016 年,Wu 等[192]在雙小行星系統90 Antiope 的雙同步橢球模型下,將雙星的表面選擇為龐加萊截面,通過穩定和不穩定流形拼接實現雙星之間的軌道轉移,如圖34 所示.Celik 和Sánchez[193]采用從小行星表面著陸點開始的反向積分方法研究了主星和次星的彈道著陸軌跡,該方法具有很強的靈活性,可以在不同約束下快速規劃出滿足約束條件的著陸軌道.

圖34 雙小行星系統90 Antiope 成員之間的轉移軌道[192]Fig.34 The spatial transfer orbits between the two bodies of the binary asteroid system 90 Antiope[192]

5 總結與展望

雙小行星系統的獨特構型和復雜引力場環境為探測任務設計帶來了極大挑戰,首當其沖的就是強攝動耦合環境中的軌道動力學特性分析與任務軌道設計.經過近三十年的研究,人們對雙小行星系統附近的軌道動力學已有了較為全面的認識,但仍有一些問題需要進一步的研究和解決.

首先,雙小行星系統的不規則引力場建模是探測器軌道動力學研究的基礎.目前已經形成了比較成熟的描述雙星不規則引力場的方法,如質點群法、球諧函數法、多面體方法,以及基于特定簡單幾何體的引力場近似方法.在這些方法中,計算簡單和高精度往往是一對矛盾,未來可針對特定場景將不同建模方法結合,形成兼顧計算效率和精度的混合方法.另外,在探測器達到目標前,小行星準確的形狀和質量信息未知,難以對其引力場進行精確建模,需要依靠探測器逐步抵近和測量,不斷提高引力場建模精度,才能確定最終的軌道方案.因此,未來需要研究通過星載計算機的有限資源,實現對雙小行星系統高效、精確甚至實時的引力場測量和計算.在雙星系統相互引力勢方面,目前已經形成了基于不同原理的建模方法,可以滿足雙星相對運動的理論和數值研究,未來可針對如何提高計算效率進行進一步的改進和提高.在雙星相對運動方面,目前已較為全面地揭示了平衡構型及其附近的動力學特性,由于大部分雙星系統處于平衡構型附近,因此目前可以滿足探測器軌道動力學研究的需要,未來可針對更精細或更長期的動力學演化特性開展研究.

基于限制性三體問題的雙小行星系統附近軌道動力學研究目前主要包括動力學建模、周期軌道求解與分析、轉移軌道設計和周期軌道維持等內容.動力學建模方面,小行星引力場模型的精度從質點/球模型、橢球模型、球諧函數模型到多面體模型精度不斷提高,逐漸逼近小行星的真實引力場.周期軌道方面,學者們對平動點及平動點周期軌道、共振周期軌道、圍繞某顆星或整個雙星系統的周期軌道進行了大量研究,主要關注了軌道族計算,以及穩定性和分岔等動力學特性.轉移軌道方面,目前已在不同動力學模型下通過流形拼接等方法實現了不同區域間的轉移,但大部分研究未考慮實際任務的需求和約束.周期軌道維持方面的研究主要包括連續推力和脈沖推力兩種形式,但目前所嘗試的維持策略種類較少,未來可進一步借鑒經典限制性三體問題中的新型維持策略,并根據雙小行星系統的特點加以改進,得到性能更優的軌道維持策略.

限制性三體問題的周期軌道較不穩定,形成條件苛刻,幾何構型和分布較為固定,一定程度上限制了其在長期、近距離探測中的應用潛力.而基于受攝二體問題的環繞雙小行星單顆星的軌道穩定性相對較好,形成條件較為寬松,分布更為廣泛,并易于實現對小行星的全球覆蓋,在長期、近距離探測和穩定停泊方面具有很大的應用價值,因而在最近受到了關注.已有的研究主要包括環繞主星的半解析軌道動力學建模和在此基礎上進行的軌道穩定性分析等.目前,該方向的研究剛剛起步,很多問題尚未解決,未來需要系統性地發展.例如,未來可研究小行星非對稱非球形引力項、太陽光壓力、太陽引力和次星軌道偏心率等對軌道穩定性的影響;還可將環繞單顆星的內問題拓展到環繞整個雙星系統的外問題,推動對雙小行星系統附近軌道動力學的理解,為探測任務提供更多的軌道類型選擇.

進入21 世紀的第三個十年,雙小行星探測將進入實踐階段.針對目前已經立項研制的雙小行星探測任務,學者們在探測目標、任務要求和實際約束等因素的綜合考量下,進行了面向實際任務的軌道動力學分析與軌道設計.周期軌道能夠提供多種類型的觀測軌道,可滿足不同的任務需求,結合軌道維持策略,能夠保證軌道的安全性和可靠性,具有一定的應用潛力.另外,針對低成本微納衛星變軌機動能力的限制,研究了基于不變流形的低能轉移和著陸/起飛軌道.實際上,目前處于工程實施階段的探測任務出于安全性的考慮,對主探測器都擬采用飛掠或往復飛掠軌道,極大限制了探測效率和探測精度.因此,未來需要基于雙小行星附近軌道動力學研究的最新成果,在復雜力學環境和多種任務約束下,提出探測效率、探測精度和安全性等綜合性能更優的軌道設計方法和控制技術.

綜上所述,目前對雙小行星探測軌道動力學的研究已經形成了比較成熟的理論體系,但某些研究方向剛剛開展,尚待系統性地發展,在比較成熟的研究方向,也仍有一些理論和技術問題尚待解決或完善.此外,需要將探測任務設計和軌道動力學的最新研究成果更好融合,從而提升探測任務的綜合效能.隨著雙小行星探測任務的不斷推進,在檢驗現有理論、方法和技術的同時,還將涌現出更多的科學和技術問題,將不斷推動軌道動力學理論、軌道設計方法和軌道控制技術的發展.