饒殷睿，殷建豐，韓 潮

（北京航空航天大學(xué) 宇航學(xué)院，北京 100191）

0 引言

隨著航天技術(shù)的不斷發(fā)展，衛(wèi)星編隊飛行因其誘人的應(yīng)用前景，逐漸成為航天領(lǐng)域的研究熱點(diǎn)。編隊飛行衛(wèi)星具高可靠性和編隊可重構(gòu)的優(yōu)勢，主要用于對地觀測、精確定位、天文觀測以及電子偵察等領(lǐng)域。這些應(yīng)用對精度的要求很高，因此研究高精度的編隊飛行衛(wèi)星系統(tǒng)相對運(yùn)動模型意義重要。文獻(xiàn)［1］建立了一種系統(tǒng)、準(zhǔn)確、穩(wěn)定、可控，可通過相對測量信息高精度地自主確定的相對運(yùn)動模型，采用球面幾何方法，嚴(yán)格定義了相對軌道要素，推導(dǎo)了開普勒軌道條件下高精度橢圓軌道近距離相對運(yùn)動方程［2］。本文稱其為經(jīng)典相對軌道要素描述的相對運(yùn)動模型。但衛(wèi)星在實(shí)際運(yùn)行中會受到地球非球形引力、第三體引力、太陽光壓、大氣阻力等多種攝動力的作用［3］。在攝動力的影響下，編隊飛行衛(wèi)星的軌道根數(shù)及相對軌道會發(fā)生變化，雖然這些攝動力與地球中心引力相比非常小，但長期作用仍可使衛(wèi)星軌道偏離衛(wèi)星應(yīng)用任務(wù)的要求［4］。因此，在需要高精度軌道預(yù)報的衛(wèi)星編隊飛行研究中，這些攝動力不能被忽略。編隊飛行衛(wèi)星的攝動分析方法主要有動力學(xué)方法和運(yùn)動學(xué)方法兩種。動力學(xué)方法主要是求解Hill方程，即將攝動的影響作為力函數(shù)添加到Hill方程右端，使其適于攝動條件下衛(wèi)星相對運(yùn)動的研究［5］。但這需求解非線性微分方程，難度較大，且僅在某些特殊情況下才能獲得解析解［6］。運(yùn)動學(xué)方法是基于運(yùn)動學(xué)模型的攝動分析方法，通過研究軌道根數(shù)在各種攝動力影響下的變化分析編隊飛行衛(wèi)星的攝動問題［7］。對此，國內(nèi)外已有相關(guān)研究［8－9］。本文在經(jīng)典相對軌道要素法的基礎(chǔ)上，考慮地球非球形引力、大氣阻力及三體引力等攝動力的影響，為相對軌道要素添加相應(yīng)的攝動項，以推導(dǎo)適于橢圓和近圓攝動軌道的預(yù)報算法，拓展相對軌道要素基本理論，對攝動條件下的衛(wèi)星高精度近距離相對運(yùn)動模型建立進(jìn)行了研究，并討論了該模型的有效性和精度。

1 相對軌道要素定義

文獻(xiàn)［1］推導(dǎo)了開普勒軌道條件下相對軌道要素的相關(guān)公式，其定義為

式中：μ為地心引力常數(shù)；a為半長軸；i為軌道傾角；e為偏心率；ω為近地點(diǎn)幅角；Ω為升交點(diǎn)赤經(jīng)；M為平近點(diǎn)角；下標(biāo)a，b分別表示參考星和伴隨星。

對衛(wèi)星近距離相對運(yùn)動，在小角度相對運(yùn)動假設(shè)條件下，由一階線性展開可導(dǎo)出從絕對軌道要素｛a，e，Ω，i，ω，M｝與 相 對 軌 道 要 素 ｛D，Δex，Δey，Δix，Δiy，ΔM′｝間相互轉(zhuǎn)換的直接表達(dá)式為

式中：為平均角速度，且

2 攝動力模型

根據(jù)文獻(xiàn)［10］的攝動力模型，引入擬平均根數(shù)定義，令ξ＝ecosω，η＝－esinω，λ＝ω＋M，以a，i，Ω，ξ，η，λ作為新的絕對軌道六要素，以消除e＝0和通約奇點(diǎn)。此處所有參數(shù)采用的單位制均為標(biāo)準(zhǔn)單位制，即長度和質(zhì)量單位各為參考橢球體赤道半徑ae和地球質(zhì)量M，且使引力常數(shù)G＝1，相應(yīng)的時間單位T＝ （（ae）3／（GM））1／2＝806.811 6341s。此處：GM＝3.986 005×1014m3／s2。

2.1 非球形引力攝動

因地球的密度分布并不均勻，其形狀亦不是球形，且相當(dāng)不規(guī)則。對人造地球衛(wèi)星的運(yùn)動來說，此形狀和密度分布的非球形部分是不可忽視的一種攝動源（地球非球形引力攝動）。對近地軌道衛(wèi)星，主要帶諧項J2項攝動是影響其運(yùn)行軌道的最主要攝動源。衛(wèi)星在中心引力場和J2項攝動作用下，以σ0表示初始時刻衛(wèi)星密切軌道根數(shù)，其瞬時軌道根數(shù)的變化為

在以上三種攝動項中，對近地衛(wèi)星軌道影響最大的為一階長期項，與其相比，短周期項和一階長周期項的影響可忽略。一階長期項可表示為

2.2 大氣阻力攝動

由于人造地球衛(wèi)星是在近地空間飛行，特別是低軌道衛(wèi)星高度僅數(shù)百公里，甚至低于200km，大氣的阻尼影響較顯著，使衛(wèi)星橢圓軌道不斷變小變圓，因此大氣也是影響衛(wèi)星運(yùn)動的一種主要攝動因素，特別是對那些面質(zhì)比較大的低軌衛(wèi)星［11］。為便于研究，作如下假設(shè)：大氣層受地球扁率的影響，即考慮扁球大氣，但不考慮日下點(diǎn)大氣密度的突變；大氣層旋轉(zhuǎn)角速度與地球自轉(zhuǎn)角速度相同；衛(wèi)星橫截面積S不變［12－13］。則大氣阻力攝動的長期攝動項系數(shù)可表示為

此處：CD＝2.2；C≈0.1；S／m≈109；ne為地球自轉(zhuǎn)角速度；μ′≈0.1；z0＝a0e0／Hp0；ρp0，Hp0分別為衛(wèi)星軌道初始近地點(diǎn)處的密度和標(biāo)高；σ0為初始近地點(diǎn)處的軌道根數(shù)；In（z）（n＝1，…，6）為第一類虛變量的貝塞爾函數(shù)［14－15］。為表達(dá)簡潔，式（21）～（25）中In（z）均簡寫為In。

2.3 三體引力攝動

日、月和大行星的引力作用對人造地球衛(wèi)星尤其是遠(yuǎn)地衛(wèi)星運(yùn)動的影響，是天體力學(xué)中的一種典型的第三體攝動問題。雖然這類攝動力也是一種保守力，但由于日、月和地球、衛(wèi)星相距都不太遠(yuǎn)，問題顯得較復(fù)雜。對高精度要求，日、月和地球不能簡單視作質(zhì)點(diǎn)，而應(yīng)精確計算日、月的位置。本文主要考慮日月攝動引起的衛(wèi)星軌道變化的三個攝動項中的二階長期項σc（t－t0），有

對月球，有

式中：

3 攝動條件下基于相對軌道要素的橢圓基準(zhǔn)軌道的相對運(yùn)動描述

3.1 相對軌道要素變換

根據(jù)前文提到的擬平均根數(shù)法，令ξi＝eicosωi，ηi＝－eisinωi，λi＝ωi＋Mi（i＝a，b），并以a，i，Ω，ξ，η，λ作為新的絕對軌道六要素，則相對軌道要素可表示為

3.2 考慮攝動影響的相對軌道要素計算

將上述各攝動力對衛(wèi)星軌道根數(shù)的影響表達(dá)式代入相對軌道要素的定義式中，并將其化為標(biāo)準(zhǔn)單位，可得

式中：

3.3 橢圓基準(zhǔn)軌道相對運(yùn)動描述

以上述攝動條件下的相對軌道要素為基本參數(shù)，建立衛(wèi)星近距離且考慮攝動力作用時相對距離和相對速度的代數(shù)關(guān)系，結(jié)果適用于攝動條件下的橢圓基準(zhǔn)軌道和近圓基準(zhǔn)軌道。

根據(jù)文獻(xiàn)［1］中的結(jié)論，可得伴隨衛(wèi)星在基準(zhǔn)星質(zhì)心軌道坐標(biāo)系中的近距離相對運(yùn)動方程和速度方程為

式中：

此處：θ為真近點(diǎn)角；u為緯度幅角；r為軌道矢徑的模長；

注：上述各式中參數(shù)，a，e，ω均為參考星的瞬時絕對軌道要素；D，Δex，Δey，Δix，Δiy，ΔM′均為瞬時相對軌道要素；r，u均為與時間相關(guān)的變量。

將攝動模型中各攝動力影響下的絕對軌道要素（式（17）～（23）、（29）、（30））及相對軌道要素變換中給出的帶攝動項的相對軌道要素（式（40）～（45））的表達(dá)式代入式（52）～（58），可得攝動條件下伴隨星相對參考星的位置速度方程。

對速度方程，需說明的是由于帶攝動項的相對軌道要素為時變的，在對位置方程進(jìn)行求導(dǎo)時也應(yīng)對其進(jìn)行求導(dǎo)。但考慮相對軌道要素對時間的導(dǎo)數(shù)項對速度的影響較小，且為使模型盡可能簡單，故將其略去。

帶攝動項的相對軌道要素的橢圓基準(zhǔn)軌道相對運(yùn)動方程的形式與經(jīng)典相對軌道要素并未復(fù)雜太多，但適用范圍卻有很大拓展。本文給出的帶攝動項的相對運(yùn)動模型，可用于非球形引力、大氣阻力、三體引力等攝動力影響下橢圓軌道衛(wèi)星的相對軌道預(yù)報。

4 數(shù)值算例

本章算例中的精確解通過高精度軌道預(yù)報器精確計算出兩顆衛(wèi)星在地心赤道坐標(biāo)系中的位置和速度，進(jìn)而計算兩星的相對位置和相對速度，最后投影到軌道坐標(biāo)系得到［16］。

算例1：非球形引力攝動在近地橢圓基準(zhǔn)軌道衛(wèi)星相對運(yùn)動中的影響分析

選擇高精度軌道預(yù)報器的攝動類型為非球形引力攝動，編隊飛行特征尺度1km，基準(zhǔn)軌道偏心率0.01，典型的初始軌道參數(shù)見表1。

考慮與不考慮非球形引力攝動相對運(yùn)動分析結(jié)果如圖1所示。衛(wèi)星運(yùn)行10圈的誤差值統(tǒng)計結(jié)果見表2。

由圖1、表2可知：在考慮非球形引力攝動影響時，對近地橢圓基準(zhǔn)軌道而言，無攝動項的橢圓相對運(yùn)動模型誤差非常大，本文給出的帶非球形引力一階長期攝動項的橢圓相對運(yùn)動模型較好地適應(yīng)了非球形引力攝動對橢圓基準(zhǔn)軌道的影響，相對軌道預(yù)報精度誤差保持在十米級的標(biāo)準(zhǔn)，相對誤差具有百分之幾的量級精度。

說明：在考慮非球形引力攝動的影響時，在相對軌道要素中并未加入短周期項，主要是因為短周期項對衛(wèi)星相對軌道的影響量級僅約1%，故可忽略。短周期項影響下衛(wèi)星相對位置和相對速度誤差如圖2所示。由圖可知：與一階長期項相比，短周期項的影響很小，故在為經(jīng)典相對軌道要素添加非球形引力攝動項時，只添加一階長期項部分既可簡化瞬時相對軌道要素計算的復(fù)雜性，又不降低太多精度，是合理的。

表1 算例1中衛(wèi)星A、B的初始軌道參數(shù)Tab.1 Orbit parameters of satellite A and B in case 1

圖1 非球形引力影響下衛(wèi)星相對位置和相對速度誤差Fig.1 Errors of relative position and relative velocity under effect of nonspherical gravitational perturbation

表2 算例1中衛(wèi)星運(yùn)行10圈后相對位置和相對速度誤差Tab.2 Errors of relative position and relative velocity after 10periods in case 1

算例2：大氣阻力攝動在近地軌道衛(wèi)星相對運(yùn)動中的影響分析

圖2 短周期項影響下衛(wèi)星相對位置和相對速度誤差Fig.2 Errors of relative position and relative velocity under effect of short periodic term

選擇高精度軌道預(yù)報器的攝動類型為非球形引力攝動＋大氣阻力攝動，編隊飛行特征尺度為1km，基準(zhǔn)軌道偏心率為0.001，典型的初始軌道參數(shù)見表3。說明：由于大氣阻力對近地軌道衛(wèi)星影響較大，且在圓軌道的情況下會更明顯，故算例中的基準(zhǔn)軌道取為近圓軌道。

考慮和未考慮大氣阻力攝動時相對運(yùn)動分析的結(jié)果如圖3所示。衛(wèi)星運(yùn)行10圈的誤差統(tǒng)計結(jié)果見表4。

由圖3和表4可知：考慮大氣阻力攝動影響時，對近地軌道來說，無大氣阻力攝動項的相對運(yùn)動模型誤差仍很大，本文給出的帶大氣阻力長期攝動項的相對運(yùn)動模型可較好地反映大氣阻力對近地軌道衛(wèi)星近距離相對運(yùn)動的影響，相對軌道預(yù)報精度誤差保持在1m級的標(biāo)準(zhǔn)，相對誤差有百分之幾的量級精度。

表3 算例2中A、B衛(wèi)星的初始軌道參數(shù)Tab.3 Orbit parameters of satellite A and B in case 2

圖3 大氣阻力影響下衛(wèi)星相對位置和相對速度誤差Fig.3 Errors of relative position and relative velocity under effect of atmospheric drag perturbation

表4 算例2中衛(wèi)星運(yùn)行10圈后相對位置和相對速度誤差Tab.4 Errors of relative position and relative velocity after 10periods in case 2

算例3：三體引力攝動在遠(yuǎn)地橢圓基準(zhǔn)軌道衛(wèi)星相對運(yùn)動中的影響分析

選擇高精度軌道預(yù)報器的攝動類型為非球形引力攝動＋日月引力攝動，編隊飛行特征尺度為1km，基準(zhǔn)軌道偏心率為0.01，典型的初始軌道參數(shù)見表5。

圖4 三體引力影響下衛(wèi)星相對位置和相對速度誤差Fig.4 Errors of relative position and relative velocity under effect of third body gravitational perturbation

考慮和未考慮三體引力攝動時相對運(yùn)動分析結(jié)果如圖4所示。衛(wèi)星運(yùn)行10圈的誤差統(tǒng)計結(jié)果見表6。由圖4和表6可知：考慮三體引力攝動影響時，對遠(yuǎn)地橢圓基準(zhǔn)軌道來說，無三體引力攝動項的橢圓相對運(yùn)動模型誤差較大，本文給出的帶三體引力長期攝動項的橢圓相對運(yùn)動模型在一定程度上減小了因日月引力攝動造成的軌道預(yù)報誤差，相對軌道預(yù)報精度誤差保持在米級的標(biāo)準(zhǔn)，相對誤差具有千分之幾的量級精度。

算例4：相對軌道要素在攝動條件下的特性分析

由前三個算例已發(fā)現(xiàn)，非球形引力攝動、大氣阻力攝動以及三體引力攝動對衛(wèi)星近距離相對運(yùn)動軌道的影響很大，不可忽略。經(jīng)典相對軌道要素與帶攝動項的相對軌道要素分別描述的相對運(yùn)動模型的誤差較大的原因是：在攝動條件下，相對軌道要素隨時間不斷變化，而經(jīng)典相對軌道要素是常量，其預(yù)報的相對運(yùn)動軌道與精確解間的誤差較大。取與算例1相同的基準(zhǔn)軌道參數(shù)，考慮非球形引力、大氣阻力、三體引力等三種攝動力的影響，計算各時刻相對軌道要素的值，并求得其與高精度軌道預(yù)報器計算得到的瞬時相對軌道要素的精確值的相對誤差，再同經(jīng)典相對軌道要素與精確值的相對誤差進(jìn)行對比，結(jié)果如圖5所示。

表5 算例3中衛(wèi)星A、B的初始軌道參數(shù)Tabl.5 Orbit parameters of satellite A and B in case 3

表6 算例3中衛(wèi)星運(yùn)行10圈后相對位置和相對速度誤差Tab.6 Errors of relative position and relative velocity after 10periods in case 3

由圖5可知：加入攝動項后，瞬時相對軌道要素與精確值的相對誤差明顯減小至約0.5%，而相對軌道要素的值的精度直接關(guān)系最后求得的衛(wèi)星相對運(yùn)動的位置和速度的精度。因此，在攝動條件下，為經(jīng)典相對軌道要素添加攝動項非常有必要且有效。

說明：由于初始給定的相對軌道要素D，ΔM′均為0，故只計算了Δex，Δey，Δix，Δiy的相對誤差。由式（40）～（51）可知，攝動項中只有大氣阻力的長期攝動項對Δiy（t）有影響，故與其他幾項相對軌道要素攝動項相比，添加與否的影響較小，圖5中經(jīng)典相對軌道要素中Δiy（t）與精確值的相對誤差和帶攝動項的相對軌道要素中Δiy（t）與精確值的相對誤差接近（相對誤差曲線基本重合）。

算例5：相對運(yùn)動模型精度分析

設(shè)典型初始軌道參數(shù)見表7，初始相對軌道要素按比例增加，并且衛(wèi)星運(yùn)行的圈數(shù)從10圈到100圈依次增加。表中：S為相對運(yùn)動尺度，從1km到10km以1km為步長遞增。仿真結(jié)果如圖6～8所示。

圖5 相對軌道要素誤差分析Fig.5 Errors of relative orbit elements

表7 算例5中A、B衛(wèi)星的初始軌道參數(shù)Tab.7 Orbit parameters of satellite A and B in case 5

圖6 兩種相對運(yùn)動模型的相對軌道要素誤差Fig.6 Relative elements errors of two relative models

圖7 經(jīng)典相對軌道要素運(yùn)動模型位置速度誤差Fig.7 Relative position and velocity errors of classical relative orbit elements model

比較圖6～8可知：典型相對軌道要素描述的衛(wèi)星近距離相對運(yùn)動模型的相對軌道要素和相對位置速度的誤差，隨圈數(shù)和相對運(yùn)動的尺度增加而急劇增大，而帶攝動項的相對軌道要素描述的運(yùn)動模型的誤差，隨圈數(shù)以及相對運(yùn)動的尺度增加的趨勢不明顯，出現(xiàn)誤差主要是由模型線性化導(dǎo)致。在低軌道上，帶攝動項的相對軌道要素運(yùn)動模型的精度明顯高于經(jīng)典相對軌道要素運(yùn)動模型。這進(jìn)一步表明本文給出的帶攝動項的相對軌道要素運(yùn)動模型有更強(qiáng)的軌道類型的適應(yīng)性和預(yù)報精度。

圖8 帶攝動的相對軌道要素運(yùn)動模型位置速度誤差Fig.8 Relative position and velocity errors of relative orbit elements with perturbation model

5 結(jié)束語

本文對攝動條件下橢圓軌道衛(wèi)星相對運(yùn)動進(jìn)行了研究。研究發(fā)現(xiàn)：本文基于文獻(xiàn)［1］給出的開普勒軌道的相對軌道要素，采用擬平均根數(shù)法，消除奇點(diǎn)情況，為其添加非球形引力、大氣阻力、三體引力部分?jǐn)z動項，由此得到帶攝動項的相對軌道要素，其形式簡單，便于拓展，適用范圍廣，且易于與精確解進(jìn)行誤差分析。算例表明，其與精確解的相對誤差保持在約0.5%，因此是正確而有效的。各種攝動力對衛(wèi)星相對運(yùn)動中的相對位置和相對速度均有很大的影響，不能忽略。文中給出的帶攝動項的相對軌道要素法提高了預(yù)報精度，其預(yù)報衛(wèi)星相對運(yùn)動軌道的誤差較采用經(jīng)典的相對軌道要素法最多可減少近90%。本文推導(dǎo)的帶攝動項相對軌道要素描述的衛(wèi)星近距離相對運(yùn)動模型，形式簡單，計算速度快，與通過積分方法得到的精確解相比，誤差保持在十米量級，提供的模型精確有效，為衛(wèi)星編隊飛行的軌道預(yù)報提供了一種新的參考模型。研究結(jié)果統(tǒng)一了開普勒軌道和攝動軌道的近圓和橢圓基準(zhǔn)軌道衛(wèi)星相對運(yùn)動描述形式，這為軌道攝動條件下橢圓軌道衛(wèi)星編隊飛行的研究提供了便利和精度保證。

［1］ 韓 潮，殷建豐.基于相對軌道要素的橢圓軌道衛(wèi)星相對運(yùn)動研究［J］.航空學(xué)報，2011，32（12）：2244－2258.

［2］ TAFF L G.Computational spherical astronomy［M］.New York：Wiley，1981.

［3］ XU Guo－chang.Orbits［M］.Berlin Heidelberg：Spinger－Verlag，2008.

［4］ 孟 鑫，李俊峰，高云峰.衛(wèi)星編隊飛行中相對軌道的J2攝動分析［J］.力學(xué)學(xué)報，2006，38（1）：89－96.

［5］ HILL G W.Researches on the lunar theory［J］.A－merican Journal of Mathematics，1978，1（3）：245－260.

［6］ 張玉錕，戴金海.衛(wèi)星編隊飛行的地球扁率攝動分析［J］.宇航學(xué)報，2002，23（3）：72－76.

［7］ 張 娟，和興鎖，鄧峰巖，等.編隊飛行星座的地球扁率攝動和大氣阻力攝動分析［J］.應(yīng)用力學(xué)學(xué)報，2006，23（3）：496－450.

［8］ 楊 宇，韓 潮.編隊飛行衛(wèi)星群描述及攝動分析［J］.中國空間科學(xué)技術(shù)，2002，22（2）：15－23.

［9］ ALFRIRND K T，SCHAUB H，GIM D W.Gravitational perturbations nonlinearity and circular orbit assumption effects on formation flying control strategies：Annual AAS Rocky Mountain Guidance and Control Conference［C］.2000.

［10］ 劉 林.航天器軌道理論［M］.北京：國防工業(yè)出版社，2000.

［11］ 劉 林.人造地球衛(wèi)星軌道力學(xué)［M］.北京：高等教育出版社，1992.

［12］ KING－HELE D G.Analysis of the orbit of 1970－114Fin its last 20days［J］.Planet，Space Sci，1976，24（1）：1－16.

［13］ KING－HELE D G.Theory of satellite orbits in an atmosphere［M］.London：Butterworths，1964.

［14］ US Government.COESA，US standard atmosphere，1976［M］.Washington DC：US Government Printing Office，1976.

［15］ WATSON G N.Theory of Bessel functions［M］.London and New York：Combridge Univ Press，1994.

［16］ 肖業(yè)倫.航天器飛行動力學(xué)原理［M］.北京：宇航出版社，1995.