雷博持 李明濤
1.中國空間技術研究院西安分院,西安710100
2.中國科學院空間科學與應用研究中心,北京100190
立方星這一概念最早由美國加州理工大學和斯坦福大學在1999年提出,立方星是標準化的微小衛星,尺寸10 cm×10 cm×10 cm,質量約1 kg,為一個基本單元1 U。以此為基礎可以有2 U、3 U等立方星。立方星具有很多特點,不僅成本低,開發時間也較短,易于實現標準化、模塊化,易于技術更新[1-2]。
近年來,利用多衛星系統的空間任務越來越多,借助標準化的微小衛星來實現分布式空間系統極具發展前景。分布式空間任務的出現,拓展了立方星在科學領域的應用。目前已有越來越多的立方星被送上太空。國外的立方星應用主要分為3個方面:以空間科學和教學為核心目標,如荷蘭代爾夫特理工大學的 Delfi-c3;以載荷任務為目標,如 QuakeSat搭載了地震電磁探測載荷;也有以新技術驗證為目標的立方星任務,如加拿大的先進納星太空實驗室CANX衛星系列[3]。國內很多高校及科研機構也在積極開展立方星的相關研究。
體積小、重量輕的自身特點決定了立方星大多沒有在軌機動能力,考慮無主動控制時立方星軌道的自然演化就成為一個重要研究方向。太陽同步軌道(SSO)具有星上光照條件好、過同一緯度的當地時間相同等特點,已被廣泛應用于各種空間任務。因此,本文選取了SSO立方星任務進行軌道演化研究,考慮地球非球形引力項J2和大氣阻力攝動的影響,以平均半長軸隨時間的變化作為攝動效應的指標。以軌道分析軟件為仿真平臺,先進行了攝動力簡化以及立方星面質比對軌道演化的影響分析,以確立數值仿真的前提條件;接著研究了軌道高度、太陽活動指數和大氣模型這3個方面對SSO立方星軌道演化的影響。最終得出結論。
研究地球非球形攝動及大氣阻力攝動對軌道演化的影響時,涉及到地心慣性坐標系和地心地固坐標系,如圖1所示。Oxiyizi(標記為 Si)為地心慣性坐標系,Oxeyeze(標記為Se)為地心地固坐標系,α為Greenwich赤經。
地心慣性坐標系Si:坐標系原點在地心,Z軸垂直于地球赤道平面;XY軸在地球赤道平面內,X軸指向春分點位置;Y軸滿足右手螺旋定則。動力學方程一般建立在慣性坐標系中。

圖1 地心慣性坐標系和地心地固坐標系
地心地固坐標系Se:坐標系原點在地心,Z軸垂直于地球赤道平面;XY軸在地球赤道平面內,X軸指向 Greenwich子午圈與赤道交點;Y軸滿足右手螺旋定則。地球非球形攝動一般在地心地固坐標系中表示。
在慣性坐標系下,衛星動力學模型可用式(1)表示。

式中:r為地心距矢量,μ是地心引力常數,r為地心距大小,f為地球中心引力之外的所有攝動力的合力和包括地球非球形引力、大氣阻力、太陽光壓、第三體引力攝動等。
將式(1)在慣性坐標系下展開并進行積分,可以得到任意時刻衛星在慣性坐標系下的運動狀態。
對于低軌道來說,衛星所受到的攝動力主要是地球非球形引力攝動、大氣阻力攝動、太陽光壓攝動以及第三體引力攝動。設f0表示地球中心引力加速度的大小,fe表示其他類型的攝動加速度的大小。與地球中心引力相比,用η=fe/f0來表征各攝動力的量級估算值,如表1所示。其中大氣阻力攝動與第三體引力攝動力均為10-7,但第三體引力攝動要小于大氣阻力攝動,且軌道高度越低,大氣阻力攝動的作用越大。由此可見,對于低軌衛星,J2攝動最重要,其次是大氣阻力攝動[4]。我們可以對攝動力進行簡化,僅考慮J2攝動和大氣阻力攝動的影響。

表1 各種攝動力的量級
在兩體問題中,認為地球是一個均勻球體,衛星運行在地球中心引力場中。但實際上地球形狀是不規則性的,質量分布也不均勻,常將地球引力勢函數表示成如下形式:

式中:r為到觀測點的地心距,φ為地心經度,λ為地心緯度,RE是赤道半徑,Pn,Pnm是正則勒讓德多項式,Cnm,Snm,Jn為引力場系數。
對引力勢函數(2)求梯度,可以得到衛星在地心地固坐標系下的加速度,再通過坐標轉換轉到地心慣性坐標系,即可以得到慣性系下衛星的動力學方程。J2項的影響要比其他項大100倍以上,通常僅考慮J2項攝動的影響。
大氣阻力是非保守力,在其影響下,衛星軌道會發生衰減。大氣阻力攝動對衛星所產生的加速度如下式:

式中:fD為大氣阻力攝動加速度,Vr為衛星相對大氣的相對速度矢量,Vr為衛星相對大氣的相對速度的大小,v和va分別為衛星和大氣相對地心慣性坐標系的速度矢量,ρ為大氣密度,CD為阻力系數。S/m為面質比參數。
另外,分析大氣阻力攝動時,必須考慮以下幾方面因素:
1)CD,阻力系數一般與衛星表面材料有關,常取2.2;
2)S/m,面質比與有效迎風面積有關,即要考慮衛星姿態和形狀,通常選用一個等效面質比來進行計算;
3)va,大氣轉速比較復雜,約為地球自轉速度的0.8~1.4倍,通常認為兩者相等;
4)ρ,大氣密度是極其復雜的問題,與軌道高度、溫度、太陽活動指數等密切相關,變化幅值非常大。已有的大氣模型也很多,如指數模型,Jacchia 77模型、NRLMSISE2000模型等,其中 NRLMSISE2000大氣模型是由美國海軍研究實驗室于2000年在MSISE-90模型的基礎上發展而出,是最新更新的一個大氣模型,與其他模型相比具有很多優勢[5]。
基于2.1節和2.2節的分析,可以得到同時考慮J2項攝動和大氣阻力攝動時,衛星在慣性系下的運動方程如下式:

給定衛星初始時刻的運動狀態后,對式(4)進行數值積分,就可以得到在J2攝動和大氣阻力攝動同時影響下,衛星軌道的演化進程[6]。
在軌道仿真軟件STK中建立一個SSO軌道,將其偏心率、近地點幅角、真近點角設置為0,降交點地方時設置為18時。改變軌道高度、設定軌道預報器就可以得到滿足文章要求的各種仿真條件。
首先分析了面質比對軌道演化的影響并進行仿真,以確立出基本的仿真初始條件。然后分別仿真分析了軌道高度、太陽活動指數以及不同大氣模型對軌道演化的影響。由于瞬時軌道根數波動性比較大,無法表現衛星受攝運動的本質,文中用平均軌道根數下的半長軸作為攝動效應的衡量指標。
由2.2節可知,在研究大氣阻力攝動時,必須考慮衛星面質比的影響,本文研究對象是1 U立方星,邊長為10 cm,質量約為1 kg,所以面質比的變化范圍為:0.01~0.014。在STK中,考慮J2項和大氣模型為 NRLMSISE2000時的大氣攝動,在0.01~0.014之間以0.0005為步長改變面質比,仿真不同面質比對500km SSO立方星軌道演化的影響,利用所得9個數據點繪圖,如圖2所示。橫坐標是面質比,縱坐標代表半長軸的衰減量。

圖2 SSO立方星軌道演化與面質比的關系
從圖2可以看出,面質比對軌道演化有很大影響,仿真分析時必須合理選擇面質比;軌道演化情況大致與面質比成線性關系。文中以下章節的仿真中,面質比默認選定為0.01m2/kg。
以NRLMSISE2000大氣模型為例,研究大氣對不同高度SSO立方星軌道演化的影響,仿真時間180d。NRLMSISE2000大氣模型中,太陽活動指數F10.7按默認值為150。在300~700km范圍內以100km為步長改變軌道高度,所得到的半長軸演化情況匯總如表2所示,從仿真結果中選出400km,600km高度軌道的演化圖,分別如圖3和4所示。由此可以得到180d內軌道衰減的情況。圖3中,在大氣阻力作用下,145d左右時軌道高度直線下降,衛星隕落,壽命終止。圖4中,仿真結束時,軌道衰減了約3.6km。

圖3 400km SSO立方星軌道演化

圖4 600km SSO立方星軌道演化

表2 軌道高度對軌道演化的影響(F10.7=150)
從表2可以看出:1)大氣阻力攝動對軌道演化(半長軸衰減)影響很大,且軌道高度越低影響越明顯;2)高度300km的軌道,壽命約23d。高度為400km時,軌道壽命約150d。而對于500km高的軌道,半長軸衰減了16km。由此可見,軌道的演化情況與軌道高度間存在非線性關系。
太陽活動指數F10.7指太陽發出的波長為10.7cm的電磁輻射強度,其變化沒有嚴格規律,一般由觀測得到。
F10.7的長期預測具有不準確性,以2015年6月1日為例,STK軌道分析軟件8.1.1根據2007年5月24日的觀測數據,預測的F10.7參數為77.7;STK軌道分析軟件9.2.1根據2009年4月29日的觀測數據,預測的F10.7參數為101.4;NOAA根據2012年11月26日的觀測數據,預測的 F10.7參數為117.2[7]。由此可見,實現以年為單位的 F10.7參數預報,誤差很大。應該采用盡可能新的預測數據。
由于F10.7的不確定性,導致大氣密度相應地也存在不確定性。因此,在考慮大氣阻力攝動時,有必要研究太陽活動指數F10.7對軌道演化的影響。
根據太陽活動指數隨年份變化的趨勢,選擇太陽活動高、低年時的 F10.7參數分別為 200和70[8]。然后分別分析太陽活動高、低年時,大氣對不同高度SSO演化情況的影響。
當軌道高度為300~1000km時,可以分別得到太陽活動高、低年時軌道半長軸的演化情況,匯總如表3。

表3 F10.7參數對軌道演化的影響
選取太陽活動低、高年時500km SSO的演化情況,分別如圖5和6所示。由圖5和6可以看出,對于500km的SSO軌道,太陽活動低、高年時的軌道演化情況差距很大,仿真結束時,軌道半長軸分別衰減了2km和37km。

圖5 太陽活動低年時的軌道演化情況

圖6 太陽活動高年時的軌道演化
從表3可以看出,對于同樣高度的SSO立方星,太陽活動高、低年時的軌道演化差別很大;太陽活動指數越高,可以忽略大氣阻力攝動的最低軌道高度也越高。一般認為180天內半長軸衰減量小于1 km時,就可以忽略大氣阻力的影響。
下面以500 km SSO軌道為例,研究軌道高度不變時,半長軸衰減量隨F10.7的變化,以得到F10.7對軌道衰減的影響曲線圖。在70~200之間改變F10.7,步長取10。利用所得數據點繪圖,如圖7所示。由圖7可知,對同一高度的軌道,太陽活動指數越高,軌道衰減就越快;軌道衰減量與F10.7呈現出非線性關系。

圖7 半長軸衰減量與F10.7之間的關系
由2.2節可知,分析大氣阻力攝動時,必須考慮大氣模型的影響。目前已經有多個大氣模型,然而這些模型的統計精度只有15% 左右。
以高度為500 km的典型 SSO軌道為例,選用不同的大氣模型,對軌道半長軸隨時間的演化進行分析,仿真初始條件與上文相同。分別選用了8個大氣模型:Cira 72模型,Exponential-Earth模型,Harris-Priester模型,Jacchia-Roberts模型,Jacchia-1971模型,MSIS 1986模型,US Standard Atmosphere模型以及NRLMSISE 2000模型。各模型中涉及到太陽活動指數F10.7的,都按照默認設置為150。
對于指數大氣模型,仿真時要輸入參考高度和標準高度,按式(5)計算[9]。

式中:H為標準高度,H0為參考高度,r為衛星質心的地心距。μ≈0.1(常取μ<0.2),H0=37.4 km,r0=H0+6371 km。因此,對于500 km軌道,計算可得參考高度是37.4 km,標準高度是60.887 km。
分別設置不同的軌道預報器,以改變大氣模型。仿真180d,得到不同大氣模型下半長軸的演化情況,如表4所示。

表4 不同大氣模型對軌道演化的影響(180d)
從表4可以看出,在不同大氣模型的影響下,軌道演化情況不盡相同,大多數模型的影響效果的量級約為16 km;與其他模型相比,指數大氣模型影響下的半長軸演化情況懸殊,反映出指數大氣模型的精度不夠高。
立方星是近年來納衛星領域的研究熱點。本文對SSO立方星任務的軌道演化情況進行了仿真分析,得出如下結論:
1)對于300~700 km軌道高度的立方體衛星,大氣阻力攝動對軌道壽命有明顯影響,導致半長軸不斷衰減,軌道越低影響越大,且影響結果與軌道壽命成非線性關系;
2)太陽活動對大氣阻力攝動有顯著影響,太陽活動高、低年時的軌道演化情況差距很大。在分析大氣阻力攝動時,必須考慮F10.7參數的影響,且影響效果與F10.7間存在非線性關系;
3)選取不同大氣模型時,軌道演化情況不完全相同。大多數大氣模型之間的軌道演化情況差距并不大,指數模型與其他模型相比存在一定偏差。
綜上所述,考慮軌道高度、太陽活動及大氣模型的耦合影響,對單顆SSO立方星任務軌道的選定具有參考意義。另外,單星的軌道演化是研究立方星編隊的基礎,星間無碰撞的安全飛行高度的確定、星間相對距離的演化情況及編隊構形控制方法可以作為下一步的研究方向。
[1] Sundaramoorthy P P,Gill E,Verhoeven C J M,Bouwmeester J.Two CubeSats with Micro-Propulsion in the QB50 Satellite Network[C]//Proceeding of the 24th Annual AIAA/USU Conference on Small Satellites,SC10-III-3.Washington:AIAA,2010:1-5.
[2] 林來興.立方體星的技術發展和應用前景[J].航天器工程,2013,23(3):90-98.(Lin Laixing.Technology development and application prospects of cubeSat[J].Spacecraft Engineering,2013,22(3):90-97.)
[3] 李軍予,伍保峰,張曉敏.立方體納衛星的發展及其啟示[J].航天器工程,2012,21(3):80-87.(Li Junyu,Wu Baofeng,Zhang Xiaomin.Development of cubeSat and its enlightenment[J].Spacecraft Engineering,2012,21(3):80-87.)
[4] 王融,熊智,等.基于受攝軌道模型的小衛星軌道攝動分析研究[J].航天控制,2007,25(3):66-70.(Wang Rong,Xiong Zhi,et al.Analysis and research of micro satellite orbit perturbation based on the perturbative orbit model[J].Aerospace Control,2007,25(3):66-70.)
[5] 盧明,李智,陳冒銀.NRLMSISE-00大氣模型的分析和驗證[J].裝備指揮技術學院學報,2010,21(4):58-61.(Lu Ming,Li Zhi,Chen Maoyin.Analysis and verification of the NRLMSISE-00 atmospheric model[J].Journal of the Academy of Equipment Command &Technology,2010,21(4):58-61.)
[6] 黃勇,李小將,王志恒,李兆銘.J2項和大氣阻力攝動作用下衛星編隊構型的演化分析[J].航天控制,2013,31(3):62-66.(Huang Yong,Li Xiaojiang,Wang Zhiheng,Li Zhaoming.The evolvement analysis of satellite formation under J2 term and atmosphere drag perturbation[J].Aerospace Control ,2013,31(3):62-66.)
[7] NOAA.Solar Cycle F10.7 cm Radio Flux Progression[EB/OL].[2013-05-20]http://www.noaa.gov.
[8] 劉暾,趙鈞.空間飛行器動力學[M].哈爾濱:哈爾濱工業大學出版社,2003:253-255.(Liu Tun,Zhao Jun.Dynamics of Spacecraft[M].Harbin:Press of Harbin Institute Technology,2003.)
[9] 劉林.航天器軌道理論[M].北京:國防工業出版社,2000.(Liu Ling.Orbit Theory of Spacecraft[M].Beijing:Defense Industry Press,2000.)