基于變尺度尋優(yōu)混合求解算法的航空發(fā)動機狀態(tài)變量模型建立方法

楊征山，仇小杰，王鐵，黃金泉

（1.中航工業(yè)航空動力控制系統(tǒng)研究所，江蘇無錫 214063；2.南京航空航天大學能源與動力學院，南京 210016；3.總參陸航部駐上海地區(qū)軍事代表室，上海 200233）

0 引言

隨著航空發(fā)動機技術的發(fā)展，對發(fā)動機控制技術的研究提出了更高、更復雜的要求。在控制技術領域，包括綜合飛行和推進控制、機載自適應模型、降低裕度的邏輯、先進多變量控制技術等方面[1-5]，開展了一系列研究。先進控制算法的研究基礎，即如何建立精確的航空發(fā)動機狀態(tài)變量模型成為研究熱點[6-9]。

目前求取航空發(fā)動機狀態(tài)變量模型的算法主要包括小擾動法、穩(wěn)態(tài)終值響應法和擬合方法[10-14]。其中，小擾動方法最簡單，且易于計算，但結(jié)果精度不高；穩(wěn)態(tài)終值響應法在處理線性模型和非線性模型的終值一致問題上效果較好，但不能保證動態(tài)過程的一致性；運用擬合法求取的結(jié)果精度較高，但需要理想的初猜值，否則會出現(xiàn)擬合過程發(fā)散的情況，同時精度隨著矩陣維數(shù)的增加而降低，耗時也隨之增加。

本文將上述3種算法加以綜合，最終獲得狀態(tài)變量線性模型。

1 狀態(tài)變量模型

設航空發(fā)動機非線性模型為

式中：x∈Rn、u∈Rr、y∈Rm，分別為發(fā)動機狀態(tài)變量、控制變量、輸出量。

當發(fā)動機處于某穩(wěn)態(tài)點（x0,u0,y0）工作時，發(fā)動機模型為

式（4）即為狀態(tài)變量線性模型（SVM）。

選取航空發(fā)動機風扇轉(zhuǎn)速nF和壓氣機轉(zhuǎn)速nC為狀態(tài)變量，發(fā)動機供油量Wfm和尾噴管面積A8為控制量，風扇出口總溫Tt22、總壓Pt22、壓氣機出口總溫Tt3、燃燒室出口總壓Pt4、高壓渦輪出口總溫Tt41、總壓Pt41為模型輸出量，則式（4）可以表示為

為解決全包線狀態(tài)下系數(shù)矩陣求解繁瑣、應用不便和模型中由于物理變量的數(shù)量級不同導致系數(shù)矩陣畸形的問題，將所有參數(shù)進行相似變換并歸一化

最終經(jīng)過相似歸一化的狀態(tài)變量線性模型為

2 傳統(tǒng)系數(shù)矩陣求解方法

2.1 小擾動法

小擾動法的原理是在航空發(fā)動機某穩(wěn)態(tài)點處加一擾動，然后用非線性部件級模型進行動態(tài)計算，當動態(tài)模型的流量達到連續(xù)準平衡收斂條件后，求出狀態(tài)變量的導數(shù)和輸出變量的增量，計算其與擾動量的比值，求出各系數(shù)矩陣。

以求取A、C矩陣為例，分別對狀態(tài)向量風扇轉(zhuǎn)速nF和壓氣機轉(zhuǎn)速nC進行擾動，在擾動某一狀態(tài)向量時，其他狀態(tài)向量和控制量保持不變，即當擾動風扇轉(zhuǎn)速nF時，令ΔnCC=0、ΔWfmC=0、ΔA8C=0，則

同理擾動壓氣機轉(zhuǎn)速nC，計算得到a12、a22、c32…、c82的值，最終求得A、C矩陣。采用相同方法求B、D

矩陣，區(qū)別是擾動的量由狀態(tài)量變?yōu)榭刂屏俊?/p>

2.2 穩(wěn)態(tài)終值響應法

系統(tǒng)的特性由狀態(tài)變量線性模型的系數(shù)矩陣決定。其中模型的狀態(tài)響應特性由A、C矩陣決定，模型的最終穩(wěn)態(tài)響應由B、D矩陣決定。通過設定A、C矩陣值，根據(jù)不同輸入階躍情況下模型的穩(wěn)態(tài)輸出變化量，應用數(shù)學解析方法計算出B、D矩陣值，具體計算公式為

2.3 擬合法

擬合法的求解步驟如下。

（1）狀態(tài)變量線性模型如式（7）描述，分別求出模型在控制量Wfm和A8作小階躍時的動態(tài)響應解析式；

（2）利用前文建立的非線性部件級模型同樣對控制量Wfm和A8作小階躍，求出部件級模型的動態(tài)響應解析式；

（3）利用線性模型動態(tài)響應與部件級非線性模型動態(tài)響應一致的原則，擬合出狀態(tài)變量模型。

3 基于變尺度法的混合求解方法

本文將上述3種算法加以綜合，利用基于變尺度法的擬合法得到最優(yōu)的A、C矩陣，使得最終的狀態(tài)變量線性模型與非線性部件級模型的動態(tài)特性相吻合，A、C矩陣的初猜值由小擾動算法求得，而B、D矩陣則通過穩(wěn)態(tài)響應終值法求出并保證狀態(tài)變量線性模型與非線性部件級模型的穩(wěn)態(tài)穩(wěn)定輸出值完全一致，最終獲得狀態(tài)變量線性模型。

近年來，變尺度優(yōu)化方法在求解無約束極值問題中取得了良好效果，其主要特點是既避免了計算2階導數(shù)矩陣以及求逆過程，收斂速度又比梯度法的快，對高維問題具有更顯著的優(yōu)越性。

無約束極值問題可以表述為

在采用牛頓法求解無約束問題過程中，其優(yōu)化搜索方向為-[?2（fxk）]-1?2（fxk），為了不計算2階導數(shù)?2（fxk）以及逆矩陣，變尺度方法運用擬牛頓法構(gòu)造1個H矩陣來逼近2階導數(shù)矩陣的逆陣，變尺度方法的計算步驟如下。

（1）設置初始值x0以及梯度允許誤差ε＆gt;0；

（2）若||?（fx0）||≤ε，則x0為近似極小值，停止迭代，若否，繼續(xù)下1步；

λ得到λ0，從而得到下個近似點

（4）迭代到第k步，得到近似點xk，計算出?（fxk），若||?（fxk）||≤ε，則xk即為最優(yōu)值，停止迭代，若否，計算

式中：ΔGk-1=?f（xk）-?f（xk-1），得到和λk，計算出 xk+1=xk+λkPk；

（5）若xk+1滿足精度要求，則xk+1為最優(yōu)解，若否，轉(zhuǎn)到第（4）步。

基于變尺度法的發(fā)動機線性模型混合求解方法計算步驟如下。

（1）對發(fā)動機非線性部件級模型的控制量Wfm和A8分別作小階躍[Δu1（t），Δu2（t）]

計算部件級模型的動態(tài)響應偏差向量為Δy(t)=[ΔnFC（t），ΔnCC（t），ΔTt22C（t），ΔPt22C（t），ΔTt3C（t），ΔPt4C（t），ΔTt41C（t），ΔPt41C（t），ΔFnC（t），ΔSmC（t）]。

式中：t=0，T，…，nT，T為采樣時間；

（2）利用小擾動法求取A、C矩陣的初猜值A0、C0；

（3）利用穩(wěn)態(tài)終值響應法中的式（9）、（10）求出B0、D0，得到線性模型；

（5）將（1）和（3）中的動態(tài)響應數(shù)據(jù)Δy（t）和Δy?（t）作為樣本，采用變尺度法，優(yōu)化如式（14）所示的目標函數(shù)

即計算與部件級非線性模型響應偏差平方和最小的線性模型系數(shù)矩陣值。通過優(yōu)化，若滿足精度要求，得到線性模型系數(shù)矩陣（A，B，C，D）；若不滿足精度要求，重復第（4）、（5）步。

限于篇幅，這里僅給出地面標準大氣條件下，Wfm=1.8kg/s和A8=0.25m2工作點處的線性模型系數(shù)矩陣求解過程。

（1）對發(fā)動機非線性部件級模型的控制量Wfm和A8分別作小階躍[Δu1（t），Δu2（t）]

計算部件級模型的動態(tài)響應偏差向量為Δy（t）；

（2）利用小擾動法求取A、C矩陣的初猜值A0、C0

（3）由式（9）、（10）求出B0、D0

得到線性模型系數(shù)矩陣初猜值A0，B0，C0，D0；

（4）對此線性模型的控制量Wfm和A8分別作小階躍[Δu1（t），Δu2（t）]

（5）利用變尺度算法優(yōu)化如式（14）的目標函數(shù)，得到第（1）步優(yōu)化所得線性模型系數(shù)矩陣（A1，B1，C1，D1）

不滿足精度要求，重復第（4）、（5）步操作，直至滿足精度要求，得到線性模型的系數(shù)矩陣（A1，B1，C1，D1）

4 狀態(tài)變量模型仿真

模擬發(fā)動機非線性部件級模型[15]和狀態(tài)變量線性模型在H=0、Ma=0、Wfm=1.8kg/s和A8=0.25m2工作點上，分別作2個控制量的小階躍仿真（其中供油量Wfm階躍0.5%，尾噴管面積A8正階躍1.5%），比較2個模型的輸出響應吻合程度，檢驗基于變尺度法的混合求解方法計算出的狀態(tài)變量線性模型的精確度，仿真結(jié)果如圖1、2所示。

圖1 燃油小階躍線性模型響應與非線性模型響應

圖2 尾噴管面積小階躍線性模型響應與非線性模型響應

從仿真結(jié)果中可見，本文提出的基于變尺度法的混合求解方法計算的線性模型與非線性部件級模型不僅在動態(tài)過程響應中吻合良好，而且能夠保證發(fā)動機最終穩(wěn)定值的一致性，具有較高的精度。

5 結(jié)束語

本文研究了發(fā)動機狀態(tài)變量線性模型的求解方法，分析了傳統(tǒng)小擾動法、穩(wěn)態(tài)終值響應法和擬合法的優(yōu)缺點，提出了基于變尺度法的混合求解方法。該方法能避免小擾動方法精度不高、穩(wěn)態(tài)終值響應法動態(tài)過程不一致以及傳統(tǒng)擬合法隨著需要擬合的矩陣維數(shù)增加而精度下降、擬合時間增加的缺點。仿真結(jié)果表明，本文方法不僅在動態(tài)過程響應中吻合良好，具有較高的穩(wěn)態(tài)精度，而且能夠保證發(fā)動機最終穩(wěn)定值的一致性，具有良好的應用前景。

[1] Skira C A, Agnello M. Control system for the next century’s fighter engines [J]. Transactions of the ASME Journal of Engineering for Gas Turbines and Power, 1992，114（4）：749-754.

[2] Paseal G. A linear matrix inequality approach to H ∞control[J].International Journal of Robust and Nonlinear Colltrol，1994，4（4）：421-448.

[3] Turner M C, Walker D J, Alford G. Design and ground-base simulation of an H ∞ limited authority flight control system for the Westland Lynx helicopter [J]. Aerospace science and Tech nology，2001，5（3）：221-234.

[4] Conners T R. Thrust stands evaluation of engine performance improvement algorithms in an F-15 airplane[R]. AIAA-92-3747.

[5] Orme J, Schkolnik G. Flight assessment of the onboard propulsion system model for the performance seeking control algorithm on an F-15 aircraft [R]. AIAA-95-2361.

[6] Leibov R. Aircraft turbine engine linear with uncertain engine values [J]. IEEE Transaction on Automatic Control，2002，47（8）：1367-1369.

[7]王磊，郭迎清，陸軍，等.用改進子空間辨識法建立航空發(fā)動機模型[J].計算機仿真，2012，29（12）：75-79.

WANG Lei, GUO Yingqing, LU Jun, et al. Establishment of aeroengine model with improved subspace model identification method [J]. Computer Simulation，2012，29 （12）：75-79. (in Chinese)

[8]李秋紅，孫健國.基于遺傳算法的航空發(fā)動機狀態(tài)變量模型建立方法[J].航空動力學報，2006，21（2）：427-431.

LI Qiuhong，SUN Jianguo. Aeroengine state variable modeling based on the genetic algorithm[J]. Journal of Aerospace Power，2006，21（2）：427-431. (in Chinese)

[9] Sugiyama N. Derivation of ABCD system matrices from nonlinear dynamic simulation of jet engine[R]. AIAA-92-3391.

[10]鄭鐵軍，王曦，李秀芹，等.建立航空發(fā)動機狀態(tài)空間模型的修正方法[J].推進技術，2005，26（1）：46-49.

ZHENG Tiejun, WANG Xi, LI Xiuqin, et al. Modified method of establishing the state space model of aeroengine[J]. Journal of Propulsion Technology，2005，26（1）：46-49. (in Chinese)

[11]苗卓廣，謝壽生，吳勇，等.基于改進粒子群算法的航空發(fā)動機狀態(tài)變量建模[J].推進技術，2012，33（1）：73-77.

MIAO Zhuoguang, XIE Shousheng, WU Yong, et al. Aeroengine state variable modeling based on the improved particle swarm optimization [J]. Journal of Propulsion Technology，2012，33（1）：73-77. (in Chinese)

[12]陸軍，郭迎清，陳小磊.線性擬合法建立航空發(fā)動機狀態(tài)變量模型[J].航空動力學報，2011，26（5）：1172-1177.

LU Jun，GUO Yingqing，CHEN Xiaolei. Establishment of aeroengine state variable model based on linear fitting method [J]. Journal of Aerospace Power，2011，26（5）：1172-1177. (in Chinese)

[13]魯峰，黃金泉，佘云峰.航空發(fā)動機狀態(tài)變量模型的QPSO尋優(yōu)混合求解法[J].推進技術，2011，32（5）：722-727.

LU Feng, HUANG Jinquan, SHE Yunfeng. State space modeling based on QPSO hybrid method for aeroengines[J]. Journal of Propulsion Technology，2011，32 （5）：722-727. (in Chinese)

[14]張海波，楊小龍，林一暉.一種求取發(fā)動機狀態(tài)變量模型的改進擬合法[J].航空動力學報，2011，26（8）：1907-1913.

ZHANG Haibo，YANG Xiaolong，LIN Yihui. An improved method of identification for aeroengine’s state variable model [J]. Journal of Aerospace Power，2011，26（8）：1907-1913. (in Chinese)

[15]周文祥，黃金泉，竇建平，等.面向?qū)ο蟮臏u扇發(fā)動機及控制系統(tǒng)仿真平臺[J].航空動力學報，2007，22（1）：119-125.

ZHOU Wenxiang，HUANG Jinquan，DOU Jianping，et al.Object-oriented simulation platform for turbofan engine and its control system [J]. Journal of Aerospace Power,2007,22（1）：119-125. (in Chinese)