火星大氣進入滑模自抗擾制導(dǎo)方法*

李 翔,朱東方,胥 彪, 聶欽博

(1. 南京航空航天大學(xué) 航天學(xué)院·南京·210016；2. 上海航天控制技術(shù)研究所·上海·201109)

0 引 言

火星是距離地球最近的類地行星,對火星的探索研究對人類了解宇宙及開發(fā)外太空資源具有重要意義?；鹦翘綔y代表了行星際太空探索和技術(shù)的頂峰,是未來深空探測領(lǐng)域的最大熱點,也是各個航天大國行星探測的首個目標(biāo)。美國和歐盟都宣布計劃2030年左右實現(xiàn)載人登火,而我國也明確指出在2020年左右完成對火星無人探測的環(huán)繞、著陸與巡視三個任務(wù)指標(biāo)[1]。

火星探測器的安全著陸是探測任務(wù)完成最為重要的保障,然而技術(shù)難度極高,在以往的火星著陸任務(wù)中成功次數(shù)不足一半。因此探測器進入、下降、著陸(Entry, Decent, Landing,EDL)過程是未來火星探測任務(wù)的關(guān)鍵技術(shù)之一,而大氣進入段的環(huán)境最為惡劣,受到火星大氣模型不確定、氣動環(huán)境復(fù)雜、強非線性、強耦合等問題的影響,對著陸精度的影響很大。為了提高探測器在有較高科學(xué)價值的特定區(qū)域精確著陸的能力,可以利用文獻[2]中總結(jié)的自主導(dǎo)航方法的優(yōu)勢,在進入段對火星探測器進行制導(dǎo)策略的設(shè)計。進入制導(dǎo)包括縱向制導(dǎo)和橫向制導(dǎo),縱向制導(dǎo)是更為重要的方面,目前這個階段制導(dǎo)律設(shè)計主要有兩種方法:標(biāo)稱軌跡制導(dǎo)法和預(yù)測校正制導(dǎo)法[3]。標(biāo)稱軌跡法就是事先按照要求設(shè)計好一條標(biāo)稱軌跡,存儲在器載計算機中,然后根據(jù)實時跟蹤誤差設(shè)計制導(dǎo)律來跟蹤這條軌跡。預(yù)測校正法是在每個制導(dǎo)周期內(nèi)根據(jù)當(dāng)前狀態(tài)預(yù)測出終端狀態(tài)值,與期望值對比產(chǎn)生誤差信號,最后通過控制器生成控制信號實現(xiàn)對飛行軌跡的控制。這兩種方法各有利弊,前者算法簡單易于實現(xiàn),但其對初始誤差和外部擾動敏感而造成該方法的制導(dǎo)精度較低,需要高性能的控制器去改善缺點。后者較前者對外界干擾具有更好的自適應(yīng)能力,但其預(yù)測過程依賴于精確的動力學(xué)模型。盡管預(yù)測校正法的制導(dǎo)精度更高,但由于模型精度不高以及器載計算機性能很難滿足計算要求,導(dǎo)致這種方法很難在工程實際中得到有效應(yīng)用。因此標(biāo)稱軌跡法仍然是近年來火星探測任務(wù)設(shè)計中優(yōu)先考慮的方法,所以進一步提高其自適應(yīng)能力和魯棒性尤為重要。阻力加速度跟蹤方法是一種很成熟的標(biāo)稱軌跡制導(dǎo)法,在航天飛機和好奇者號任務(wù)中都得到成功應(yīng)用。航天飛機中采用了傳統(tǒng)的比例、積分、微分控制方法[4],雖然技術(shù)比較成熟,但是在一些線性假設(shè)條件下求得的,且增益系數(shù)的整定比較麻煩。文獻[5]運用反饋線性化(Feedback Linearization, FL)方法設(shè)計了標(biāo)稱軌跡跟蹤控制,證明了該方法的優(yōu)越性。文獻[6]設(shè)計了對低升阻比航天器適用的魯棒制導(dǎo)律。文獻[7]運用了非線性預(yù)測控制的方法設(shè)計跟蹤控制器,提高了制導(dǎo)精度,但增大了計算負(fù)擔(dān)。文獻[8]對于火星大氣進入段軌跡跟蹤控制,設(shè)計了基于連續(xù)有限時間滑?？刂破?具有較高的魯棒性。文獻[9]利用滑?？刂品椒▉硐薪鐢_動的影響,并利用基于徑向基函數(shù)神經(jīng)網(wǎng)絡(luò)方法用來逼近未知擾動。文獻[10]為了減小大氣密度和空氣動力學(xué)系數(shù)的有界誤差的影響,使用了基于命令生成跟蹤器的模型參考自適應(yīng)方法。文獻[11]結(jié)合軌跡在線更新策略和動態(tài)逆方法,給出了一種再入飛行器閉環(huán)穩(wěn)定的魯棒跟蹤方法。文獻[12]針對再入飛行器設(shè)計了具有前向補償?shù)闹茖?dǎo)控制一體化系統(tǒng),具有良好的制導(dǎo)性能和魯棒性。

預(yù)測-校正制導(dǎo)算法需要精確的動力學(xué)模型來預(yù)測最終的狀態(tài)變量,而由于火星大氣環(huán)境的復(fù)雜多變,存在各種不確定性,難以得到精確的模型,限制了這種方法的應(yīng)用。故本文采用經(jīng)典的基于阻力加速度跟蹤的標(biāo)稱軌跡制導(dǎo)方法。反饋線性化方法在解決非線性系統(tǒng)控制律設(shè)計方面表現(xiàn)出了一定的優(yōu)越性,但不確定因素很大程度上影響了其在進入制導(dǎo)段的跟蹤效果。雖然非線性預(yù)測控制和基于神經(jīng)網(wǎng)絡(luò)的控制方法能有效提高制導(dǎo)精度,但其計算量較大,器載計算機處理能力有限,所以在實際應(yīng)用中存在一定的困難。由于滑動模態(tài)控制(Sliding Mode Control, SMC)對系統(tǒng)參數(shù)變化及擾動不靈敏,具有較強的魯棒性,而且控制器設(shè)計簡單易于實現(xiàn),可靠性較高,因此,本文基于反饋線性化的思想,設(shè)計了一種將滑?？刂婆c擴張狀態(tài)觀測器(Extended State Observer, ESO)結(jié)合的方法來得出跟蹤制導(dǎo)律,可利用觀測器估計出系統(tǒng)模型不確定性導(dǎo)致的總體誤差,并對控制系統(tǒng)實時進行補償,以得到更加精確有效的控制量,減小跟蹤誤差。

1 系統(tǒng)建模

1.1 動力學(xué)建模

本文將火星探測器看作一個質(zhì)點,在整個大氣進入過程中,只受到氣動力和火星引力的作用,且探測器的質(zhì)量恒定不變。本文采用雙錐體探測器外部構(gòu)型[13],使質(zhì)心偏離對稱中心軸線來獲得升力,通過改變升力的大小或調(diào)整升力的方向來控制探測器的飛行軌跡。探測器的受力示意圖如圖1所示,O為探測器質(zhì)心,R為作用在探測器上的氣動力,可分解為氣動升力L和氣動阻力D。G為火星的引力。飛行過程中配平攻角為15°,并假設(shè)無側(cè)滑角。σ為滾轉(zhuǎn)角,即探測器關(guān)于速度矢量的轉(zhuǎn)動角,Lx(Lsinσ)為升力沿著橫向的分量,Ly(Lcosσ)為縱向分量。通過改變探測器滾轉(zhuǎn)角的大小可以調(diào)整縱向分量的大小,從而改變其縱向運動。滾轉(zhuǎn)角的符號并不影響縱向分量的大小,故可通過改變滾轉(zhuǎn)角的符號對探測器的橫向運動進行控制。

(a)側(cè)面

(b)背面圖1 探測器受力示意圖Fig.1 Force acting on the probe

為了簡化運動方程計算,將火星看作一個標(biāo)準(zhǔn)球體,且不考慮其自轉(zhuǎn)及表面風(fēng)力的影響?；谌杂啥鹊幕鹦谴髿膺M入段動力學(xué)模型如下式

(1)

式中，r為探測器質(zhì)心到火星質(zhì)心的距離,V為探測器速度,φ為緯度,θ為經(jīng)度。γ是飛行路徑角,即速度和當(dāng)?shù)厮矫娴膴A角。ψ是飛行方向角,即速度在水平面上的投影和正東方向之間的夾角。g為重力加速度,為了便于計算,L和D重新定義為氣動升力加速度和阻力加速度,可由式(2)計算得到,

(2)

其中m=2802kg為探測器質(zhì)量,S=15.9m2為探測器參考面積。升力系數(shù)CL和阻力系數(shù)CD為氣動系數(shù),升阻比CL/CD小于0.25。本文采用簡化的指數(shù)大氣密度模型,

(3)

其中ρ為高度h=r-r0處的大氣密度,r0=3396.2km為火星半徑,hs為火星高度參數(shù),取9354m,ρ0為火星表面的標(biāo)準(zhǔn)大氣密度,取0.0158kg/m3。

1.2 動力學(xué)分析

因為在進入大氣過程中,探測器處于無動力飛行狀態(tài),滿足能量守恒定律,其能量可以由下式表示

(4)

(5)

其中Ei為初始能量,Ef為終端能量。標(biāo)準(zhǔn)能量取值范圍為0到1。

由航程的定義知,當(dāng)以時間為變量時,探測器航程只與速度有關(guān)。而由能量變化率可知,當(dāng)以能量作為自變量時,探測器的航程由阻力加速度的大小唯一確定,所以只要能夠跟蹤上探測器的阻力加速度,就能跟蹤其軌跡路線,這就將軌跡跟蹤的問題轉(zhuǎn)化為了阻力加速度跟蹤的問題。故本文以阻力加速度作為跟蹤變量,以標(biāo)準(zhǔn)能量作為自變量進行探測器運動狀態(tài)的分析。

2 制導(dǎo)律設(shè)計

進入制導(dǎo)策略包括縱向制導(dǎo)和橫向制導(dǎo)。由于假設(shè)攻角為配平攻角,故兩個通道均是將滾轉(zhuǎn)角作為唯一的控制變量。為了實現(xiàn)縱向制導(dǎo)和橫向制導(dǎo)的解耦,分別采用控制滾轉(zhuǎn)角大小和符號的策略。滾轉(zhuǎn)角的大小直接與探測器升力在縱平面的投影相關(guān),而其符號則會改變探測器的航向,所以縱向制導(dǎo)策略就是通過控制滾轉(zhuǎn)角的大小來調(diào)整縱程,橫向制導(dǎo)策略則是通過改變滾轉(zhuǎn)角的符號來修正探測器的橫程偏差。

2.1 縱向制導(dǎo)律設(shè)計

根據(jù)前文可知,縱向制導(dǎo)采用跟蹤阻力加速度的方法,因控制變量為滾轉(zhuǎn)角,首先要得到阻力加速度和滾轉(zhuǎn)角的顯性關(guān)系。由式(1)和式(2)計算出阻力加速度的一階導(dǎo),如下式所示

(6)

繼續(xù)對上式求導(dǎo),得阻力加速度二階導(dǎo)和滾轉(zhuǎn)角的關(guān)系式

(7)

其中u=cosσ,a和b是各狀態(tài)量的非線性函數(shù)形式,如下式

(8)

(9)

由式(7)可知,對阻力加速度求二階導(dǎo)后,表達(dá)式中顯含控制量u,所以系統(tǒng)相對階數(shù)為2,可用下式表示

(10)

其中D的一階和二階導(dǎo)數(shù)為狀態(tài)變量,輸入變量為u,輸出變量為D。

考慮到在實際大氣進入時存在模型誤差,記氣動系數(shù)誤差為ΔCL、ΔCD,大氣密度誤差為Δρ,則實際阻力加速度和升力加速度可由下式得到

(11)

(12)

將Ds,Ls代入表達(dá)式(8)和(9)中得非線性函數(shù)a和b的實際值as=a+Δa,bs=b+Δb。則式(10)可寫成

=a+Δa+(b+Δb)u=a+bu+Δa+Δbu

(13)

由于文獻[14]中只將as作為不確定量,沒有考慮第二項bsu的影響,所以本文更進一步地令受控系統(tǒng)的總不確定量為Δa+Δbu=w(t),則控制系統(tǒng)可寫成

(14)

(15)

對系統(tǒng)(15)利用直接反饋線性化理論,其基本思路就是選擇虛擬控制量,從而抵消原系統(tǒng)中的非線性因素, 使系統(tǒng)實現(xiàn)線性化。把式(14)右邊用一個函數(shù)N代表,即

a+bu+w(t)=N(t)

(16)

那么相對于輸入量N,非線性系統(tǒng)(15)就轉(zhuǎn)化成了一個線性受控系統(tǒng)

(17)

N就稱為原系統(tǒng)的虛擬控制輸入量,可按照線性系統(tǒng)控制律設(shè)計方法求出。

控制系統(tǒng)的目的是調(diào)整控制量u使輸出變量阻力加速度Ds能夠跟蹤上期望的阻力加速度Dr,定義跟蹤誤差

(18)

(19)

對于系統(tǒng)(19),定義滑模面

s=c1e1+c2e2

(20)

利用二次型性能指標(biāo)最優(yōu)法[15]選取滑模面系數(shù)c1、c2。首先將方程寫成如下狀態(tài)空間形式

(21)

可知系統(tǒng)滿足可控性。式(21)的二次型最優(yōu)性能指標(biāo)為

(22)

eTQe=e1Q11e1+e2Q21e1+e1Q12e2+e2Q22e2

(23)

(24)

則性能指標(biāo)變?yōu)?/p>

(25)

且有如下方程

(26)

其中

(27)

故可得

(A12P+Q12)e1+Q22e2=0

(28)

由此可取滑模面系數(shù)c1=A12P+Q12、c2=Q22。選取合適的Q陣,得c1=0.05、c2=1。

為改善滑?？刂葡到y(tǒng)到達(dá)段的品質(zhì),采用如下趨近律

(29)

選取合適的參數(shù)α可保證當(dāng)系統(tǒng)狀態(tài)遠(yuǎn)離滑動模態(tài)即滑模面s時,能以較大速度趨近于滑模面,系統(tǒng)狀態(tài)趨近于滑模面時,能保持較小的控制增益,有效降低抖振。這里其值取為0.85。參數(shù)k取為0.038。

下面證明控制律的穩(wěn)定性。取李雅普諾夫函數(shù)Vs=s2/2,顯然在原點鄰域內(nèi)Vs>0,且

(30)

由式(20)和式(29)可得出,

=-k|s|αsgn(s)

(31)

結(jié)合式(17)可得虛擬控制輸入量N的表達(dá)式

(32)

代入式(16)可得控制律

(33)

式中a和b由式(8)和(9)計算得到。

對于式(33)中的未知不確定量w(t),本文通過設(shè)置一個擴張狀態(tài)觀測器對其數(shù)值大小進行估計,實現(xiàn)對控制系統(tǒng)的動態(tài)補償。將不確定量看作一個未知狀態(tài)變量x3,其導(dǎo)數(shù)為h(t),則控制系統(tǒng)(15)變換為

(34)

由于傳統(tǒng)的利用非線性函數(shù)設(shè)計的ESO參數(shù)較多,且其整定大多依靠經(jīng)驗,故本文采用只有一個可調(diào)參數(shù)的線性擴張狀態(tài)觀測器(Linear Extended State Observer, LESO),利于實際工程應(yīng)用[16]。其三階形式如下

(35)

其中z1和z2為對狀態(tài)變量x1和x2的估計量,z3為對不確定量w(t)的估計量。eo1為觀測誤差,l1,l2,l3為觀測器增益。下面給出增益系數(shù)的整定方法,令式(35)減去式(34)得,

(36)

令狀態(tài)誤差向量eo=[eo1,eo2,eo3]T,則有

(37)

如果矩陣Ao的特征值全部在左半復(fù)平面,且h(t)是有界的,則可使系統(tǒng)(36)均對原點穩(wěn)定,LESO就能跟蹤式(34)中的擴張狀態(tài)。這里不確定量w(t)及其導(dǎo)數(shù)h(t)都是有界的。將Ao的特征值全部設(shè)為-ωo,則其特征方程與期望的特征多項式相等,即

λo(s)=|sI-Ao|=(s+ωo)3

(38)

最后將利用LESO得到的對未知不確定量的估計值z3代入控制輸入u中,即引入相應(yīng)補償,抑制模型參數(shù)不確定變化導(dǎo)致的影響,最終得到縱向制導(dǎo)律為

(39)

2.2 橫向制導(dǎo)律設(shè)計

由第一章動力學(xué)模型可知,通過改變滾轉(zhuǎn)角符號,Lx(Lsinσ)的符號也隨之反轉(zhuǎn),探測器的升力指向縱平面的另一側(cè),這樣探測器的航向角速度便會反向。本文橫向制導(dǎo)律的設(shè)計就是根據(jù)這個特點,設(shè)計一個與探測器橫向運動有關(guān)的邊界區(qū)間,使探測器在此區(qū)間范圍內(nèi)飛行,當(dāng)橫程超出邊界時,滾轉(zhuǎn)角反號,調(diào)整橫向運動,使得探測器到達(dá)開傘點的橫向距離偏差盡可能小。本文設(shè)計的橫向制導(dǎo)律如下

(40)

(41)

3 仿真結(jié)果及分析

本節(jié)采用第一節(jié)給出的結(jié)構(gòu)參數(shù)及動力學(xué)模型(1)進行仿真計算。仿真的初始條件和開傘點參數(shù)參考文獻[18]。具體數(shù)據(jù)如下表1和表2所示。

表1 仿真初始參數(shù)

表2 終端參數(shù)

根據(jù)開傘之后著陸段的初始狀態(tài)要求,當(dāng)探測器速度達(dá)到400m/s或飛行高度低于8km時,探測器降落傘打開,此時大氣進入段制導(dǎo)過程結(jié)束。

將本文基于反饋線性化理論的滑?？刂婆cLESO結(jié)合的方法與文獻[18]給出的反饋線性化方法進行仿真對比,兩種方法的橫向制導(dǎo)均采用本文設(shè)計的制導(dǎo)律。首先在無模型不確定性情況下分別用兩種方法進行仿真。如圖2所示為阻力加速度跟蹤曲線,圖3為探測器經(jīng)緯度曲線,根據(jù)具體結(jié)果可計算得反饋線性化方法的最終開傘點位置偏差為770m,本文方法的偏差為451m,這表明在不存在模型誤差的情況下,兩種方法均能實現(xiàn)對參考軌跡的良好跟蹤。

圖2 阻力加速度曲線Fig.2 Drag acceleration profile

圖3 經(jīng)緯度曲線 Fig.3 Latitude and longitude profile

由于在實際火星大氣進入段存在著模型不確定性,對制導(dǎo)精度有重要影響,所以通過人為施加誤差,將兩種方法進行仿真對比,來驗證本文設(shè)計的縱向制導(dǎo)律有更好的魯棒性。本文假設(shè)的誤差如下表3所示。

表3 模型誤差

在存在模型誤差即氣動系數(shù)誤差及大氣密度誤差的情況下,利用上述兩種方法對標(biāo)準(zhǔn)軌跡進行跟蹤。由于事先已經(jīng)知道施加的參數(shù)誤差大小,故可通過模型計算出總體誤差量的大小,與實際對不確定量w(t)的估計值進行對比來評估觀測器的觀測效果。

首先考慮氣動正拉偏的情況,從阻力加速度跟蹤結(jié)果圖4看出,在整個進入階段,本文方法設(shè)計的縱向制導(dǎo)律對阻力加速度的跟蹤效果優(yōu)于反饋線性化方法。仿真結(jié)果圖5和圖6分別是滾轉(zhuǎn)角變化曲線和漏斗邊界,從圖中能夠看出,在橫 向制導(dǎo)律的作用下,探測器的滾轉(zhuǎn)角在大氣進入過程中進行了數(shù)次變號,調(diào)整橫向運動方向,使得探測器在漏斗區(qū)間內(nèi)飛行,能夠保證最終橫程誤差在很小范圍內(nèi)。探測器飛行高度變化曲線如圖7所示,兩種方法的最終開傘點位置的高度均大于8km,滿足終端條件約束。圖8所示的經(jīng)緯度曲線表明,相比于反饋線性化方法,利用SMC結(jié)合LESO方法的探測器最終開傘點經(jīng)緯度位置離期望點更近,精度明顯更高。圖9為LESO對未知不確定量的估計曲線,與實際誤差量進行對比,能夠看出觀測器變量z3的估計誤差較小,觀測器表現(xiàn)出了良好的性能。

圖4 阻力加速度曲線Fig.4 Drag acceleration profile

圖5 滾轉(zhuǎn)角曲線Fig.5 Bank angle profile

圖6 漏斗邊界Fig.6 Funnel boundary

圖7 高度曲線Fig.7 Altitude profile

圖8 經(jīng)緯度曲線 Fig.8 Latitude and longitude profile

圖9 LESO估計曲線Fig.9 LESO estimation profile

在氣動負(fù)拉偏情況下的阻力加速度跟蹤效果如圖10所示。可以看出,相比反饋線性化方法,本文方法的跟蹤效果同樣更好。

圖10 阻力加速度曲線Fig.10 Drag acceleration profile

存在模型誤差情況下,開傘點位置的具體仿真結(jié)果對比如表4所示。

表4 仿真結(jié)果

以上仿真結(jié)果的對比表明,在存在模型誤差的情況下,反饋線性化方法對參數(shù)變化比較敏感,制導(dǎo)精度不高。本文設(shè)計的方法在縱向制導(dǎo)方面具有更優(yōu)越的抗不確定性干擾能力,結(jié)合橫向制導(dǎo)律的作用,表現(xiàn)出很好的軌跡跟蹤效果,使得開傘點誤差保持在較小的范圍。

4 結(jié) 論

本文針對火星探測器大氣進入段,利用直接反饋線性化理論對跟蹤模型進行線性化,設(shè)計了一種滑動模態(tài)控制結(jié)合線性擴張狀態(tài)觀測器的縱向跟蹤制導(dǎo)律,利用滑模控制方法魯棒性較高的特點以及通過對不確定量進行實時補償來提高制導(dǎo)精度。此外,還給出了橫向制導(dǎo)律以減小探測器的橫向偏差。在存在氣動參數(shù)和大氣密度誤差的情況下進行仿真,并與反饋線性化方法進行比較,仿真結(jié)果表明,觀測器實現(xiàn)了對不確定量的良好估計并實時對滑?？刂坡蛇M行精確補償,增強了制導(dǎo)律對模型參數(shù)不確定性的魯棒性,使得開傘點位置精度更高。故所設(shè)計的制導(dǎo)方法有效地改善了探測器大氣進入階段的制導(dǎo)性能,從而減小了開傘點偏差,并為探測器最終安全著陸到預(yù)期位置提供了保障。