牟春陽，李世中

（中北大學(xué) 機(jī)電工程學(xué)院，太原 030051）

0 引言

發(fā)動機(jī)數(shù)學(xué)模型是發(fā)動機(jī)性能等用數(shù)學(xué)表述的一種形式。由于發(fā)動機(jī)控制系統(tǒng)設(shè)計(jì)、控制系統(tǒng)數(shù)字仿真和半物理仿真等都要基于發(fā)動機(jī)非線性部件級模型展開[1～4]，所以發(fā)動機(jī)建模是發(fā)動機(jī)控制系統(tǒng)研究的重要部分，只有所建立的數(shù)學(xué)模型能夠在全包線、全狀態(tài)下模擬真實(shí)發(fā)動機(jī)的工作情況，才能保證基于部件級模型的研究成果的有效性。為此早在20世紀(jì)70年代發(fā)動機(jī)數(shù)學(xué)模型的研究就引起的廣泛關(guān)注[5～7]。

建立發(fā)動機(jī)數(shù)學(xué)模型的方法有解析法和試驗(yàn)法[8]。目前國內(nèi)的發(fā)動機(jī)數(shù)學(xué)模型大多數(shù)為通過求解共同工作方程建立的部件級氣動熱力學(xué)模型或小偏差線性化模型，但這種方法有它的不足之處，由于發(fā)動機(jī)具有強(qiáng)非線性，只有初猜值在真實(shí)解附近，該方法才能保證迭代收斂。這在實(shí)際計(jì)算中保證全飛行包線內(nèi)都收斂是很困難的。為解決部件級建模中解非線性方程組時(shí)的迭代問題，本文基于容腔動力學(xué)建立渦輪發(fā)動機(jī)的數(shù)學(xué)模型，認(rèn)為動態(tài)時(shí)流量不再平衡，壓力隨流量變化而變化，并通過該數(shù)學(xué)模型各截面參數(shù)與試驗(yàn)數(shù)據(jù)比對，驗(yàn)證模型的準(zhǔn)確性。

1 容腔動力學(xué)

一臺噴氣發(fā)動機(jī)包含眾多的小容腔，每一個(gè)小容腔可以存儲熱能和一定質(zhì)量的空氣/燃?xì)?。由于空氣和燃?xì)馐沁B續(xù)介質(zhì)，一個(gè)小容腔內(nèi)的能量和質(zhì)量動力學(xué)行為均是分布式系統(tǒng)，而小容腔的大小決定了用集中參數(shù)系統(tǒng)來近似容腔動力學(xué)的逼近程度。把空氣/燃?xì)馔ǖ赖男∪萸粍澐值脑叫?，則逼近越準(zhǔn)確[8]。對于分排渦扇發(fā)動機(jī)，其容腔劃分如圖1灰色邊框區(qū)域。

圖1 分排渦扇發(fā)動機(jī)空氣和燃?xì)馊萸?/p>

容腔內(nèi)的質(zhì)量存儲效應(yīng)導(dǎo)致容腔內(nèi)的壓力和溫度發(fā)生變化。這些變化是容腔內(nèi)空氣或燃?xì)獾膲毫蜏囟却嬖跓崃W(xué)聯(lián)系的結(jié)果。考慮如圖2所示的單位空氣容腔。

圖2 空氣/燃?xì)鈫挝蝗莘e

控制面內(nèi)空氣（或燃?xì)猓┵|(zhì)量的變化率為：

容腔內(nèi)的空氣以溫度T、壓力p和密度ρ來表征。假設(shè)空氣在標(biāo)稱工作點(diǎn)附近為理想氣體，容腔內(nèi)壓力的微小變化量為：

式中，V是控制面包圍的體積，為常量。上述方程對時(shí)間微分，得到：

因而容腔內(nèi)壓力的變化率與溫度和質(zhì)量的變化率相關(guān)。然而，由于上述方程右邊第一項(xiàng)要遠(yuǎn)小于第二項(xiàng)，可以舍去第一項(xiàng)，這樣，上式中壓力項(xiàng)就近似正比于容腔內(nèi)質(zhì)量的變化率。

即：

2 分排渦扇發(fā)動機(jī)動態(tài)模型建模方法

分排渦扇發(fā)動機(jī)工作過程中，各部件之間滿足以下共同工作條件：

壓氣機(jī)與高壓渦輪之間空氣流量連續(xù)；

高壓渦輪與低壓渦輪之間燃?xì)饬髁窟B續(xù)；

低壓渦輪與噴管之間空氣流量連續(xù)；

風(fēng)扇進(jìn)口與發(fā)動機(jī)出口總流量連續(xù)；

高壓轉(zhuǎn)子功率平衡；

低壓轉(zhuǎn)子功率平衡。

根據(jù)容腔動力學(xué)描述，將每個(gè)流量連續(xù)共同工作條件對應(yīng)一個(gè)容腔模塊，計(jì)算時(shí)根據(jù)容腔進(jìn)出口截面流量差更新模型中相對應(yīng)的壓力狀態(tài)；另外，轉(zhuǎn)子的功率平衡約束是一組描述發(fā)動機(jī)轉(zhuǎn)速動態(tài)變化的微分方程，根據(jù)轉(zhuǎn)子功率平衡方程，計(jì)算過程中可采用歐拉積分的方法對高低壓轉(zhuǎn)子轉(zhuǎn)速進(jìn)行更新。

分排渦扇發(fā)動機(jī)動態(tài)模型可以表示成如下形式：

式(6)中，模型狀態(tài)x1，輸入變量u和輸出變量y分別為：

在給定的狀態(tài)、輸入條件下，通過發(fā)動機(jī)部件計(jì)算可以求得對應(yīng)的發(fā)動機(jī)整機(jī)及各個(gè)部件性能參數(shù)。部件計(jì)算后，根據(jù)部件容腔模塊更新發(fā)動機(jī)動態(tài)模型中的壓力狀態(tài)參數(shù)，根據(jù)功率平衡方程更新高低壓轉(zhuǎn)子轉(zhuǎn)速。

分排渦扇發(fā)動機(jī)動態(tài)模型中將風(fēng)扇后作為1號容腔模塊，由1號容腔更新風(fēng)扇后總壓Pt2；在壓氣機(jī)與燃燒室之間加入2號容腔，由2號容腔更新壓氣機(jī)后總壓Pt3；在高低壓渦輪之間加入3號容腔，由3號容腔更新高壓渦輪后總壓Pt5；在低壓渦輪后與尾噴管之間加入4號容腔，由4號容腔更新低壓渦輪出口總壓Pt6。部件容腔的放置位置如圖3中所示。

圖3 分排渦扇發(fā)動機(jī)模型部件容腔位置示意圖

容腔計(jì)算：

1）風(fēng)扇容腔

對風(fēng)扇容腔效應(yīng)進(jìn)行計(jì)算：

其中，Tt2為風(fēng)扇后總溫，VF為風(fēng)扇處容腔體積，容腔入口流量WaFAN,in=Wa2，Wa2為風(fēng)扇空氣流量，容腔出口流量WaFAN,out=Wa3+Wa82，Wa3為壓氣機(jī)空氣流量，Wa82為外涵道空氣流量。

2）壓氣機(jī)容腔

其中，Tt3為壓氣機(jī)后總溫，VC為壓氣機(jī)處容腔體積，Wa3,in為壓氣機(jī)空氣流量Wa3，Wa3,out=Wg5-Wf1為高壓渦輪燃?xì)饬髁繙p去主燃燒室燃油流量。

3）高壓渦輪容腔

其中，Tt5為高壓渦輪后總溫，Vht為高壓渦輪處容腔體積，Wa5,in為高壓渦輪燃?xì)饬髁縒g5，Wa5,out為低壓渦輪燃?xì)饬髁縒g6。

4）低壓渦輪容腔

其中，Tt6為低壓渦輪后總溫，Vlt為低壓渦輪處容腔體積，Wa6,in為低壓渦輪部件空氣流量的計(jì)算結(jié)果Wg6，Wa6,out為尾噴管出口燃?xì)饬髁縒g81。

另外，考慮兩個(gè)轉(zhuǎn)子的動力學(xué)方程。

根據(jù)動能定理，高、低壓渦輪輸出功率與壓氣機(jī)、增壓級、風(fēng)扇部件的消耗功率之間滿足：

上式中，η為渦輪軸的機(jī)械效率，J為轉(zhuǎn)子轉(zhuǎn)動慣量。

轉(zhuǎn)子轉(zhuǎn)速n與轉(zhuǎn)動角速度ω之間的換算關(guān)系為：

將式(16)帶入式(14)、式(15)中，整理可得描述發(fā)動機(jī)轉(zhuǎn)子功率平衡的微分方程：

3 仿真驗(yàn)證

為保證動態(tài)模型程序的實(shí)時(shí)性，并使軟件在不同環(huán)境下具有良好的兼容性，采用標(biāo)準(zhǔn)C語言進(jìn)行代碼編寫。根據(jù)發(fā)動機(jī)動態(tài)數(shù)學(xué)模型功能計(jì)算需求，按照“高內(nèi)聚、低耦合”的模塊劃分原則，采用面向?qū)ο笏枷氪罱òl(fā)動機(jī)動態(tài)數(shù)學(xué)模型。發(fā)動機(jī)模型程序主要由16個(gè)子功能模塊組成，分別為氣動函數(shù)計(jì)算模塊、矩陣運(yùn)算模塊、插值計(jì)算模塊、發(fā)動機(jī)尺寸模塊、發(fā)動機(jī)轉(zhuǎn)子模塊、發(fā)動機(jī)輸入模塊、進(jìn)氣道計(jì)算模塊、風(fēng)扇計(jì)算模塊、壓氣機(jī)計(jì)算模塊、燃燒室計(jì)算模塊、高壓渦輪計(jì)算模塊、低壓渦輪計(jì)算模塊、外涵道計(jì)算模塊、噴管計(jì)算模塊、容腔計(jì)算模塊和發(fā)動機(jī)上層組織模塊。

分排渦扇發(fā)動機(jī)動態(tài)模型程序模塊之間的組織關(guān)系如圖4所示。

圖4 分排渦扇發(fā)動機(jī)模型程序模塊組織關(guān)系

將模型封裝為MATLAB/Simulink中S函數(shù)，模型輸入為飛行高度H、馬赫數(shù)Ma和燃油流量wf，模型輸出為各截面特征參數(shù)，模型方框圖如圖5所示。

圖5 封裝為S函數(shù)的Simulink圖

模型的輸入為0高度0馬赫數(shù)，燃油輸入值為試驗(yàn)數(shù)據(jù)，如圖6所示。

圖6 燃油流量

試驗(yàn)數(shù)據(jù)與模型輸出對比如圖7～圖9所示。

圖7 高壓轉(zhuǎn)子轉(zhuǎn)速對比曲線

圖8 低壓轉(zhuǎn)子轉(zhuǎn)速對比曲線

圖9 壓氣機(jī)出口總壓對比曲線

綜上，模型仿真數(shù)據(jù)與試驗(yàn)數(shù)據(jù)存在一定的差異（其中壓氣機(jī)出口總壓相差較大，約6%左右，高、低壓轉(zhuǎn)子轉(zhuǎn)速相差較小，約1%左右），其差異存在的主要原因是部件特性存在一定的誤差。

4 結(jié)論

分排渦扇發(fā)動機(jī)動態(tài)數(shù)學(xué)模型是進(jìn)行發(fā)動機(jī)性能仿真分析和進(jìn)行發(fā)動機(jī)控制系統(tǒng)設(shè)計(jì)、分析和數(shù)值驗(yàn)證的重要基礎(chǔ)。本文詳細(xì)介紹了基于容腔動力學(xué)約束的非線性動態(tài)數(shù)學(xué)模型建模方法，該方法可以有效避免發(fā)動機(jī)性能計(jì)算中的迭代過程，具有更好的計(jì)算穩(wěn)定性。文中通過仿真實(shí)例驗(yàn)證了模型的收斂性。該模型是發(fā)動機(jī)實(shí)時(shí)動態(tài)仿真和基于在線實(shí)時(shí)模型的先進(jìn)控制系統(tǒng)設(shè)計(jì)的基礎(chǔ)。

[1]Koenig R W,Fisbbaco L H.Geneng-A program for calculating design and off-design performance for turbojet and turbofan engines[R].NASA TND-6552,1972.

[2]Fisbbacb L H, Koenig R W.Geneng II-A program for calculating design and off-design performance of two- and three-spool turbofans with as many as three nozzles[R].NADA TND-6553,1972.

[3]陶金偉.航空發(fā)動機(jī)組態(tài)建模仿真技術(shù)研究[D].南京:南京航空航天大學(xué),2009.12.

[4]鄒先權(quán).渦軸發(fā)動機(jī)自適應(yīng)模型的建立與魯棒控制[D].南京:南京航空航天大學(xué),2008.12.

[5]于龍江,樸英.一種簡化算法的航空發(fā)動機(jī)全狀態(tài)數(shù)學(xué)模型[J].航空動力學(xué)報(bào),2008,23(3):510-515.

[6]楊剛,孫健國,黃向華,等.一種不需要迭代的發(fā)動機(jī)輔助變量建模方法[J].航空動力學(xué)報(bào),2003,18(2):289-294.

[7]夏超,王繼強(qiáng),商國軍,等.基于Matlab/Simulink的航空發(fā)動機(jī)部件級建模與分析[J].航空發(fā)動機(jī),2012,38(4):31-33.

[8]張新國,譯.飛機(jī)發(fā)動機(jī)控制——設(shè)計(jì)、系統(tǒng)分析和健康監(jiān)視[M].北京:航空工業(yè)出版社,2011.