劉俊堯，賈建芳，李文強，宋錦武

(1. 中北大學 電氣與控制工程學院，山西 太原 030051； 2. 中北大學 信息與通信工程學院，山西 太原 030051；3. 廣州廣電計量檢測股份有限公司，廣東 廣州 510656)

0 引 言

航天器軌道轉(zhuǎn)移優(yōu)化技術(shù)是空間技術(shù)的重要組成部分. 目前對于軌道轉(zhuǎn)移優(yōu)化的研究大多基于脈沖變軌[1]. 在真實工作狀態(tài)下進行軌道轉(zhuǎn)移時，由于發(fā)動機推力制約，脈沖式變軌并不符合軌道轉(zhuǎn)移的實際情況[2]，特別是發(fā)動機推力較小的航天器在遠距離軌道轉(zhuǎn)移時，所需時間較長，存在很長的推力弧段，此時脈沖式變軌假設不再成立[3]. 與脈沖式軌道轉(zhuǎn)移技術(shù)相比，基于有限推力的軌道轉(zhuǎn)移技術(shù)能夠真實反映航天器的工作狀態(tài)，在提高研究結(jié)果真實性、 精確性方面有很大的優(yōu)越性.

基于有限推力的航天器軌道轉(zhuǎn)移優(yōu)化問題的本質(zhì)是泛函空間最優(yōu)控制問題[4]，根據(jù)優(yōu)化過程中對動力學模型有限維近似處理方式不同，大致分為間接法、 直接法和混合法三種. 間接法[5-7]的本質(zhì)是以極小值原理以及非線性規(guī)劃作為理論基礎，通過相應的邊界條件將軌道轉(zhuǎn)移優(yōu)化問題轉(zhuǎn)化為兩點邊值求解問題，雖能得到很光滑的解, 但由于對初始猜測異常要求很高，求解時困難較大； 直接法[8-10]的一般思路是對待優(yōu)化軌道進行離散處理，利用優(yōu)化算法在每一個軌道段內(nèi)求解若干個待優(yōu)化參數(shù)，轉(zhuǎn)化后待優(yōu)化參數(shù)較多，需要在很大范圍內(nèi)進行參數(shù)搜索，計算性能隨著節(jié)點數(shù)增加會有較大程度下降，算法魯棒性較差； 混合法[11-13]研究以優(yōu)化算法的實現(xiàn)為主要方向，對問題本身具有很大依賴性，求解時需要針對具體問題對算法做大量修改.

智能優(yōu)化算法[14-17]是一種啟發(fā)式優(yōu)化算法，一般是針對具體問題設計相關(guān)的算法，智能優(yōu)化算法技術(shù)性、 應用性強，計算速度快且理論要求不高，得到國內(nèi)外學者廣泛關(guān)注并且發(fā)展迅速.

針對有限推力下遠程導引變軌問題，本文采用改進遺傳算法求解遠程軌道能量最優(yōu)的多脈沖變軌策略，將求解得到的多脈沖變軌策略在變軌點轉(zhuǎn)換為推力弧段上的有限推力，通過發(fā)動機推力計算出航天器加速度，對加速度進行積分求得航天器瞬時速度增量，將速度增量不斷累加完成變軌，實現(xiàn)有限推力下遠程導引變軌，并對有限推力下遠程導引變軌策略速度增量方向進行修正，修正后變軌策略精度有較大程度提高，與目標軌道位置誤差大幅減小，能夠反映航天器在真實工作狀態(tài)下的軌道轉(zhuǎn)移情況.

1 遠程導引變軌原理

1.1 特殊點變軌

根據(jù)特殊點變軌特點，變軌方案如圖1 所示.

圖1 特殊點變軌方法Fig.1 Special point change method

特殊點變軌大致分為四步：

1) 調(diào)整近地點高度： 初始軌道遠地點處施加跡向脈沖量Δv1；

2) 調(diào)整軌道傾角： 緯度幅角為0或π處施加法向脈沖量Δv2，將異面軌道調(diào)整為共面軌道；

3) 調(diào)整遠地點高度： 轉(zhuǎn)移軌道近地點處施加跡向脈沖量Δv3；

4) 調(diào)整軌道偏心率： 在目標軌道遠地點處施加跡向脈沖量Δv4，實現(xiàn)橢圓軌道向圓軌道轉(zhuǎn)換.

1.2 脈沖量計算

二體條件下，基于軌道動力學方程，根據(jù)變軌前的軌道參數(shù)和變軌后軌道參數(shù)，可以求得每次變軌所需要脈沖量，每次施加脈沖量計算如下：

第1次變軌

(1)

第2次變軌

(2)

(3)

第3次變軌

(4)

第4次變軌

(5)

式中：n表示軌道角速度；e表示軌道偏心率；hp表示近地點高度；ha表示遠地點高度； Δi表示軌道傾角修正量； ΔΩ表示升交點赤經(jīng)修正量；r表示航天器地心距；a表示軌道半長軸.

1.3 脈沖量轉(zhuǎn)換為有限推力方法

根據(jù)初始軌道參數(shù)和目標軌道參數(shù)求解得到脈沖量，將脈沖量根據(jù)航天器發(fā)動機推力、 比沖在發(fā)動機工作時間段內(nèi)轉(zhuǎn)換為有限推力，將長時間變軌過程轉(zhuǎn)換為短時間內(nèi)的連續(xù)變軌，轉(zhuǎn)換方法如圖2 所示.

圖2 脈沖量向有限推力轉(zhuǎn)換Fig.2 Pulse thrust to finite thrust conversion

圖2 中，ti表示連續(xù)變軌過程中的時間節(jié)點，Δvi表示速度增量.

在軌道轉(zhuǎn)換過程中，需要確定發(fā)動機開始工作時間、 終止工作時間以及發(fā)動機產(chǎn)生的加速度.

發(fā)動機工作時間為

(6)

式中：tb為發(fā)動機開始工作時刻；tf為發(fā)動機工作終止時刻；M為發(fā)動機消耗燃料的質(zhì)量；F為發(fā)動機提供的推力. 發(fā)動機燃料消耗量M由航天器質(zhì)量、 發(fā)動機比沖以及速度增量決定，計算公式為

(7)

式中：m為航天器質(zhì)量； Δv為速度增量；isp為發(fā)動機比沖.

1.4 有限推力轉(zhuǎn)換修正

在脈沖變軌向有限推力變軌轉(zhuǎn)換過程中，在軌道轉(zhuǎn)移弧段對發(fā)動機推力產(chǎn)生的航天器加速度進行積分，求得節(jié)點處的速度增量. 在qsw本地軌道坐標系中，X′軸是從地球中心指向衛(wèi)星的單位矢量；Z′軸是垂直衛(wèi)星軌道平面且方向向上的單位矢量；Y′軸是由(x,y,z)形成規(guī)則確定的單位矢量，如圖3 所示.

圖3 本地軌道坐標系Fig.3 Local orbital coordinate system

圖3 中，O-XYZ是ECI坐標系，O′-X′Y′Z′是qsw本地軌道坐標系. 在qsw本地軌道坐標系中，軌道為圓軌時，Y′軸方向為軌道切線方向，與航天器速度方向一致； 軌道為近圓軌道時，Y′軸方向與橢圓軌道切線方向存在少量偏差，與航天器速度方向基本保持一致.

在脈沖變軌向有限推力變軌轉(zhuǎn)換過程中，將qsw本地軌道坐標系Y′軸方向近似為速度增量方向，通過確定節(jié)點處航天器速度矢量，計算速度增量方向.

2 遺傳算法

遺傳算法是一種高效、 并行的全局搜索方法. 它的求解過程是基于概率搜索技術(shù)，以目標函數(shù)值作為搜索目標，同時在多個搜索點對目標進行搜索. 標準遺傳算法全局最優(yōu)搜索能力強，但是在實踐應用中，也存在局部搜索能力差的問題，無法保證算法完全收斂. 針對標準遺傳算法存在的問題，對遺傳算法的選擇和交叉操作進行改進，在對群體根據(jù)適應度值高低進行排序的基礎上，用個體最優(yōu)與其他偶數(shù)位的所有個體進行交叉，每次交叉產(chǎn)生兩個新的個體，在交叉過后，對新產(chǎn)生的群體進行多點變異產(chǎn)生子群體，再計算其適應度值，然后與父群體合并，并且根據(jù)適應度值進行排序，取前NP個個體為新群體，進行下一次遺傳操作. 改進遺傳算法運算流程如圖4 所示.

二體條件下，式(1)中的n1、e1、hp由初始軌道根數(shù)確定； 式(3)中的n3、e3由初始軌道根數(shù)及Δv1確定，h′a為目標軌道遠地點高度，ha為初始軌道遠地點高度； 式(5)中的n4、a4、e4由目標軌道根數(shù)確定. 選取第1次變軌時轉(zhuǎn)移軌道近地點高度h′p為唯一優(yōu)化變量，參數(shù)設置如下： 基因數(shù)目D=1； 染色體數(shù)目NP=100； 最大遺傳代數(shù)G=1 000； 交叉概率Pc=0.8； 變異概率Pm=0.1. 構(gòu)建目標函數(shù)

fmin=Δv=∑(|Δv1|+|Δv2|+|Δv3|).

(8)

圖4 改進遺傳算法流程圖Fig.4 Improved genetic algorithm flow chart

改進遺傳算法求解h′p最優(yōu)解時，先計算Δv1、 Δv2、 Δv3的適應度并按升序排列，選擇適應度最高的Δv1、 Δv2、 Δv3染色體個體與其他偶數(shù)位的所有染色體個體進行交叉并進行多點變異產(chǎn)生子代染色體Δv′1、 Δv′2、 Δv′3，計算子代染色體適應度并按升序排列，將Δv′1、 Δv′2、 Δv′3染色體與Δv1、 Δv2、 Δv3染色體合并，計算合并后的適應度并按升序排列，選取前NP個適應度最高值作為此代最優(yōu)個體.

3 仿真算例與分析

3.1 軌道根數(shù)設置

初始軌道根數(shù)和目標軌道根數(shù)設置如表1. 表1 中，a表示軌道半長軸；e表示偏心率；i表示軌道傾角；Ω表示升交點赤經(jīng)；ω表示近地點幅角；anv表示真近點角.

由于初始軌道與目標軌道共面，無需改變軌道傾角，僅需要施加3次脈沖量，實現(xiàn)遠程導引變軌.

表1 軌道根數(shù)參數(shù)

3.2 仿真結(jié)果與分析

以第一次變軌時轉(zhuǎn)移軌道近地點高度h′p為優(yōu)化變量，采用改進遺傳算法求解總速度增量的最小值. 結(jié)果顯示，當h′p為6 982.039 km時，總速度增量最小值為66.699 m/s，其中第1次變軌的速度增量Δv1為28.855 m/s，第2次變軌的速度增量Δv2為33.122 m/s，第3次變軌的速度增量Δv3為 4.722 m/s，航天器遠程導引變軌序列如表2 所示.

表2 脈沖變軌方案

將初始得到的脈沖變軌方案在變軌點處轉(zhuǎn)換為推力弧段上的有限推力，設置航天器初始質(zhì)量為1 000 kg，發(fā)動機比沖為220 s，分別施加50 N, 100 N, 200 N的推力，轉(zhuǎn)換結(jié)果如表3 所示.

表3 50 N推力下變軌方案

表4 100 N推力下變軌方案

由表3～表5 可知，有限推力變軌方案消耗的燃料會多于脈沖變軌方案； 推力越小，軌道上推力弧段越長，變軌結(jié)果誤差越大. 由表5可知，推力為200 N時，航天器在經(jīng)過3次變軌后，軌道半長軸為6 998.295 km，與目標軌道半長軸誤差為1.715 km，相對位置誤差為3.89 km. 轉(zhuǎn)換后的連續(xù)推力變軌方案速度不斷迭加，在不同節(jié)點處不斷變化，但速度增量方向一直保持脈沖變軌處方向不變，對轉(zhuǎn)換結(jié)果造成較大影響，并且變軌持續(xù)時間，誤差越大. 對200 N推力下變軌方案進行修正，結(jié)果如表6 所示.

表5 200 N推力下變軌方案

表6 有限推力變軌方案修正結(jié)果

由表6 可知，有限推力變軌轉(zhuǎn)換修正方案精度更高，與目標軌道半長軸誤差為0.957 km，相對位值誤差為1.33 km. 經(jīng)過3次變軌，軌道根數(shù)與目標軌道根數(shù)對比如表7 所示.

表7 軌道根數(shù)對比

表7 中，真近點角差異是由于推力弧段的存在，導致變軌結(jié)束時刻有差異. 由表7 可知，脈沖變軌方案轉(zhuǎn)換為有限推力方案，僅有軌道半長軸、 偏心率、 真近點角存在少量差異，軌道傾角、 升交點赤經(jīng)、 近地點幅角保持一致； 有限推力修正方案精度更高，可以反映航天器在真實飛行環(huán)境中的變軌情況. 航天器有限推力方案變軌仿真結(jié)果如圖5 所示.

由圖5 可知，變軌后軌道高度提升，軌道仍為圓軌； 由圖6 可知，航天器軌道半長軸在機動時間段內(nèi)不斷增大，非機動時間段內(nèi)基本保持不變，經(jīng)過3次變軌，軌道半長軸與目標軌道半長軸高度一致.

圖5 有限推力方案3D仿真圖Fig.5 Finite thrust scheme 3D simulation diagram

與脈沖變軌相比，有限推力變軌所需燃料更多，是由于有限推力變軌弧段存在引力損失； 不進行推力方向修正的變軌方案的誤差要高于進行推力方向修正的變軌方案，是因為通過修正推力方向可找到時間節(jié)點處推力的方位角和俯仰角.

4 結(jié) 論

本文對有限推力下遠程導引變軌進行了研究. 通過改進遺傳算法的選擇和交叉操作，求解得到燃料最優(yōu)的多脈沖變軌方案，在脈沖點處將脈沖量轉(zhuǎn)換為推力弧段上的有限推力，實現(xiàn)多脈沖機動遠程導引變軌向有限推力下遠程導引變軌的轉(zhuǎn)換，并對有限推力變軌方案速度增量方向進行修正，修正結(jié)果精度有較大程度提高，證明方法確實有效. 基于有限推力的遠程導引變軌方案可以反映航天器在真實工作狀態(tài)下環(huán)境中變軌的過程，并取得較高的精度.