基于混沌多項(xiàng)式的指令魯棒優(yōu)化及在飛行控制中的應(yīng)用

曹 瑞，沈海東，劉燕斌 ，陸宇平

(1.南京航空航天大學(xué)自動(dòng)化學(xué)院，江蘇南京 211106;2.南京航空航天大學(xué)航天學(xué)院，江蘇南京 210016)

1 引言

在實(shí)踐中，通常會(huì)存在由于測(cè)量誤差、缺乏相關(guān)的操縱知識(shí)、或者無(wú)法達(dá)到某些參數(shù)的測(cè)量條件使得參數(shù)無(wú)法測(cè)量的情況.特別是對(duì)于行星著陸器和高超聲速飛行器，由于環(huán)境未知以及系統(tǒng)間耦合無(wú)法完全表述等原因，導(dǎo)致系統(tǒng)不完全可知.如果不考慮系統(tǒng)存在的不確定性進(jìn)行軌跡和控制設(shè)計(jì)可能對(duì)任務(wù)產(chǎn)生不利后果.因此，必須在計(jì)算過(guò)程中對(duì)不確定性進(jìn)行適當(dāng)?shù)奶幚?，以便能夠充分了解參?shù)、初始條件或邊界條件中誤差或不確定帶來(lái)的影響.針對(duì)隨機(jī)系統(tǒng)，為了獲得可靠的預(yù)測(cè)結(jié)果或使預(yù)測(cè)誤差最小，人們付出大量的努力研究不確定系統(tǒng)分析和預(yù)測(cè)算法.對(duì)于非線性動(dòng)力系統(tǒng)，蒙特卡羅(Monte-Carlo，MC)[1]、馬爾科夫鏈[2]、局部線性化、無(wú)跡卡爾曼濾波[3]、Fokker-Planck-Kolmogorov框架[4]、混沌多項(xiàng)式和細(xì)胞間映射[5]等方法可以用于量化由隨機(jī)變量引起的狀態(tài)不確定性.

其中，混沌多項(xiàng)式(polynomial chaos，PC)是最有前途的不確定量化方法之一.與MC方法相比，PC 方法具有計(jì)算量小、精度高等優(yōu)點(diǎn).MC方法因其相對(duì)簡(jiǎn)單而廣受歡迎，通常是分析此類(lèi)隨機(jī)問(wèn)題優(yōu)先考慮的方法.但是MC方法為獲得高精度計(jì)算結(jié)果，所需要的重復(fù)過(guò)程會(huì)帶來(lái)非常大的計(jì)算負(fù)擔(dān).此外，MC方法在計(jì)算上是不可伸縮的，并且可能出現(xiàn)統(tǒng)計(jì)上的不一致.PC方法提供的從隨機(jī)問(wèn)題到確定性問(wèn)題的轉(zhuǎn)換方案，由于其準(zhǔn)確性和計(jì)算效率，在隨機(jī)系統(tǒng)軌跡優(yōu)化的背景下十分具有吸引力.Norbert Wiener[6]提出了Hermite齊次混沌理論，其基于隨機(jī)變量的概率分布特性，結(jié)合正交多項(xiàng)式將不確定方程轉(zhuǎn)化為一組增廣的確定方程.之后，Ghanem和Spanos[7]對(duì)該理論進(jìn)行改進(jìn)，克服了應(yīng)用于非高斯隨機(jī)系統(tǒng)時(shí)先驗(yàn)的收斂性問(wèn)題，并將其應(yīng)用于許多實(shí)際系統(tǒng)的不確定性量化.研究表明，只有根據(jù)隨機(jī)變量的概率分布選擇適當(dāng)?shù)恼欢囗?xiàng)式，才能實(shí)現(xiàn)完全收斂[8].選擇合適的正交基之后，需要使用Galerkin投影[9]或隨機(jī)配置[10]進(jìn)一步解決原隨機(jī)系統(tǒng)的不確定量化問(wèn)題.Galerkin投影法可以采用比隨機(jī)配置方法更少的方程數(shù)目，實(shí)現(xiàn)對(duì)隨機(jī)模型的準(zhǔn)確預(yù)測(cè)[11]，因此被廣泛使用.Galerkin投影法利用投影原理，最小化基于混沌多項(xiàng)式的估計(jì)系統(tǒng)和原隨機(jī)系統(tǒng)之間的誤差，從而產(chǎn)生了一組增廣的確定方程.之后，可以使用歐拉方程或龍格-庫(kù)塔方程對(duì)增廣系統(tǒng)進(jìn)行求解，進(jìn)而預(yù)測(cè)不確定系統(tǒng)的運(yùn)動(dòng)學(xué)特性.然而，不確定量化的實(shí)現(xiàn)過(guò)程并不像Galerkin投影那樣簡(jiǎn)單明了.這對(duì)于復(fù)雜的非線性應(yīng)用來(lái)說(shuō)仍是一個(gè)具有挑戰(zhàn)性的過(guò)程[12].混沌多項(xiàng)式方法最早由Wiener[6]引入，他使用Hermite多項(xiàng)式對(duì)具有高斯分布變量的隨機(jī)過(guò)程進(jìn)行模擬.根據(jù)Cameron和Martin[13]研究得知，基于混沌多項(xiàng)式的擴(kuò)展在L2意義下可以收斂于任何具有有限方差的隨機(jī)過(guò)程，這適用于大多數(shù)的物理系統(tǒng).Liu和Zhang[14]使用來(lái)自Askey方案[15]的正交多項(xiàng)式將Martin的結(jié)果推廣到各種連續(xù)和離散的分布中，并證明L2在相應(yīng)的Hilbert泛函空間中收斂.通常將其稱(chēng)為廣義PC(generalized PC，gPC)方法.gPC方法已經(jīng)應(yīng)用于各個(gè)學(xué)科之中，包括飛行器軌跡優(yōu)化、翼型設(shè)計(jì)[16-17]、代理模型[18]和參數(shù)敏感性研究[19]等.

雖然gPC方法廣泛應(yīng)用于系統(tǒng)分析和軌跡優(yōu)化之中，但很少將其應(yīng)用于飛行器控制指令優(yōu)化設(shè)計(jì).在控制環(huán)境下，gPC 應(yīng)用的發(fā)展主要?dú)w功于Nagy 和Braatz，他們?cè)诠I(yè)應(yīng)用中證明了參數(shù)不確定對(duì)非線性系統(tǒng)的影響[20].但是，工業(yè)中的應(yīng)用是一個(gè)驗(yàn)證技術(shù)，并沒(méi)有對(duì)穩(wěn)定性進(jìn)行分析.文獻(xiàn)[21]將gPC應(yīng)用于簡(jiǎn)單的懸停和飛行雙線性系統(tǒng)中，主要用于分析系統(tǒng)在不確定下的穩(wěn)定性.2010 年，Prabhakar等人[22]首次使用gPC演示了一種用于高超聲速飛行器的新型不確定性傳播框架.之后，Dutta和Bhattacharya[23]對(duì)隨機(jī)過(guò)程的先驗(yàn)概率密度函數(shù)(probability density function，PDF)進(jìn)行貝葉斯估計(jì)，從而對(duì)這項(xiàng)工作進(jìn)行了擴(kuò)展.與線性估計(jì)器相比，貝葉斯估計(jì)器和gPC的組合在非線性高超聲速飛行器不確定性分析問(wèn)題上取得了很好的表現(xiàn).在gPC隨機(jī)框架中，目標(biāo)函數(shù)如何設(shè)計(jì)十分重要，通常采用隨機(jī)參數(shù)的期望或方差作為目標(biāo)函數(shù).Bhattacharya和Fisher[24]在2011年對(duì)此進(jìn)行了研究，他們證明了gPC隨機(jī)框架中的目標(biāo)函數(shù)等價(jià)于標(biāo)準(zhǔn)的二次目標(biāo)函數(shù).Xiong等人[25]也證明了這一點(diǎn)，并且他們將這些優(yōu)化原則應(yīng)用于超敏問(wèn)題和范德波爾問(wèn)題.上述例子均是對(duì)單變量隨機(jī)系統(tǒng)(參數(shù)或初始條件不確定)的分析或設(shè)計(jì)，雖然Xiong等人[26]對(duì)系統(tǒng)存在多元、混合分布隨機(jī)變量的問(wèn)題進(jìn)行了討論，但不在優(yōu)化的范圍內(nèi).

因此，本文基于gPC方法對(duì)多維、混合分布隨機(jī)系統(tǒng)分析和飛行器控制指令優(yōu)化設(shè)計(jì)應(yīng)用方面進(jìn)行研究.通過(guò)gPC方法將隨機(jī)系統(tǒng)轉(zhuǎn)換為增廣的確定系統(tǒng)，并設(shè)計(jì)合適的目標(biāo)函數(shù)，生成一個(gè)增廣的最優(yōu)控制問(wèn)題(augmented optimal control problem，AOCP).通 過(guò)求解該最優(yōu)控制問(wèn)題，從而獲得具有魯棒特性的優(yōu)化控制指令.優(yōu)化求解方法通常包括直接和間接兩種方法.間接方法利用龐特里亞金原理，從最小化哈密頓能量的角度來(lái)闡述問(wèn)題.直接方法將最優(yōu)控制問(wèn)題轉(zhuǎn)化為有限維非線性規(guī)劃問(wèn)題[27].與間接法相比，直接法易于使用，因此被廣泛采用.本文利用gPC方法將多維、混合分布隨機(jī)系統(tǒng)分析，不確定性預(yù)測(cè)分析與優(yōu)化控制相結(jié)合，在考慮結(jié)構(gòu)不確定性的情況下獲得魯棒性更強(qiáng)的優(yōu)化控制方案，擴(kuò)展了前人在該領(lǐng)域的工作.

本文首先介紹混沌多項(xiàng)式擴(kuò)展方法，描述如何通過(guò)混沌多項(xiàng)式對(duì)隨機(jī)信息進(jìn)行數(shù)學(xué)表述;并利用Galerkin投影法將隨機(jī)變量的混沌多項(xiàng)式引入常微分方程中.然后，將隨機(jī)變量的均值和方差考慮至優(yōu)化問(wèn)題的成本函數(shù)中，并利用偽譜方法對(duì)控制指令進(jìn)行魯棒優(yōu)化.最后，通過(guò)仿真驗(yàn)證了本文所提方法在控制指令設(shè)計(jì)上的有效性，并對(duì)存在一維隨機(jī)變量和多維、混合分布隨機(jī)變量的系統(tǒng)進(jìn)行分析.仿真結(jié)果表明該方法計(jì)算效率高、精度高，并且生成的控制指令對(duì)存在參數(shù)或初始條件不確定的系統(tǒng)具有很好的魯棒性.

2 基于混沌多項(xiàng)式的控制指令優(yōu)化方法

在本節(jié)中，基于混沌多項(xiàng)式擴(kuò)展法和Galerkin投影法，對(duì)隨機(jī)動(dòng)力學(xué)系統(tǒng)轉(zhuǎn)換為高維空間中的等效確定動(dòng)力學(xué)系統(tǒng)的過(guò)程進(jìn)行介紹.并將上述方法應(yīng)用于最優(yōu)控制，生成考慮系統(tǒng)不確定性的魯棒優(yōu)化控制方案.

2.1 隨機(jī)系統(tǒng)的混沌多項(xiàng)式擴(kuò)展方法

2.1.1 混沌多項(xiàng)式

其中:ω為隨機(jī)事件，φi[Δ(ω)]表示隨機(jī)變量Δ(ω)的第i個(gè)gPC基函數(shù).函數(shù)φ是L2(Ω，F(xiàn)，)上的正交基，滿(mǎn)足式(2)所示的關(guān)系式:

其中:δij=0，i/=j;δij=1，i=j;E[·]表示概率測(cè)度d(ω)=ρ[Δ(ω)]dω的期望，其中ρ[Δ(ω)]為概率密度函數(shù).之后，采用Δ代表Δ(ω).對(duì)于隨機(jī)變量Δ，可以根據(jù)其概率分布選擇合適的正交基函數(shù)φ，使其權(quán)函數(shù)的形式與概率密度函數(shù)ρ(Δ)相同.這些正交多項(xiàng)式屬于Askey多項(xiàng)式[15];在希爾伯特空間中，由它們構(gòu)成了完備基.表1總結(jié)了Δ概率分布函數(shù)與相對(duì)應(yīng)多項(xiàng)式之間的關(guān)系[28].

表1 基于隨機(jī)變量PDF的正交多項(xiàng)式選取Table 1 Orthogonal polynomial selection based on PDF of stochastic variable

2.1.2 基于Galerkin投影引入隨機(jī)信息

本節(jié)描述了一個(gè)形式為˙x=f(x，Δ)的隨機(jī)系統(tǒng)采用PC框架轉(zhuǎn)換為高維空間中等效確定系統(tǒng)的過(guò)程，其中:x∈Rn，Δ ∈Rd，n為系統(tǒng)維數(shù)，d為隨機(jī)變量數(shù)目.假設(shè)隨機(jī)微分方程的解為x(t，Δ)，則對(duì)于n階系統(tǒng)，x ∈Rn中每個(gè)分量可近似表達(dá)為式(3)所示形式:

通過(guò)執(zhí)行Galerkin投影，將誤差el設(shè)置為零，可以獲得最優(yōu)的動(dòng)力學(xué)系統(tǒng)近似方程

其中:DΔ為Δ的定義域;ρ為隨機(jī)變量Δ的概率密度分布函數(shù).結(jié)合式(4)-(5)獲得的n(P +1)維確定性常微分方程，然后通過(guò)數(shù)值求解得到隨機(jī)系統(tǒng)的近似響應(yīng).

基于上述流程實(shí)現(xiàn)了將Rn中的隨機(jī)動(dòng)力學(xué)方程轉(zhuǎn)化為Rn(P+1)中的確定動(dòng)力學(xué)方程.式(5)中正交多項(xiàng)式φ的階數(shù)P由隨機(jī)變量Δ的維數(shù)d和運(yùn)動(dòng)方程階數(shù)r決定，滿(mǎn)足如下條件:

該表達(dá)式給出了具有d個(gè)隨機(jī)變量，最高階為r的隨機(jī)系統(tǒng)進(jìn)行不確定量化所需的混沌正交多項(xiàng)式項(xiàng)數(shù).下面以一階線性系統(tǒng)(6)為例，對(duì)這一理論的應(yīng)用進(jìn)行展示:

其中a為均勻分布于[0，1]的隨機(jī)變量.

根據(jù)表1可以看出，隨機(jī)變量a概率分布所對(duì)應(yīng)的混沌正交基函數(shù)是勒讓德(Legendre)多項(xiàng)式.由于勒讓德多項(xiàng)式的定義域?yàn)閇-1，1]，因此，需要將隨機(jī)變量a轉(zhuǎn)換為a=-0.5+0.5Δ的形式，且Δ ∈[-1，1];則隨機(jī)變量a的混沌多項(xiàng)展開(kāi)式如式(7)所示:

其中:φi為隨機(jī)變量a對(duì)應(yīng)的第i階正交基函數(shù);ai為第i階正交基函數(shù)對(duì)應(yīng)的系數(shù).

在一階線性系統(tǒng)中，式(3)和式(4)的具體表達(dá)形式分別如式(8)和式(9)所示:

其中:φj為的第j階正交基函數(shù);xj為第j階正交基函數(shù)對(duì)應(yīng)的系數(shù);φk為的第k階正交基函數(shù);為第k階正交基函數(shù)對(duì)應(yīng)的系數(shù).

將e(t，Δ)投影到φk上并將其設(shè)為0，可以得到式(10):

根據(jù)式(10)可以得到

對(duì)式(11)進(jìn)一步展開(kāi)可以獲得

由于φ是正交函數(shù)，所以式(12)可以進(jìn)一步簡(jiǎn)化為如下形式:

對(duì)式(13)進(jìn)行改寫(xiě)，將其轉(zhuǎn)換為如下表達(dá)式:

矩陣Apc∈R(P+1)×(P+1)由式(13)確定，并且Xpc=[x0x1··· xP]T為不確定量化后增廣矩陣多對(duì)應(yīng)的狀態(tài)量.

上述過(guò)程完成了隨機(jī)方程x(t，Δ)∈R向確定方程Xpc(t)∈R(P+1)的轉(zhuǎn)化，進(jìn)而可以采用常微分方程的標(biāo)準(zhǔn)數(shù)值算法進(jìn)行求解.文獻(xiàn)[29]推導(dǎo)了線性動(dòng)力學(xué)x(t，Δ)∈R的廣義表達(dá)式，并證明了這種變換在動(dòng)力學(xué)中也保持線性關(guān)系.

2.1.3 多維不確定性描述

考慮系統(tǒng)存在多維、混合分布隨機(jī)變量的情況，如果它們是相互獨(dú)立的，則PDF為每個(gè)隨機(jī)變量概率密度分布函數(shù)的乘積.因此，即使在混合分布情況下，也可以很好地捕捉系統(tǒng)特性.下面對(duì)多維、混合分布隨機(jī)系統(tǒng)的混沌多項(xiàng)展開(kāi)式(polynomial chaos expansion，PCE)建立過(guò)程進(jìn)行簡(jiǎn)單介紹.假設(shè)隨機(jī)變量初始的PCE系數(shù)a由均值μa和方差δa組成.對(duì)于僅包含參數(shù)不確定性的問(wèn)題，隨機(jī)參數(shù)的均值和標(biāo)準(zhǔn)差分別置于第1項(xiàng)和第2項(xiàng)，其他值均為零.并且，其PCE系數(shù)的維數(shù)為P +1，則a=[μaδa0···0]P+1.但是，在多維的情況下，情況略有不同(下面將詳細(xì)描述).若采用符號(hào)eijk表示式(12)中的三重變量積分，采用γk表示歸一化因子(Δ)，則不確定量化后的增廣系統(tǒng)可以采用式(15)進(jìn)行表示:

其中:

下面以具有兩個(gè)不確定量的三階線性系統(tǒng)為例，介紹多維、混合分布隨機(jī)系統(tǒng)的不確定量化過(guò)程.假設(shè)兩個(gè)隨機(jī)變量Δ1和Δ2的概率密度函數(shù)分別為高斯分布和均勻分布.然后，針對(duì)變量Δ1和Δ2分別應(yīng)用埃爾米特(Hermite)和勒讓德正交多項(xiàng)式，計(jì)算對(duì)應(yīng)的三重積分eijk和歸一化因子γk，其具體表達(dá)形式如下所示:

其中:s表示正交條件，滿(mǎn)足式(20)所示形式;N+表示整數(shù)集.

除了采用式(18)和式(19)所示的表達(dá)式求解eijk，也可以采用積分方法進(jìn)行求解.若僅考慮一個(gè)隨機(jī)變量，則多項(xiàng)式的階數(shù)由下標(biāo)i，j，k給出.若考慮多維隨機(jī)變量，多項(xiàng)式基函數(shù)可以通過(guò)分級(jí)順序索引方法[4]進(jìn)行建立.該方法將生成維數(shù)為[(P +1)×d]的索引數(shù)組，例如:考慮一個(gè)具有兩個(gè)不確定量(d=2)的三階系統(tǒng)(n=3)混沌多項(xiàng)式擴(kuò)展問(wèn)題.表2展示了使用上述方法所建立的多維、混合分布隨機(jī)系統(tǒng)的多項(xiàng)式正交基，其中:φ的下標(biāo)對(duì)應(yīng)基函數(shù)的階數(shù);Δ的下標(biāo)與隨機(jī)變量相對(duì)應(yīng).

表2中多級(jí)引索因子α決定了隨機(jī)變量Δ1和Δ2均值μ(Δ)和方差δ(Δ)在初始PCE系數(shù)中的位置.根據(jù)表2可得隨機(jī)變量Δ1和Δ2初始PCE系數(shù)分別為

對(duì)于隨機(jī)變量采用埃爾米特多項(xiàng)式或勒讓德多項(xiàng)式進(jìn)行展開(kāi).隨機(jī)系統(tǒng)歸一化因子γk的表達(dá)式分別如式(21)和式(22)所示:

其中γi為第i個(gè)隨機(jī)變量的歸一化因子.

表2 n=3，d=2的多級(jí)引索方法Table 2 Multi-index method for n=3，d=2

根據(jù)式(21)-(22)可知，若隨機(jī)系統(tǒng)中多維不確定變量相互獨(dú)立，則該系統(tǒng)的歸一化因子γk為各不確定變量對(duì)應(yīng)歸一化因子的乘積.然后結(jié)合方程(13)，可以獲得不確定量化后的增廣矩陣A.考慮式(23)所示的連續(xù)時(shí)間線性系統(tǒng)，其存在一個(gè)受參數(shù)b影響的控制輸入u.

LTI:

假設(shè)變量x為隨機(jī)變量(上述線性系統(tǒng)中的參數(shù)a和初始條件x0均為不確定參數(shù))，控制變量是確定值，即b是常數(shù).則對(duì)隨機(jī)系統(tǒng)(23)進(jìn)行不確定量化，獲得的增廣線性系統(tǒng)可由式(24)給出:

其中:狀態(tài)向量定義為X=[x0x1··· xp1]T，控制向量定義為u=ns為隨機(jī)系統(tǒng)維數(shù);p2=nc(P +1)，nc為控制輸入維數(shù).狀態(tài)轉(zhuǎn)換矩陣第k行，第j列元素的表達(dá)式如式(25)所示:

其中ai定義同式(7)一致.假設(shè)參數(shù)b=1，即系統(tǒng)中控制參數(shù)確定，則相應(yīng)增廣系統(tǒng)的控制矩陣為B=[1 0··· 0]T，其維數(shù)為其中PCE系數(shù)由式(24)中狀態(tài)向量X的各分量表示.若系統(tǒng)初始條件不確定，則式(24)中的狀態(tài)初始向量為X0=[μxδx0··· 0]T;若僅參數(shù)a存在不確定性，初始條件確定，則X0中僅存在μx這一非零項(xiàng).下面將介紹如何將獲得的增廣矩陣用于控制優(yōu)化設(shè)計(jì)，以充分利用已知的不確定性信息.

2.2 隨機(jī)系統(tǒng)最優(yōu)控制指令設(shè)計(jì)

對(duì)于一個(gè)隨機(jī)系統(tǒng)，為獲得考慮系統(tǒng)不確定性的最優(yōu)控制指令u*，必須對(duì)原成本函數(shù)進(jìn)行修正.對(duì)于不確定系統(tǒng)控制優(yōu)化，需要將不確定參數(shù)的信息(均值和方差)考慮到成本函數(shù)中.不確定參數(shù)的信息可以使用混沌多項(xiàng)式和Galerkin投影法產(chǎn)生增廣系統(tǒng)的狀態(tài)進(jìn)行描述.在PC框架中，均值和標(biāo)準(zhǔn)差定義為

根據(jù)上述定義可知，不確定參數(shù)的均值總是對(duì)應(yīng)于增廣系統(tǒng)的第1個(gè)狀態(tài)x0，而方差是歸一化因子與增廣系統(tǒng)其他狀態(tài)平方乘積的和.因此，考慮不確定性的增廣最優(yōu)控制問(wèn)題可以通過(guò)下式進(jìn)行描述:

其中:J0表示與隨機(jī)變量無(wú)關(guān)的性能指標(biāo);Φ[μ(x)，σ(x)]表示終端優(yōu)化形式，目的是使隨機(jī)變量在軌跡末端的變化最小;表示積分優(yōu)化形式，目的是最小化整個(gè)任務(wù)段隨機(jī)變量的均值或方差. k1和k2分別為J0和隨機(jī)變量統(tǒng)計(jì)信息的權(quán)重值.

例如，若希望最小化最終期望值，則可以將式(28)所表示的隨機(jī)系統(tǒng)期望項(xiàng)寫(xiě)入成本函數(shù)中.

其中:W表示由歸一化因子組成的權(quán)重矩陣;X表示增廣矩陣的狀態(tài)量.

若期望最小化最終方差，則可以將式(30)所表示的隨機(jī)系統(tǒng)方差項(xiàng)寫(xiě)入成本函數(shù)中.

由于不確定參數(shù)的均值總是對(duì)應(yīng)于增廣系統(tǒng)的第1個(gè)狀態(tài)，因此

故，可以將式(30)進(jìn)一步簡(jiǎn)化為

因此，可以對(duì)式(28)(32)進(jìn)行積分，以積分形式考慮在成本函數(shù)中或直接設(shè)置為終端條件.然后求解式(27)所示的優(yōu)化問(wèn)題獲得隨機(jī)系統(tǒng)的魯棒優(yōu)化控制指令信號(hào).下面將概述隨機(jī)系統(tǒng)控制指令魯棒優(yōu)化設(shè)計(jì)過(guò)程.最后通過(guò)仿真，介紹增廣系統(tǒng)的設(shè)計(jì)過(guò)程，并驗(yàn)證基于混沌多項(xiàng)式的控制指令魯棒優(yōu)化設(shè)計(jì)的有效性.

提出的隨機(jī)系統(tǒng)控制指令魯棒優(yōu)化方案可歸納為以下步驟，其實(shí)現(xiàn)流程如圖1所示，其中各隨機(jī)變量相應(yīng)的混沌正交基的初始階數(shù)由P0=(d+r)!/(d!r!)-1確定.

圖1 基于混沌多項(xiàng)式的控制指令魯棒優(yōu)化方案Fig.1 Control command robust optimization strategy based on polynomial chaos

步驟1為系統(tǒng)中各隨機(jī)變量設(shè)置相應(yīng)的混沌多項(xiàng)展開(kāi)式，并根據(jù)第2.1.3節(jié)所示內(nèi)容建立多級(jí)引索，進(jìn)而獲取多維、混合分布隨機(jī)系統(tǒng)的多項(xiàng)式正交基;

步驟2根據(jù)Galerkin投影法將隨機(jī)動(dòng)力學(xué)方程Xsto轉(zhuǎn)換為等價(jià)的增廣確定性模型XAOCP;

步驟3根據(jù)Galerkin投影法將隨機(jī)約束條件Ssto轉(zhuǎn)換為增廣的確定性等價(jià)約束SAOCP;

步驟4根據(jù)最小期望/協(xié)方差對(duì)隨機(jī)成本函數(shù)進(jìn)行量化，根據(jù)式(27)所示成本函數(shù)形式將隨機(jī)系統(tǒng)成本函數(shù)Jsto轉(zhuǎn)為增廣確定系統(tǒng)的成本函數(shù)JAOCP;

步驟5判斷增廣系統(tǒng)的運(yùn)動(dòng)學(xué)模型XAOCP以及約束條件SAOCP是否滿(mǎn)足可以近似描述隨機(jī)系統(tǒng)動(dòng)力學(xué)方程不確定邊界|Xsto|boundary和約束條件不確定邊界|Ssto|boundary的要求.若滿(mǎn)足，則進(jìn)入步驟6;若不滿(mǎn)足，則增加混沌正交基的階數(shù)并返回步驟1;

步驟6通過(guò)偽譜法對(duì)增廣系統(tǒng)的優(yōu)化問(wèn)題(27)進(jìn)行求解;

步驟7判斷成本函數(shù)中隨機(jī)量的統(tǒng)計(jì)信息Js是否小于設(shè)定值δ(人為設(shè)定).若滿(mǎn)足，則獲得具有魯棒性的最優(yōu)控制指令;否則對(duì)成本函數(shù)JAOCP中隨機(jī)信息統(tǒng)計(jì)量對(duì)應(yīng)的權(quán)重值k2進(jìn)行調(diào)整，并返回步驟6.

3 隨機(jī)系統(tǒng)魯棒優(yōu)化方法仿真示例

在此考慮一個(gè)簡(jiǎn)單非線性系統(tǒng)，不考慮系統(tǒng)存在不確定時(shí)，其標(biāo)稱(chēng)系統(tǒng)的優(yōu)化控制問(wèn)題描述如下所示:

考慮系統(tǒng)存在兩個(gè)隨機(jī)變量:參數(shù)a概率分布滿(mǎn)足正態(tài)分布a～N(2.5，0.1);初始狀態(tài)量x0的概率分布滿(mǎn)足均勻分布x0～U(0.9，1.1).根據(jù)式(33)可以看出，此非線性系統(tǒng)的控制輸入也受到隨機(jī)變量的影響.

首先考慮僅存在參數(shù)a不確定的情況，參數(shù)a的初始PCE系數(shù)為ai=[2.5 0.1 0]，增廣系統(tǒng)的初始條件為X0=[1 0 0]T，則增廣系統(tǒng)如式(34)所示.該增廣系統(tǒng)最優(yōu)控制問(wèn)題的目標(biāo)函數(shù)設(shè)置為式(35)所示結(jié)構(gòu)

在標(biāo)稱(chēng)系統(tǒng)下求解優(yōu)化控制問(wèn)題(33)所得的控制指令記為CJ;在隨機(jī)系統(tǒng)下基于混沌多項(xiàng)式魯棒優(yōu)化方法生成的控制指令記為CJAOCP，其控制輸入和對(duì)應(yīng)的增廣系統(tǒng)狀態(tài)分別如圖2(a)-2(b)所示.

根據(jù)圖2(b)可知，在控制指令CJAOCP作用下增廣系統(tǒng)狀態(tài)隨時(shí)間逐漸收斂于零.

將CJ和CJAOCP分別作用于隨機(jī)系統(tǒng)進(jìn)行1000次蒙特卡羅分析，以驗(yàn)證基于混沌多項(xiàng)式的控制指令魯棒優(yōu)化方案的有效性.圖3展示的MC仿真結(jié)果證明了本文所提方法可以有效利用已知的不確定信息，生成對(duì)不確定性有抑制作用的控制指令.

圖2 控制量和PCE系數(shù)變化曲線Fig.2 Change curve of control and PCE coefficients

圖3 一維不確定性系統(tǒng)MC分析Fig.3 MC analysis for one dimensional uncertainty system

然后考慮參數(shù)a和初始狀態(tài)x0不確定的情況，根據(jù)第2.1.3節(jié)內(nèi)容對(duì)隨機(jī)系統(tǒng)進(jìn)行不確定量化，以獲得式(36)所示的增廣矩陣;其中，參數(shù)a的初始PCE系數(shù)、增廣系統(tǒng)初始狀態(tài)和歸一化因子分別為

圖4(a)展示了二維不確定情況、一維不確定情況以及標(biāo)稱(chēng)情況下的控制曲線;圖4(b)展示了二維不確定情況下的增廣系統(tǒng)狀態(tài)變化圖.圖5展示了CJ和二維不確定系統(tǒng)優(yōu)化控制指令CJAOCP分別作用于隨機(jī)系統(tǒng)進(jìn)行1000次蒙特卡羅仿真的結(jié)果.

圖4 控制量和PCE系數(shù)變化曲線Fig.4 Change curve of control and PCE coefficients

圖5 二維不確定性系統(tǒng)MC分析Fig.5 MC analysis for two dimensional uncertainty system

根據(jù)圖3和圖5的仿真結(jié)果可知:雖然隨機(jī)初始條件會(huì)對(duì)系統(tǒng)變化帶來(lái)非常明顯影響，但是基于混沌多項(xiàng)式的魯棒優(yōu)化方法生成的控制指令仍然具有良好的不確定抑制效果.本節(jié)以非線性系統(tǒng)為例，詳細(xì)介紹基于混沌多項(xiàng)式的魯棒優(yōu)化方法對(duì)隨機(jī)系統(tǒng)控制指令的優(yōu)化設(shè)計(jì)過(guò)程，并仿真驗(yàn)證所提方案的有效性.下節(jié)將介紹本文所提方法在飛行控制中的應(yīng)用.

4 飛行控制中的分析應(yīng)用

4.1 隨機(jī)系統(tǒng)運(yùn)動(dòng)學(xué)特性預(yù)測(cè)

本節(jié)考慮F-16飛行器模型存在參數(shù)不確定性，分析本文所提方法獲取的增廣系統(tǒng)對(duì)隨機(jī)系統(tǒng)運(yùn)動(dòng)學(xué)特性的預(yù)測(cè)性能.為簡(jiǎn)單起見(jiàn)，假設(shè)系統(tǒng)參數(shù)變化依賴(lài)于單個(gè)隨機(jī)變量Δ.在這個(gè)例子中，考慮了F-16飛行器短周期近似模型，該模型由文獻(xiàn)[30]給出.

其中狀態(tài)向量為x=[α q xe]T，α是迎角，q為俯仰角速率，xe為所測(cè)量的升降舵狀態(tài).控制輸入u=δec是升降舵指令，單位為度.飛行器模型中各矩陣數(shù)值如下所示:

假設(shè)括號(hào)內(nèi)的參數(shù)是不確定的，并且其概率分布滿(mǎn)足均勻分布特點(diǎn)，與其標(biāo)稱(chēng)值有20%的偏差.這種不確定性是由于阻尼項(xiàng)Cxq不確定性所引起的;Cxq在大迎角條件下很難建立準(zhǔn)確模型.文獻(xiàn)[31]設(shè)計(jì)了一種基于俯仰角速率q反饋的頻域控制(控制律在標(biāo)稱(chēng)系統(tǒng)上進(jìn)行設(shè)計(jì))，控制器的表達(dá)形式如下所示:

式(38)表示一種俯仰角速率反饋控制器.將上述表達(dá)式轉(zhuǎn)換為狀態(tài)空間形式(Ac，Bc，Cc)并對(duì)原系統(tǒng)進(jìn)行增廣，得到閉環(huán)系統(tǒng)

圖6和圖7展示了本文所提方法獲取的增廣系統(tǒng)對(duì)隨機(jī)系統(tǒng)(39)動(dòng)力學(xué)特性的預(yù)測(cè)能力.圖6中的圓表示混沌增廣系統(tǒng)的特征值，實(shí)點(diǎn)代表對(duì)隨機(jī)系統(tǒng)(39)進(jìn)行蒙特卡羅仿真所得到的系統(tǒng)特征值.

圖7顯示了上述參數(shù)中存在±20%不確定性時(shí)系統(tǒng)的俯仰角速率響應(yīng).由PC所預(yù)測(cè)的俯仰角速率曲線和通過(guò)蒙特卡羅仿真獲得的俯仰角速率曲線分別由虛線和暗實(shí)線表示.

圖6 ±20%參數(shù)不確定性的短周期閉環(huán)特征值分布Fig.6 Closed loop eigenvalue distributions for a short-period mode for±20%parameter uncertainty

圖7 PC和MC系統(tǒng)對(duì)±20%參數(shù)不確定性的響應(yīng)Fig.7 PC and MC system response to±20%parameter uncertainty

從圖6可以得知，混沌增廣系統(tǒng)捕捉到了隨機(jī)系統(tǒng)的特征值分布.根據(jù)圖7可以得知，PC系統(tǒng)預(yù)測(cè)的邊界與MC仿真曲線邊界十分一致.此仿真結(jié)果表明:獲得的增廣系統(tǒng)可以用于預(yù)測(cè)隨機(jī)系統(tǒng)運(yùn)動(dòng)學(xué)特性以及分析隨機(jī)系統(tǒng)穩(wěn)定性，并且這種方法比MC方法更有效，計(jì)算效率更高.

4.2 控制指令的魯棒優(yōu)化設(shè)計(jì)

F-16飛行器系統(tǒng)的不確定參數(shù)滿(mǎn)足均勻分布的特點(diǎn)，與其標(biāo)稱(chēng)值存在±20%的偏差，其不確定量化后的增廣矩陣A由表示.然后，根據(jù)第2.2節(jié)提出的控制指令優(yōu)化方案對(duì)增廣后系統(tǒng)的控制指令進(jìn)行魯棒優(yōu)化設(shè)計(jì)，并與文獻(xiàn)[31]所提魯棒控制器和式(38)所示控制器作用到原不確定系統(tǒng)中進(jìn)行仿真，并對(duì)仿真結(jié)果進(jìn)行對(duì)比.

由于是對(duì)飛行控制系統(tǒng)進(jìn)行控制器設(shè)計(jì)，因此需要考慮控制輸入的飽和約束.

在此設(shè)計(jì)除去配平舵面需求外，滿(mǎn)足額外控制需求的舵面約束為δec∈[-8°，8°].根據(jù)圖7可知控制器(38)對(duì)不確定抑制能力弱，因此在圖8中僅展示了不確定系統(tǒng)在多胞魯棒控制器[31]、基于混沌多項(xiàng)式的魯棒優(yōu)化控制指令下，仿真獲得的控制輸入和根據(jù)式(27)計(jì)算所得性能指標(biāo)對(duì)比結(jié)果.

最后，將多胞魯棒控制器[31]、基于混沌多項(xiàng)式的魯棒優(yōu)化控制指令和控制器(38)作用到隨機(jī)系統(tǒng)(37)中進(jìn)行蒙特卡羅仿真，其仿真對(duì)比結(jié)果如圖9所示.然后考慮存在參數(shù)和初始狀態(tài)不確定情況下，分析基于上述控制指令的控制效果.假設(shè)俯仰角速率的隨機(jī)初始狀態(tài)滿(mǎn)足正態(tài)分布q0～N(2，0.25)，在參數(shù)和初始狀態(tài)不確定條件下的蒙特卡羅仿真和性能指標(biāo)對(duì)比結(jié)果如圖10所示.

圖8 兩種控制算法下控制輸入與性能指標(biāo)Fig.8 Control input and performance indicator of two control algorithms

圖9 單維不確定系統(tǒng)MC分析Fig.9 MC analysis for single-dimensional uncertainty system

圖10 性能指標(biāo)以及多維、混合分布不確定系統(tǒng)的MC分析Fig.10 Performance indicator and MC analysis for multi-dimensional，mixed-distribution uncertainty system

根據(jù)圖8-10的仿真結(jié)果，可以看出本文提出的方法可以充分利用系統(tǒng)中不確定參數(shù)信息，其優(yōu)化產(chǎn)生的控制指令比魯棒控制律在不確定性抑制方面具有更好的性能.根據(jù)圖8可以得出本文提出的魯棒優(yōu)化方法與文獻(xiàn)[31]中的魯棒控制器相比，可以在有效抑制系統(tǒng)不確定性的基礎(chǔ)上降低所需的控制能量.根據(jù)圖10可以得知，本文提出的方法對(duì)存在多維、混合分布隨機(jī)變量的系統(tǒng)同樣具有很好的不確定抑制能力，可以處理系統(tǒng)存在初始狀態(tài)不確定的問(wèn)題.

通過(guò)對(duì)F-16飛行器的仿真分析可以得知，該方法可以有效應(yīng)用于飛行控制中的動(dòng)力學(xué)預(yù)測(cè)以及具有魯棒特性的控制指令生成.

5 結(jié)論

本文提出一種基于混沌多項(xiàng)式的魯棒優(yōu)化方法，將混沌多項(xiàng)式方法應(yīng)用于飛行器控制設(shè)計(jì)中.本文利用混沌多項(xiàng)式和Galerkin投影法將隨機(jī)動(dòng)力學(xué)系統(tǒng)轉(zhuǎn)換為高維空間中的等效確定性動(dòng)力學(xué)系統(tǒng).之后，將隨機(jī)變量的均值和標(biāo)準(zhǔn)差考慮至優(yōu)化問(wèn)題的成本函數(shù)中，并根據(jù)偽譜法對(duì)控制指令進(jìn)行魯棒優(yōu)化.最后，通過(guò)仿真驗(yàn)證了本文所提方法可以充分利用系統(tǒng)中不確定參數(shù)的信息，生成具有良好魯棒性能的控制指令.并對(duì)存在一維隨機(jī)變量和多維、混合分布隨機(jī)變量的系統(tǒng)進(jìn)行分析.分析結(jié)果表明，該方法可以有效應(yīng)用于飛行控制中的動(dòng)力學(xué)預(yù)測(cè)，并且生成的控制指令對(duì)存在初始條件和參數(shù)不確定的系統(tǒng)具有良好的魯棒性.