微小推力下小衛(wèi)星的軌道遞推與推力標(biāo)定

張軍華，陳建林，孫沖，袁建平

1. 西北工業(yè)大學(xué) 航天學(xué)院航天飛行動(dòng)力學(xué)技術(shù)重點(diǎn)實(shí)驗(yàn)室，西安 710072

2. 西北工業(yè)大學(xué) 民航學(xué)院，西安 710072

1 引言

隨著商業(yè)航天的發(fā)展，低軌小衛(wèi)星平臺(tái)受到越來(lái)越廣泛地關(guān)注。小衛(wèi)星相對(duì)于大衛(wèi)星來(lái)說(shuō)，由于體積小、質(zhì)量輕、成本低，具有更靈活的應(yīng)用潛力，因此近年來(lái)小衛(wèi)星廣泛應(yīng)用于近地航天任務(wù)，已經(jīng)成為航天工程的研究熱點(diǎn)之一[1-2]。

通常，小衛(wèi)星受限于質(zhì)量約束，無(wú)法像大衛(wèi)星平臺(tái)一樣攜帶大量的燃料，因此高比沖微小推力推進(jìn)器小衛(wèi)星的常用推進(jìn)方式，主要包括電推進(jìn)、太陽(yáng)帆推進(jìn)、電帆推進(jìn)、冷氣推進(jìn)等[3-10]。隨著航天技術(shù)的發(fā)展，當(dāng)前的微小衛(wèi)星任務(wù)對(duì)軌道預(yù)測(cè)和估計(jì)的精度要求不斷提高，這不僅對(duì)衛(wèi)星軌道動(dòng)力學(xué)建模精度提出了較高的要求，還對(duì)軌道預(yù)測(cè)和估計(jì)方法提出了精度要求。用于軌道的預(yù)測(cè)和估計(jì)的數(shù)值積分方法，如Dormand-Prince方法、Gauss-Jackson方法、Adams-Bashforth-Moulton方法等往往在誤差控制、穩(wěn)定性等方面存在問(wèn)題[11-15]。而采用簡(jiǎn)化的動(dòng)力學(xué)模型往往導(dǎo)致軌道預(yù)測(cè)的精度無(wú)法滿足實(shí)際需求[16-19]。

進(jìn)一步，由于微小推力的量級(jí)小，難以嚴(yán)格在地面完成高精度的推力標(biāo)定，需要衛(wèi)星在軌對(duì)其推力進(jìn)行標(biāo)定。常見(jiàn)的推力標(biāo)定包括參數(shù)標(biāo)定、飛輪標(biāo)定、矢量標(biāo)定等[20-23]，這些標(biāo)定方法采用衛(wèi)星載荷對(duì)標(biāo)定所需參數(shù)進(jìn)行估計(jì)。但小衛(wèi)星搭載載荷有限，上述標(biāo)定方法無(wú)法適用于小衛(wèi)星的推力在軌標(biāo)定。另一方面，微小推進(jìn)器存在推力損失等不確定性因素，如何保證推力在軌標(biāo)定的精度，是個(gè)亟待解決的問(wèn)題。

本文針對(duì)小衛(wèi)星的高精度軌道預(yù)測(cè)和微小推力在軌標(biāo)定問(wèn)題，對(duì)主要的攝動(dòng)力進(jìn)行了高精度建模。利用變步長(zhǎng)龍格庫(kù)塔積分RK78對(duì)衛(wèi)星軌道進(jìn)行了外推，通過(guò)與STK仿真結(jié)果對(duì)比驗(yàn)證了模型的精確性；使用無(wú)跡卡爾曼濾波器(UKF)同時(shí)估計(jì)小衛(wèi)星軌道并標(biāo)定衛(wèi)星推力，證明了在軌估計(jì)算法的可行性。本文對(duì)軌道遞推和推力標(biāo)定算法的研究，可以應(yīng)用于小衛(wèi)星的實(shí)時(shí)在軌運(yùn)算中，實(shí)時(shí)掌握小衛(wèi)星自身的軌道狀態(tài)和推力大小，同時(shí)算法的速度較快，有利于小衛(wèi)星在線運(yùn)算任務(wù)的實(shí)現(xiàn)。

2 微小推力下小衛(wèi)星動(dòng)力學(xué)模型

2.1 問(wèn)題描述

小衛(wèi)星在近地軌道上運(yùn)動(dòng)時(shí)，除了受到地球中心引力場(chǎng)的引力外，還受到地球非球形引力、太陽(yáng)光壓和第三體引力等攝動(dòng)影響，小衛(wèi)星采用微小推力進(jìn)行軌道控制。對(duì)微小推力下小衛(wèi)星的在軌運(yùn)行情況進(jìn)行預(yù)測(cè)，不僅需要精確的軌道運(yùn)動(dòng)模型，而且對(duì)軌道預(yù)測(cè)和估計(jì)的精度也有較高的要求。此外，由于微小推力量級(jí)小，且小衛(wèi)星任務(wù)過(guò)程中存在不確定性，需要將微小推力進(jìn)行在軌標(biāo)定。

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

首先建立三個(gè)坐標(biāo)系[24]，如圖1所示。

圖1 坐標(biāo)系示意圖Fig.1 Illustration of coordinate frames

地球慣性坐標(biāo)系EXYZ(JDate)，坐標(biāo)系原點(diǎn)在地心E，X軸指向當(dāng)前時(shí)刻春分點(diǎn)，Z軸垂直于地球赤道平面向上，Y軸滿足右手法則。

地球慣性坐標(biāo)系J2 000，坐標(biāo)系原點(diǎn)在地心E，X'軸指向J2 000時(shí)刻的平春分點(diǎn)，Z軸垂直于地球赤道平面向上，Y'軸滿足右手法則。

衛(wèi)星軌道坐標(biāo)系Sxyz，坐標(biāo)系原點(diǎn)在衛(wèi)星質(zhì)心S，x軸沿衛(wèi)星矢徑方向，從地心指向衛(wèi)星，y軸垂直衛(wèi)星矢徑并指向速度方向，z軸垂直于軌道平面，遵循右手法則。

除此3個(gè)坐標(biāo)系外，建立地球固聯(lián)坐標(biāo)系(ECEF)，其定義為：坐標(biāo)系原點(diǎn)在地心，X''軸從地心指向赤道與國(guó)際參考子午線的交點(diǎn)，Z''軸沿著地球自轉(zhuǎn)軸從地心指向北極點(diǎn)，根據(jù)右手準(zhǔn)則Y''軸從地心指向赤道與90°東經(jīng)子午線的交點(diǎn)。

另外，假設(shè)衛(wèi)星本身采用三軸對(duì)地穩(wěn)定，因此衛(wèi)星的本體坐標(biāo)系與軌道坐標(biāo)系是相重合的。

衛(wèi)星受到微小推力加速度aF的作用，且該微小作用力在本體坐標(biāo)系中的方向不變，如圖2所示。因此，aF在軌道坐標(biāo)系下可以由θ和φ兩個(gè)方位角表示，aF=[aFcosφcosθ,aFcosφsinθ,aFsinφ]T。

圖2 微小推力在軌道坐標(biāo)系下的示意Fig.2 Illustration of micro-thrust in frame Sxyz

衛(wèi)星受到微小推力、地球扁率(J2、J22等)、太陽(yáng)光壓以及第三體(日月)引力等攝動(dòng)，則衛(wèi)星在地心慣性坐標(biāo)系下的動(dòng)力學(xué)方程可以描述為：

(1)

(2)

式中：ax,nm、ay,nm和az,nm分別表示地球形狀攝動(dòng)加速度在地球固聯(lián)坐標(biāo)系中三個(gè)方向上的分量：

(3)

式中：R⊕表示地球半徑；地勢(shì)系數(shù)Cnm,Snm由EGM96S引力模型給出；Vnm和Wnm分別為：

(4)

(5)

在地球慣性系下的太陽(yáng)光壓和太陽(yáng)、月球等第三體引力引起的攝動(dòng)力為：

(6)

(7)

式中：r⊕,rSun,rMoon,r分別表示地球慣性系下衛(wèi)星到太陽(yáng)的位置矢量；太陽(yáng)的位置矢量；月球的位置矢量和衛(wèi)星的位置矢量；G表示萬(wàn)有引力常數(shù)；AU表示地球到太陽(yáng)的距離；P表示單位AU距離下的太陽(yáng)光壓力；ν表示地球的陰影函數(shù)；M表示太陽(yáng)或者月球的質(zhì)量；Cr表示輻射壓力系數(shù)；m表示衛(wèi)星的質(zhì)量；A表示衛(wèi)星的受照橫截面積。

因此，其他攝動(dòng)總加速度可以表示為：

(8)

式中：aFx、aFy、aFz、aJx、aJy、aJz、adx、ady、adz分別表示微小推力加速度，地球形狀攝動(dòng)加速度，其他攝動(dòng)加速度在慣性坐標(biāo)系下的三個(gè)坐標(biāo)軸的分量。

3 衛(wèi)星軌道遞推算法

衛(wèi)星常見(jiàn)的星上定軌方式是基于GNSS信號(hào)，但是由于GNSS輸出的軌道信息一般在位置和速度上存在誤差，直接應(yīng)用GNSS信號(hào)進(jìn)行軌道確定，會(huì)造成衛(wèi)星軌道定位不精確。為了解決這一問(wèn)題，提出一種變步長(zhǎng)龍格庫(kù)塔算法RK78進(jìn)行軌道外推，具體的算法模型如下：

Runge-Kutta方法是一種間接采用泰勒級(jí)數(shù)展開(kāi)而求解常微分方程初值問(wèn)題的數(shù)值方法。該方法用某些點(diǎn)上函數(shù)值的不同組合來(lái)近似初值問(wèn)題的解析解，從而避免函數(shù)高階偏導(dǎo)數(shù)的計(jì)算，并獲得高階的方法[26]。將Runge-Kutta公式與初值問(wèn)題解的泰勒展開(kāi)式進(jìn)行比較，使之有某些項(xiàng)相同，便可以確定組合系數(shù)，從而建立所需的方法。

由m個(gè)函數(shù)值線性組合構(gòu)成的顯式Runge-Kutta公式為：

(9)

式中：K1和Ki分別表示時(shí)間段開(kāi)始時(shí)的斜率和每一時(shí)間段中點(diǎn)的斜率；ci,αi,βij均為待定系數(shù)；h為步長(zhǎng)。可以看出，方程(9)含有m個(gè)函數(shù)值，所以稱為m階Runge-Kutta公式。本文使用的RK78算法表示采用七階-八階變步長(zhǎng)Runge-Kutta算法，它用七階方法提供候選解，八階方法控制誤差，是一種自適應(yīng)步長(zhǎng)的常微分方程數(shù)值解法，其整體截?cái)嗾`差為O(h8) 。其變步長(zhǎng)的原理如下：

對(duì)于一步h1對(duì)應(yīng)其誤差估計(jì)Δ1，而h0對(duì)應(yīng)其誤差估計(jì)Δ0，兩者關(guān)系為：

(10)

針對(duì)兩個(gè)誤差估計(jì)值，如果Δ1>Δ0，可以根據(jù)方程(10)確定所要減小步長(zhǎng)的大小；如果Δ1<Δ0，可以根據(jù)方程(10)確定能夠在安全范圍內(nèi)增加步長(zhǎng)的大小。

針對(duì)微小推力軌道外推問(wèn)題，根據(jù)上述Runge-Kutta公式，使用動(dòng)力學(xué)模型(8)帶入到七階Runge-Kutta公式(9)中。在此基礎(chǔ)上，設(shè)置合適的變步長(zhǎng)，利用八階Runge-Kutta公式設(shè)置截?cái)嗾`差，通過(guò)積分和迭代，對(duì)衛(wèi)星軌道進(jìn)行外推。

4 基于UKF的推力在軌標(biāo)定算法

微小推力由于在軌保持時(shí)間較長(zhǎng)，為了保證微小推力在軌工作的有效性，需要對(duì)推力進(jìn)行在軌標(biāo)定。本文采用UKF對(duì)微小推力進(jìn)行高精度推力在軌標(biāo)定。

4.1 UKF的基本原理

無(wú)跡卡爾曼濾波(Unscented Kalman Filter，UKF)，是無(wú)跡變換(Unscented Transform，UT)與標(biāo)準(zhǔn)卡爾曼濾波體系的結(jié)合，通過(guò)無(wú)跡變換使非線性系統(tǒng)方程適用于線性假設(shè)下的標(biāo)準(zhǔn)卡爾曼體系，該算法的核心思想是采用UT變換，利用一組Sigma采樣點(diǎn)來(lái)描述隨機(jī)變量的高斯分布，然后通過(guò)非線性函數(shù)的傳遞，再利用加權(quán)統(tǒng)計(jì)線性回歸技術(shù)來(lái)近似非線性函數(shù)的后驗(yàn)均值和方差[27]。

UT變換的特點(diǎn)是對(duì)非線性函數(shù)的概率密度分布進(jìn)行近似，而不是對(duì)非線性函數(shù)進(jìn)行近似，即使系統(tǒng)模型復(fù)雜，也不增加算法實(shí)現(xiàn)的難度；而且所得到的非線性函數(shù)的統(tǒng)計(jì)量的準(zhǔn)確性可以達(dá)到三階；除此之外，它不需要計(jì)算雅可比矩陣，可以處理不可導(dǎo)非線性函數(shù)。

算法模型如下：

設(shè)非線性系統(tǒng)：

(11)

式中：Xk為nx維的系統(tǒng)狀態(tài)向量；Zk為nz維的系統(tǒng)觀測(cè)向量；vk為系統(tǒng)噪聲；協(xié)方差矩陣為Pv；nk為觀測(cè)噪聲；協(xié)方差矩陣為Pn；式中：vk、nk都是高斯白噪聲，且互不相關(guān)。算法具體計(jì)算步驟如下：

1)選定濾波初值。

(12)

2)計(jì)算k-1時(shí)刻的2n+1個(gè)σ(Sigma)樣本點(diǎn)，構(gòu)造(2n+1)維向量χi為σ向量：

(13)

λ=α2(n+k)-n

(14)

3)狀態(tài)更新，計(jì)算k時(shí)刻的一步預(yù)測(cè)模型值。

(15)

4)計(jì)算k時(shí)刻的一步預(yù)測(cè)增廣樣本點(diǎn)。

(16)

5)觀測(cè)更新，計(jì)算k時(shí)刻的觀測(cè)變量的一步預(yù)測(cè)均值、方差和協(xié)方差。

(17)

6)計(jì)算增益矩陣。

(18)

7)計(jì)算濾波值。

(19)

4.2 在軌標(biāo)定算法

利用UKF對(duì)微小推力的大小及方向進(jìn)行在軌標(biāo)定的過(guò)程如下：

(1)動(dòng)力學(xué)模型改進(jìn)

將笛卡爾動(dòng)力學(xué)模型進(jìn)行狀態(tài)擴(kuò)維，即包含帶估計(jì)的參數(shù)。

(20)

(2)基本方程

利用上一部分的軌道遞推算法，對(duì)一個(gè)測(cè)量周期內(nèi)的航天器狀態(tài)進(jìn)行預(yù)測(cè)，寫(xiě)出預(yù)測(cè)方程；然后根據(jù)UKF基本原理中的公式分別寫(xiě)出狀態(tài)估值方程、濾波增益方程、一步預(yù)測(cè)均方誤差方程和估計(jì)均方誤差方程等。

(3)更新

1)設(shè)置位置誤差和速度誤差，設(shè)計(jì)GNSS數(shù)值生成模擬器。

2)融合GNSS的測(cè)量值，在UKF濾波器內(nèi)對(duì)狀態(tài)預(yù)測(cè)進(jìn)行更新。實(shí)時(shí)估計(jì)狀態(tài)x、m、aF、φ、θ。

5 仿真校驗(yàn)

5.1 軌道遞推仿真結(jié)果

針對(duì)微小推力軌道外推問(wèn)題，使用變步長(zhǎng)龍格庫(kù)塔RK78算法進(jìn)行衛(wèi)星軌道預(yù)測(cè)。初始參數(shù)設(shè)置如下：

1)設(shè)置初始軌道為圓軌道，參數(shù)軌道半徑為755 3 km，軌道傾角86.5°，偏心率0，升交點(diǎn)赤經(jīng)40°，近地點(diǎn)幅角60°，真近點(diǎn)角30°，地球引力常數(shù)3.986×1014m3/s2。通過(guò)軌道要素變換，可以求解衛(wèi)星的初始位置和速度，如表1所示。

表1 初始條件

2)設(shè)置衛(wèi)星質(zhì)量為194 kg，微小推力12 mN，平均分配到三軸推力[6.928 2,6.928 2,6.928 2]TmN，推力加速度為[3.571 24,3.571 24,3.571 24]T×10-5m/s2。此外，考慮到微小推力較小，衛(wèi)星推進(jìn)器的質(zhì)量流量忽略不計(jì)。

3)設(shè)置太陽(yáng)光壓攝動(dòng)參數(shù)：輻射壓力系數(shù)為1.21，衛(wèi)星表面積為3.88 m2，地球到太陽(yáng)的距離AU=149 597 870 000 m，單位AU距離下的太陽(yáng)光壓力P=4.56×10-6N/m2。

4)設(shè)置太陽(yáng)引力攝動(dòng)太陽(yáng)引力常數(shù)為：1.327 124 38×1020m3/s2，月球引力常數(shù)為：4.902 793 455×1012m3/s2。

利用C++軟件編寫(xiě)整個(gè)推演過(guò)程的程序，和STK軟件模擬衛(wèi)星軌道推演24 h，比較不同情況下所編寫(xiě)的C++程序和STK的結(jié)果對(duì)比。根據(jù)對(duì)比結(jié)果，有如下結(jié)論：

1)由表2～表6可以看出，在完整考慮4種主要軌道攝動(dòng)和微小推力作用下，本文所設(shè)計(jì)的C++程序和STK模擬結(jié)果具有較高的一致性，24 h軌道遞推的位置誤差均在22 m以下，速度誤差小于0.012 m/s，幾乎可以忽略不計(jì)。

2)由表2和表3對(duì)比可以看出，在僅考慮微小推力情況時(shí)(不考慮其他軌道攝動(dòng))，本文所設(shè)計(jì)的C++程序和STK模擬的結(jié)果具有高度一致性，兩種方法得到的24 h軌道遞推的位置差異小于0.3 m，速度差異小于0.001 m/s。而由表4和表2對(duì)比可以看出，當(dāng)考慮4種主要軌道攝動(dòng)時(shí)，兩種方法得到的24 h軌道遞推的位置差異增加至30 m量級(jí)，速度差異增加至0.01 m/s量級(jí)。

表2 無(wú)攝動(dòng)+無(wú)推力結(jié)果對(duì)比

表3 無(wú)攝動(dòng)+12 mN推力結(jié)果對(duì)比

表4 考慮4種主要軌道攝動(dòng)+無(wú)推力結(jié)果對(duì)比

3)由表3與表6對(duì)比可以看出，地球形狀攝動(dòng)力是引起本文所設(shè)計(jì)的C++程序和STK模擬結(jié)果差異的最大因素，位置差異為22 m左右，主要原因可能是采用的地球引力模型之間存在輕微差異，其他三種主要攝動(dòng)力并不會(huì)引起另種方法之間的較大位置差異。

表5 考慮4種主要軌道攝動(dòng)+12 mN推力結(jié)果對(duì)比

表6 僅考慮地球扁率攝動(dòng)+12 mN推力結(jié)果對(duì)比

總之，盡管本文所設(shè)計(jì)的C++程序?qū)Φ厍蛐螤顢z動(dòng)力的處理與STK內(nèi)置的地球形狀攝動(dòng)力有所差異，但本文所設(shè)計(jì)的C++程序仍可實(shí)現(xiàn)在微小推力作用下的高精度軌道預(yù)測(cè)。

為了驗(yàn)證其他影響因素對(duì)所設(shè)計(jì)的C++程序的影響，本文分別針對(duì)軌道高度、推力大小和推演時(shí)間對(duì)程序的結(jié)果進(jìn)行分析。

圖3和圖4分別給出了不同軌道高度下(200 km～1 500 km)C++程序和STK模擬之間的位置誤差和速度誤差，其中除軌道高度外，其余初始參數(shù)均與本節(jié)第一部分一致，程序考慮12 mN推力且無(wú)攝動(dòng)下的結(jié)果。通過(guò)結(jié)果可以看出，不同軌道高度下，程序的位置誤差均在0.33 m以下，速度誤差均在4×10-4m/s以下，這說(shuō)明本文的程序在不同軌道高度下均可以保證高精度軌道遞推。

圖3 不同軌道高度下位置誤差Fig.3 Position errors at different orbital altitudes

圖4 不同軌道高度下速度誤差Fig.4 Velocity errors at different orbital altitudes

圖5和圖6分別給出了不同推力大小下(0 mN-20 mN)C++程序和STK模擬之間的位置誤差和速度誤差，其中除推力大小外，其余初始參數(shù)均與本節(jié)第一部分一致，程序考慮存在推力且無(wú)攝動(dòng)下的結(jié)果。通過(guò)結(jié)果可以看出，不同推力大小下，程序的位置誤差均在0.6 m以下，速度誤差均在6×10-4m/s以下，這說(shuō)明本文的程序在不同推力大小下均可以保證高精度軌道遞推。

圖5 不同推力大小下位置誤差Fig.5 Position errors with different thrusts

圖6 不同推力大小下速度誤差Fig.6 Velocity errors with different thrusts

圖7和圖8分別給出了不同推演時(shí)間下(1～14 d)C++程序和STK模擬之間的位置誤差和速度誤差，其中除推演時(shí)間外，其余初始參數(shù)均與本節(jié)第一部分一致，程序考慮存在4種攝動(dòng)和12 mN推力下的結(jié)果。通過(guò)結(jié)果可以看出，隨著推演時(shí)間的增加位置誤差和速度誤差都會(huì)變大，14 d內(nèi)的位置誤差在140 m以下，位置誤差在0.35 m/s以下，這說(shuō)明本文的程序隨著推演時(shí)間的增大，誤差會(huì)累積增大，但還能保證正常的精度要求。

圖7 不同推演時(shí)間下位置誤差Fig.7 Position errors with different propagation time

圖8 不同推演時(shí)間下速度誤差Fig.8 Velocity errors with different propagation time

5.2 在軌標(biāo)定仿真結(jié)果

根據(jù)上文的算法，利用C/C++對(duì)微小推力大小和其本體坐標(biāo)系下的方向在3 h內(nèi)進(jìn)行在軌標(biāo)定的仿真驗(yàn)證。

針對(duì)本文中的微小推力在軌標(biāo)定問(wèn)題，使用無(wú)跡卡爾曼濾波器(UKF)進(jìn)行衛(wèi)星軌道估計(jì)和推力標(biāo)定，通過(guò)預(yù)測(cè)和更新，完成對(duì)微小推力的標(biāo)定。衛(wèi)星軌道遞推的初始參數(shù)與5.1節(jié)第一部分初始參數(shù)一致。

利用所編寫(xiě)的C++程序在3 h時(shí)間間隔內(nèi)對(duì)衛(wèi)星微小推力進(jìn)行推力在軌標(biāo)定，初始位置各分量誤差為1 000 m，速度各分量誤差為0.05 m/s，微小推力各方向初始誤差為0.000 01 m/s2，測(cè)量精度為1 m，濾波器更新時(shí)間間隔為10 s。表7展示了一次軌道標(biāo)定結(jié)果，可以看出位置誤差、速度誤差和推力誤差均收斂，推力精度占比收斂至0.429 4%。

表7 推力在軌標(biāo)定3 h結(jié)果

如圖9～圖12以圖示的形式給出了UKF推力在軌標(biāo)定的結(jié)果，圖9展示了3 h內(nèi)的位置估計(jì)誤差，圖10為速度誤差，圖11為推力誤差，圖12為推力精度占比，由上圖可以看出，各誤差和占比均收斂。

圖9 推力在軌標(biāo)定的位置誤差Fig.9 Position errors of on-orbit thrust calibration

圖10 推力在軌標(biāo)定的速度誤差Fig.10 Velocity errors of on-orbit thrust calibration

圖11 推力在軌標(biāo)定的推力誤差Fig.11 Thrust errors of on-orbit thrust calibration

圖12 推力在軌標(biāo)定的推力精度占比Fig.12 Thrust accuracy ratio of on-orbit thrust calibration

此外，通過(guò)改變C++程序中的位置誤差、速度誤差和位置誤差，該程序仍能保證各誤差和占比均收斂。

6 結(jié)論

對(duì)于微小推力下的小衛(wèi)星軌道遞推和推力在軌標(biāo)定問(wèn)題，本文提出了基于變步長(zhǎng)龍格庫(kù)塔的軌道遞推算法，利用高精度的軌道動(dòng)力學(xué)模型，能夠?qū)π⌒l(wèi)星軌道進(jìn)行高精度的推演。在此基礎(chǔ)上，設(shè)計(jì)基于UKF推力標(biāo)定算法，對(duì)小衛(wèi)星施加的微小推力進(jìn)行在軌標(biāo)定。仿真實(shí)例中對(duì)軌道遞推算法與STK進(jìn)行比較，位置誤差均在22 m以下，速度誤差均在0.02 m/s以下；同時(shí)推力在軌標(biāo)定各誤差均收斂，且精度得到保證，仿真數(shù)據(jù)說(shuō)明了本文所設(shè)計(jì)算法的可行性。本文設(shè)計(jì)的算法基于C++程序，相較于其他程序算法具有更快的計(jì)算速率，有利于實(shí)際任務(wù)中實(shí)現(xiàn)星上的在線計(jì)算。在下一步的研究中，致力于消除地球形狀攝動(dòng)力對(duì)算法帶來(lái)的誤差。