張 冉,李小娟,韓 潮,張亞航,于俊慧
(1.北京空間飛行器總體設計部·北京·100094; 2.北京航空航天大學 宇航學院·北京·100191)
電推進系統以其比沖大、推力小的技術特點在深空探測、近地軌道阻力補償、軌道位置保持等領域得到應用[1]。近年來,電推進器在GEO軌道轉移中嶄露頭角。歐洲先進通信衛星ARTEMIS(2001)和美國先進極高頻軍事通信衛星AEHF-1(2010)原計劃以化學推進和電推進結合的混合方式進行GTO-GEO軌道的轉移,但因不同的故障導致化學推進失效,電推進系統成為GEO軌道轉移的主推力器,超時工作挽救了任務。波音公司基于702SP平臺發展了全電推進通信衛星,迅速投入商業市場;截至目前,已有五顆702SP平臺全電推進衛星發射入軌,單星質量在1950 kg~2300 kg之間,每顆衛星配備了4臺XIPS-25離子推力器,總推力大小為165 mN,GEO入軌時間為6個月;其中兩次發射采用獵鷹九號“一箭雙星”的發射方式,極大降低了發射成本。空中客車防務與航天基于E3000EOR平臺也發展了全電推通信衛星,單星質量3550 kg~5300 kg之間,每顆衛星配備一臺SPT140D穩態等離子推力器,推力大小約290 mN,GTO-GEO軌道轉移時間為4~6個月;目前已有三顆E3000EOR平臺的通信衛星發射入軌。全電推進通信衛星采用比沖更高的電推進系統進行GTO-GEO的軌道轉移,可以大幅減少燃料質量,有效降低衛星總重或提高有效載荷的質量比,延長任務的壽命,降低發射成本[2-4]。因此,在GTO-GEO軌道轉移中采用電推進技術已經成為中小型通信衛星發展的趨勢,研究電推進GEO軌道轉移優化和制導方法具有重要的工程意義。
全電推進GEO衛星的軌道轉移設計主要分為軌跡優化和制導控制兩個方面。求解航天器軌跡優化最優控制問題的數值方法分為直接法、間接法和結合二者的混合法[5-6]。直接法通過離散變量和參數將最優控制問題轉換為參數優化問題,然后通過非線性規劃算法獲得最優解,直接法的計算量比較大。間接法通過變分法和龐特里亞金極大值原理構造兩點邊值問題(Two-Point-Boundary-Value-Problem,TPBVP),然后通過打靶法求解靶函數;間接法的困難在于靶函數的收斂域很小,靶函數對于初值特別敏感。一些學者提出同倫法[7-9]和非線性濾波參數估計算法[10]等解決初值敏感的數值問題。混合法則結合了直接法和間接法各自的優點,通常假定推力結構,由間接法的必要條件保證控制的最優性。在對控制結構適當假設和簡化的前提下,一些學者獲得了軌跡優化的解析求解方法。Kechichian[11]提出了保持偏心率為零的前提下使半長軸和軌道傾角最快改變的分段常值控制律;田百義等[12]在小推力GEO入軌問題中假設平面內外控制解耦,平面內保持遠地點高度不變,平面外在升降交點采取控制的解析策略。高揚[13]提出三種分別控制軌道半長軸、偏心率和軌道傾角的切向、慣性和偏航控制律。解析法的控制律只能達到次優的控制結果,與最優轉移相比有一定差別。在小推力軌道轉移制導方法研究方面,許多文獻提出了跟蹤最優軌跡的控制思路[14-15],通過設計一組控制參數,修正實際軌道與參考軌道之間的偏差;另外,基于李雅普諾夫函數理論的制導律也得到利用,主要問題是李雅普諾夫函數的定義方式和增益的計算[16-17]。
本文主要研究全電推進GEO衛星的入軌轉移制導策略問題,提出了一種跟蹤最優參考軌跡的分段常值制導策略,該策略形式簡單,易于星上執行,并且控制結果與最優控制律相差不大,可為我國全電推進通信衛星平臺的研究提供參考。
本文采用經典軌道要素的高斯攝動方程作為軌道轉移優化的控制模型,形式如下:
(1)
式中,μ為地球引力常數,a為軌道半長軸,e為軌道偏心率,i為軌道傾角,Ω為升交點赤經,ω為近地點幅角,θ為真近點角;p為軌道半通徑,h為軌道角動量大小,r為衛星地心距,v為衛星速度大小;控制加速度a在衛星速度方向、垂直速度方向和動量矩方向的三個分量分別為ft,fn,fh。a是除二體加速度外的其他外力加速度矢量,包含發動機推力和攝動加速度矢量兩部分;本文攝動加速度僅考慮了地球非球形引力場J2項的影響。
設衛星質量為m,全電推進發動機具有上限推力大小Tmax和常值比沖Isp,g0為赤道重力加速度,式(1)改寫為矩陣的形式:
(2)

對于全電推進GEO衛星的燃料最優軌道轉移問題,性能指標為
(3)
轉移時間tf是相應時間最短問題轉移時間tfmin的倍數,用系數Ctf表示,Ctf>1。
利用變分法和極大值原理,小推力燃料最優GEO軌道轉移問題可以表述為下面的兩點邊值問題:
(4)
其中,λ是與軌道要素對應的協態變量,λm是與質量對應的協態變量,λθ,λΩ,λω分別是與真近點角、升交點赤經、近地點幅角對應的協態變量,Ctf是飛行時間系數;下標“0”表示初始時刻的狀態,下標“f”表示變軌終止時刻的狀態。
上述TPBVP問題需要利用打靶法求解協態變量的初值λ(t0)、λm(t0)。對于電推進GEO軌道的轉移問題,軌道轉移過程包含幾百圈軌道,轉移時間長達半年左右,靶函數對于初值非常敏感。為了解決數值求解困難,本文采用了文獻[10]提出的UKF非線性濾波參數估計算法。
本文仿真中采取的初始軌道是GTO轉移軌道,軌道要素如表1所示。飛行時間系數Ctf設為1.5,使用UKF濾波參數估計算法求解不同推力值Tmax情況下的燃料最優GEO軌道轉移問題。為避免結果陷入局部極值,本文采取逐步縮小推力幅值Tmax的延拓思路,每一步計算時以上一步計算結果為初值,并根據控制結構篩查全局最優的極值。

表1 軌道轉移初始與目標軌道要素
圖1(a)顯示了燃料最優GEO軌道轉移問題不同推力幅值Tmax條件下軌道半長軸、偏心率和軌道傾角隨時間的變化歷程,為便于比較,橫軸的時間做了去單位化處理(t=t/tf)。Tmax較大時,軌道要素的變化呈現出“階梯”形,這是由Bang-bang控制結構造成的。隨著Tmax的減小,軌道要素的變化變得平滑,局部放大后仍存在“階梯”特性。不同Tmax條件下,軌道要素呈現出相似的變化趨勢,可以考慮用一條平滑曲線去擬合軌道要素的變化規律。
第1節獲得了燃料最優GEO轉移的最優軌跡和最優控制曲線,這些軌道和控制信息可以作為實際工程中飛行制導的參考,然而直接存儲上述軌道與控制信息有以下幾個缺點:①軌道圈數多,轉移時間長,存儲數據量大;②最優控制律的表達式與協態變量相關,而協態變量沒有其他物理意義;③狀態變量和控制變量與時間關聯,攝動力和控制偏差會造成相應狀態在指定時間達到預定數值,使存儲信息失效;④軌道要素隨時間的小幅度振蕩影響制導律的設計。
為解決以上問題,本文利用擬合方式獲得簡單解析形式的曲線,更適合在星上計算機存儲和利用。根據最優軌跡的形狀,本文選擇簡單的多次多項式對最優軌跡進行擬合,數值試驗表明,采用四次及以上的多項式對最優軌跡的擬合效果比較好。
圖1(b)為使用四次多項式對推力幅值Tmax為2.5 N、1.0 N、0.5 N和0.2 N條件下的燃料最優GEO軌道轉移軌跡進行擬合時獲得的參考軌道,四條擬合曲線非常接近。Ctf固定時,燃料最優軌道轉移軌跡的軌道要素變化趨勢與推力大小關系不大,在設計小推力燃料最優轉移問題的參考軌道時,僅需計算較大推力條件下的軌道轉移工況,通過曲線擬合方式可以快速得到小推力條件下的參考軌道。

(a)最優軌跡

(b)擬合軌跡
第一節的最優軌跡對應的最優控制是一組開環的控制,這些數據中控制力的大小和方向等數據量非常大,不適合在星上計算機存儲;控制力的方向與當前協態變量和狀態變量以及時間相關,協態變量本身沒有具體的物理含義,這給制導律的設計帶來了很多不便。在衛星的實際變軌過程中,不確定的擾動因素、未建模的攝動力、導航與控制偏差會導致衛星偏離參考軌道,必須設計含閉環反饋的制導控制律。在獲得解析化的參考軌道之后,多圈長時的軌道轉移問題可以分解為多個單圈轉移設計問題,根據參考軌道的形式,將整個軌道改變量以一種近乎優化的方式分配到每一個軌道周期內,本節基于此基礎設計單圈內的軌道制導律。
電推進控制加速度[ft,fn,fh]T與引力加速度相比是小量,使用偏近點角E代替真近點角θ作為第六個軌道要素,dE/dt可以近似為:
(5)
式中,n為軌道角速度。將式(5)代入式(1),可得到以偏近點角為獨立變量的微分方程:
(6)
在推力加速度[ft,fn,fh]T作用下,衛星軌道偏近點角由E1到E2的過程中軌道要素的變化量具有解析表達式,軌道要素X的改變量具有以下形式:

(7)
根據式(7),根據不定積分的公式,可以獲得TX,NX,HX的表達式,具體形式可以參考文獻[18]。
衛星在偏近點角由E1到E2的過程中滿推力工作時,質量的改變量為:
(8)
圖2顯示了燃料最優GEO軌道轉移的最優控制的曲線和設計的常值推力控制結構,其中ut,un,uh是推力加速度單位矢量在速度方向、垂直速度方向和動量矩方向上的三個分量。細線表示的是最優控制曲線,黑粗實線是根據最優控制曲線形狀設計的分段常值控制結構。

(a)軌道控制前中期

(b)軌道控制末期圖2 最優和分段常值推力控制結構Fig.2 Optimal and piecewise constant thrust control structure
在軌道轉移的前中期,平面外有兩個推力弧段:[E1,E2]和[E2,E3]段,平面外推力方向分別用β1和β2描述;平面內的推力方向角用α1和α2描述。常值推力結構見表2所示,在軌道轉移的前中期,每個軌道周期內,在三段常值推力的作用下,軌道要素的變化量為:

(9)
軌道轉移的末期,平面外有四個推力弧段,弧段用五個設計變量描述:E1,E2,E3,E4,E5;平面內根據E3和π的大小關系,分為五個推力弧段;為減少變量數目,分別用四個方向角α1,α2,β1,β2描述推力矢量的方向,軌道轉移末期的常值推力結構如表3所示。在軌道轉移末期,在五段常值推力作用下,軌道要素的變化量為:
ΔX=
(10)

表3 軌道末期的常值推力弧段
在考慮J2攝動時,單圈內的控制結構如圖3所示。在軌道轉移的前期和中期,控制序列包含三個連續的推力弧段和一個漂移弧段,設計變量為
σ=[E1,E2,E3,α1,α2,β1,β2]
(11)
性能指標為:
J=(E3-esinE3)-(E1-esinE1)
(12)
約束條件為:

0≤α1,α2,β1≤π/2,-π/2≤β2≤π/2
0 (13) (a)軌道控制前中期 (b)軌道控制末期圖3 單圈內的控制序列示意圖Fig.3 Schematic diagram of control sequence in a single loop 在軌道轉移的末期,單圈內控制序列包括五個推力弧段和兩個漂移弧段,設計變量為: σ=[E1,E2,E3,E4,E5,α1,α2,β1,β2] (14) 性能指標為: J=(E4-esinE4)-(E2-esinE2)+(E1-esinE1)-(E5-esinE5) (15) 約束條件為: (16) 其中ΔX為根據式(10)所計算的軌道要素改變量。 以上完成了控制單圈內數學規劃問題的數學建模,這是一個具有簡單邊界約束的優化問題,該類問題有較成熟的優化算法;控制律充分利用了最優控制的變化規律,給出了較好的初值;每一圈的優化結果可以作為下一圈的設計變量初值,提高了優化求解速率。 選取表1所示的初末軌道要素,最大推力值為Tmax=0.2N,在最優軌跡計算時考慮J2項攝動,最優軌跡的飛行時間為265.53天,燃料消耗123.84 kg;采用四次多項式擬合生成參考軌道。圖4顯示了本文提出的分段常值推力跟蹤制導策略的控制結果,半長軸、偏心率和軌道傾角的跟蹤效果非常好。在考慮J2項攝動的情況下,變軌結束后軌道半長軸控制誤差為0.83 km,偏心率誤差為0.00145,軌道傾角控制誤差為0.0025°。由近地點幅角的變化規律可以看出,跟蹤控制并不會抑制近地點幅角ω的漂移。分段常值跟蹤控制的飛行時間為266.96天,燃料消耗132.72 kg。本文提出的分段常值跟蹤制導方法與最優軌跡相比,飛行時間延長0.54%,燃料消耗增加了7.17%。燃料消耗的主要原因是簡化的分段常值的推力結構控制效率低于連續變化的最優控制結構。如果能根據最優控制結構的形狀,找出一組簡化的數學解析式來模擬推力值與偏近點角的關系,燃料消耗能夠進一步優化。 圖4 GTO-GEO轉移分段常值制導軌跡、擬合軌跡和最優軌跡Fig.4 GTO-GEO transfer piecewise constant guidance trajectory, fitting trajectory and optimal trajectory 本文針對全電推進GEO衛星的入軌轉移制導策略進行了研究分析,提出了一種跟蹤參考軌道的制導策略,并結合算例分析了制導策略的效果。根據最優軌跡的變化規律,提出采用多次多項式擬合最優軌跡生成參考軌道的方法,參考軌道形式簡單,能夠表征最優轉移軌道的變化特點,易于在星上計算機存儲和使用。跟蹤參考軌道的閉環反饋控制策略將多圈軌道轉移的控制量分配到每一個軌道單圈內,保證了制導的燃料近優性。在軌道單圈內,本文提出一種分段常值推力控制結構,這種控制結構形式簡單,在此假設下軌道要素改變量有解析解,簡化了制導策略設計;并且常值推力結構易于工程實施。仿真結果顯示,在考慮J2攝動條件下,分段常值的制導策略比最優控制解多消耗7.17%的燃料,對參考軌道的跟蹤效果好,控制精度高。上述設計可為我國全電推進GEO通信衛星的入軌控制提供參考。為進一步提高跟蹤制導策略的效率,后續可研究一類模擬最優控制結構的解析連續變化的曲線,并進一步考慮更多的工程約束條件。




4 結 論