基于改進(jìn)粒子群算法的月地轉(zhuǎn)移軌道優(yōu)化

李遠(yuǎn)平，蔡遠(yuǎn)利，

(西安交通大學(xué)電子與信息學(xué)部，陜西，西安 710049)

近年來(lái)，國(guó)際上再次掀起月球探測(cè)高潮，美國(guó)、俄羅斯、中國(guó)、歐空局和日本都制定了各自的載人登月和月球基地計(jì)劃[1]。我國(guó)原計(jì)劃于2019年底實(shí)施“嫦娥五號(hào)”月球采樣返回任務(wù)，軌道器與返回器組合體從月球停泊軌道加速后，進(jìn)入月地轉(zhuǎn)移軌道，在地面高度約5000 km處軌道器分離，返回器采用跳躍式再入，著陸場(chǎng)選擇為中緯度的內(nèi)蒙古四子王旗地區(qū)。月地轉(zhuǎn)移軌道優(yōu)化設(shè)計(jì)作為其中的技術(shù)難題之一，需要滿足月球停泊軌道、月地轉(zhuǎn)移和再入地球大氣層等方面的復(fù)雜約束條件，盡量減小返回窗口內(nèi)的逃逸速度增量，降低返回燃料消耗，同時(shí)獲得最優(yōu)的再入點(diǎn)參數(shù)，為跳躍式再入控制創(chuàng)造良好的起始條件。

月球返回軌道的設(shè)計(jì)方案可分為兩種：一是月球自由返回軌道設(shè)計(jì)[2]，采用多段拼接方法[3]，適用于短期的月球探測(cè)返回任務(wù)；二是月地快速轉(zhuǎn)移軌道設(shè)計(jì)，從月球停泊軌道加速后進(jìn)入月地轉(zhuǎn)移軌道[4]，適用于各類登月返回任務(wù)。文獻(xiàn)[5]在圓形限制性三體模型下研究了地月空間對(duì)稱自由返回軌道的設(shè)計(jì)方法及其軌道特性。文獻(xiàn)[6]分析了地月自由返回時(shí)近地距的各種影響因素，包括出發(fā)時(shí)刻、轉(zhuǎn)移時(shí)間、傾角組合等，提出了一種可使返回時(shí)近地距逐漸向設(shè)定的目標(biāo)值靠攏的單向搜索方法。文獻(xiàn)[7]分析了探測(cè)器低能量返回軌道的存在性與軌道動(dòng)力學(xué)特性，提出了一種混合自適應(yīng)遺傳算法，用于尋找月球探測(cè)器返回地球所需的最低能量以及對(duì)應(yīng)的轉(zhuǎn)移軌道。文獻(xiàn)[8―9]對(duì)從月球返回地球的軌道進(jìn)行了分析，給出了快速月地轉(zhuǎn)移軌道設(shè)計(jì)方法。文獻(xiàn)[10]采用了自起步解析算法，構(gòu)造出多次點(diǎn)火的軌道機(jī)動(dòng)序列，可實(shí)現(xiàn)從任意近月停泊軌道到地球再入狀態(tài)的轉(zhuǎn)移，但在某些情況下需要消耗較多的燃料。

對(duì)于復(fù)雜約束條件下的優(yōu)化設(shè)計(jì)問(wèn)題，采用群智能優(yōu)化算法具有明顯的優(yōu)勢(shì)：算法結(jié)構(gòu)簡(jiǎn)單、對(duì)初始值要求不高、不需要繁瑣的公式推導(dǎo)。常用的群智能優(yōu)化算法包括遺傳算法[11]、模擬退火算法[12]、粒子群算法[13]、蟻群算法[14]、蜂群算法[15]、螢火蟲算法[16]和離子運(yùn)動(dòng)算法[17]。

本文基于月地快速轉(zhuǎn)移軌道的思想，建立了多段拼接的解析法軌道計(jì)算模型，考慮了月地轉(zhuǎn)移軌道和再入彈道參數(shù)等諸多約束條件，以月球停泊軌道的逃逸速度增量最小為優(yōu)化目標(biāo)，提出了一種改進(jìn)粒子群算法，求解月心真赤道坐標(biāo)系內(nèi)的逃逸點(diǎn)軌道參數(shù)，搜索最優(yōu)的逃逸點(diǎn)時(shí)刻、升交點(diǎn)經(jīng)度和緯度幅角。

1 月地返回窗口選擇

月球停泊軌道通常采用圓軌道，軌道高度hMP和軌道傾角iMP(相對(duì)月球真赤道面)由月球軌道交會(huì)對(duì)接任務(wù)確定，月地返回窗口選擇主要是確定逃逸點(diǎn)時(shí)刻tMP的取值范圍，而逃逸速度增量ΔvMP、升交點(diǎn)經(jīng)度?MP和緯度幅角uMP由下一步的月地轉(zhuǎn)移軌道優(yōu)化過(guò)程來(lái)確定。

利用噴氣推進(jìn)實(shí)驗(yàn)室JPL的DE/LE405歷表來(lái)準(zhǔn)確計(jì)算月球軌道位置及天平動(dòng)參數(shù)，瞄準(zhǔn)著陸場(chǎng)中心位置，考慮月地轉(zhuǎn)移時(shí)間、地心段軌道傾角、返回當(dāng)?shù)貢r(shí)間和返回航程等約束條件，在一年或者更長(zhǎng)的時(shí)間段內(nèi)采用幾何解析法進(jìn)行覆蓋式搜索。

2 月地轉(zhuǎn)移軌道優(yōu)化方法

2.1 軌道計(jì)算模型

在全局優(yōu)化設(shè)計(jì)時(shí)，為了提高計(jì)算效率，采用雙二體引力模型來(lái)描述月地轉(zhuǎn)移軌道的運(yùn)動(dòng)規(guī)律。

1) 月球停泊軌道逃逸計(jì)算

先利用逃逸點(diǎn)時(shí)刻的月球停泊軌道參數(shù)，計(jì)算出探測(cè)器在月心真赤道坐標(biāo)系內(nèi)的逃逸點(diǎn)位置矢量rMP和速度矢量vMP，將其轉(zhuǎn)換到J2000.0月心天球慣性坐標(biāo)系，再考慮逃逸速度增量ΔvMP后，得到月地轉(zhuǎn)移軌道的起始點(diǎn)位置矢量rM0和速度矢量vM0。

2) 月心段轉(zhuǎn)移軌道計(jì)算

探測(cè)器從月球停泊軌道逃逸后，進(jìn)入一條月心雙曲線軌道。在月心天球慣性坐標(biāo)系內(nèi)計(jì)算探測(cè)器飛經(jīng)月球引力影響球出口點(diǎn)的歷元時(shí)刻tM1、位置矢量rM1和速度矢量vM1，將其轉(zhuǎn)換到J2000.0地心天球慣性坐標(biāo)系，得到相應(yīng)的位置矢量rE1和速度矢量vE1。

3) 地心段轉(zhuǎn)移軌道計(jì)算

探測(cè)器飛出月球引力影響球后，進(jìn)入一條地心橢圓軌道。在地心天球慣性坐標(biāo)系內(nèi)計(jì)算開普勒軌道參數(shù)并推算至再入點(diǎn)，得到月地轉(zhuǎn)移時(shí)間TET、地心軌道傾角iE、升交點(diǎn)赤經(jīng)?E、近地點(diǎn)幅角ωE和近地點(diǎn)高度hEP，以及再入點(diǎn)歷元時(shí)刻tR、當(dāng)天時(shí)間τR、位置矢量rR和速度矢量vR。

4) 再入彈道參數(shù)計(jì)算

計(jì)算再入相對(duì)速度VR、再入角θR、再入點(diǎn)地理經(jīng)度λR、再入點(diǎn)地理緯度φR、再入彈道面傾角iR和再入段航程LR。

2.2 優(yōu)化目標(biāo)函數(shù)

通常多目標(biāo)優(yōu)化問(wèn)題可描述如下：

式中：x為決策向量，x∈Rm；J為目標(biāo)函數(shù)向量，J∈Rn；fj(x)為目標(biāo)函數(shù)；gj(x)為約束條件。

對(duì)于月地轉(zhuǎn)移軌道優(yōu)化問(wèn)題，目標(biāo)函數(shù)可簡(jiǎn)化為：

式中，Δv1、Δv2分別為月球逃逸速度增量的下界和上界。

約束條件如下：

式中：hEP0為參考近地點(diǎn)高度；TET0為標(biāo)稱月地轉(zhuǎn)移時(shí)間；ΔTET為月地轉(zhuǎn)移時(shí)間偏差允許值；τR0為標(biāo)稱再入點(diǎn)當(dāng)天時(shí)間；ΔRτ為再入點(diǎn)當(dāng)天時(shí)間偏差允許值；iE0為標(biāo)稱地心軌道傾角；ΔiE為地心軌道傾角偏差允許值；ΔiR為再入彈道面傾角偏差允許值；ΩE0為標(biāo)稱地心軌道升交點(diǎn)赤經(jīng)；ΔEΩ為地心軌道升交點(diǎn)赤經(jīng)偏差允許值；ωE0為標(biāo)稱地心軌道近地點(diǎn)幅角；ΔEω為地心軌道近地點(diǎn)幅角偏差允許值；θR0為標(biāo)稱再入角；ΔRθ為再入角偏差允許值；LR0為標(biāo)稱再入航程；ΔLR為再入航程偏差允許值。

為了求解上述多約束條件下的優(yōu)化問(wèn)題，定義如下罰函數(shù)：

上述問(wèn)題轉(zhuǎn)化為：

式中：χ為收斂域條件標(biāo)識(shí)，0為不滿足，1為滿足。

2.3 優(yōu)化求解方法

由于搜索空間大，約束條件多，目標(biāo)函數(shù)復(fù)雜，本文提出了一種改進(jìn)粒子群算法。

2.3.1 標(biāo)準(zhǔn)粒子群算法

粒子群算法(Particle Swarm Optimization, PSO)的基本概念源于對(duì)于鳥群捕食行為的簡(jiǎn)化社會(huì)模型的模擬，在解空間追隨最優(yōu)的粒子進(jìn)行搜索[18]。

在PSO算法中，每個(gè)粒子在飛行過(guò)程所經(jīng)歷過(guò)的最好位置，就是粒子本身找到的最優(yōu)解，稱之為個(gè)體極值(pi)。整個(gè)粒子群所經(jīng)歷過(guò)的最好位置，就是整個(gè)群體找到的最優(yōu)解，稱之為全局極值(pg)。實(shí)際操作中，通過(guò)適應(yīng)度值(fitness value)來(lái)評(píng)價(jià)粒子的“好壞”程度。

設(shè)粒子的群體規(guī)模為N，第j次迭代時(shí)第i個(gè)粒子的位置表示為，所經(jīng)歷過(guò)的“最好”位置記為，其速度表示為，群體中所有粒子所經(jīng)歷過(guò)的“最好”位置表示為。在標(biāo)準(zhǔn)PSO算法中，第i個(gè)粒子的速度和位置更新規(guī)則如下：

式中：j為迭代次數(shù)；c1、c2均為學(xué)習(xí)因子，取值為常數(shù)；r1、r2均為[0，1]區(qū)間上的隨機(jī)數(shù)；w為慣性權(quán)重(Inertia Weight)。

設(shè)適應(yīng)度函數(shù)為F(X)，第i個(gè)粒子個(gè)體極值的更新公式為：

所有粒子的全局極值更新公式為：

標(biāo)準(zhǔn)PSO算法在應(yīng)用時(shí)，存在收斂速度慢和早熟的缺陷[19]。

2.3.2 改進(jìn)粒子群算法

為了克服標(biāo)準(zhǔn)PSO算法的缺點(diǎn)，針對(duì)上述月地轉(zhuǎn)移軌道優(yōu)化問(wèn)題，提出一種改進(jìn)粒子群算法(Modified Particle Swarm Optimization, MPSO)，利用適應(yīng)度函數(shù)F(X)來(lái)更新慣性權(quán)重w，并對(duì)粒子的速度加以約束。

對(duì)于月地轉(zhuǎn)移軌道優(yōu)化問(wèn)題，可將適應(yīng)度函數(shù)F(X)的最大化變換為目標(biāo)函數(shù)J的最小化。

對(duì)不同的粒子位置，采用不同的慣性權(quán)重

式中，JL為標(biāo)函數(shù)的上限值，取值為10。

隨著目標(biāo)函數(shù)J的最小化，慣性權(quán)重→0，加快了粒子位置的收斂速度。

定義粒子的位置向量X為：

其搜索空間如下：

式中：[t1,t2]為月球停泊軌道逃逸點(diǎn)時(shí)間窗口；[Δv1,Δv2]為月球停泊軌道逃逸速度增量區(qū)間，[850,860]m/s；[?1,?2]為月球停泊軌道升交點(diǎn)經(jīng)度區(qū)間，[0, 2π]；[u1,u2]為月球停泊軌道緯度幅角區(qū)間，[0,2π]。

粒子位置的初值為各維參數(shù)區(qū)間內(nèi)的隨機(jī)數(shù)組合，即：

式中，r為[0，1]區(qū)間上的隨機(jī)數(shù)。

粒子速度的初值為：

粒子速度的最大約束值為：

為進(jìn)一步提高M(jìn)PSO算法的收斂能力，利用月球停泊軌道逃逸點(diǎn)時(shí)刻tMP與再入點(diǎn)當(dāng)天時(shí)間τR的相關(guān)性，當(dāng)約束條件g2(x)≤0滿足時(shí)，對(duì)第一維位置參數(shù)引入下面的隨機(jī)反饋控制：

2.3.3 算法收斂性分析

利用式(24)進(jìn)行位置更新的粒子，最終會(huì)收斂到一個(gè)穩(wěn)定的位置[20]：

將式(35)代入式(24)中，可得：

式中，c=c1r1+c2r2。

設(shè)yj=pc-Xj，算法的速度和位置更新公式可變換為：

式(37)中，系數(shù)矩陣的特征方程為：

求解式(38)可得：

算法的收斂條件為max{|λ1|, |λ2|}＜1，代入式(39)可得：

算法的收斂區(qū)域見圖1。

圖1 改進(jìn)粒子群算法的收斂區(qū)域Fig.1 The convergence region of MPSO algorithm

對(duì)于MPSO算法，由于w→0，當(dāng)滿足c1+c2≤2時(shí)，可保證收斂性。

3 算例及分析

利用本文所提出的MPSO算法進(jìn)行月地轉(zhuǎn)移軌道優(yōu)化設(shè)計(jì)，所使用的計(jì)算機(jī)配置為2.53 GHz單核CPU、4G內(nèi)存，編程環(huán)境為MATLAB 7.11.0。

設(shè)定月球停泊軌道高度hMP為210 km、軌道傾角iMP為45°，著陸場(chǎng)中心點(diǎn)地理經(jīng)度λL=111.43°、地心緯度ΦL=42.16°，標(biāo)稱月地轉(zhuǎn)移時(shí)間TET0=3.0 d、偏差允許值ΔTET=0.5 d，標(biāo)稱地心軌道傾角iE0=45.0°、偏差允許值ΔiE=3.0°，標(biāo)稱升交點(diǎn)赤經(jīng)?E0=156.6°、偏差允許值ΔEΩ=10.0°，標(biāo)稱近地點(diǎn)幅角ωE0=10.8°、偏差允許值ΔEω=10.0°，參考近地點(diǎn)高度hEP0=45.0 km，標(biāo)稱再入點(diǎn)當(dāng)天時(shí)間τR0=22:00:00(UTC)、偏差允許值ΔRτ=0.5 h，標(biāo)稱再入角θR0=-6.0°、偏差允許值ΔRθ=0.5°，標(biāo)稱再入航程LR0=6500 km、偏差允許值ΔLR=2000 km，再入彈道面傾角偏差允許值ΔiR=0.5°。

月地返回時(shí)間窗口為2017-01-30 03:00:00～2017-01-31 03:00:00(UTC)，粒子的群體規(guī)模數(shù)N為20。采用標(biāo)準(zhǔn)PSO算法，慣性權(quán)重w=0.73，學(xué)習(xí)因子c1=c2=0.75，搜索步數(shù)設(shè)置為300，運(yùn)行一次計(jì)算程序耗時(shí)為54 s，部分粒子的目標(biāo)函數(shù)值一度陷入局部極值，經(jīng)過(guò)228步才收斂，但最終結(jié)果仍不滿足全部約束條件，目標(biāo)函數(shù)值搜索狀態(tài)見圖2。

采用MPSO改進(jìn)算法，搜索步數(shù)設(shè)置為100，運(yùn)行一次計(jì)算程序耗時(shí)為18 s。每次運(yùn)行可獲得一個(gè)優(yōu)化結(jié)果，10次運(yùn)行結(jié)果見表1。

綜合MPSO改進(jìn)算法的運(yùn)行結(jié)果，最小收斂步數(shù)為21，最大收斂步數(shù)為68，平均收斂步數(shù)為38，月地轉(zhuǎn)移軌道優(yōu)化結(jié)果可歸納為兩組，包含所有局部最優(yōu)解，滿足全部約束條件，見表2。兩組優(yōu)化結(jié)果的目標(biāo)函數(shù)值相差僅為0.008，選取逃逸速度增量最小的第1組結(jié)果作為全局最優(yōu)解。

圖2 目標(biāo)函數(shù)值搜索狀態(tài)(標(biāo)準(zhǔn)PSO算法)Fig.2 The search status of objective function(PSO)

一次軌道優(yōu)化過(guò)程的目標(biāo)函數(shù)值搜索狀態(tài)見圖3，逃逸速度增量搜索狀態(tài)見圖4。其中，經(jīng)過(guò)25步搜索后，全部20個(gè)粒子收斂一致。

為驗(yàn)證上述優(yōu)化結(jié)果的正確性，繪制出目標(biāo)函數(shù)的全局等高線和曲面見圖5、圖6。

從圖5和圖6中，可以明顯看到存在兩組局部最優(yōu)解，使目標(biāo)函數(shù)達(dá)到局部極小值，與MPSO算法的優(yōu)化結(jié)果一致，表明該算法對(duì)解決搜索空間大、約束條件多的月地轉(zhuǎn)移軌道優(yōu)化問(wèn)題是正確有效的。

表1 10次MPSO算法運(yùn)行結(jié)果Table 1 Ten computational results of MPSO algorithm

表2 月地轉(zhuǎn)移軌道優(yōu)化結(jié)果Table 2 The optimized result of moon-to-earth transfer trajectory

圖3 目標(biāo)函數(shù)值搜索狀態(tài)(MPSO算法)Fig.3 The search status of objective function(MPSO)

圖4 逃逸速度增量搜索狀態(tài)(MPSO算法)Fig.4 The search status of escape velocity increment(MPSO)

圖5 目標(biāo)函數(shù)全局等高線Fig.5 The contour map of objective function

圖6 目標(biāo)函數(shù)全局曲面Fig.6 The curved surface of objective function

4 結(jié)論

本文結(jié)合我國(guó)月球探測(cè)三期工程的實(shí)際需求，針對(duì)返回中緯度著陸場(chǎng)，考慮月球停泊軌道、月地轉(zhuǎn)移和跳躍式再入等方面的復(fù)雜約束條件，提出了一種改進(jìn)的粒子群算法，實(shí)現(xiàn)了返回窗口內(nèi)逃逸速度增量的最小化，降低了返回燃料消耗，有效解決了多約束條件下的月地轉(zhuǎn)移軌道優(yōu)化問(wèn)題：

(1) 通過(guò)罰函數(shù)法將多約束條件下的月地轉(zhuǎn)移軌道優(yōu)化問(wèn)題轉(zhuǎn)化為無(wú)約束優(yōu)化問(wèn)題，在收斂域內(nèi)外使用不同的目標(biāo)函數(shù)，便于采用粒子群算法進(jìn)行優(yōu)化求解。

(2) 用標(biāo)準(zhǔn)PSO算法求解月地轉(zhuǎn)移軌道優(yōu)化問(wèn)題時(shí)收斂較慢且最終結(jié)果仍不滿足全部約束條件，采用MPSO改進(jìn)算法后能夠較快收斂且滿足全部約束條件，表明改進(jìn)算法對(duì)解決復(fù)雜優(yōu)化問(wèn)題的有效性。

(3) 在預(yù)定月地返回時(shí)間窗口內(nèi)，月地轉(zhuǎn)移軌道優(yōu)化的多次運(yùn)行結(jié)果可歸納為兩組，逃逸速度增量最小的那組結(jié)果是全局最優(yōu)解。

(4) 通過(guò)目標(biāo)函數(shù)的全局等高線和曲面圖分析，驗(yàn)證了優(yōu)化結(jié)果的正確性。

本文提出的優(yōu)化方法，具有求解效率高、收斂性好、適應(yīng)性強(qiáng)等優(yōu)點(diǎn)，能夠有效避免陷入局部極值，可適用于解空間有界的多參數(shù)優(yōu)化問(wèn)題。