一種考慮狀態(tài)靈敏度的火星進(jìn)入軌跡優(yōu)化方法

彭玉明 尤偉 馬彥涵

（上海衛(wèi)星工程研究所，上海 200240）

1 引言

大氣進(jìn)入過程的運動方程形式復(fù)雜，進(jìn)入軌跡對控制變量高度敏感且進(jìn)入過程的非線性約束較強，因此進(jìn)入軌跡優(yōu)化問題一直是研究的熱點。目前火星探測器大都沿用“海盜號”的氣動外形，依靠滾轉(zhuǎn)角的變化實現(xiàn)軌跡控制，控制能力較弱，因此在設(shè)計標(biāo)稱軌跡時還需要考慮控制系統(tǒng)的跟蹤能力。傳統(tǒng)的軌跡優(yōu)化方法首先根據(jù)約束條件確定一條軌跡，然后再設(shè)計控制系統(tǒng)，若控制系統(tǒng)無法跟蹤標(biāo)稱軌跡，則重新規(guī)劃，如此反復(fù)直到滿足要求。在對軌跡進(jìn)行優(yōu)化設(shè)計時，完全不考慮控制系統(tǒng)的跟蹤能力，設(shè)計過程完全依靠經(jīng)驗費時費力，且不一定能取得滿意的效果。因此本文提出一種考慮控制飽和的進(jìn)入軌跡優(yōu)化算法，把末端狀態(tài)的靈敏度加入到目標(biāo)函數(shù)中，并與末端高度最高、燃料最省通過罰函數(shù)法結(jié)合到一起形成最終的優(yōu)化目標(biāo)[1-3]。

2 Simpson 法

直接配點法的核心思想是通過離散化把最優(yōu)控制問題轉(zhuǎn)化成非線性規(guī)劃問題，直接配點法是將整個時間過程劃分成N 段，每一段的兩端稱為節(jié)點，然后用多項式逼近每一段的狀態(tài)和控制變量[4-5]。根據(jù)多項式階次的不同，直接配點法又可以分為低階的梯形法、Simpson 法和高階的四階、五階方法，本文采用三階Simpson 法。

在時間[t0,tf]內(nèi)將連續(xù)時間分成N 段（tf為未端時刻），每一個子區(qū)間為[ti,ti+1]，記Ti為相鄰2個節(jié)點的時間間隔，s=（t -ti）Ti（s為歸一化變量；t為時間），在每個子區(qū)間上狀態(tài)變量x 可以用三次Hermite多項式表示

式中 c0、c1、c2、c3為多項式系數(shù)。

其邊界條件為

式中 x1和x2分別表示相鄰兩個節(jié)點的狀態(tài)。

把邊界條件代入式(1)中得

解上述線性方程組得

在子區(qū)間的中點處，即s=0.5 時，將式(3)代入式(1)得

式中 xmi為中點處的狀態(tài)量；xi為第i個節(jié)點的狀態(tài)量；xi+1為第i+1個節(jié)點的狀態(tài)量；fi，fi1+分別表示函數(shù) f（x,u,t）在第i個子區(qū)間兩端點處的函數(shù)值，（其中，u為控制量）即

式中iu為第i個節(jié)點的控制量。

式中 ui1+為第i+1個節(jié)點的控制量。

當(dāng)由公式(5)近似得到的中點處的導(dǎo)數(shù)與動力學(xué)方程在中點處的值足夠接近時，就認(rèn)為多項式很好的逼近動力學(xué)微分方程，如圖1所示。

圖1 Simpson 法Fig.1 Simpson method

原來的動力學(xué)方程由N×m個Hermite-Simpson 積分形式的約束代替

式中 umi為中點處的控制變量，由兩端線性插值得到；tmi為中點對應(yīng)的時刻。

定義決策變量Z為

式中 x為各節(jié)點狀態(tài)變量構(gòu)成的狀態(tài)向量；mu為中點處控制量構(gòu)成的控制向量。

定義Δi為各節(jié)點處狀態(tài)量逼近誤差

這樣動力學(xué)微分方程約束轉(zhuǎn)變成N×m個等式約束，問題轉(zhuǎn)化為有約束的非線性規(guī)劃問題，搜索決策變量最優(yōu)值，使得Δi趨于零。

3 考慮狀態(tài)靈敏度的進(jìn)入段最優(yōu)軌跡

3.1 優(yōu)化目標(biāo)

由于火星大氣特別稀薄，降落傘減速效果有限，為了使探測器能夠安全著陸，希望開傘高度越高越好，所以目標(biāo)函數(shù)中應(yīng)包含末端高度；另一方面，探測器上攜帶的燃料有限，過多的消耗燃料是不可取的，因此燃料最省也是優(yōu)化的目標(biāo)之一。然而，如果簡單的以末端高度最高或燃料最省為性能指標(biāo)進(jìn)行優(yōu)化，得到的控制變量的形式很可能是最大–最小–最大的形式，雖然控制量留取一定余量，但是在存在較大不確定性擾動的情況下，會導(dǎo)致控制飽和，更加不利的是，像“海盜號”這一類小升阻比的探測器，軌跡控制能力較弱，加之進(jìn)入段惡劣的氣動加熱環(huán)境以及大量不確定性因素，探測器很可能會因為控制飽和而無法跟蹤標(biāo)稱軌跡，因此在選取目標(biāo)函數(shù)時還要考慮控制系統(tǒng)能否跟蹤上標(biāo)稱軌跡[6-7]。末端高度、燃料消耗以及跟蹤控制性能等幾個指標(biāo)是相互對立的，因此采用罰函數(shù)法把3個性能指標(biāo)統(tǒng)一到一起，本文選取目標(biāo)函數(shù)J為

式中 cS,cu為加權(quán)系數(shù)；h（tf）為探測器末端高度；Js為末端狀態(tài)對進(jìn)入軌跡的靈敏度函數(shù)；Ju為進(jìn)入過程中關(guān)于燃料消耗的目標(biāo)函數(shù)

式中 σ為滾轉(zhuǎn)角。

3.2 約束條件

進(jìn)入軌跡除了動力學(xué)方程的約束外，還考慮了初始條件約束、終端條件約束、過程約束和控制約束。初始狀態(tài)約束包括初始高度0r 、速度0v 、飛行路徑角0γ和航程 s0，即

探測器最后要到達(dá)指定開傘區(qū)域，并且滿足開傘條件，由于末端高度是目標(biāo)函數(shù)之一，因此這里不做約束，末端條件約束為

式中fs為末端狀態(tài)靈敏度；fv為末端速度。末端時刻 tf自由無約束。

進(jìn)入過程中，為了結(jié)構(gòu)和設(shè)備的安全，需要考慮過載約束[8]。鑒于本文的研究對象為彈道升力式飛行器，其軸向和法向均有可能產(chǎn)生較大過載，因此本文過載約束na取總過載約束

式中 L為升力加速度；D為阻力加速度；na為過載；nmax為最大過載。

對動壓的限制應(yīng)從其峰值和末端時刻來考慮，進(jìn)入過程中動壓峰值q 應(yīng)不超過給定的最大值qmax。

式中 ρ為大氣密度；V為飛行速度。

為了安全打開降落傘，末端時刻動壓應(yīng)不超過給定的最大值

式中 c為常數(shù)；b 是與探測器半徑相關(guān)的常數(shù)[9]。

由于只考慮縱向運動，所以假設(shè)σ ∈ [ 0,π ]，為了進(jìn)行側(cè)向制導(dǎo)，控制量需要留取一定余量，控制約束為

式中 σmin為最小滾轉(zhuǎn)角；σmax為最大滾轉(zhuǎn)角。同時由于探測器上姿態(tài)執(zhí)行機(jī)構(gòu)的限制，滾轉(zhuǎn)角角速度和角加速度還必須滿足以下條件

式中 σ˙min為最小滾轉(zhuǎn)角速率；為最大滾轉(zhuǎn)角速率；為最小滾轉(zhuǎn)角加速度；為最大滾轉(zhuǎn)角加速度。

3.3 狀態(tài)靈敏度計算

控制系統(tǒng)的抗擾動性能體現(xiàn)在末端狀態(tài)對進(jìn)入軌跡的靈敏度上，靈敏度越小，說明末端狀態(tài)對進(jìn)入軌跡的變化越不敏感，也就是說抗擾動性能越強，反之亦然。因此，若是把所有時間點上狀態(tài)擾動對末端狀態(tài)的影響用靈敏度函數(shù)表示出來，然后加入到目標(biāo)函數(shù)中，這樣在設(shè)計最優(yōu)軌跡時就考慮了控制系統(tǒng)的抗擾動性能。

根據(jù)靈敏度的定義，任意軌跡X（t t0,x0）相對初始時刻狀態(tài) x0的靈敏度矩陣為

為簡化表達(dá)把S（t t0,x0）表示成 S（t,t0），它滿足以下微分方程

則反饋控制律變?yōu)?/p>

靈敏度矩陣微分方程變?yōu)?/p>

根據(jù)靈敏度矩陣的定義，容易證明上述靈敏度矩陣具有如下性質(zhì)

式中 I7×7為7 維單位矩陣。

在進(jìn)入過程中，主要關(guān)注狀態(tài)擾動對末端狀態(tài)的影響，因此定義ξ（x（tf）,tf）函數(shù)

根據(jù)定義ξ（x（tf）,tf）相對t 時刻的狀態(tài)擾動的靈敏度矩陣可以表示為

所以ξ（x（tf）,tf）相對t 時刻的狀態(tài)擾動的靈敏度矩陣等于S(tf,t)。很顯然只有矩陣中的每一個元素同時最小才能使靈敏度最小[10]。通過積分得到末端狀態(tài)對所有時刻狀態(tài)擾動的靈敏度，因此靈敏度函數(shù)JS為

式中Si,j（tf,t ）表示靈敏度矩陣中第(i,j)個元素；ci為加權(quán)系數(shù)。

4 仿真分析

4.1 仿真參數(shù)

探測器進(jìn)入火星大氣層的初始約束條件為：[r0,v0,γ0,s0]=[3 522 000,6 000,-0.2,0]，末端約束條件v（tf）≤350m/s ；s（tf）=935km ；最大過載值nmax=8gn；最大動壓值q=10 0 00Pa，末端動壓 qf≤750Pa；最大加熱率=70kW/m2；控制變量約束為，11.5°≤ σ≤168.5°，-20(°) s≤ σ˙≤20(°) s，-5 (°) s2≤≤5(°) s2。

不考慮擾動對末端飛行路徑角的影響，靈敏度函數(shù)的加權(quán)系數(shù)為 c1=1,c2=1,c3=0,c4=0.01；反饋增益為k=-[0.01 0.005 50 0.001]。

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

圖2、3 是加權(quán)系數(shù) cs和cu取不同值時所得到的最優(yōu)控制曲線。

圖2 進(jìn)入段最優(yōu)控制軌跡Fig.2 Entry optimal control trajectory

圖3 cu=10 時進(jìn)入段最優(yōu)控制軌跡Fig.3 Entry optimal control trajectory while cu=10

從圖中可以看出，當(dāng)cs=0，即目標(biāo)函數(shù)中不包括靈敏度，最優(yōu)控制軌跡大部分時間處于控制約束邊緣，當(dāng)存在擾動時，控制系統(tǒng)很可能會因為控制飽和而無法跟蹤標(biāo)稱軌跡。當(dāng)uc=10，即考慮燃料消耗以后，控制軌跡變得平滑，但是仍有大部分時間處于飽和狀態(tài)。隨著 cs的不斷增大，最優(yōu)控制軌跡離約束邊界越來越遠(yuǎn)，抗擾動性能逐漸增強，但代價是末端高度越來越低，因此實際應(yīng)用時需要根據(jù)具體情況采取不同方案。表1 給出了cs和cu取不同值時的末端高度。最優(yōu)軌跡確定以后采用蒙特卡洛仿真方法進(jìn)一步驗證考慮控制飽和的火星進(jìn)入軌跡優(yōu)化方法的有效性。具體仿真參數(shù)見表2。

表1 不同參數(shù)情況下的末端高度Tab.1 Terminal altitude for different parameters

表2 進(jìn)入段蒙特卡洛仿真參數(shù)Tab.2 Entry Monte Carlo simulation parameters

圖4 是采用線性反饋控制律跟蹤最優(yōu)軌跡得到的實際控制軌跡，如前文所料，cs越大，控制量飽和的時間越少，當(dāng)cs=0時，即不考慮狀態(tài)靈敏度，在200s 以后控制量一直處于飽和狀態(tài)，而考慮了狀態(tài)靈敏度以后，在150s 以后基本不會出現(xiàn)控制量飽和。圖5 是 cs取不同值時末端狀態(tài)誤差分布。從圖中可以看到僅考慮末端高度時，末端誤差高度、速度、航程誤差分別約為2000m，100m，10km，隨著 cs的逐漸增大，目標(biāo)函數(shù)中狀態(tài)靈敏度的權(quán)重越來越大，末端狀態(tài)誤差逐漸減小，當(dāng)cs=0.005時末端高度、速度和航程誤差分別減小到100m，5m/s，1km，進(jìn)入著陸精度顯著提高，這說明考慮控制飽和的進(jìn)入軌跡優(yōu)化算法是有效的。

圖4 反饋控制Fig.4 Feedback control

圖5 末端狀態(tài)誤差分布Fig.5 Terminal state error distribution

5 結(jié)束語

本文采用直接配點，對火星大氣進(jìn)入軌跡進(jìn)行優(yōu)化設(shè)計，并考慮了狀態(tài)靈敏度對進(jìn)入軌跡的影響，仿真算例表明，在目標(biāo)函數(shù)中加入狀態(tài)靈敏度以后，進(jìn)入軌跡不再靠近邊界，控制飽和的現(xiàn)象明顯減少，進(jìn)而減小了末端狀態(tài)估計誤差。

References)

[1]Leavitt J,Mease K.Feasible Trajectory Generation for Atmospheric Entry Guidance[J].Journal of Guidance,Control,and Dynamics,2007,30(2)：473-481.

[2]高濱.火星探測器著陸技術(shù)[J].航天返回與遙感,2009,30(1)：1-9,20.GAO Bin.Mars Exploration Entry,Descent and Landing Technologies[J].Spacecraft Recovery & Remote Sensing,2009,30(1)：1-9,20.(in Chinese)

[3]于瑩瀟, 田佳林.美國新型載人火星探測技術(shù)方案[J].航天返回與遙感,2009,30(3)：1-7,15.YU Yingxiao,TIAN Jialin.USA New Mars Exploration Technology Concept[J].Spacecraft Recovery & Remote Sensing,2009,30(3)：1-7,15.(in Chinese)

[4]Shen H J,Seywald H,Powell R W.Desensitizing the Minimum-fuel Power Descent for Mars Pinpoint Landing[J].Journal of Guidance Control and Dynamics,2010,33(1)：108-115.

[5]Seywald H,Kumar R R.Desensitized Optimal Trajectories[R].AMA report No 03-16,2003.

[6]包進(jìn)進(jìn), 榮偉.火星探測器進(jìn)入階段穩(wěn)定性分析[J].航天返回與遙感,2011,32(4)：6-13.BAO Jinjin,RONG Wei.Mars Explorer Stability Analysis in the Entry Phase[J].Spacecraft Recovery & Remote Sensing,2011,32(4)：6-13.(in Chinese)

[7]Seywald H.Desensitized Optimal Trajectories with Control constraints[R].AMA report No.03-28,2003.

[8]Shen H J,Seywald H,Powell R W.Desensitizing the Pin-point Landing Trajectory on Mars[C].AAS/AIAA Astrodynamics Specialist Conference,Honolulu,HI,2008.

[9]Tang S,Conway B A.Optimization of Low-thrust Interplanetary Trajectories Using Collocation and Nonlinear Programming[J].Journal of Guidance Control and Dynamics,1995,18(3)：599-604.

[10]陳曉, 張偉, 彭玉明.基于器間測量的火星進(jìn)入過程實時高精度導(dǎo)航[J].航天返回與遙感,2012,33(6)：17-23.CHEN Xiao,ZHANG Wei,PENG Yuming.Mars Entry Real-time Navigation Based on Orbiter Tracking Data[J].Spacecraft Recovery & Remote Sensing,2012,33(6)：17-23.(in Chinese)