楊洪偉 陳楊 寶音賀西 李俊峰

（清華大學(xué)航天航空學(xué)院，北京 100084）

1 引言

在深空探測任務(wù)中，往往需要通過行星引力輔助來達到減少燃料消耗的目的。在利用行星引力輔助的深空探測任務(wù)分析與初步軌道設(shè)計中，通常采用等效脈沖引力輔助模型（下文簡稱“簡化模型”），將行星引力輔助效果看作是一個瞬時的速度脈沖［1］，通過圓錐曲線拼接技術(shù)［2］設(shè)計軌道。但是在實際情況下，引力輔助是連續(xù)的過程：航天器在行星影響球內(nèi)飛行需要一定的時間，飛行的過程中會受到太陽和其他行星的引力等攝動，并且這段時間內(nèi)行星的位置速度會發(fā)生變化。這些因素會導(dǎo)致引力輔助結(jié)束時刻航天器實際的狀態(tài)量與簡化模型中不同。而引力輔助結(jié)束后航天器還要長時間的飛行，初始時刻狀態(tài)量很小的不同就能引起終端時刻很大的響應(yīng)，因此在實際的工程應(yīng)用中必須采用精確的引力輔助模型設(shè)計深空探測器飛行軌道。此外，現(xiàn)在大多數(shù)研究都是基于簡化模型，通過精確模型和簡化模型計算結(jié)果對比可以驗證簡化模型的準確性和可信度。

針對精確動力學(xué)模型中的引力輔助軌道設(shè)計問題，文獻［3］給出了在簡化引力輔助模型基礎(chǔ)上引入行星攝動設(shè)計精確軌道的方法，但是其方法存在影響球內(nèi)外引力攝動不連續(xù)問題并且過程復(fù)雜。為了將簡化模型的設(shè)計結(jié)果迭代至精確動力學(xué)模型中，需要描述航天器在引力輔助行星附近的軌道和飛行時間。而以行星為中心的B 平面參數(shù)［4］就能夠準確反映航天器相對行星的雙曲線軌道的大小和空間取向，并且在打靶修正時B平面參數(shù)與修正脈沖有良好的線性關(guān)系。文獻［5－6］采用B平面參數(shù)研究了火星探測的精確軌道設(shè)計。本文將脈沖等效模型求得的引力輔助參數(shù)化為目標B平面約束，使用微分修正法改進地心逃逸速度，最終得到了精確模型中的引力輔助軌道。

2 引力輔助簡化模型下的軌道計算

在用圓錐曲線拼接法近似設(shè)計軌道的時候，采用引力輔助等效脈沖模型，即忽略行星影響球半徑，且認為引力輔助使航天器獲得瞬時速度脈沖。文獻［7］在設(shè)計木星探測軌道時提出了引力輔助簡化模型下軌道計算的有效方法。引力輔助前后，航天器的日心位置與引力輔助行星位置rPl相同。

航天器引力輔助時相對于引力輔助行星的軌道為雙曲線軌道，雙曲線進入剩余速度和離開剩余速度大小相等

引力輔助對航天器產(chǎn)生的等效脈沖ΔvGA為

雙曲線剩余速度的偏轉(zhuǎn)角為δ根據(jù)下式確定

式中 μPl為引力輔助行星引力常數(shù)；rp為引力輔助半徑。

圖1 引力輔助坐標系Fig.1 Gravity assist reference frame

為了求得引力輔助給航天器提供的速度增量和引力輔助后航天器的速度，建立行星心引力輔助坐標系P－ξηζ（見圖1）。坐標系原點建立在引力輔助行星的中心，ξ軸沿進入雙曲線剩余速度方向，ζ軸垂直于由引力輔助前航天器的日心速度v－與引力輔助行星的日心速度vPl所決定的平面，η 軸與ξ、ζ軸構(gòu)成右手坐標系。三坐標系的單位矢量分別記為i、j和k。

由雙曲線軌道的性質(zhì)可知，離開引力輔助行星的雙曲線剩余速度在以進入雙曲線剩余速度為軸、頂角為δ的圓錐面上，且大小等于v∞。根據(jù)幾何關(guān)系可以得到

3 精確模型引力輔助軌道設(shè)計

3.1 目標B平面參數(shù)

圖2 B平面示意［5］Fig.2 B plane

航天器經(jīng)過行星引力輔助時的雙曲線軌道可以通過B平面來描述。B平面參數(shù)在進行軌道微分修正的時候具有良好的收斂性。B 平面定義為過引力輔助行星中心并垂直于接近引力輔助行星的雙曲線漸近線的平面［4］（見圖2）。

S為射入漸近線方向的單位矢量，即

為了方便起見，參考方向單位矢量N 取黃道平面法線方向。T 矢量和R 矢量定義為

B 矢量由引力輔助行星中心指向B 平面與射入漸近線的交點，其大小等于雙曲線軌道的半虛軸。B平面參數(shù)BT和BR定義如下：

根據(jù)引力輔助等效脈沖模型確定的雙曲線軌道進入和離開雙曲線剩余速度即可計算微分修正時的B平面參數(shù)。首先，雙曲線軌道平面法向量n為

雙曲線軌道半長軸為

B 矢量的大小等于雙曲線軌道的半虛軸大小b

因此，目標B 矢量可以表示為

從而可以求得目標約束B平面參數(shù)BTtar和BRtar。

3.2 實際B平面參數(shù)

微分修正時可以根據(jù)進入行星影響球時的位置r、速度v求出實際S 矢量，進而求出實際B平面參數(shù)

3.3 時間參數(shù)

在精確模型引力輔助軌道設(shè)計時，通常將脈沖等效模型搜索得到的引力輔助時間tGA認為是航天器實際飛行到達近行星點的時刻。但是由于模型不同，在精確動力學(xué)模型中積分至tGA時刻求得航天器飛行軌道，航天器相對行星的雙曲線軌道一般并沒有到達近行星點或已經(jīng)超過近行星點。航天器的當前狀態(tài)距離近行星點的時間可以根據(jù)其位置、速度信息計算，公式為

式中 H為雙曲角。

3.4 微分修正

建立描述航天器深空飛行的精確動力學(xué)模型，考慮太陽引力、各大行星的引力和太陽光壓力。引力輔助微分修正的設(shè)計變量選為地心逃逸速度vp，約束變量有目標B 平面參數(shù)BTtar和BRtar，另外，為了保證脈沖搜索得到的引力輔助時間為近行星點時間，需要設(shè)定航天器積分結(jié)束時距離近星點時間約束τtar為0。約束變量的響應(yīng)和設(shè)計變量的擾動所滿足的函數(shù)關(guān)系線性化后有

式中 K為約束變量對設(shè)計變量的敏度矩陣。由于精確動力學(xué)模型較為復(fù)雜，敏度矩陣K 通常根據(jù)數(shù)值方法求得。那么可以得到為了消除末端誤差從而需要對設(shè)計變量施加的擾動

在滿足引力輔助約束后，由于存在模型誤差，航天器還無法如二體問題設(shè)計結(jié)果那樣飛到目標天體，因此，還需要進行第二次以任務(wù)目標星體為約束的軌道修正。若探測目標為小行星，末端約束應(yīng)為航天器和小行星的日心位置誤差；若探測目標為火星等大行星，則末端約束應(yīng)設(shè)定為目標火心軌道傾角、近火點或火星中心B 平面參數(shù)。以第一次引力輔助修正結(jié)果為初值，可以較快得到滿足精度的設(shè)計結(jié)果。

4 算例及分析

本節(jié)以地球出發(fā)利用火星甩擺交會8號小行星花神星的軌道為例，說明本文設(shè)計方法的有效性。8號花神星是1顆主帶小行星，也是與太陽距離最近的較大的主帶小行星?；ㄉ裥堑能壍酪厝绫?所示。

表1 J2000日心黃道坐標系下花神星的軌道要素（MJD 55600）Tab.1 Orbit elements of Flora in J2000heliocentric ecliptic reference frame（MJD 55600）

4.1 等效脈沖模型軌道初步設(shè)計

在進行初步軌道設(shè)計時，基于日心二體模型，引力輔助采用等效脈沖模型。設(shè)計變量選為地球出發(fā)時間t0、火星甩擺時間tM、小行星到達時間tf以及火星引力輔助變量rp和ψ，優(yōu)化搜索算法為粒子群算法（PSO）。火星甩擺時的最低高度設(shè)為300km，出發(fā)時間選為2020年1月1日（MJD 58849）至2025年12月31日（MJD 61040）。軌道初步設(shè)計結(jié)果如表2所示。

表2 簡化模型引力輔助初步設(shè)計結(jié)果（J2000日心黃道坐標系）Tab.2 Preliminary design result by the simplified model（J2000heliocentric ecliptic reference frame）

4.2 精確軌道迭代設(shè)計

選擇地心停泊軌道高度為300km，赤道傾角為22°。計算得到的地心逃逸軌道參數(shù)如表3所示。

表3 地心軌道設(shè)計結(jié)果（J2000地心黃道坐標系）Tab.3 Design result of earth－centred orbit（J2000earth－centred ecliptic reference frame）

以逃逸雙曲線近地點速度為設(shè)計變量，兩次微分修正結(jié)果如表4、表5所示。

表4 第一次微分修正結(jié)果（J2000火心黃道坐標系）Tab.4 Result of first differential correction（J2000 Mars－centred ecliptic reference frame）

表5 第二次微分修正結(jié)果（J2000日心黃道坐標系）Tab.5 Result of second differential correction（J2000 Mars－centred ecliptic reference frame）

精確模型中航天器的日心飛行軌道如圖3所示。航天器到達火星引力輔助時相對火星的雙曲線軌道如圖4所示。

圖3 日心軌道（J2000日心黃道系）Fig.3 Sun－centred trajectory（J2000heliocentric ecliptic reference frame）

圖4 火星甩擺時相對火星的軌道（J2000火心黃道系）Fig.4 Orbit during Mars gravity assist（J2000 Mars－centred ecliptic reference frame）

4.3 結(jié)果分析

由表4和表5可以看到，兩次微分修正的迭代次數(shù)較小，說明本文設(shè)計含精確引力輔助軌道方法的高效性。其次，從圖3可知本文的設(shè)計方法能夠使得航天器在精確動力學(xué)模型下成功實現(xiàn)引力輔助，并在末端時刻精確到達目標位置。此外，根據(jù)軌道設(shè)計結(jié)果，簡化模型脈沖搜索、第一次微分修正和第二次微分修正所得到的設(shè)計變量，即地心逃逸近地點速度大小，有所不同，但差別均非常?。ㄒ姳?～5）。這說明簡化等效脈沖模型精度較高，使用其對引力輔助過程進行近似有較高的合理性和可信度。

5 結(jié)束語

本文給出了精確動力學(xué)模型下，通過微分修正設(shè)計含引力輔助軌道的一種有效方法。精確模型相比于傳統(tǒng)的簡化模型，在動力學(xué)方程中引入了攝動，并將引力輔助看作連續(xù)的過程，符合工程實際情況。另外，通過算例分析可知，精確模型與簡化模型計算參數(shù)相差不大，說明使用引力輔助等效脈沖模型的簡化模型精度較好，且以簡化模型結(jié)果為基礎(chǔ)，通過微分修正幾次迭代即可得到滿足目標約束的精確引力輔助模型下的結(jié)果。

在本文的基礎(chǔ)上，可以進一步研究用B 平面設(shè)計多次精確引力輔助探測軌道方法。此外，精確引力輔助模型下航天器小推力軌跡優(yōu)化也是下一步研究內(nèi)容。

［1］SIMS J A.Delta－V Gravity－Assist Trajectory Design：Theory and Practice ［D］.Dissertation，School of Aeronautics and Astronautics，Purdue University，MI，1997.

［2］JV BREAKWELL，LM PERKO.Matched asymptotic expansions，patched conics and the computation of interplanetary trajectories［C］.AIAA／AAS Astrodynamics Specialist Conference，AIAA Monterey，CA，1965：685－689.

［3］STEPHEN S BAYLISS.Precision targeting for multiple swing by interplanetary trajectories［D］.Measurement Systems Laboratory，Massachusetts Institute of Technology，Cambridge，Massachusetts，1970.

［4］KIZNER W.A method of describing miss distances for lunar and interplanetary trajectories［J］.Planetary and Space Science，1961（7）：125－131.

［5］陳楊，趙國強，寶音賀西，等.精確動力性模型下的火星探測軌道設(shè)計［J］.中國空間科學(xué)技術(shù)，2011，31（1）：8－15.CHEN YANG，ZHAO GUOQIANG，BAOYIN HEXI，et al.Orbit design for Mars exploration by the accurate dynamic model［J］.Chinese Space Science and Technology，2011，31（1）：8－15.

［6］趙國強，寶音賀西，李俊峰.基于B 平面的火星探測直接轉(zhuǎn)移軌道設(shè)計方法 ［J］.中國空間科學(xué)技術(shù)，2012，32（1）：1－7.ZHAO GUOQIANG，BAOYIN HEXI，LI JUNFENG.Direct transfer trajectory design for Mars Exploration Using B－plane［J］.Chinese Space Science and Technology，2012，32（1）：1－7.

［7］陳楊，寶音賀西，李俊峰.木星探測軌道分析與設(shè)計 ［J］.天文學(xué)報，2012，53（2）：106－118.CHEN YANG，BAOYIN HEXI，LI JUNFENG.Jupiter exploration mission analysis and trajectory ［J］.ACTA ASTRONOMICA SINICA，2012，53（2）：106－118.