基于遺傳算法的渦噴發(fā)動機身份證模型建立

陳 煜，黃金泉，羅啟君，魯 峰

（1.南京航空航天大學(xué)能源與動力學(xué)院，南京210016；2.中航工業(yè)貴州黎陽航空發(fā)動機（集團）有限公司，貴州平壩561102）

0 引言

航空發(fā)動機模型計算精度問題受到工程建模人員廣泛關(guān)注。模型精度除了受所建氣動熱力學(xué)模型與航空發(fā)動機真實工作過程吻合程度的影響，更受制于航空發(fā)動機部件特性[1]和其他部件參數(shù)的準確性影響[2]。由于同系列的其他發(fā)動機在制造與安裝上存在公差，部件的效率和流量隨使用時間變化發(fā)生蛻化，導(dǎo)致發(fā)動機部件特性存在個體差異。因此有必要開發(fā)發(fā)動機模型的功能模塊，真實反映發(fā)動機部件的特性，解決發(fā)動機模型與真實發(fā)動機之間失配性問題，建立發(fā)動機身份證模型，即個體模型。

20世紀90代初以來，國內(nèi)外學(xué)者提出一系列航空發(fā)動機部件特性和模型參數(shù)修正方法[3-4]。Stamatis等提出1種根據(jù)發(fā)動機實際測量數(shù)據(jù)校正發(fā)動機部件特性的方法，將發(fā)動機模型的殘差向量以及待修正參數(shù)誤差的加權(quán)平方和作為代價函數(shù)，采用非線性廣義最小殘差法，解共同工作方程[5-6]；Kong Changduk等提出基于遺傳算法生成部件特性圖的方法，將壓氣機部件特性圖擬合為3次多項式，并根據(jù)遺傳算法尋優(yōu)求解多項式的常系數(shù)[7]；Joachim Kurzke介紹了基于模型的燃氣渦輪部件特性修正方法，研究大氣環(huán)境條件對發(fā)動機流道內(nèi)變量的影響，給出詳細的校準方法[8]；段守付等提出用加權(quán)函數(shù)法和卡爾曼濾波方法對模型和特性圖進行修正[9]；陳策等提出1種將BP神經(jīng)網(wǎng)絡(luò)應(yīng)用于渦扇發(fā)動機部件特性的方法，插值和外推出較為準確的部件特性數(shù)據(jù)[10]；吳虎等對發(fā)動機部件特性修正方法進行發(fā)展，以試車數(shù)據(jù)與模型計算數(shù)據(jù)之間誤差最小為尋優(yōu)目標，采用單純形法修正通用部件特性[11]。上述修正方法在原理上存在相似性，在求解非線性方程組或者進行非線性尋優(yōu)的方法不同。

本文采用具有高效、并行、全局搜索能力的遺傳算法修正發(fā)動機部件特性，建立渦噴發(fā)動機身份證模型。不改變模型的核心計算部分，通過遺傳算法尋優(yōu)求出高/低壓壓氣機部件特性的調(diào)節(jié)因子，按比例縮放通用部件特性，以提高發(fā)動機模型精度。該方法通用性強，計算速度快。對于出廠試車后的發(fā)動機，通過試車數(shù)據(jù)匹配發(fā)動機計算模型，減小人工調(diào)整部件特性圖的工作量；在發(fā)動機使用壽命期內(nèi)由于性能正常蛻化引起發(fā)動機模型失配時，該方法自動調(diào)整發(fā)動機模型輸出以匹配發(fā)動機測量數(shù)據(jù)。

1 渦噴發(fā)動機模型

建立雙轉(zhuǎn)子渦噴發(fā)動機模型。該發(fā)動機具有3級低壓壓氣機、5級高壓壓氣機、單級高/低壓渦輪轉(zhuǎn)子和收斂型尾噴口，其結(jié)構(gòu)如圖1所示。

圖1 雙轉(zhuǎn)子渦噴發(fā)動機模型

部件級模型的輸入?yún)?shù)為大氣壓力Pamb、大氣溫度TAT、高度H、馬赫數(shù)Ma、主燃油流量Wfuel，輸出參數(shù)為發(fā)動機轉(zhuǎn)速、各截面氣動參數(shù)（總溫、總壓等）及性能參數(shù)。按部件法建模思想，分別建立進氣道、低壓壓氣機、高壓壓氣機、燃燒室、高壓渦輪、低壓渦輪和噴管的氣動熱力學(xué)模型，如圖2所示。從以下3方面建立發(fā)動機部件級模型[12]：（1）沿發(fā)動機流道的各部件氣動熱力計算；（2）建立反映發(fā)動機動態(tài)、穩(wěn)態(tài)工作過程中各部件間共同工作關(guān)系（流量平衡、轉(zhuǎn)速平衡、功率平衡）的控制方程；（3）用Newton-Raphson法求解控制方程。

圖2 模型流路計算

建立渦噴發(fā)動機各部件氣動熱力學(xué)非線性方程組[13]含有6個獨立變量，分別為低壓軸轉(zhuǎn)速nL，低壓壓氣機增壓比πLPC、高壓壓氣機增壓比πHPC、高壓軸轉(zhuǎn)速nH、高壓渦輪落壓比πHPT、低壓渦輪落壓比πLPT。這也是發(fā)動機穩(wěn)態(tài)仿真的初猜值，對應(yīng)的6個共同工作方程如下。

低壓壓氣機流量與高壓壓氣機流量平衡

高壓壓氣機流量與高壓渦輪流量平衡

高壓渦輪流量與低壓渦輪流量平衡

低壓渦輪流量與內(nèi)涵噴管流量平衡

穩(wěn)態(tài)求解時，功率平衡方程

式中：Wfuel為燃油流量；W為氣體流量；下標數(shù)字為各截面標識；Wbleed為高壓壓氣機引氣流量；WcoolHPT、WcoolLPT分別為高、低壓轉(zhuǎn)子的冷卻流量；ηHPT、ηLPT分別為高、低壓轉(zhuǎn)子的機械效率；Eext為從高壓軸輸出的功率。

對于動態(tài)過程的仿真，式（5）、（6）轉(zhuǎn)化為轉(zhuǎn)子動力學(xué)方程

式中：JH、JL分別為高、低壓軸轉(zhuǎn)動慣量。

在動態(tài)仿真前，先進行穩(wěn)態(tài)計算作為動態(tài)仿真的初始條件，然后根據(jù)下一時刻模型輸入?yún)?shù)，采用改進的歐拉法求解上述2個轉(zhuǎn)子動力學(xué)方程，解出nL和nH后，采用Newton-Raphson法求解發(fā)動機部件組成的非線性方程組。

2 遺傳算法修正渦噴發(fā)動機模型

采用遺傳算法修正渦噴發(fā)動機模型的原理為：選取發(fā)動機模型測量參數(shù)與試車數(shù)據(jù)的偏差作為被調(diào)節(jié)的目標函數(shù)，分別設(shè)定高、低壓壓氣機部件特性圖的調(diào)節(jié)因子，通過遺傳算法尋優(yōu)求出使目標函數(shù)最小的調(diào)節(jié)因子，按比例縮放壓氣機部件特性，以達到修正模型精度目的。

2.1 調(diào)節(jié)因子的選取

發(fā)動機高、低渦輪部件相對換算轉(zhuǎn)速與折合流量變化范圍很窄，共同工作點只在設(shè)計點附近很小的區(qū)域內(nèi)變化。而壓氣機的折合流量與相對換算轉(zhuǎn)速隨物理轉(zhuǎn)速變化范圍很大，且具有較好的線性和單調(diào)性的特點。因此，保持發(fā)動機渦輪部件特性不變，通過設(shè)立低壓壓氣機的流量調(diào)節(jié)因子λLpcMas和效率調(diào)節(jié)因子λLpcEff、高壓壓氣機的流量調(diào)節(jié)因子λHpcMas和效率調(diào)節(jié)因子λHpcEff，調(diào)節(jié)高、低壓壓氣機部件特性以匹配渦輪部件特性，達到修正發(fā)動機模型的目的。

2.2 修正的目標函數(shù)

修正的目標函數(shù)即為遺傳算法的適應(yīng)度函數(shù)。本文選取3個測量參數(shù)（nL、nH和低壓渦輪后溫度T6）作為模型修正的參考。設(shè)定發(fā)動機測量參數(shù)的偏差[14]為

式中：ai為加權(quán)因子，本文令其皆為1。

2.3 基于遺傳算法尋優(yōu)調(diào)節(jié)因子

利用遺傳算法模擬自然選擇和遺傳中發(fā)生的復(fù)制、交叉和變異等現(xiàn)象。從初始種群出發(fā)，通過隨機選擇、交叉、變異操作，產(chǎn)生更適合環(huán)境的個體。在每代遺傳算法中，根據(jù)個體在問題域中的適應(yīng)度值和從自然遺傳學(xué)中借鑒來的再造方法進行個體選擇，產(chǎn)生1個新的近似解。種群進化后得到問題的最優(yōu)解。

本文采用二進制編碼，輪盤賭選擇機制，均勻交叉、基本變異。控制參數(shù)包括：群體大小M=20，交叉概率PC=0.8，變異概率Pm=0.02，最大進化代數(shù)G=20。遺傳算法參數(shù)尋優(yōu)時調(diào)用發(fā)動機模型用于適應(yīng)度函數(shù)值計算，自動尋優(yōu)的計算結(jié)構(gòu)如圖3所示。

利用遺傳算法修正發(fā)動機高、低壓壓氣機流量特性和效率特性的具體步驟如下：

圖3 基于遺傳算法的渦噴發(fā)動機模型修正

（1）為減小計算量，按試車數(shù)據(jù)提取具有代表意義的數(shù)據(jù)，按燃油量從大到小排序；

（2）按數(shù)據(jù)表順序，將大氣條件、飛行狀態(tài)條件、燃油量輸入發(fā)動機模型后輸出仿真計算結(jié)果，該結(jié)果與發(fā)動機測量參數(shù)共同構(gòu)成目標函數(shù)；

（3）判斷是否達到最大遺傳代數(shù)。若未達到則采用遺傳算法改變高、低壓壓氣機特性圖修正系數(shù)，并回到（2）；若達到則記錄基于低壓壓氣機相對換算轉(zhuǎn)速NLLPC，cor的λLpcMas、λLpcEff、λHpcMas和λHpcEff。判斷數(shù)據(jù)表中的數(shù)據(jù)驗算是否完畢，若執(zhí)行完成，生成總修正系數(shù)表，否則將輸入的標識移至下一行，回到（2）。總修正系數(shù)表格式見表1。

表1 以LPC相對換算轉(zhuǎn)速為基準的修正系數(shù)

3 渦噴發(fā)動機仿真分析

3.1 修正個體差異模型

因存在個體差異，需要建立身份證模型來準確進行發(fā)動機仿真。采用地面試車數(shù)據(jù)修正模型進行仿真，結(jié)果如圖4~6所示。

圖4 nL轉(zhuǎn)速比較

圖5 nH轉(zhuǎn)速比較

圖6 T6比較

將原模型和基于遺傳算法修正后模型的計算結(jié)果與試車數(shù)據(jù)比較，修正后模型的計算結(jié)果與試車數(shù)據(jù)接近，總體精度有所提高。通過改變目標函數(shù)中選取的發(fā)動機測量參數(shù)對應(yīng)的加權(quán)因子大小，改善計算精度。目標函數(shù)中選擇較多的發(fā)動機測量參數(shù)，有利于提高發(fā)動機模型的總體仿真精度。

3.2 修正性能蛻化后模型

在發(fā)動機工作一段時間后，各部件性能發(fā)生一定量的蛻化，若不調(diào)整發(fā)動機計算模型，仿真計算結(jié)果將產(chǎn)生偏差。通過定期采集的試車數(shù)據(jù)，可以自動修正模型輸出結(jié)果。本節(jié)采用表2列出的蛻化量進行仿真驗證[15]，結(jié)果如圖7~9所示。

表2 發(fā)動機3000循環(huán)各部件性能蛻化量

圖7 發(fā)動機性能蛻化后nL轉(zhuǎn)速比較

圖8 發(fā)動機性能蛻化后nH轉(zhuǎn)速比較

基于遺傳算法修正后，模型的輸出結(jié)果得到改善，總體的仿真精度有所提高。渦輪部件蛻化帶來模型輸出偏差時，通過修正壓氣機部件特性進行補償。以圖7為例，由圖可知：在低轉(zhuǎn)速下，原模型計算結(jié)果與試車數(shù)據(jù)偏差較大，修正后，模型修正了低壓軸轉(zhuǎn)速；在高轉(zhuǎn)速下，原模型計算結(jié)果接近于試車數(shù)據(jù)，修正后，模型改善了仿真精度。

圖9 發(fā)動機性能蛻化后T6比較

3.3 修正后模型動態(tài)仿真

在發(fā)動機蛻化后的修正模型中，帶入修正系數(shù)表，發(fā)動機動態(tài)仿真如圖10所示。仿真結(jié)果表明：帶入遺傳算法求解的修正系數(shù)表進行動態(tài)仿真，模型運算順暢，Newton-Raphson求解模型非線性方程組時收斂。尚無試車動態(tài)數(shù)據(jù)可供對比，因此不能定量比較動態(tài)模型的計算精度。

圖10 發(fā)動機動態(tài)計算仿真

4 結(jié)論

（1）高、低壓壓氣機部件特性流量與效率調(diào)節(jié)因子可作為被修正參數(shù)，修正發(fā)動機部件特性。

（2）建立渦噴發(fā)動機修正模型，采用遺傳算法尋優(yōu)的方式可使模型輸出與發(fā)動機試車數(shù)據(jù)的偏差最小，實現(xiàn)發(fā)動機穩(wěn)態(tài)仿真。

（3）以低壓壓氣機相對換算轉(zhuǎn)速為基準，對修正系數(shù)表插值計算，實時求出部件特性圖調(diào)節(jié)因子。通過調(diào)節(jié)因子修正部件特性圖，可使發(fā)動機模型動態(tài)仿真計算收斂流暢。

[1]Drummond C K，F(xiàn)ollen G J，Putt C W．Gas turbine system simulation：an objected-oriented approach[R]．NASA-TM-106044.

[2]范作民，孫春林，白杰.航空發(fā)動機故障診斷導(dǎo)論[M].北京：科學(xué)出版社，2004：1-2.FAN Zuomin，SUN Chunlin，BAI Jie.Introduction of aeroengine fault diagnosis[M].Beijing：Science Press，2004：1-2.（in Chinese）

[3]Doel D L.Interpretation of weight-least squares gas path analysis results[J].Transactions of ASME，2005，（125）：624-633.

[4]Carl JD，Susan M K.Generation of linear dynamic models from a digital nonlinear simulation[R].NASA-TP-1388.

[5]Stamatis A，Mathioudakis K，Papailiou K D.Adaptive simulation of gas turbine performance[J].Journal of Engineering for Gas Turbines and Power，1990，112（2）：168-175.

[6]Stamatis A，Mathioudakis K，Berbs G，et al.Jet engine fault detection with discrete operation gas path analysis[J].Journal of Propulsion，1991，7（6）：1043-1048.

[7]Changduk Kong，Seonghee Kho，Jayoung Ki.Component map generation of a gas turbine using genetic algorithms[J].Journal of Engineering for Gas Turbines and Power，2006，128（1）：92-96.

[8]Joachim K.Model based gas turbine parameter corrections[R].ASME 2003-GT-38234.

[9]段守付，樊思齊，盧燕.航空發(fā)動機自適應(yīng)建模技術(shù)研究[J].航空動力學(xué)報，1999，14（4）：440-442.DUAN Shoufu，F(xiàn)AN Siqi，LU Yan.Adaptive modelling technique for aeroengine[J].Journal of Aerospace Power，1999，14（4）：440-442.（in Chinese）

[10]陳策，李軍.渦扇發(fā)動機部件特性的BP網(wǎng)絡(luò)研究[J].航空動力學(xué)報，2004，19（1）：61-64.CHEN Ce，LI Jun.A study of turbofan component characteristics based on back-propagation network[J].Journal of Aerospace Power，2004，19（1）：61-64.（in Chinese）

[11]吳虎，肖洪，蔣建軍.渦扇發(fā)動機部件特性自適應(yīng)模擬[J].推進技術(shù)，2005，26（5）：430-433.WU Hu，XIAO Hong，JIANG Jianjun.Adaptive simulation of turbofan engine component performance[J].Journal of Propulsion Technology，2005，26（5）：430-433.（in Chinese）

[12]孫健國，黃金泉，葉志峰，等.現(xiàn)代航空動力裝置控制[M].北京：航空工業(yè)出版社，2001：125-130.SUN Jianguo，HUANG Jinquan，YE Zhifeng，et al.The control of modern aircraft power plant[M].Beijing：Aviation Indurstry Press，2001：125-130.（in Chinese）

[13]周文祥，黃金泉，竇建平，等.面向?qū)ο蟮臏u扇發(fā)動機及控制系統(tǒng)仿真平臺[J].航空動力學(xué)報，2007，22（1）：119-125.ZHOU Wenxiang，HUANG Jinquan，DOU Jianping，et al.Objectedoriented simulation platform for turbine engine and its control system[J].Journal of Aerospace Power，2007，22（1）：119-125.（in Chinese）

[14]肖洪，劉振俠，廉筱純.兩種渦扇發(fā)動機部件特性自適應(yīng)模型對比[J].中國民航大學(xué)學(xué)報，2008，26（3）：17-19.XIAO Hong，LIU Zhenxia，LIAN Xiaochun.Two simulation models of turbofan component performance[J].Journal of Civil Aviation University of China，2008，26（3）：17-19.（in Chinese）

[15]Sanjay G.Propulsion control and health management research at NASA Glenn research center[R].NASA-TM-2002-211590.