張 磊于登云張 熇
(1 北京空間飛行器總體設(shè)計(jì)部,北京 100094) (2 中國(guó)航天科技集團(tuán)公司,北京 100048)
月球采樣返回是月球探測(cè)任務(wù)中的重要內(nèi)容之一[1],俄羅斯和美國(guó)都已完成了若干次采樣返回任務(wù)。在這些任務(wù)中,探測(cè)器都是由月-地轉(zhuǎn)移軌道直接再入大氣。作為一種已被任務(wù)驗(yàn)證且被普遍采用的飛行方案,直接再入大氣的月球返回軌道對(duì)我國(guó)即將開展的月球采樣返回任務(wù)具有一定的參考價(jià)值。
直接再入大氣的月球返回軌道,以再入大氣點(diǎn)為界面分為月-地轉(zhuǎn)移軌道和再入大氣軌跡兩部分。本文采用兩級(jí)修正策略使月-地轉(zhuǎn)移軌道滿足再入大氣界面約束(再入點(diǎn)地心距、慣性再入角、傾角),進(jìn)而通過(guò)再入大氣仿真及落點(diǎn)匹配,使落點(diǎn)地理位置滿足約束要求。此外,本文還對(duì)用于月球返回軌道設(shè)計(jì)的軌道動(dòng)力學(xué)模型進(jìn)行了攝動(dòng)項(xiàng)影響分析。
在地心赤道慣性坐標(biāo)下,月-地轉(zhuǎn)移軌道設(shè)計(jì)中的軌道動(dòng)力學(xué)模型可由下式表示
式(1)中, μe為地球引力常數(shù), rep為探測(cè)器在地心赤道慣性系下的位置矢量, aN為N 體引力攝動(dòng)加速度, aNSE為地球非球形攝動(dòng)加速度, aNSM為月球非球形攝動(dòng)加速度, aD為大氣阻力攝動(dòng)加速度, aR為太陽(yáng)光壓攝動(dòng)加速度。若完全考慮所有攝動(dòng)項(xiàng),則會(huì)造成軌道計(jì)算效率降低,因此,需要通過(guò)攝動(dòng)項(xiàng)影響分析,得到兼顧計(jì)算精度和計(jì)算效率的軌道動(dòng)力學(xué)模型。表1 所示為各攝動(dòng)項(xiàng)被忽略時(shí)造成的再入點(diǎn)偏差。
表1 忽略某一攝動(dòng)項(xiàng)時(shí)對(duì)月-地轉(zhuǎn)移軌道影響的再入點(diǎn)偏差Table 1 Reentry point deviations resulting from neglecting a certain perturbation terms
可以看出:太陽(yáng)引力攝動(dòng)對(duì)月-地轉(zhuǎn)移軌道的影響最大,也最為嚴(yán)重,以至于探測(cè)器不能進(jìn)入大氣;其次是月球非球形攝動(dòng)的影響,其對(duì)到達(dá)再入點(diǎn)時(shí)刻、再入點(diǎn)位置、慣性再入角造成一定的偏差;地球非球形攝動(dòng)的影響、太陽(yáng)光壓攝動(dòng)、大氣阻力攝動(dòng)對(duì)月-地轉(zhuǎn)移軌道的影響非常微小。本文采用如下軌道動(dòng)力學(xué)模型,該模型考慮日、月引力攝動(dòng),忽略其它攝動(dòng)項(xiàng)。式(2)中, rmp表示探測(cè)器相對(duì)月球的位置矢量,rsp表示探測(cè)器相對(duì)太陽(yáng)的位置矢量,月球、太陽(yáng)相對(duì)地球的位置矢量rem 、res 由NASA 噴氣推進(jìn)實(shí)驗(yàn)室(JPL)DE405 星歷得到, rmp=rep-rem,rsp=rep-res。盡管采用該模型計(jì)算的月地轉(zhuǎn)移軌道存在偏差(主要由月球非球形攝動(dòng)引起),但還應(yīng)考慮到采用該模型計(jì)算的月-地轉(zhuǎn)移軌道初始狀態(tài)與采用式(1)模型計(jì)算所得的差異非常小,如表2 所示。
表2 不同模型置入?yún)?shù)的偏差Table 2 Insertion point parameters deviations of different models
可見,升交點(diǎn)赤經(jīng)Ωm、近月點(diǎn)幅角ωm的偏差均遠(yuǎn)小于0.1°,置入速度vpm偏差遠(yuǎn)小于0.1m/s,這樣的偏差小于發(fā)動(dòng)機(jī)執(zhí)行誤差甚至導(dǎo)航誤差,從工程設(shè)計(jì)的角度看,采用式(2)模型進(jìn)行標(biāo)稱月-地軌道設(shè)計(jì)即可,而不必采用更為復(fù)雜的軌道動(dòng)力學(xué)模型。
在球形旋轉(zhuǎn)地球、無(wú)風(fēng)條件下,建立無(wú)量綱化的再入飛行動(dòng)力學(xué)方程[2-3]??紤]彈道式和彈道-升力式兩種再入方式。彈道式再入,動(dòng)力學(xué)方程中不包含氣動(dòng)升力項(xiàng)。彈道-升力式再入,僅考慮其縱向運(yùn)動(dòng)過(guò)程,動(dòng)力學(xué)方程中不包含sinσ項(xiàng), σ為速度傾側(cè)角。
彈道-升力式再入通過(guò)σ調(diào)整其縱向航程, σ由再入制導(dǎo)律給出[4]。這里對(duì)問題進(jìn)行簡(jiǎn)化,考慮再入地心航程角為一給定值,可由此確定一個(gè)常值σ。為此,引入下面的方程
式中, sgo為需要完成的地心航程角, τ為時(shí)間變量,Vks為航跡速度, γ為航跡傾斜角,是航跡速度矢量與當(dāng)?shù)厮矫骈g的角度,偏向上方時(shí)為正, rs為地心距。由于sgo(0)為給定值,因此, sgo(τf)是σ的函數(shù),記作
采用割線法求解滿足sgof(σ)=0 的速度傾側(cè)角,如式(5)所示
月-地轉(zhuǎn)移軌道初步設(shè)計(jì)是在雙二體模型下進(jìn)行的,由初步設(shè)計(jì)得到的設(shè)計(jì)變量初值,直接用于上述軌道力學(xué)模型進(jìn)行軌道傳播會(huì)產(chǎn)生較大的偏差,因此,必須對(duì)設(shè)計(jì)初值進(jìn)行修正,這也是標(biāo)稱月-地轉(zhuǎn)移軌道設(shè)計(jì)的基礎(chǔ)。初值修正是一個(gè)正向參數(shù)搜索的過(guò)程,對(duì)初值進(jìn)行傳播得到終端狀態(tài),計(jì)算終端偏差并根據(jù)該偏差改正初值,重復(fù)該過(guò)程使終端偏差最小。可用于初值修正的數(shù)值方法和非數(shù)值方法很多,本文采用在轉(zhuǎn)移軌道設(shè)計(jì)常用的微分校正法[5-8]。
初值修正的設(shè)計(jì)變量p 選為:轉(zhuǎn)移軌道置入速度vpm、升交點(diǎn)赤經(jīng)Ωm、近月點(diǎn)幅角ωm;目標(biāo)變量q選為再入界面終端參數(shù):再入點(diǎn)地心距ren、慣性再入角γen、傾角ie。在初值修正中,若直接采用終端參數(shù)作為目標(biāo)變量,迭代過(guò)程的收斂性不易保證。為此,本文引入橢圓軌道B 平面參數(shù),并提出初值的兩級(jí)修正策略。
對(duì)月-地轉(zhuǎn)移軌道,考慮節(jié)省燃耗,地心段軌道通常是橢圓軌道,因此,雙曲線軌道B 平面參數(shù)[9-10]在這里就不再適用,但可以參考其定義引入橢圓軌道B 平面參數(shù)[8]。橢圓軌道B 平面參數(shù)Bp記為
式(6)中, BT、BR為B 矢量的兩個(gè)標(biāo)量, b =TF P 為探測(cè)器從月-地轉(zhuǎn)移軌道置入點(diǎn)至近地點(diǎn)的飛行時(shí)間。B 平面參數(shù)與軌道根數(shù)間的關(guān)系由下式表示
式(7)、(8)中, ie為橢圓軌道傾角, ωe為橢圓軌道近地點(diǎn)幅角。
初值的兩級(jí)修正策略為:
1)第一級(jí)采用B 平面參數(shù)作為目標(biāo)變量,積分終止條件為近地點(diǎn),可適當(dāng)放寬迭代終止條件;
2)第二級(jí)使用第一級(jí)修正所得的設(shè)計(jì)參數(shù)作為迭代初值,并采用終端參數(shù)作為目標(biāo)變量,積分終止條件為飛行時(shí)間。
前文敘述了通過(guò)兩級(jí)修正獲得標(biāo)稱軌道初始狀態(tài),該初始狀態(tài)可以保證月-地轉(zhuǎn)移軌道滿足轉(zhuǎn)移時(shí)間、再入高度、再入角、傾角這些約束條件,但不能保證落點(diǎn)地理位置滿足要求,因此,還需要進(jìn)行落點(diǎn)匹配。大致過(guò)程如下:
1)對(duì)設(shè)計(jì)初值p0進(jìn)行兩級(jí)修正,第一級(jí)修正的結(jié)果為, 第二級(jí)修正的結(jié)果為;
3)由ren、ven生成再入動(dòng)力學(xué)方程狀態(tài)變量初值x0, 積分再入動(dòng)力學(xué)方程,得到落點(diǎn)地理位置λu、φu和地心航程角θnu;
4)比較φu與要求的落點(diǎn)緯度, 判斷是否滿足要求,若不滿足,則調(diào)整傾角ie , 更新橢圓軌道B平面參數(shù)qB和終端參數(shù)qE, 將第二級(jí)修正的結(jié)果作為設(shè)計(jì)初值p0, 返回步驟1),若滿足,則進(jìn)行下一步;
5)比較λu與要求的落點(diǎn)緯度, 判斷是否滿足要求,若不滿足,則調(diào)整置入時(shí)刻ti, 將第二級(jí)修正的結(jié)果作為設(shè)計(jì)初值p0, 返回步驟1),若滿足,則輸出設(shè)計(jì)結(jié)果。
圖1 所示為結(jié)合落點(diǎn)仿真的標(biāo)稱月地轉(zhuǎn)移軌道設(shè)計(jì)流程。
調(diào)整軌道傾角ie的目標(biāo)是為了使落點(diǎn)緯度φu滿足約束要求, 由再入點(diǎn)狀態(tài)ren、ven可以生成一組瞬時(shí)軌道根數(shù)ae、ee、ie、Ωe、ωe, 結(jié) 合由再入仿真得到的地心航程角θ, 可方便的計(jì)算落點(diǎn)位置單位矢量
式 中, 坐 標(biāo) 轉(zhuǎn) 換 矩 陣 Toei= RZ(ωe+ fen+θnu)RX(ie)RZ(Ωe), fen為再入點(diǎn)真近點(diǎn)角。則落點(diǎn)緯度由下式計(jì)算
因此,可以將落點(diǎn)緯度φu看作傾角ie的函數(shù)
保持ae、ee、Ωe、ωe不變,通過(guò)調(diào)整傾角ie, 可使落點(diǎn)緯度滿足要求,本文采用New ton-Raphson 迭代法
圖1 標(biāo)稱月地轉(zhuǎn)移軌道設(shè)計(jì)流程圖Fig.1 Flow chart of nominal moon-to-earth transfer orbit design
流程圖的外層是通過(guò)調(diào)整月-地轉(zhuǎn)移軌道置入時(shí)刻進(jìn)行落點(diǎn)經(jīng)度匹配,如下式所示
經(jīng)度匹配實(shí)質(zhì)上是對(duì)升級(jí)點(diǎn)赤經(jīng)的調(diào)整,從式(7)、(8)可以看出,B 平面參數(shù)與升交點(diǎn)赤經(jīng)無(wú)關(guān),因此,經(jīng)度匹配對(duì)置入時(shí)刻的調(diào)整,并不會(huì)影響用于緯度匹配的B 平面參數(shù)。
相應(yīng)的橢圓軌道B 平面參數(shù)按以下步驟更新,由再入點(diǎn)狀態(tài)ren、ven計(jì)算B 平面參數(shù)qB
式(17)、(18)中, ωe是由ren、ven計(jì)算的瞬時(shí)近地點(diǎn)幅角。
B 矢量的模b 可由下式計(jì)算
由月-地轉(zhuǎn)移軌道特性可知,月-地轉(zhuǎn)移軌道近地點(diǎn)地心距rp與慣性再入角余弦在一定范圍呈線性關(guān)系,如式(20)所示,且與其它因素幾乎無(wú)關(guān)。
因此,在慣性再入角確定的條件下,近地點(diǎn)地心距rp可認(rèn)為是確定的。由于對(duì)軌道傾角ie的調(diào)整非常微小,所以由ren、ven確定的一組瞬時(shí)軌道根數(shù),除ie外,其它都認(rèn)為是不變的,這樣在計(jì)算新B 平面參數(shù)時(shí), b、ωe保持不變,只需將新計(jì)算的軌道傾角i*e帶入即可,而過(guò)近地點(diǎn)的飛行時(shí)間TP F 也保持不變。這樣調(diào)整后的B 平面參數(shù)與新的終端參數(shù)具有很好的對(duì)應(yīng)性,表現(xiàn)為第一級(jí)修正后的初值用于第二級(jí)修正時(shí),迭代2、3 次即可以預(yù)定精度收斂,而不調(diào)整B 平面參數(shù)或直接以終端參數(shù)作為目標(biāo)變量時(shí),往往收斂速度慢,甚至不能以預(yù)定精度收斂。
設(shè)置約束條件如下
1)月-地轉(zhuǎn)移軌道置入日期:2017年10月;
2)月球端約束:環(huán)月停泊軌道高度100km,環(huán)月停泊軌道傾角40°;
3)地球端約束:再入點(diǎn)地心距6 500km;
4)地-月轉(zhuǎn)移軌道飛行時(shí)間:72h;
5)落點(diǎn)地理位置:東經(jīng)112°、北緯42°。
考慮彈道式再入和彈道升力再入兩種地球大氣再入方式,參數(shù)如下
1)彈道式再入:降段再入,軌道傾角60°左右,再入角-12°;
2)彈道-升力式再入:升段再入,軌道傾角43°左右,再入角-6°,按地心航程角70°計(jì)算常值σ。
探測(cè)器參數(shù)
1)彈道式再入:彈道系數(shù)BC=60;
2)彈道-升力式再入:質(zhì)量m=2 800kg,迎風(fēng)面面積S =5m2,升力系數(shù)CL=0.443,阻力系數(shù)CD=1.11。
采用月-地轉(zhuǎn)移軌道初步設(shè)計(jì)結(jié)果作為設(shè)計(jì)初值。如表3 所示。軌道動(dòng)力學(xué)模型采用地-月-日-器限制性地體模型,不考慮其它攝動(dòng)項(xiàng)。月球、太陽(yáng)位置采用JPL 的DE405 星歷計(jì)算得到,軌道積分及再入仿真均采用RKF7(8)積分器,計(jì)算結(jié)果如表4、表5 所示。月地返回軌道的星下點(diǎn)軌跡如圖2、圖3 所示,再入軌跡如圖4、圖5 所示。
表3月-地轉(zhuǎn)移軌道設(shè)計(jì)初值Table 3 Preliminary design results of moon-to-earth transfer orbits
表4 標(biāo)稱月-地轉(zhuǎn)移軌道軌道設(shè)計(jì)結(jié)果Table 4 Design results of nominal moon-to-earth transfer orbits
表5 再入初始狀態(tài)Table 5 Initial state parameters of atmospheric reentry
圖2 彈道式再入的月球返回軌道星下點(diǎn)軌跡Fig.2 Subsatellite point track of moon-return orbit with ballistic reentry
圖3 彈道-升力式再入的月球返回軌道星下點(diǎn)軌跡Fig.3 Subsatellite point track of moon-return orbit with ballistic-lift reentry
圖4 彈道式再入彈道Fig.4 Trajectory of ballistic reentry
圖5 彈道-升力式再入彈道Fig.5 Trajectory of ballistic-lift reentry
直接大氣再入的月球返回軌道,要求月-地轉(zhuǎn)移軌道滿足大氣再入界面約束,并使探測(cè)器到達(dá)預(yù)定落點(diǎn),本文將月-地轉(zhuǎn)移軌道設(shè)計(jì)同大氣再入仿真相結(jié)合,根據(jù)落點(diǎn)偏差調(diào)整月-地轉(zhuǎn)移軌道設(shè)計(jì)參數(shù),使得月球返回軌道滿足要求。兩組分別采用彈道式再入和彈道-升力式再入方式的月球返回軌道算例表明,該方法是有效的。
)
[1]韓鴻碩, 陳杰.21 世紀(jì)國(guó)外深空探測(cè)發(fā)展計(jì)劃及進(jìn)展[J].航天器工程, 2008, 17(3):1-22
[2]李惠峰, 凌源.基于預(yù)測(cè)校正法的低升阻比飛船再入軌跡設(shè)計(jì)[J].北京航空航天大學(xué)學(xué)報(bào),2009, 35(6):705-708
[3]雍恩米, 唐國(guó)金, 陳磊.基于Gauss 偽譜方法的高超音速飛行器再入軌跡快速優(yōu)化[J].宇航學(xué)報(bào), 2008, 29(6):1766-1772
[4]王希季 主編.航天器進(jìn)入與返回技術(shù)[M].北京:中國(guó)宇航出版社, 1991
[5]楊維廉.發(fā)射極月衛(wèi)星的轉(zhuǎn)移軌道研究[J].航天器工程,1997, 6(3):19-33
[6]楊維廉.擊中月球的轉(zhuǎn)移軌道研究[J].飛行力學(xué),1998, 16(4):20-25
[7]周文艷, 楊維廉.月球探測(cè)器轉(zhuǎn)移軌道的特性分析[J].空間科學(xué)學(xué)報(bào),2004, 24(5):354-359
[8]張磊, 于登云, 張熇.繞月自由返回軌道的設(shè)計(jì)與分析[J].航天器工程, 2010, 19(2):128-135
[9]Kizner W.A method of describing miss distances for lunar and interplanetary trajectories[J].Ballistic Missile and Space Technology,1961(3)
[10]李立濤, 楊滌, 崔祜濤.奔月軌道的一種求解方法[J].宇航學(xué)報(bào),2003, 24(2):150-155