雷博持 李明濤

1.中國空間技術(shù)研究院西安分院，西安710100

2.中國科學(xué)院空間科學(xué)與應(yīng)用研究中心，北京100190

立方星這一概念最早由美國加州理工大學(xué)和斯坦福大學(xué)在1999年提出，立方星是標(biāo)準(zhǔn)化的微小衛(wèi)星，尺寸10 cm×10 cm×10 cm，質(zhì)量約1 kg，為一個(gè)基本單元1 U。以此為基礎(chǔ)可以有2 U、3 U等立方星。立方星具有很多特點(diǎn)，不僅成本低，開發(fā)時(shí)間也較短，易于實(shí)現(xiàn)標(biāo)準(zhǔn)化、模塊化，易于技術(shù)更新［1-2］。

近年來，利用多衛(wèi)星系統(tǒng)的空間任務(wù)越來越多，借助標(biāo)準(zhǔn)化的微小衛(wèi)星來實(shí)現(xiàn)分布式空間系統(tǒng)極具發(fā)展前景。分布式空間任務(wù)的出現(xiàn)，拓展了立方星在科學(xué)領(lǐng)域的應(yīng)用。目前已有越來越多的立方星被送上太空。國外的立方星應(yīng)用主要分為3個(gè)方面:以空間科學(xué)和教學(xué)為核心目標(biāo)，如荷蘭代爾夫特理工大學(xué)的 Delfi-c3;以載荷任務(wù)為目標(biāo)，如 QuakeSat搭載了地震電磁探測載荷;也有以新技術(shù)驗(yàn)證為目標(biāo)的立方星任務(wù)，如加拿大的先進(jìn)納星太空實(shí)驗(yàn)室CANX衛(wèi)星系列［3］。國內(nèi)很多高校及科研機(jī)構(gòu)也在積極開展立方星的相關(guān)研究。

體積小、重量輕的自身特點(diǎn)決定了立方星大多沒有在軌機(jī)動(dòng)能力，考慮無主動(dòng)控制時(shí)立方星軌道的自然演化就成為一個(gè)重要研究方向。太陽同步軌道(SSO)具有星上光照條件好、過同一緯度的當(dāng)?shù)貢r(shí)間相同等特點(diǎn)，已被廣泛應(yīng)用于各種空間任務(wù)。因此，本文選取了SSO立方星任務(wù)進(jìn)行軌道演化研究，考慮地球非球形引力項(xiàng)J2和大氣阻力攝動(dòng)的影響，以平均半長軸隨時(shí)間的變化作為攝動(dòng)效應(yīng)的指標(biāo)。以軌道分析軟件為仿真平臺，先進(jìn)行了攝動(dòng)力簡化以及立方星面質(zhì)比對軌道演化的影響分析，以確立數(shù)值仿真的前提條件;接著研究了軌道高度、太陽活動(dòng)指數(shù)和大氣模型這3個(gè)方面對SSO立方星軌道演化的影響。最終得出結(jié)論。

1 衛(wèi)星軌道動(dòng)力學(xué)模型

1.1 坐標(biāo)系定義

研究地球非球形攝動(dòng)及大氣阻力攝動(dòng)對軌道演化的影響時(shí)，涉及到地心慣性坐標(biāo)系和地心地固坐標(biāo)系，如圖1所示。Oxiyizi(標(biāo)記為 Si)為地心慣性坐標(biāo)系，Oxeyeze(標(biāo)記為Se)為地心地固坐標(biāo)系，α為Greenwich赤經(jīng)。

地心慣性坐標(biāo)系Si:坐標(biāo)系原點(diǎn)在地心，Z軸垂直于地球赤道平面;XY軸在地球赤道平面內(nèi)，X軸指向春分點(diǎn)位置;Y軸滿足右手螺旋定則。動(dòng)力學(xué)方程一般建立在慣性坐標(biāo)系中。

圖1 地心慣性坐標(biāo)系和地心地固坐標(biāo)系

地心地固坐標(biāo)系Se:坐標(biāo)系原點(diǎn)在地心，Z軸垂直于地球赤道平面;XY軸在地球赤道平面內(nèi)，X軸指向 Greenwich子午圈與赤道交點(diǎn);Y軸滿足右手螺旋定則。地球非球形攝動(dòng)一般在地心地固坐標(biāo)系中表示。

1.2 動(dòng)力學(xué)模型

在慣性坐標(biāo)系下，衛(wèi)星動(dòng)力學(xué)模型可用式(1)表示。

式中:r為地心距矢量，μ是地心引力常數(shù)，r為地心距大小，f為地球中心引力之外的所有攝動(dòng)力的合力和包括地球非球形引力、大氣阻力、太陽光壓、第三體引力攝動(dòng)等。

將式(1)在慣性坐標(biāo)系下展開并進(jìn)行積分，可以得到任意時(shí)刻衛(wèi)星在慣性坐標(biāo)系下的運(yùn)動(dòng)狀態(tài)。

2 地球非球形攝動(dòng)和大氣阻力攝動(dòng)分析

對于低軌道來說，衛(wèi)星所受到的攝動(dòng)力主要是地球非球形引力攝動(dòng)、大氣阻力攝動(dòng)、太陽光壓攝動(dòng)以及第三體引力攝動(dòng)。設(shè)f0表示地球中心引力加速度的大小，fe表示其他類型的攝動(dòng)加速度的大小。與地球中心引力相比，用η=fe/f0來表征各攝動(dòng)力的量級估算值，如表1所示。其中大氣阻力攝動(dòng)與第三體引力攝動(dòng)力均為10-7，但第三體引力攝動(dòng)要小于大氣阻力攝動(dòng)，且軌道高度越低，大氣阻力攝動(dòng)的作用越大。由此可見，對于低軌衛(wèi)星，J2攝動(dòng)最重要，其次是大氣阻力攝動(dòng)［4］。我們可以對攝動(dòng)力進(jìn)行簡化，僅考慮J2攝動(dòng)和大氣阻力攝動(dòng)的影響。

表1 各種攝動(dòng)力的量級

2.1 J2項(xiàng)攝動(dòng)

在兩體問題中，認(rèn)為地球是一個(gè)均勻球體，衛(wèi)星運(yùn)行在地球中心引力場中。但實(shí)際上地球形狀是不規(guī)則性的，質(zhì)量分布也不均勻，常將地球引力勢函數(shù)表示成如下形式:

式中:r為到觀測點(diǎn)的地心距，φ為地心經(jīng)度，λ為地心緯度，RE是赤道半徑，Pn，Pnm是正則勒讓德多項(xiàng)式，Cnm，Snm，Jn為引力場系數(shù)。

對引力勢函數(shù)(2)求梯度，可以得到衛(wèi)星在地心地固坐標(biāo)系下的加速度，再通過坐標(biāo)轉(zhuǎn)換轉(zhuǎn)到地心慣性坐標(biāo)系，即可以得到慣性系下衛(wèi)星的動(dòng)力學(xué)方程。J2項(xiàng)的影響要比其他項(xiàng)大100倍以上，通常僅考慮J2項(xiàng)攝動(dòng)的影響。

2.2 大氣阻力攝動(dòng)

大氣阻力是非保守力，在其影響下，衛(wèi)星軌道會(huì)發(fā)生衰減。大氣阻力攝動(dòng)對衛(wèi)星所產(chǎn)生的加速度如下式:

式中:fD為大氣阻力攝動(dòng)加速度，Vr為衛(wèi)星相對大氣的相對速度矢量，Vr為衛(wèi)星相對大氣的相對速度的大小，v和va分別為衛(wèi)星和大氣相對地心慣性坐標(biāo)系的速度矢量，ρ為大氣密度，CD為阻力系數(shù)。S/m為面質(zhì)比參數(shù)。

另外，分析大氣阻力攝動(dòng)時(shí)，必須考慮以下幾方面因素:

1)CD，阻力系數(shù)一般與衛(wèi)星表面材料有關(guān)，常取2.2;

2)S/m，面質(zhì)比與有效迎風(fēng)面積有關(guān)，即要考慮衛(wèi)星姿態(tài)和形狀，通常選用一個(gè)等效面質(zhì)比來進(jìn)行計(jì)算;

3)va，大氣轉(zhuǎn)速比較復(fù)雜，約為地球自轉(zhuǎn)速度的0.8～1.4倍，通常認(rèn)為兩者相等;

4)ρ，大氣密度是極其復(fù)雜的問題，與軌道高度、溫度、太陽活動(dòng)指數(shù)等密切相關(guān)，變化幅值非常大。已有的大氣模型也很多，如指數(shù)模型，Jacchia 77模型、NRLMSISE2000模型等，其中 NRLMSISE2000大氣模型是由美國海軍研究實(shí)驗(yàn)室于2000年在MSISE-90模型的基礎(chǔ)上發(fā)展而出，是最新更新的一個(gè)大氣模型，與其他模型相比具有很多優(yōu)勢［5］。

2.3 J2項(xiàng)攝動(dòng)和大氣攝動(dòng)影響下的衛(wèi)星運(yùn)動(dòng)方程

基于2.1節(jié)和2.2節(jié)的分析，可以得到同時(shí)考慮J2項(xiàng)攝動(dòng)和大氣阻力攝動(dòng)時(shí)，衛(wèi)星在慣性系下的運(yùn)動(dòng)方程如下式:

給定衛(wèi)星初始時(shí)刻的運(yùn)動(dòng)狀態(tài)后，對式(4)進(jìn)行數(shù)值積分，就可以得到在J2攝動(dòng)和大氣阻力攝動(dòng)同時(shí)影響下，衛(wèi)星軌道的演化進(jìn)程［6］。

3 仿真分析

在軌道仿真軟件STK中建立一個(gè)SSO軌道，將其偏心率、近地點(diǎn)幅角、真近點(diǎn)角設(shè)置為0，降交點(diǎn)地方時(shí)設(shè)置為18時(shí)。改變軌道高度、設(shè)定軌道預(yù)報(bào)器就可以得到滿足文章要求的各種仿真條件。

首先分析了面質(zhì)比對軌道演化的影響并進(jìn)行仿真，以確立出基本的仿真初始條件。然后分別仿真分析了軌道高度、太陽活動(dòng)指數(shù)以及不同大氣模型對軌道演化的影響。由于瞬時(shí)軌道根數(shù)波動(dòng)性比較大，無法表現(xiàn)衛(wèi)星受攝運(yùn)動(dòng)的本質(zhì)，文中用平均軌道根數(shù)下的半長軸作為攝動(dòng)效應(yīng)的衡量指標(biāo)。

3.1 面質(zhì)比的選定

由2.2節(jié)可知，在研究大氣阻力攝動(dòng)時(shí)，必須考慮衛(wèi)星面質(zhì)比的影響，本文研究對象是1 U立方星，邊長為10 cm，質(zhì)量約為1 kg，所以面質(zhì)比的變化范圍為:0.01～0.014。在STK中，考慮J2項(xiàng)和大氣模型為 NRLMSISE2000時(shí)的大氣攝動(dòng)，在0.01～0.014之間以0.0005為步長改變面質(zhì)比，仿真不同面質(zhì)比對500km SSO立方星軌道演化的影響，利用所得9個(gè)數(shù)據(jù)點(diǎn)繪圖，如圖2所示。橫坐標(biāo)是面質(zhì)比，縱坐標(biāo)代表半長軸的衰減量。

圖2 SSO立方星軌道演化與面質(zhì)比的關(guān)系

從圖2可以看出，面質(zhì)比對軌道演化有很大影響，仿真分析時(shí)必須合理選擇面質(zhì)比;軌道演化情況大致與面質(zhì)比成線性關(guān)系。文中以下章節(jié)的仿真中，面質(zhì)比默認(rèn)選定為0.01m2/kg。

3.2 軌道高度對軌道演化的影響

以NRLMSISE2000大氣模型為例，研究大氣對不同高度SSO立方星軌道演化的影響，仿真時(shí)間180d。NRLMSISE2000大氣模型中，太陽活動(dòng)指數(shù)F10.7按默認(rèn)值為150。在300～700km范圍內(nèi)以100km為步長改變軌道高度，所得到的半長軸演化情況匯總?cè)绫?所示，從仿真結(jié)果中選出400km，600km高度軌道的演化圖，分別如圖3和4所示。由此可以得到180d內(nèi)軌道衰減的情況。圖3中，在大氣阻力作用下，145d左右時(shí)軌道高度直線下降，衛(wèi)星隕落，壽命終止。圖4中，仿真結(jié)束時(shí)，軌道衰減了約3.6km。

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

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

表2 軌道高度對軌道演化的影響(F10.7=150)

從表2可以看出:1)大氣阻力攝動(dòng)對軌道演化(半長軸衰減)影響很大，且軌道高度越低影響越明顯;2)高度300km的軌道，壽命約23d。高度為400km時(shí)，軌道壽命約150d。而對于500km高的軌道，半長軸衰減了16km。由此可見，軌道的演化情況與軌道高度間存在非線性關(guān)系。

3.3 太陽活動(dòng)對軌道演化的影響

太陽活動(dòng)指數(shù)F10.7指太陽發(fā)出的波長為10.7cm的電磁輻射強(qiáng)度，其變化沒有嚴(yán)格規(guī)律，一般由觀測得到。

F10.7的長期預(yù)測具有不準(zhǔn)確性，以2015年6月1日為例，STK軌道分析軟件8.1.1根據(jù)2007年5月24日的觀測數(shù)據(jù)，預(yù)測的F10.7參數(shù)為77.7;STK軌道分析軟件9.2.1根據(jù)2009年4月29日的觀測數(shù)據(jù)，預(yù)測的F10.7參數(shù)為101.4;NOAA根據(jù)2012年11月26日的觀測數(shù)據(jù)，預(yù)測的 F10.7參數(shù)為117.2［7］。由此可見，實(shí)現(xiàn)以年為單位的 F10.7參數(shù)預(yù)報(bào)，誤差很大。應(yīng)該采用盡可能新的預(yù)測數(shù)據(jù)。

由于F10.7的不確定性，導(dǎo)致大氣密度相應(yīng)地也存在不確定性。因此，在考慮大氣阻力攝動(dòng)時(shí)，有必要研究太陽活動(dòng)指數(shù)F10.7對軌道演化的影響。

根據(jù)太陽活動(dòng)指數(shù)隨年份變化的趨勢，選擇太陽活動(dòng)高、低年時(shí)的 F10.7參數(shù)分別為 200和70［8］。然后分別分析太陽活動(dòng)高、低年時(shí)，大氣對不同高度SSO演化情況的影響。

當(dāng)軌道高度為300～1000km時(shí)，可以分別得到太陽活動(dòng)高、低年時(shí)軌道半長軸的演化情況，匯總?cè)绫?。

表3 F10.7參數(shù)對軌道演化的影響

選取太陽活動(dòng)低、高年時(shí)500km SSO的演化情況，分別如圖5和6所示。由圖5和6可以看出，對于500km的SSO軌道，太陽活動(dòng)低、高年時(shí)的軌道演化情況差距很大，仿真結(jié)束時(shí)，軌道半長軸分別衰減了2km和37km。

圖5 太陽活動(dòng)低年時(shí)的軌道演化情況

圖6 太陽活動(dòng)高年時(shí)的軌道演化

從表3可以看出，對于同樣高度的SSO立方星，太陽活動(dòng)高、低年時(shí)的軌道演化差別很大;太陽活動(dòng)指數(shù)越高，可以忽略大氣阻力攝動(dòng)的最低軌道高度也越高。一般認(rèn)為180天內(nèi)半長軸衰減量小于1 km時(shí)，就可以忽略大氣阻力的影響。

下面以500 km SSO軌道為例，研究軌道高度不變時(shí)，半長軸衰減量隨F10.7的變化，以得到F10.7對軌道衰減的影響曲線圖。在70～200之間改變F10.7，步長取10。利用所得數(shù)據(jù)點(diǎn)繪圖，如圖7所示。由圖7可知，對同一高度的軌道，太陽活動(dòng)指數(shù)越高，軌道衰減就越快;軌道衰減量與F10.7呈現(xiàn)出非線性關(guān)系。

圖7 半長軸衰減量與F10.7之間的關(guān)系

3.4 大氣模型對軌道演化的影響

由2.2節(jié)可知，分析大氣阻力攝動(dòng)時(shí)，必須考慮大氣模型的影響。目前已經(jīng)有多個(gè)大氣模型，然而這些模型的統(tǒng)計(jì)精度只有15% 左右。

以高度為500 km的典型 SSO軌道為例，選用不同的大氣模型，對軌道半長軸隨時(shí)間的演化進(jìn)行分析，仿真初始條件與上文相同。分別選用了8個(gè)大氣模型:Cira 72模型，Exponential-Earth模型，Harris-Priester模型，Jacchia-Roberts模型，Jacchia-1971模型，MSIS 1986模型，US Standard Atmosphere模型以及NRLMSISE 2000模型。各模型中涉及到太陽活動(dòng)指數(shù)F10.7的，都按照默認(rèn)設(shè)置為150。

對于指數(shù)大氣模型，仿真時(shí)要輸入?yún)⒖几叨群蜆?biāo)準(zhǔn)高度，按式(5)計(jì)算［9］。

式中:H為標(biāo)準(zhǔn)高度，H0為參考高度，r為衛(wèi)星質(zhì)心的地心距。μ≈0.1(常取μ＜0.2)，H0=37.4 km，r0=H0+6371 km。因此，對于500 km軌道，計(jì)算可得參考高度是37.4 km，標(biāo)準(zhǔn)高度是60.887 km。

分別設(shè)置不同的軌道預(yù)報(bào)器，以改變大氣模型。仿真180d，得到不同大氣模型下半長軸的演化情況，如表4所示。

表4 不同大氣模型對軌道演化的影響(180d)

從表4可以看出，在不同大氣模型的影響下，軌道演化情況不盡相同，大多數(shù)模型的影響效果的量級約為16 km;與其他模型相比，指數(shù)大氣模型影響下的半長軸演化情況懸殊，反映出指數(shù)大氣模型的精度不夠高。

4 結(jié)論

立方星是近年來納衛(wèi)星領(lǐng)域的研究熱點(diǎn)。本文對SSO立方星任務(wù)的軌道演化情況進(jìn)行了仿真分析，得出如下結(jié)論:

1)對于300～700 km軌道高度的立方體衛(wèi)星，大氣阻力攝動(dòng)對軌道壽命有明顯影響，導(dǎo)致半長軸不斷衰減，軌道越低影響越大，且影響結(jié)果與軌道壽命成非線性關(guān)系;

2)太陽活動(dòng)對大氣阻力攝動(dòng)有顯著影響，太陽活動(dòng)高、低年時(shí)的軌道演化情況差距很大。在分析大氣阻力攝動(dòng)時(shí)，必須考慮F10.7參數(shù)的影響，且影響效果與F10.7間存在非線性關(guān)系;

3)選取不同大氣模型時(shí)，軌道演化情況不完全相同。大多數(shù)大氣模型之間的軌道演化情況差距并不大，指數(shù)模型與其他模型相比存在一定偏差。

綜上所述，考慮軌道高度、太陽活動(dòng)及大氣模型的耦合影響，對單顆SSO立方星任務(wù)軌道的選定具有參考意義。另外，單星的軌道演化是研究立方星編隊(duì)的基礎(chǔ)，星間無碰撞的安全飛行高度的確定、星間相對距離的演化情況及編隊(duì)構(gòu)形控制方法可以作為下一步的研究方向。

［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］ 林來興.立方體星的技術(shù)發(fā)展和應(yīng)用前景［J］.航天器工程，2013，23(3):90-98.(Lin Laixing.Technology development and application prospects of cubeSat［J］.Spacecraft Engineering，2013，22(3):90-97.)

［3］ 李軍予，伍保峰，張曉敏.立方體納衛(wèi)星的發(fā)展及其啟示［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］ 王融，熊智，等.基于受攝軌道模型的小衛(wèi)星軌道攝動(dòng)分析研究［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大氣模型的分析和驗(yàn)證［J］.裝備指揮技術(shù)學(xué)院學(xué)報(bào)，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項(xiàng)和大氣阻力攝動(dòng)作用下衛(wèi)星編隊(duì)構(gòu)型的演化分析［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］ 劉暾，趙鈞.空間飛行器動(dòng)力學(xué)［M］.哈爾濱:哈爾濱工業(yè)大學(xué)出版社，2003:253-255.(Liu Tun，Zhao Jun.Dynamics of Spacecraft［M］.Harbin:Press of Harbin Institute Technology，2003.)

［9］ 劉林.航天器軌道理論［M］.北京:國防工業(yè)出版社，2000.(Liu Ling.Orbit Theory of Spacecraft［M］.Beijing:Defense Industry Press，2000.)