彭玉明 尤偉 馬彥涵

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

1 引言

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

2 Simpson 法

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

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

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

其邊界條件為

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

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

解上述線性方程組得

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

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

式中iu為第i個(gè)節(jié)點(diǎn)的控制量。

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

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

圖1 Simpson 法Fig.1 Simpson method

原來(lái)的動(dòng)力學(xué)方程由N×m個(gè)Hermite-Simpson 積分形式的約束代替

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

定義決策變量Z為

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

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

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

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

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

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

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

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

3.2 約束條件

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

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

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

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

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

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

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

為了安全打開(kāi)降落傘，末端時(shí)刻動(dòng)壓應(yīng)不超過(guò)給定的最大值

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

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

4 仿真分析

4.1 仿真參數(shù)

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

不考慮擾動(dòng)對(duì)末端飛行路徑角的影響，靈敏度函數(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取不同值時(shí)所得到的最優(yōu)控制曲線。

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

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

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

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

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

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

圖4 反饋控制Fig.4 Feedback control

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

5 結(jié)束語(yǔ)

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

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]高濱.火星探測(cè)器著陸技術(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]于瑩瀟, 田佳林.美國(guó)新型載人火星探測(cè)技術(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), 榮偉.火星探測(cè)器進(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]陳曉, 張偉, 彭玉明.基于器間測(cè)量的火星進(jìn)入過(guò)程實(shí)時(shí)高精度導(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)