祁 瑞,徐世杰

(北京航空航天大學宇航學院，北京 100191)

相比限制性二體模型，限制性三體模型因其對空間引力資源的更深入利用，在深空探測軌道設(shè)計方面表現(xiàn)出巨大的應用優(yōu)勢.特別是平面圓型限制性三體問題(PCR3BP，planar circular restricted 3-body problem)共線平動點附近周期軌道的不變流形，以其漸近性質(zhì)和運動分界特點被廣泛用于平動點附近編隊飛行、低能奔月以及木星多顆衛(wèi)星探測等相關(guān)任務研究，并取得了一系列有價值的理論與工程實踐結(jié)果[1-2].然而，這些研究只適用于自治系統(tǒng)PCR3BP，對于橢圓限制性三體問題(ER3BP，elliptic restricted 3-body problem )，由于不存在合適的參考系使系統(tǒng)擺脫對時間的顯性依賴，相關(guān)研究難以推廣.

對于不變流形的計算，傳統(tǒng)的方法是用不變特征向量方向上的一條短線段來代替局部不變流形，對該線段上的點進行迭代來確定全局不變流形.然而有三方面的原因會導致數(shù)值計算出的不變流形位置不精確，時間跨度很大時尤其如此.首先，局部不變流形是利用直線段來近似的；其次，點的數(shù)值迭代有舍入誤差；最后，進行高次迭代時，曲線上的點會漸漸散開.Parker和Chua[3]給出了一個較好的算法，通過沿不變流形改變點的數(shù)量，使這些點均勻分散在流形上.對于ER3BP中的時間周期不變流形，由瞬時不變集出發(fā)采用傳統(tǒng)方法進行計算是很自然的選擇.然而，Shadden等人[4]對雙漩渦流場的研究結(jié)果表明，該時間周期系統(tǒng)的不變流形并不與瞬時不動點相連接.這說明，不變流形的傳統(tǒng)算法在ER3BP系統(tǒng)中不再適用.

在流體力學領(lǐng)域，學者采用拉格朗日擬序結(jié)構(gòu)(LCS，Lagrangian coherent structures)作為時變流場中的不變流形來研究動力系統(tǒng)的相空間.在較早期的文獻中，LCS的定義比較模糊，給問題分析和算法設(shè)計造成不便.2001年，Haller[5]提出將LCS定義為有限時間Lyapunov指數(shù)(FTLE，finite-time Lyapunov exponent )域中的脊.有限時間Lyapunov指數(shù)是經(jīng)典Lyapunov指數(shù)的變形，用于度量系統(tǒng)對初值的敏感依賴性，最早由Lorenz引入研究大氣模型中的混沌現(xiàn)象.近年來，Gawlik[6]，Qi[7]等人將LCS引入天文動力學的研究并取得成功.

本文首先以單擺系統(tǒng)為例闡述了LCS與穩(wěn)定、不穩(wěn)定流形之間的關(guān)系，進而將LCS理論應用于天體力學的研究中，以質(zhì)量比為0.1、偏心率為0.2的虛擬系統(tǒng)為例，研究了ER3BP中的時間周期不變流形，通過數(shù)值方法研究了其運動分界面本質(zhì)和軌道不變特性.

1 相關(guān)的基礎(chǔ)理論

1.1 橢圓限制性三體問題模型

假設(shè)天體M1和M2圍繞其公共質(zhì)心作橢圓運動，航天器在M1和M2的引力合力作用下運動，且航天器質(zhì)量足夠小以至于不會影響兩個主天體的橢圓運動.

引入無因次量綱：長度單位DU為M2真近點角f達90°時的兩主天體間距；質(zhì)量單位MU為兩主天體質(zhì)量之和；時間單位TU為兩主天體軌道周期與2π的比(本文中出現(xiàn)的數(shù)量，未經(jīng)特殊聲明，均為無因次量綱.).

建立M1-M2質(zhì)心旋轉(zhuǎn)脈動系：原點位于M1和M2質(zhì)心；x軸由M1指向M2；y軸由x軸在軌道平面內(nèi)逆時針旋轉(zhuǎn)90°確定；z軸垂直于軌道面完成右手系；通過3個軸同等比例的實時伸縮，保持M1和M2質(zhì)心在x軸上靜止.

為了使ER3BP下的航天器運動方程與CR3BP下的方程形式相似，以f為系統(tǒng)自變量.在上述坐標系和無因次量綱下，容易得到航天器運動微分方程[8]：

(1)

式中：

(2)

1.2 拉格朗日擬序結(jié)構(gòu)

這里采用Haller的定義，將LCS定義為有限時間Lyapunov域中的脊線[5].粗略地說，F(xiàn)TLE就是在相流作用下，鄰近的兩點在有限時間內(nèi)的平均最大分離速率，度量了系統(tǒng)對于初值的敏感依賴性.

在開區(qū)域D?Rn上考查動力系統(tǒng)：

(3)

在相流映射φ的作用下，t0時刻的點x在t0+T時刻到達φ(t0+T;t0;x).假設(shè)x受到無窮小攝動δx(0)，則該攝動在t0+T時刻的變分為：

δx(T)=φ(t0+T;t0;x+δx(0))-φ(t0+T;t0;x)

(4)

只保留線性部分，得到：

(5)

(6)

其中，λmax是矩陣Δ的最大特征值，ξ是相應于λmax的特征向量.

此時，定義有限時間Lyapunov指數(shù)為：

(7)

從而，式(6)可改寫為：

(8)

對于整個定義域D，由式(7)對每一點都賦以一標量值，這就得到了FTLE域.由動力系統(tǒng)理論可知，不變流形位于不同側(cè)的鄰近兩點會以指數(shù)速率快速分離，即不變流形對應于較大的FTLE.為此，將LCS定義為FTLE域中的脊線.分析可知，對于T>0，LCS蘊含穩(wěn)定流形；對于T<0，LCS蘊含不穩(wěn)定流形.

2 單擺相空間中的LCS

為了驗證上述分析，考查單擺系統(tǒng)：

(9)

在經(jīng)典的動力系統(tǒng)理論中，已經(jīng)得到了系統(tǒng)(9)的相圖(見圖1(a)).(±π,0)是兩個鞍點；(0,0)是中心.軌道a和b是連接兩個鞍點的異宿軌道；軌道c是環(huán)繞中心的周期解；軌道d是旋轉(zhuǎn)解.顯然，異宿軌道是單擺相空間中的運動分界面，將所有軌道劃分為周期解和旋轉(zhuǎn)解.

圖1 單擺系統(tǒng)的相空間及其FTLE域(T=10)

在考查的定義域內(nèi)取格點劃分為300×300，T=10.對每一格點計算Lyapunov指數(shù)，并在圖1(b)中以等高圖繪出FTLE域.在圖中可以清晰地看出脊線結(jié)構(gòu)，即LCS，它們是趨近于兩個鞍點的穩(wěn)定流形.類似地，當取T=-10，我們可以得到不穩(wěn)定流形.

3 時間周期不變流形的性質(zhì)

在PCR3BP模型下，Koon等人[2]的研究表明，共線平動點附近的不變流形管作為軌道的不變集，是穿越軌道與非穿越軌道的分界面.本節(jié)將通過數(shù)值方法驗證，ER3BP模型下的時間周期不變流形依然具有上述性質(zhì).

下面的數(shù)值仿真研究基于系統(tǒng)參數(shù)μ=0.1和e=0.2.

3.1 時間周期不變流形是運動分界面

圖2 U1截面處的FTLE域(T=2,f=π/2)

從上圖中可以清楚地看到一條封閉脊線，即LCS，這是時間周期不變流形在U1截面上的截交線.如圖所示，在該截線附近隨機布置試驗點陣.其中，位于不變流形截線外部的點以實心圓點表示，以其為初狀態(tài)的軌道示于圖3(a)；位于不變流形截線內(nèi)部的點以五角星符號表示，相應的軌道示于圖3(b).

從圖3(a)和圖3(b)中易見，由不變流形截線外部出發(fā)的軌道在M1-M2系L1點區(qū)域附近折返回M1，而由截線內(nèi)部出發(fā)的軌道穿越了L1平動點區(qū)域，達到了M2附近.這說明，時間周期不變流形依然是穿越軌道與非穿越軌道的分界面.

3.2 時間周期不變流形是軌道不變集

在使用流形生長法對PCR3BP中的不變流形(以穩(wěn)定流形為例)進行計算時，由Lyapunov軌道上沿穩(wěn)定流形方向的微小偏移點出發(fā)，反向積分得到一條位于穩(wěn)定流形上的軌道.實際上，這種處理本身已經(jīng)蘊含了不變流形是軌道不變集的結(jié)論.

對于時間周期不變流形，由于難以直觀地表達其整體形狀，為了證明類似的結(jié)論，本文將采用一系列Poincare截面進行采樣.

圖3 非穿越軌道和穿越軌道

該典型軌道進入M1-M2系平動點區(qū)域，并表現(xiàn)出漸近收斂于準周期軌道的趨勢.考察該軌道與截面x=0.3、x=0.4和x=0.5的截交，以截交點處的軌道能量確定出相應的Poincare截面：

(10)

(11)

(12)

在對應的真近點角下，各截面上的FTLE域等高圖及該軌道在各截面的截交點繪于圖5(a)～(c)中.

由圖5(a)～(c)可以看到，在各截面上，軌道的交點恰好處于相應的LCS上.這反映出，該條軌道恰好落在時間周期不變流形上，從而表明時間周期不變流形是軌道的不變集.

圖4 U1截面處的FTLE域(T=2,f=π/2)以及一條典型軌道

圖5 U2、U3和U4截面處的FTLE域

4 結(jié) 論

本文對于橢圓三體模型中的時間周期不變流形得到了如下結(jié)論：

(1)時間周期不變流形是穿越軌道與非穿越軌道的分界面；

(2)時間周期不變流形是軌道不變集；

本文中利用LCS研究時間周期不變流形的方法可以推廣到對時間任意依賴的動力系統(tǒng)中.

參 考 文 獻

[1]Qi R, Xu S J, Xu M. Impulsive control for formation flight about libration points[J]. Journal of Guidance, Control, and Dynamics, 2012, 35(2): 484-496

[2]Koon W S, Lo M W, Marsden J E, Ross S D. Heteroclinic connections between periodic orbits and resonance transitions in celestial mechanics[J]. Chaos, 2000, 10: 427-469

[3]Parker T S, Chua L O. Practical numerical algorithms for chaotic systems[M]. New York: Springer-Verleg, 1989

[4]Shadden S C, Lekien F, Marsden J E. Definition and properties of Lagrangian coherent structures from finite-time Lyapunov exponents in two-dimensional aperiodic flows[J]. Physica D, 2005, 212: 271-304

[5]Haller G. Distinguished material surfaces and coherent structures in three-dimensional fluid flows[J]. Physica D, 2001, 149: 248-277

[6]Gawlik E S, Du Toit P C, Campagnola S, et al. Lagrangian coherent structures in the planar elliptic restricted three-body problem[J]. Celestial Mechanics and Dynamical Astronomy, 2009, 103: 227-249

[7]Qi R, Xu S J, Zhang Y, Wang Y. Earth-to-moon low energy transfer using time-dependent invariant manifolds [C]. AIAA/AAS Astrodynamics Specialist Conference, Minneapolis, Minnesota， August 13-16, 2012

[8]Szebehely V G. Theory of orbits: the restricted problem of three bodies[M]. New York: Academic Press, 1967