采用虛擬引力場的太陽帆/電推進混合推力機動軌道設計

2018-09-07 10:30:00孫沖袁建平方群崔堯汪明曉劉易斯

西北工業大學學報 2018年4期

孫沖, 袁建平, 方群, 崔堯, 汪明曉, 劉易斯

(1．西北工業大學航天飛行動力學技術國家級重點實驗室，陜西西安710072;2．西北工業大學深圳研究生院，深圳518057;3．上海宇航系統工程研究所，上海201109; 4．國家航天局空間碎片監測與應用中心，北京100012)

連續小推力軌道設計在空間機動軌道設計方面有廣泛的應用[1-2]。與傳統脈沖推力或大推力相比，連續小推力軌道機動容易控制，機動能力強；另外其比沖小，實現軌道機動所需要的燃料很少，因而成為研究熱點。文獻[3]中，Mcinnes首次提出了使用太陽帆來改變引力場的大小，設計了懸浮軌道，并推導了實現懸浮軌道所需推力的解析解。Malcolm采用太陽帆產生切向固定推力，設計了從地球軌道到太陽極軌的轉移軌道，解決了軌道優化問題[4]。然而，由于單獨采用電推進或者太陽帆產生的推力較小，軌道轉移常常需要很長時間。針對這個問題，Mengali等在文獻[5]中提出了太陽帆和太陽帆/電推進結合的混合推進方法，并采用間接優化方法進行了軌道優化。在其文章中，太陽帆產生徑向推力，用于改變引力場的引力常數，而與此同時電推進產生固定切向推力。在假定太陽帆能夠產生50%的中心引力(μ2=0.5μ1)的條件下，混合推力能夠減少燃料的消耗并且縮短空間機動任務時間。然而該方法是作為一種間接優化方法，并沒有考慮推力約束問題，很難應用于實際空間機動軌道設計中。

本文針對連續小推力軌道設計問題，提出了一種基于虛擬引力場的混合連續小推力軌道設計方法[6-7]。該方法將虛擬引力場與平均算法結合，通過混合連續小推力形成虛擬引力場，進而采用平均算法求解在虛擬引力場下固定大小切向力作用下的機動軌道解析解。進一步，本文采用四階龍格庫塔積分方法分析解析解的精度，并提出了誤差因子參數，從而保證解析解的精度。本文采用虛擬引力場方法，將航天器機動軌道設計優化問題轉化成參數優化問題。該方法主要分為3步：①設計并優化滿足任務約束(包括邊界約束和推力約束)的虛擬引力場，以及固定切向推力的大小；②求解實現虛擬引力場的推力大小；③考慮太陽帆的動力學模型，求解采用太陽帆來實現虛擬引力場的控制律。

1 模型

1.1 虛擬中心引力場的定義

虛擬中心引力場是指引力場中心不在地心上，且其引力常數不等于地心引力常數的中心引力場。該引力場可用參數(r0,μ2)來定義,其中r0為虛擬引力場中心與地心相對位置矢量,稱為虛擬引力場位置參數;μ2為虛擬引力場引力常數[6-7]。

1.2 虛擬中心引力場下虛擬圓錐曲線軌道

根據虛擬中心引力場參數r0維數的不同,虛擬中心引力場可以分為2種,如圖1所示;

圖1 虛擬中心引力場示意圖

設定地心引力場引力常數為μ1,航天器在地心慣性坐標系下位置矢量為r1,假定在航天器上施加推力加速度為at,則有:

(1)

假設地心慣性坐標系為OXYZ,虛擬中心引力場坐標系為O′X′Y′Z。推力與地心引力場的合力形成的虛擬中心引力場參數為(r0,μ2),則在虛擬中心引力場坐標系下有:

(2)

(3)

由(2)式、(3)式可知,在虛擬中心引力場坐標系下,航天器運動方程可表示為如下形式:

(4)

實現虛擬中心引力場所需推力為

(5)

將(3)式帶入(5)式, 可得

(6)

由上述分析可知,通過調整推力,可形成一種虛擬中心引力場。在該虛擬中心引力場中,航天器運動軌跡為圓錐曲線。若該虛擬圓錐曲線軌道滿足軌道轉移約束,則可實現軌道機動。

1.3 固定切向力作用下軌道攝動方程

在切向推力作用下,平面軌道高斯攝動方程為

(7)

(8)

徑向推力以及周向推力和切向推力之間的關系為

(9)

而真近點角和偏近點角的關系為

(10)

因此,將公式(7)中自變量轉換為偏近點角,則軌道半長軸和偏心率的變化[9]為

(11)

式中，E為偏近點角,ft為連續切向推力。在給定軌道交會的初始點A和終點B的條件下,采用混合推力方法設計機動軌道就是求解引力場參數(μ2,rvg)以及切向推力ft。若選擇合適參數X=(μ2,rvg,fsep)時,在虛擬引力場中和切向推力的同時作用下,航天器在虛擬引力場可從A點轉移到目標點B點,從而滿足軌道交會的邊界約束。這就是采用基于虛擬引力場方法的混合連續小推力軌道設計的基本思路。

采用虛擬引力場的混合連續小推力軌道設計方法來設計軌道,實際上是設計在虛擬引力場中的常值切向推力轉移軌道。轉移軌道需要滿足兩點邊界值約束,以及推力約束。當初始軌道的偏心率很小,推力較小的條件下[9],在虛擬引力場中的切向推力軌道近似解析解為

(12)

當初始點和目標點已知時,上式中的參數E已知。邊界約束可以轉化為

(13)

由此可見,采用虛擬引力場下的連續推力軌道設計方法,可將所有滿足邊界約束的軌道用一組參數X=(μvg,rvg,fsep)表示,即轉移軌道參數化。這樣轉移軌道設計與優化問題就可以轉化成為參數設計與優化問題。

由上分析可知,采用虛擬引力場下的混合連續小推力軌道設計方法,需要給航天器施加推力,一部分推力用于形成虛擬引力場,另一部分推力用于形成虛擬引力場中的切向力。具體推力分配:將形成虛擬引力場所需要的連續推力與在切向力的合力分解為徑向力和切向力,然后采用太陽帆和電推進分別產生分解的徑向力和切向力。由公式(6)可知,實現虛擬引力場所需要的推力是有虛擬引力場參數決定的如下式所示,

Trvg=F(μvg,rvg,fsep)

1.4 太陽帆/電推進混合最近最優軌道設計

1.4.1 太陽帆推進模型

由公式(5)可知,實現虛擬引力場需要的力隨著軌道半徑的變化而變化。因而采用太陽帆來實現虛擬引力場,其推力大小和方向也都是實時變化的。在此采用太陽帆推進的理想模型[4]

Fsail=2PAcos2αn

(14)

式中，P=(S0/c)(r0/r)2,為太陽光壓率。

采用理想推進模型,太陽帆產生的推力表示如下

(15)

式中，a1,a2,a3均為常數。而太陽帆產生的力在軌道的徑向和切向的分量為

(16)

由此可知太陽帆推進可以實現部分徑向和周向的正推力,而剩余部分需要用電推進來實現。

1.4.2 太陽帆電推進模型

本文選擇文獻[5]中的電推進器,其產生的小推力和電推進器功率參數相關,而需要的燃料消耗也由點推進器的功率決定。點推進器的功率是可變的,如公式(17)所示

(17)

則其產生的推力加速度在徑向周向分量為:

(18)

1.4.3 推力分配

采用虛擬引力場混合小推力軌道設計方法設計軌道,需要推力來形成虛擬引力場以及在虛擬引力場中的切向力。本文采用太陽帆/電推進混合推力來實現虛擬引力場,并提供在虛擬引力場中的常值推力。由于形成虛擬引力場所需的徑向力和周向力是固定的,而太陽帆產生的推力在徑向和周向是耦合的,因此僅僅采用太陽帆無法直接提供形成虛擬引力場需要的推力。在此,本文在選擇太陽帆的同時,另外選擇安裝2個電推進器,分別用于產生徑向和周向的推力,使其與太陽帆光壓推進產生的徑向力和切向力相互匹配,共同作用產生形成引力勢場,并提供常值推力。在此假定形成引力勢場所需要的徑向和周向的推力加速度為Tr,Tc,而太陽帆產生的徑向和周向的推力加速度為assr,assc。電推進產生的徑向和周向的推力加速度為asepr,asepc,則有如下公式

(19)

由太陽帆的數學模型可知,其產生的徑向力和軸向力與太陽帆的位置和太陽帆姿態有關。采用虛擬引力場和常值推力結合的方法,其太陽帆的在每一時刻的位置是一定的,因此太陽帆產生的推力,實際上是由太陽帆的姿態角來確定的。為了能夠最大程度利用太陽帆產生的推力,從而節省燃料,需要對其姿態角進行優化,在此假定太陽帆的姿態角α和軌道機動的轉移角θ有如下的關系:

α=a2·θ2+a1·θ+a0

(20)

式中，a0,a1,a2為待優化的系數項。

當給定虛擬引力場后,軌道轉移過程中需要電推進器產生推力的大小為:

(21)

對比本文與文獻[5]可知,本文的推力分配與其不同。文獻[5]中軌道機動的徑向推力完全由太陽帆提供,而切向力由電推進提供。本文采用虛擬引力場方法設計機動軌道,太陽帆推力一部分提供徑向力一部分提供切向力,這樣可以減少電推進所需提供的推力。因而本文的推力分配方法比較靈活,可節省軌道機動所需燃料。

1.4.4 混合小推力軌道設計與優化

由上述分析可知,采用基于虛擬引力場的混合連續小推力軌道設計方法,可求出大量的滿足邊界約束的轉移軌道,如圖3所示。假定虛擬引力場參數以及切向力參數的參數范圍為

圖2 混合小推力軌道交會

圖3 混合小推力轉移軌道優化示意圖

則軌道優化問題即可轉化為參數優化問題,優化的參數為虛擬引力場參數和切向推力加速度,以及太陽帆光壓推進的推力系數,而優化的目標函數可以為軌道轉移過程中,電推進器消耗燃料質量:

J=min(mfuel)

(22)

式中，mfuel為任務結束時航天器的質量。

或是轉移軌道的任務時間:

J=min(Tf)

(23)

或是同時考慮燃料消耗和飛行時間的權重函數

J=min(k1·mfuel+k2·Tf)

(24)

式中，k1,k2為權重,k1+k2=1。

結合上述分析可知,整體軌道設計流程包括以下幾個步驟:

a) 給出初始軌道和目標軌道;

b) 給出滿足推力約束的虛擬引力場參數范圍和切向推力范圍;

c) 在b)的范圍中,選擇虛擬引力場參數。

d) 將初始軌道狀態(ri,vi)轉化為虛擬引力

場中的軌道狀態(rvg,vvg);

e) 優化切向力參數,求解滿足邊界約束條件的轉移軌道;

f) 重復a)～e),得到大量的滿足約束的轉移軌道;

g) 采用Matlab軟件中的Fmincon函數優化轉移軌道。

2 誤差分析

由第1節分析可知,采用虛擬引力場下混合連續小推力軌道設計方法,可求解出滿足推力約束和邊界值約束的解析解。然而在此需強調,公式(10)中的解析解是存在誤差的。在此為了保證該解析方法的精度,本節采用四階龍格庫塔數值積分方法求解了軌道積分數值解,求出了解析解的誤差。在此基礎上,分析誤差產生的原因并給出解析解的應用范圍。具體的仿真條件見表1。

在此,將解析解的誤差分為位置誤差和速度誤差,并分別用Δεr,Δεv表示

Δεr=ranslytical-rnumerical=f(x)

Δεv=vanalytical-vnumerical=g(x)

對比圖4圖5,以及公式(10)、公式(11)可知,解析解與數值解的誤差和飛行時間、引力場參數、以及軌道參數有關。切向推力越大、飛行時間越長,該解析解的誤差會越大;初始軌道(虛擬引力場中)的偏心率越大、半長軸越大其誤差越大。這一結論也與文獻[1]得出的結論一致。為描述解析解與數值解的誤差,在此定義一個誤差因子參數,如(25)式所示

(25)

對比公式(11)及(12)可知,解析解的誤差因子會隨著切向推力、軌道半長軸、軌道轉移時間(偏近點角)的增大而增大,當(e<0.7)時,偏心率越大,誤差參數越大。如圖5所示,采用解析方法可以得到大量的轉移軌道,同時采用四階龍格庫塔數值方法,就可以求解并分析這些解析解的速度誤差、位置誤差以及誤差因子之間的關系。

以地球-火星轉移軌道為例,采用虛擬引力場方法可以得到大量的轉移軌道,在其中任意選擇n條轉移軌道(n=157)。采用數值積分方法,可以得出這些解析解的速度誤差和位置誤差與誤差因子的關系,如圖6所示。

圖6 速度、位置誤差隨誤差因子的變化

圖6中,藍色的線條為位置誤差;紅色的線條為速度誤差;綠色的線條為誤差因子。由圖6可知,速度誤差、位置誤差和誤差因子是同步的,當誤差因子很小時,其對應的速度誤差和位置誤差也很小;相反,當誤差因子比較大時,其速度誤差和位置誤差都很大。由圖6可知，速度誤差和位置誤差都隨著誤差因子的增大而增大。當保證誤差因子小于0.002 5時,速度誤差小于0.001%,而位置誤差小于0.002%。由以上分析可知,在采用虛擬引力場常值切向力解析方法設計轉移軌道時,只需要保證誤差因子k

3 數值仿真

表1 地球火星軌道轉移虛擬引力場參數和常值推力

表2 太陽帆姿態角系數

圖7 虛擬引力場混合推力方法的轉移軌道

圖8 切向推力方法的轉移軌道(相同時間)

圖7為采用虛擬引力場混合推力方法設計的轉移軌道；圖8為在相同轉移時間條件下，采用常值切向力方法設計的轉移軌道。圖9表述轉移軌道所需的周向推力和徑向推力以及太陽帆提供的周向推力和徑向推力。圖10是轉移軌道太陽帆姿態角。圖11是為相同情況下(切向推力大小相等)得到的轉移軌道。圖12為實現虛擬引力場所需要的推力加速度。圖13和圖14表示推力秒流量和推力方向角。由表5可知，相比于常值切向推力，虛擬引力場混合推力方法能夠很大程度上節省燃料，縮短時間。假定軌道轉移時間相等，虛擬引力場混合推力方法所需要的燃料僅為常值推力方法的。采用虛擬引力場混合連續小推力方法，所需要的電推進推力大大減小了，因而虛擬引力場混合推力方法能夠擴大小推力軌道設計方法的應用范圍。另一方面，在同樣大小推力條件下，采用常值推力方法，需要的軌道轉移時間和燃料消耗是相當大的。