馬 力，李 冬，原渭蘭

（1.海軍駐沈陽地區(qū)發(fā)動機專業(yè)軍事代表室，沈陽 110015；2.海軍航空兵學(xué)院，遼寧葫蘆島 125001；3.海軍航空工程學(xué)院，山東煙臺 264001）

基于迭代修正影響系數(shù)矩陣的某型發(fā)動機單機氣路模型修正研究

馬 力1，李 冬2，原渭蘭3

（1.海軍駐沈陽地區(qū)發(fā)動機專業(yè)軍事代表室，沈陽 110015；2.海軍航空兵學(xué)院，遼寧葫蘆島 125001；3.海軍航空工程學(xué)院，山東煙臺 264001）

發(fā)動機由于制造、安裝工藝等不同，其部件特性存在差異，使得每臺發(fā)動機的特性不完全相同，因此建立的發(fā)動機模型不能很好地與實際使用的發(fā)動機相匹配。針對此情況，本文引入影響系數(shù)矩陣的方法，修正發(fā)動機氣路模型，使修正后的模型能更好的符合發(fā)動機實際使用情況。仿真結(jié)果表明：迭代修正的影響系數(shù)矩陣方法修正的模型精度優(yōu)于一次修正的結(jié)果，參數(shù)精度提高最大將近1.6%。進一步指出模型修正精度與所選的待修正部件參數(shù)對目標(biāo)性能參數(shù)的敏感性相關(guān)。

模型修正；影響系數(shù)矩陣；迭代修正；敏感性

不同發(fā)動機有不同的特性，甚至同種類型的發(fā)動機由于制造、安裝工藝等不同，其部件特性有差異，使得每臺發(fā)動機的特性也會不完全相同［1－2］。熱力學(xué)建模技術(shù)和一些發(fā)動機性能仿真計算軟件已經(jīng)能夠越來越準(zhǔn)確的建立發(fā)動機模型。但是，這些技術(shù)仍然難以在沒有足夠發(fā)動機部件特性數(shù)據(jù)（例如空氣流量、壓比、效率等）的情況下建立發(fā)動機模型，甚至在設(shè)計點也會產(chǎn)生很大的偏差［3］。特別是對發(fā)動機故障診斷人員來說，發(fā)動機特性的偏差對分析單臺發(fā)動機性能衰退有重要影響。并且準(zhǔn)確的發(fā)動機模型是開展性能分析計算的前提和基礎(chǔ)。因此，發(fā)動機單機氣路模型修正顯得十分有必要。

文獻［4］基于遺傳算法，研究了發(fā)動機氣路性能模型的修正問題。修正結(jié)果表明：遺傳算法能夠很好地修正模型，并且敏感性參數(shù)在發(fā)動機穩(wěn)態(tài)模型修正中起了關(guān)鍵的作用。Stamatis等人介紹了修正系數(shù)法（MF，Modification Factors），并用廣義最小殘差法優(yōu)化修正系數(shù)［5］。Lambiris等人介紹了加權(quán)誤差函數(shù)并使用了多種算法去優(yōu)化修正系數(shù)［6］。文獻［7］基于奇異值分解的方法，利用敏感度分析進行快速選擇并對修正系數(shù)進行優(yōu)化。

本文在上述研究基礎(chǔ)上，利用迭代修正影響系數(shù)矩陣的方法，對發(fā)動機模型進行修正，使之能更好地符合發(fā)動機實際使用情況。并與一次修正方法的修正精度進行對比，分析了待修正部件參數(shù)對目標(biāo)性能參數(shù)的敏感性對模型修正精度的影響。

1 基于影響系數(shù)矩陣的模型修正

發(fā)動機單個狀態(tài)點的性能自適應(yīng)修正是一個反向的數(shù)學(xué)問題。文獻［8］中Escher在他們的發(fā)動機氣路故障診斷方法中應(yīng)用了影響系數(shù)矩陣（ICM，Influence Coefficient Matrix）的逆矩陣模型，有效解決了基于模型的逆向問題，該方法同時適用于發(fā)動機模型的自適應(yīng)修正。

1.1 一次修正模型

在整個自適應(yīng)修正過程中，有兩類參數(shù)需要定義：

（1）待修正的部件參數(shù)——也叫獨立自適應(yīng)參數(shù)或者發(fā)動機設(shè)計點部件參數(shù)，例如空氣流量Wa2、高壓壓氣機增壓比πcH、壓氣機效率ηcH、渦輪效率ηTH等。

（2）目標(biāo)性能參數(shù)——也稱為關(guān)聯(lián)自適應(yīng)參數(shù)或者發(fā)動機性能表征參數(shù)，例如推力F、耗油率sfc、尾噴管出口總溫T*8等。

發(fā)動機待修正的部件參數(shù)與目標(biāo)性能參數(shù)之間的熱力關(guān)系可以表征為

z＝h（x）（1）

其中，目標(biāo)性能參數(shù)向量z∈RM，M是目標(biāo)性能參數(shù)個數(shù)；待修正的部件參數(shù)向量x∈RN，N是待修正的部件參數(shù)個數(shù)；h（）是多目標(biāo)判別函數(shù)?；鶞?zhǔn)點用下標(biāo)0表示，則式（1）可以在基準(zhǔn)點按Taylor級數(shù)展開，即

式（2）中，H．O．T為高階小量，其影響不大，可以忽略。因此基準(zhǔn)點的帶修正部件參數(shù)的偏離量和目標(biāo)性能參數(shù)的偏離量之間的關(guān)系式可以線性化為

式（3）中，H為影響系數(shù)矩陣?？梢赃x擇這些待修正的部件參數(shù)設(shè)計點作為初值，與之相對應(yīng)的目標(biāo)性能參數(shù)與實際值相比會存在一定的誤差。通過對影響系數(shù)矩陣H求逆，得到自適應(yīng)系數(shù)矩陣H－1，就可以用相對應(yīng)基準(zhǔn)點的目標(biāo)性能參數(shù)的誤差來修正自適應(yīng)參數(shù)的偏差。當(dāng)N＝M時，

很多時候，待修正的部件參數(shù)N和目標(biāo)性能參數(shù)M并不相等。此時，可以得到最小二乘解。

其中，

理想的線性修正效果如圖1所示。

在待修正的部件參數(shù)與初始值的偏差不是很大，并且氣動熱力模型在初始設(shè)計點近似線性的情況下，發(fā)動機部件參數(shù)從初始值變化到合適的數(shù)值時，Δx的線性預(yù)測結(jié)果將是這個偏差良好的估計值，其計算流程如圖2所示。

1.2 迭代修正模型

如果待修正的部件參數(shù)和目標(biāo)性能參數(shù)之間關(guān)系是強非線性的，用1.1節(jié)的一次修正模型計算得到的目標(biāo)性能參數(shù)與實際值相比會有較大誤差。為了改進修正模型的計算效果，可以通過Newton-Raphson算法進行反復(fù)迭代來建立修正模型，計算的收斂過程見圖3。

當(dāng)預(yù)測值十分接近實際值的時候，圖3所示的收斂過程則停止迭代。收斂準(zhǔn)則為預(yù)測參數(shù)和實際參數(shù)的標(biāo)準(zhǔn)差RMS小于閾值δ。

其中，選取δ＝0．001。

每一個目標(biāo)性能參數(shù)計算誤差εi可以用模型的預(yù)測參數(shù)值和實際的發(fā)動機參數(shù)值來表示。

模型中選擇低壓壓氣機出口空氣流量Wa22、低壓壓氣機增壓比πcL、低壓壓氣機效率ηcL、高壓壓氣機增壓比πcH、高壓壓氣機效率ηcH、低壓渦輪效率ηTL、低壓渦輪流量Wst5、高壓渦輪效率ηTH和高壓渦輪流量Wst45等9個參數(shù)作為待修正的部件參數(shù)；選擇低壓壓氣機出口溫度T2、低壓壓氣機出口壓力P2、高壓壓氣機出口溫度T3、高壓壓氣機出口壓力P3、低壓渦輪出口溫度T5、尾噴管出口溫度T8、發(fā)動機推力F和耗油率sfc等8個參數(shù)作為目標(biāo)性能參數(shù)。

1.3 迭代修正模型的New ton-Raphson解法

對于多元非線性方程組Z＝F（X），已知第k次試取值：

殘差

求第k＋1次試取值X（k＋1），力求使Z（k＋1）＝0。

使用N－R法，則有：

其中，

迭代的目標(biāo)是求出使Z（k＋1）＝0的第（k＋1）次試給值X（k＋1），對于多元方程ΔZ＝H·ΔX，若要使Z（k＋1）＝0，則應(yīng)取：

以新的試取值重新進行計算，一直到RMS＜δ。迭代修正模型的修正流程見圖4。

2 一次修正模型與迭代修正模型的求解及其比較

利用發(fā)動機穩(wěn)態(tài)模型，保持其他部件性能參數(shù)不變，每次只給一個性能參數(shù)x1、一個小增量Δx1，使x′1＝x1＋Δx1，重復(fù)發(fā)動機流程計算，求出zi＋x1（i＝1，2，…，M），由式（18）確定全部對x1的近似偏導(dǎo)數(shù)值，即系數(shù)矩陣第一行數(shù)值：

用同樣的方法可依次求得對其他性能參數(shù)xi（i＝2，3，…，M）的近似偏導(dǎo)數(shù)值，就可以得到整個系數(shù)矩陣H的所有元素值。由計算得到，本文所選取修正的發(fā)動機設(shè)計點的系數(shù)影響矩陣數(shù)值見表1。

因為發(fā)動機工作狀態(tài)不能完全定在一個點，其狀態(tài)有一定的不穩(wěn)定性，主要由以下幾個方面的原因造成：

（1）發(fā)動機通常在非標(biāo)準(zhǔn)條件下工作；

（2）噪聲干擾和測量誤差；

（3）發(fā)動機流場的不均勻分布，甚至在同一個軸向位置不同徑向和周向位置的流場也會有所偏差。

因此，選擇發(fā)動機工作狀態(tài)穩(wěn)定在設(shè)計狀態(tài)附近的時候作為測量點，考慮到制造公差等原因，認為該工作點就是發(fā)動機的設(shè)計點。如果測量點選擇不是在海平面標(biāo)準(zhǔn)條件，有必要對其推力等性能參數(shù)進行修正?；谏鲜鲈?，所測量的氣路性能參數(shù)并不能保證和真實的性能參數(shù)完全一致。假設(shè)測量噪聲分布屬于正態(tài)分布，測量的準(zhǔn)確性由標(biāo)準(zhǔn)差和偏離值所決定，如圖5所示。因此這些氣路參數(shù)的測量值只能作為真實參數(shù)的一個近似值，而通過實測參數(shù)進行自適應(yīng)修正的部件參數(shù)也同樣只能作為一個近似值。

理論上，非線性修正模型因為考慮了發(fā)動機模型的強非線性關(guān)系，所以其修正的結(jié)果要比線性模型更加接近于真實值；但是非線性修正模型需要更多次數(shù)的迭代，因此非線性修正模型所運行的時間會更長。線性修正算法是一種線性運算，算法本身結(jié)構(gòu)簡單，迭代次數(shù)較少；并且在一定偏差范圍內(nèi)，線性修正算法也能達到很高的精度。因此本文采用基于影響系數(shù)矩陣的模型修正方法。

選取發(fā)動機在設(shè)計狀態(tài)下，以實測的參數(shù)為模型修正基準(zhǔn)?？紤]到發(fā)動機的狀態(tài)不能完全穩(wěn)定在一個點，因此以獲取的實際測量參數(shù)均值作為模型修正基準(zhǔn)。其中，實測參數(shù)包括各個部件截面參數(shù)和推力、耗油率等參數(shù)。設(shè)置9個待修正部件參數(shù)，分別利用一次修正和迭代修正系數(shù)矩陣對8個相同的目標(biāo)性能參數(shù)進行修正，具體結(jié)果見表2和表3。

表3和圖6同時對比了分別經(jīng)過一次修正模型和迭代修正模型修正后的待修正部件參數(shù)的偏差量。從表2中可知迭代修正模型的結(jié)果要明顯好于一次修正模型的結(jié)果，因此迭代修正模型將首選作為基于影響系數(shù)矩陣的修正模型。圖7給出了迭代修正模型的迭代收斂過程。

在上述基于影響系數(shù)矩陣修正模型中，選擇了9個待修正的部件參數(shù)，8個目標(biāo)性能參數(shù)。但是改變待修正部件參數(shù)的個數(shù)，修正結(jié)果會有所差別，結(jié)果見表4和表5。仍然選擇全部8個目標(biāo)性能參數(shù)，但只選8個待修正的部件參數(shù)（未選ηTH）。從表4和表5中可以看出，目標(biāo)性能參數(shù)與初始值的偏差變化并不明顯，但是待修正部件參數(shù)的偏差有比較明顯的差別。

保持其他待修正部件參數(shù)不變，每次僅改變一個待修正部件參數(shù)值，使其增大1%，基于發(fā)動機穩(wěn)態(tài)模型，對所有待修正參數(shù)進行目標(biāo)參數(shù)敏感性分析。選擇不同待修正部件參數(shù)組合對自適應(yīng)修正效果有明顯影響。由計算得到，不包括πcL、ηTL的待修正參數(shù)組合的修正結(jié)果沒有明顯的變化，這兩個部件參數(shù)對發(fā)動機性能參數(shù)的影響敏感性不強；而沒有選擇敏感性強的待修正部件參數(shù)的修正模型，可能會導(dǎo)致非線性模型結(jié)果的發(fā)散。

因此在基于影響系數(shù)矩陣的修正模型中，待修正部件參數(shù)必須要將敏感性強的部件參數(shù)包括進來，如果漏選這些參數(shù)，那么修正模型很難得到正確的結(jié)果，迭代也可能趨于發(fā)散。

3 總結(jié)

實際中發(fā)動機由于生產(chǎn)、制造工藝等原因，即使同一批次生產(chǎn)的發(fā)動機，每臺發(fā)動機的特性也不盡相同。因此要對所建立的發(fā)動機模型進行修正。采用迭代修正影響系數(shù)矩陣的方法對發(fā)動機模型進行修正，取得了不錯的效果，最大誤差約為0.12%。在模型修正過程中，盡量選取敏感性強的待修正部件參數(shù)，才能保證模型修正過程不致發(fā)散，修正結(jié)果精度高。

［1］ 謝小平.某型發(fā)動機氣路模型修正研究［D］.煙臺：海軍航空工程學(xué)院飛行工程系，2011.

［2］ 李冬，王永華，顏青欽，趙凱.基于非線性約束模型的發(fā)動機部件性能衰退研究［J］.海軍航空工程學(xué)院學(xué)報，2012：27（3）：301-307.

［3］ 翟高蘭.航空發(fā)動機非線性參數(shù)估計方法［D］.南京：南京航空航天大學(xué)能源與動力學(xué)院，2008.

［4］ 李冬，王冠超，曹明川.基于遺傳算法的發(fā)動機氣路性能模型修正研究［J］.燃氣渦輪試驗與研究.2012：25（4）：46- 50.

［5］ Stamatis A.Mathioudakis K，Papailiou K D.Adaptive simulation of gas turbine performance［J］.Gas Turbine Power.1990：112（2）：168-175.

［6］ Lambiris B，Mathioudakis K，Papailiou K D.Adaptivemodeling of jet engine performance with application to condition monitoring［J］.ISABE.1991（NO 91-7058）.

［7］ Stamatis A.Mathioudakis K，Papailiou K D.Optimal measured parameter and health index selection for gas turbine performance status and fault diagnosis［J］.Gas Turbine Power.1992（114）：209-216.

［8］ Escher PC，Singh R.An objective-oriented diagnostics computer program suitable for industrial gas turbine［C］.21st（CIMAC）International Congress on Combustion Engines.Switzeland，1995，May 15-18.

Research of X Engine Single Gas Path M odel M odification based on Influence Coefficient M atrix of Iterative M odification

MA Li1，LIDong2，YUANWei-lan3
（1.Engine Military Representatives Office of Navy in Shenyang，Shenyang，110015，China；2.Naval Aviation Academy，Liaoning Huludao，125001，China；3.Naval Aeronautical Engineering University，Shandong Yantai，264001，China）

Because of the difference ofmanufacture and installation，component characteristics existdifference，whichmakes characteristics of the engine differ，therefore enginemodel doesn'tmatch factual condition well.Aiming at this condition，influence coefficient matrix was introduced，and engine gas path model wasmodified，modelmodified can corresponding to factual condition.Simulated result indicates：model precision modified by influence coefficientmatrix of iterativemodification is superior to oncemodification，its precision approximates to 1.6%maximally.And precision ofmodelmodified correlates to the sensitivity of selected component parameter to objective performance parameter.

modelmodification，influence coefficientmatrix；iterativemodification；sensitivity

V235.13

A

1009-2889（2014）03-0050-06

2013-12-18改稿日期：2014-02-10

馬力（1977-），男，遼寧沈陽人，工程碩士，工程師，主要從事發(fā)動機設(shè)計與性能分析。通訊作者：李冬.E-mail：happyli.dong＠163.com。