吳愛國，鞏志浩

哈爾濱工業(yè)大學(深圳)，機電工程與自動化學院，深圳 518055

燃料最優(yōu)是目前世界各國發(fā)射火星探測器的一項重要指標。自Howard于1962年首次提出利用氣動力輔助軌道轉(zhuǎn)移以來[1]，人們開始研究利用大氣阻力對火星探測器進行減速，使其成為環(huán)火星探測器，這種方式能夠大大減少火星探測器在捕獲過程中的燃料消耗?；鹦翘綔y中的氣動輔助技術(shù)主要有氣動捕獲和氣動剎車2種，其中氣動捕獲技術(shù)是指探測器僅在一次穿越火星大氣后即可在大氣阻力的作用下到達目標軌道，而氣動剎車則需要多次穿越大氣。相比于氣動剎車來說，氣動捕獲可以大大縮短軌道轉(zhuǎn)移的時間。氣動捕獲軌道的優(yōu)化問題是關(guān)系到氣動捕獲能否成功以及最終效果的關(guān)鍵問題，其實際上是探測器在穿越火星大氣時，以燃料最優(yōu)或能量最優(yōu)為性能指標，滿足一系列狀態(tài)約束和邊界約束的非線性優(yōu)化問題。將火星探測器由捕獲軌道以接近最小的燃料消耗和相應的進入條件轉(zhuǎn)移到指定的目標軌道，根據(jù)軌道力學和可控集理論可以求得入口條件和出口條件[2]。氣動捕獲軌道的優(yōu)化問題一般來說沒有解析解，現(xiàn)有的軌跡優(yōu)化方法主要包括間接法、直接法等[3]。

在利用間接法解決軌跡優(yōu)化問題時，首先需要依據(jù)Pontryagin極大值原理將其轉(zhuǎn)化為一個Hamiltonian兩點邊值問題，然后進行求解[4-5]。采用間接法求解軌跡優(yōu)化問題能夠保證最優(yōu)解滿足一階最優(yōu)性必要條件，但是間接法存在求解過程復雜、收斂域小、對初值依賴性高且難以估計等缺點[6]。

相對于間接法來說，直接法在軌跡優(yōu)化中的應用更為廣泛。其利用離散變換將問題轉(zhuǎn)化為非線性規(guī)劃問題，然后根據(jù)優(yōu)化目標對其進行尋優(yōu)。其中配點法是廣泛應用的間接法，它同時離散控制變量和狀態(tài)變量，并通過多項式插值對控制變量和狀態(tài)變量進行擬合，通過優(yōu)化配點處的控制變量搜尋滿足約束的最優(yōu)軌跡[7-8]。正交配點法是最常用的一種配點法，又稱作偽譜法。配點一般選擇正交多項式的根，利用全局插值多項式得到狀態(tài)的表達式，通過對插值多項式求導將微分方程轉(zhuǎn)化為代數(shù)約束，通過積分得到終端狀態(tài)[9-11]。根據(jù)所選用的插值多項式，可以將偽譜法分為Legendre偽譜法、Gauss偽譜法、Chebyschev偽譜法、Radau偽譜法[6]。

近年來，啟發(fā)式智能優(yōu)化算法憑借其魯棒性、適應性以及全局優(yōu)化性能受到越來越多的關(guān)注，已經(jīng)被應用到再入軌跡優(yōu)化和高超聲速飛行器上升軌跡的優(yōu)化等問題中，其無需使用哈密頓函數(shù)、也不用計算其導數(shù)即可獲得接近最佳的結(jié)果[12-13]。另外，也有學者將啟發(fā)式智能優(yōu)化算法和配點法相結(jié)合，使用較高階多項式來近似控制量，并將此控制變量在配點處進行離散，使用啟發(fā)式智能算法優(yōu)化正交配點的值[14-15]。

基于鴿子獨特的導航能力，文獻[16]提出了一種新型的群體智能優(yōu)化方法：鴿群算法。文獻[16]利用所提出的鴿群算法給出了空間機器人路徑規(guī)劃的優(yōu)化方法。文獻[17]引入捕食逃逸機制對鴿群算法進行改進，并將其用于無人機緊密編隊協(xié)同控制的內(nèi)環(huán)控制器的優(yōu)化設(shè)計。文獻[18]提出了基于自適應種群變異的鴿群算法，并將其應用于航天器集群的軌道優(yōu)化，提高了種群演化深度和收斂速度。文獻[19]則將鴿群算法與正交配點法相結(jié)合，用于再入航天器的軌道設(shè)計。文獻[20] 針對基本鴿群算法的尋優(yōu)缺陷利用量子行為規(guī)則對其進行了改進，并利用改進后的鴿群算法用于無人機緊密編隊控制。

本文考慮火星探測器的氣動捕獲軌道的優(yōu)化問題。以火星探測器穿越火星大氣時的動力學方程為基礎(chǔ)，分析并確定了能量最優(yōu)氣動捕獲軌道，給出了優(yōu)化的性能指標。以傾側(cè)角為控制序列，將氣動捕獲軌道優(yōu)化問題轉(zhuǎn)化為多參數(shù)優(yōu)化問題。為了解決該多參數(shù)優(yōu)化問題，本文通過引入?yún)?shù)對基本鴿群算法進行了改進，并將其用于火星探測器的氣動捕獲軌道優(yōu)化。

1 氣動捕獲軌道優(yōu)化問題建模

1.1 大氣內(nèi)飛行動力學方程的建立

假設(shè)火星大氣相對于火星靜止，并且認為飛行過程中探測器的質(zhì)量變化可忽略，則探測器在火星大氣內(nèi)飛行的動力學方程為[21]

(1)

式中：探測器的狀態(tài)由6個變量r、λ、θ、v、γ、ψ(分別為火心距、經(jīng)度、緯度、飛行速度、飛行路徑角及偏航角)確定；gM為火星重力加速度；σ為滾轉(zhuǎn)角；D為阻力加速度；L為升力加速度。D、L、gM的表達式分別為

式中：Sr為探測器的有效截面積；m為探測器的質(zhì)量；CD為阻力系數(shù)；CL為升力系數(shù)；μ為火星引力常數(shù)；ρ為火星大氣密度，不同高度的火星大氣密度不同，火星大氣密度與高度之間的關(guān)系為

ρ=ρ0exp(-h/hs)

(2)

式中：hs=8 805.7 m為火星大氣標高；ρ0=0.0147 4 kg/m3為火星大氣參考密度。又由式(1) 可知，火心距、飛行速度和飛行路徑角不改變軌道平面，并且其變化與其他3個狀態(tài)無關(guān)，所以當不考慮軌道平面的變化時，上述大氣內(nèi)飛行動力學方程可以簡化為

(3)

1.2 能量最優(yōu)氣動捕獲軌道的確定

利用火星大氣對探測器進行減速，通過調(diào)整姿態(tài)進而改變所受升力的大小來控制飛行軌跡，可以認為此過程是不消耗能量的無動力飛行過程，這里的能量最優(yōu)指的是由捕獲軌道轉(zhuǎn)移到目標軌道時所消耗的能量最少。由捕獲軌道轉(zhuǎn)移到目標軌道所需要的速度增量為[2]

(4)

(5)

總速度增量為

ΔV=ΔV1+ΔV2

(6)

則

(7)

又

(8)

(9)

式中：rca為捕獲軌道遠火點半徑；rcp為捕獲軌道近火點半徑；rta為目標軌道遠火點半徑；rtp為目標軌道近火點半徑。由式(9)可知，對于給定的rcp，ΔV在區(qū)間rca∈[0,rta]內(nèi)單調(diào)遞減，在區(qū)間rca∈[rta,+∞]內(nèi)單調(diào)遞增。因此，當rca=rta時，ΔV取得最小值，即當探測器穿越大氣后的后捕獲軌道的遠地點高度與目標軌道高度相等，也就是二者相切時，軌道轉(zhuǎn)移所需的能量消耗最少。

由動量矩守恒原理和能量守恒原理可知

rtavca=ravfcosγf

(10)

式中：vf為探測器飛出大氣時的速度；γf為探測器飛出大氣時的航跡角；ra為火星大氣的高度；vca為捕獲軌道遠火點處速度。根據(jù)能量守恒有

(11)

聯(lián)立式(10)和式(11)可得

(12)

由式(12)可得

(13)

因此，當捕獲軌道遠火點與目標軌道遠火點相切時有

(14)

式中:vta為目標軌道遠火點處的速度，當目標軌道確定時其為一常值。由式(14)可知，探測器到達目標軌道高度時所需脈沖變軌的大小與飛出角度γf密切相關(guān)?？梢?，飛出角度γf越小，在到達目標軌道高度時所需的脈沖ΔV也就越小。

氣動捕獲軌道優(yōu)化問題實質(zhì)上就是對探測器穿越火星大氣時的軌跡進行控制，使得探測器在飛出大氣時恰好到達或接近能量最優(yōu)氣動捕獲軌道。由前述分析可知，要使氣動捕獲軌道為最優(yōu)氣動捕獲軌道，則氣動捕獲軌道的遠地點應與目標軌道相切，飛出時的航跡角盡量小，且飛出速度與飛出時的航跡角滿足式(13)。取最優(yōu)控制指標為

(15)

式中：

為對應飛出航跡角滿足式(13)的飛出速度。進入條件為

(16)

(17)

式中：γfmax為保證捕獲軌道與目標軌道相切的最大飛出航跡角。探測器大氣內(nèi)飛行時的總飛行約束為

(18)

為將氣動捕獲軌道優(yōu)化問題轉(zhuǎn)化為多參數(shù)優(yōu)化問題，選取傾側(cè)角σ作為控制變量，將控制變量在配點處進行離散，并以這些離散點作為節(jié)點構(gòu)造Lagrange插值多項式來近似擬合控制變量，選取Chebyschev多項式的根作為配點，有

τl=cos(πl(wèi)/n)l=0,1,…,n

(19)

式中：n為插值多項式的階數(shù)，并且插值點落在[-1,1]之間，對t∈(t0,tf)進行時域轉(zhuǎn)換

(20)

在這里由于大氣內(nèi)飛行的時間未知，所以將飛出時間tf也作為待優(yōu)化的變量。由于控制量σ是關(guān)于時間t的連續(xù)函數(shù)，很難得到其精確表達式。為此，按照式(19)和式(20)的時間轉(zhuǎn)換方式對控制量σ進行離散化，即

σl=σ(t(τl))l=0,1,…,n

因此優(yōu)化參數(shù)為

(21)

于是，對于氣動捕獲軌道的優(yōu)化問題就轉(zhuǎn)化為對大氣內(nèi)飛行時間tf以及控制序列σ在配點處的離散值的參數(shù)優(yōu)化問題。針對這個多參數(shù)優(yōu)化問題，第2節(jié)將提出一種改進鴿群算法對上述待優(yōu)化的參數(shù)向量進行優(yōu)化。

2 鴿群算法的分析與改進

2.1 原始鴿群算法及其存在的問題

鴿群算法是受鴿子覓食的啟發(fā)提出的一種仿生學算法[16]。當鴿子距離目標較遠時，通過感知磁場在腦海中形成地圖，來不斷接近目標；隨著鴿子不斷接近目的地，其導航工具由磁場變?yōu)槟康牡馗浇牡貥?。熟悉地標的鴿子直接飛往目的地，其他鴿子則跟隨那些熟悉地標的鴿子向目的地飛行。

設(shè)某優(yōu)化問題的目標函數(shù)為J(ξ)，其中優(yōu)化變量為ξ。鴿群算法在處理優(yōu)化問題時模擬鴿子導航的過程。假設(shè)在優(yōu)化時需要N只鴿子，每只鴿子的位置Xi,i=1,2,…,N表示優(yōu)化變量ξ的值，在優(yōu)化時化還需要對應鴿子的速度Vi,i=1,2,…,N。假設(shè)優(yōu)化變量是m維的，則每只鴿子的位置和速度分別記作在采用鴿群算法優(yōu)化時，設(shè)鴿群中有N只鴿子，它們的位置和速度分別為Xi和Vi。該算法求解優(yōu)化問題時包括2部分迭代過程，分別稱作地圖迭代和地標迭代。在地圖迭代中，鴿群中的每一個體通過種群中的最優(yōu)信息來更新自身的位置和速度，更新公式為

(22)

式中：R為迭代因子；rand為[0,1]區(qū)間上的隨機值；t=1,2,…,tmax為當前迭代次數(shù)，tmax為迭代總次數(shù)；i=1,2,…,N代表第i只鴿子；Vi(t)為第i個鴿子在第t次迭代的速度；Xi(t)為第i個鴿子在第t次迭代的位置；G(t-1)為Xi(t-1) (i=1,2,…,N)中的某一個，其滿足：

當?shù)螖?shù)到達最大迭代次數(shù)的75%以后，進入地標迭代的過程，選取Xc(t)為第t次迭代鴿群的中心點，作為鴿子飛行的參考方向?qū)γ恳恢圾澴舆M行更新，更新公式為

(23)

式中：

其中：fitness[Xi(t)]為當前迭代下鴿子的適應度；N(t)為第t次迭代下種群的數(shù)量；ε為防止分母為零而附加的充分小的數(shù)。

可見，在地標迭代過程中，每次迭代舍去部分個體，既保證了算法的較優(yōu)信息又提高了迭代收斂的速率。原始鴿群算法的全局搜索能力主要體現(xiàn)在地圖迭代的過程中，動態(tài)權(quán)值e-RN體現(xiàn)了算法的探索能力，在鴿群向當前最優(yōu)值收斂的過程中保持鴿群的多樣性。然而這種方式具有局限性，在鴿群數(shù)目一定的情況下，R值的選取在很大程度上決定了全局搜索的結(jié)果。

圖1為R在不同取值下動態(tài)因子的變化情況?？梢钥闯觯擱取較大值0.3時，迭代到20次e-Rt就已衰減接近為零，之后的迭代過程主要體現(xiàn)在向最優(yōu)值靠近，失去了探索能力，算法易早熟；當R取較小值0.002 6時，在迭代完成時e-Rt仍然具有較大的取值，算法不易收斂；為了平衡二者的關(guān)系，可以通過指定迭代完成時期望的動態(tài)因子的取值自適應地選取R值。令迭代完成時期望的動態(tài)因子的值為α,有e-RT=α,則R=-(1/T)lnα, 其中T為迭代次數(shù)。

圖1 R的取值對動態(tài)權(quán)值的影響Fig.1 Influence of R on dynamic weight

通過自適應的方式選取R值，雖然能夠在一定程度上緩解易早熟和不易收斂的問題，但是卻不能從根本上很好地平衡全局搜索與局部搜索之間的關(guān)系。在實際應用中，在迭代初期，往往更注重于全局搜索能力，而在迭代后期則希望局部搜索能力更強一些。由圖1可以看出，當取α=0.01 時，在迭代后期動態(tài)因子較小，算法主要在最優(yōu)值附近進行局部搜索，然而在迭代初期，動態(tài)因子下降很快，并不能很好地進行全局搜索；當取α=0.4時，雖然使得動態(tài)因子下降速度減緩，全局搜索能力增強，但是卻削弱的算法的局部搜索能力。

2.2 改進的鴿群算法

本節(jié)在原有鴿群算法的負指數(shù)動態(tài)權(quán)值的基礎(chǔ)上，為了平衡全局搜索與局部搜索提出一種廣義的動態(tài)權(quán)值項。即將式(23)修改為

(24)

式中：

(25)

且有

0≤f(t)≤1

(26)

為了更加方便直觀地確定參數(shù)ω、b的取值， 引入另外2組參數(shù)(t1,α1)、(t2,α2)，令

(27)

表示當?shù)螖?shù)為t1時，動態(tài)權(quán)值取值為α1；當?shù)螖?shù)為t2時，動態(tài)權(quán)值取值為α2。即

(28)

求解可得

(29)

因此，可以通過指定k、(t1,α1)、(t2,α2)來唯一確定動態(tài)權(quán)值。

下面對各參數(shù)的取值情況進行討論。對式(25) 求導可得

(30)

當k≥1時，對式(25)取極限可得

(31)

此時k、(t1,α1)、(t2,α2)可任意取值都可以滿足式(26)的約束。

當k<1時，此時k的取值與(t1,α1)、(t2,α2)相互制約，由式(26)可得

(32)

即

1-e-ω b≤k

(33)

又由式(29)可得

(34)

(35)

不妨取t1=T/2、t2=T，則有

(36)

(37)

即

(38)

又

(39)

(40)

(41)

(42)

由式(40)可知此時k≥0。當取k=0時，由式(29) 有b=0，又由式(25)可得

(43)

此時改進算法退化為原始的鴿群算法，也就是說原始鴿群算法是改進后算法的一種特殊情況。

改進后的鴿群算法能夠通過調(diào)整動態(tài)權(quán)值很好地平衡全局搜索與局部搜索之間的關(guān)系，既保證了算法的收斂性，又大大加強了算法的探索能力。圖2為改進后的鴿群算法在k、(t1=N/2,α1)、(t2=N,α2=0.1)的不同取值下動態(tài)權(quán)值的變化情況，其中N=200。通過調(diào)節(jié)α1和k的值可得到不同的動態(tài)權(quán)值曲線。

可以看出，當k=0且b=0時，算法即為原始的鴿群算法，即原算法是改進后的一種特殊情況。當k≠0時，可以通過調(diào)節(jié)(t1,α1)、(t2,α2)選擇合適的動態(tài)權(quán)值曲線，在搜索前期使動態(tài)權(quán)值緩慢下降，保持搜索活力，提高全局搜索的能力；在搜索后期，使得動態(tài)權(quán)值取得較小的值，保證算法的收斂。相比于原算法，改進后的算法不僅兼容了原算法中指數(shù)形式下降的動態(tài)權(quán)值，而且對其進行了擴展，在形式選擇上提供了更大的自由度，提升了全局搜索與局部搜索之間的平衡。

圖2 動態(tài)權(quán)值隨迭代次數(shù)的變換Fig.2 Evolution of dynamic weights

3 氣動捕獲軌道優(yōu)化

采用改進的鴿群算法對氣動捕獲軌道進行優(yōu)化。鴿群中的每一只鴿子所在的位置代表所求問題的一個可能解，其初值在搜索空間內(nèi)隨機生成。選擇鴿群中鴿子的數(shù)目為m，則鴿群的位置X可表示為

(44)

式中：Xk(k=1,2,…,m)為每一只鴿子當前所得到的待優(yōu)化變量的值，包括探測器飛出大氣的時間tf，以及以傾側(cè)角為控制序列在配點處的離散值。鴿群的更新速度V與X同維，首先對每一只鴿子所優(yōu)化得到的控制量進行Lagrange插值得到控制序列，然后根據(jù)所得到的控制序列和優(yōu)化得到的飛出時間tf對動力學方程式(3)進行積分得到各個狀態(tài)的值，計算此時的目標函數(shù)值，然后進行迭代尋優(yōu)。

鴿群算法本身并不具備處理約束的能力，在進行迭代尋優(yōu)的過程中為了保證滿足過程約束和終端約束，將這些約束條件以罰函數(shù)的方式加到目標函數(shù)當中，經(jīng)過修正的目標函數(shù)為

(45)

4 仿真分析

表1 火星探測器參數(shù)Table 1 Parameters of Mars explorers

表2 火星參數(shù)Table 2 Parameters of Mars

采用改進鴿群算法，取傾側(cè)角為五階插值多項式時，改進鴿群算法的迭代收斂曲線以及優(yōu)化結(jié)果如圖3所示，圖3(a)為探測器飛行高度隨時間變化的情況，可見進入高度和飛出高度相等，都等于火星大氣高度，并且最小飛行高度大于所要求的35 000 m。圖3(b)為飛行速度隨時間變化的曲線，飛出速度與飛出航跡角基本滿足式(13)。圖3(c)為航跡角隨時間變化的曲線，圖3(d)和圖3(e) 為路徑約束，分別表示過載和熱流率隨時間的變化情況，可見大氣內(nèi)飛行過程可以滿足過程約束，圖3(f)為控制變量隨時間變化的曲線。綜上，通過改進鴿群算法，以傾側(cè)角為控制變量可以得到一條接近能量最優(yōu)的氣動捕獲軌道。

圖3 氣動捕獲軌道優(yōu)化結(jié)果Fig.3 Performance of optimized aerocapture orbit

5 結(jié) 論

1) 在對火星探測器的氣動捕獲過程分析的基礎(chǔ)上，從脈沖變軌的角度給出了氣動捕獲軌道的優(yōu)化指標。

2) 對原始鴿群算法的局限性進行分析，提出一種能夠平衡收斂速度與全局搜索關(guān)系的改進算法。

3) 將氣動捕獲軌道優(yōu)化問題轉(zhuǎn)化為多參數(shù)優(yōu)化問題，并利用提出的改進個鴿群算法對其進行優(yōu)化。