牟春陽(yáng),李世中
(中北大學(xué) 機(jī)電工程學(xué)院,太原 030051)
發(fā)動(dòng)機(jī)數(shù)學(xué)模型是發(fā)動(dòng)機(jī)性能等用數(shù)學(xué)表述的一種形式。由于發(fā)動(dòng)機(jī)控制系統(tǒng)設(shè)計(jì)、控制系統(tǒng)數(shù)字仿真和半物理仿真等都要基于發(fā)動(dòng)機(jī)非線性部件級(jí)模型展開(kāi)[1~4],所以發(fā)動(dòng)機(jī)建模是發(fā)動(dòng)機(jī)控制系統(tǒng)研究的重要部分,只有所建立的數(shù)學(xué)模型能夠在全包線、全狀態(tài)下模擬真實(shí)發(fā)動(dòng)機(jī)的工作情況,才能保證基于部件級(jí)模型的研究成果的有效性。為此早在20世紀(jì)70年代發(fā)動(dòng)機(jī)數(shù)學(xué)模型的研究就引起的廣泛關(guān)注[5~7]。
建立發(fā)動(dòng)機(jī)數(shù)學(xué)模型的方法有解析法和試驗(yàn)法[8]。目前國(guó)內(nèi)的發(fā)動(dòng)機(jī)數(shù)學(xué)模型大多數(shù)為通過(guò)求解共同工作方程建立的部件級(jí)氣動(dòng)熱力學(xué)模型或小偏差線性化模型,但這種方法有它的不足之處,由于發(fā)動(dòng)機(jī)具有強(qiáng)非線性,只有初猜值在真實(shí)解附近,該方法才能保證迭代收斂。這在實(shí)際計(jì)算中保證全飛行包線內(nèi)都收斂是很困難的。為解決部件級(jí)建模中解非線性方程組時(shí)的迭代問(wèn)題,本文基于容腔動(dòng)力學(xué)建立渦輪發(fā)動(dòng)機(jī)的數(shù)學(xué)模型,認(rèn)為動(dòng)態(tài)時(shí)流量不再平衡,壓力隨流量變化而變化,并通過(guò)該數(shù)學(xué)模型各截面參數(shù)與試驗(yàn)數(shù)據(jù)比對(duì),驗(yàn)證模型的準(zhǔn)確性。
一臺(tái)噴氣發(fā)動(dòng)機(jī)包含眾多的小容腔,每一個(gè)小容腔可以存儲(chǔ)熱能和一定質(zhì)量的空氣/燃?xì)狻S捎诳諝夂腿細(xì)馐沁B續(xù)介質(zhì),一個(gè)小容腔內(nèi)的能量和質(zhì)量動(dòng)力學(xué)行為均是分布式系統(tǒng),而小容腔的大小決定了用集中參數(shù)系統(tǒng)來(lái)近似容腔動(dòng)力學(xué)的逼近程度。把空氣/燃?xì)馔ǖ赖男∪萸粍澐值脑叫?,則逼近越準(zhǔn)確[8]。對(duì)于分排渦扇發(fā)動(dòng)機(jī),其容腔劃分如圖1灰色邊框區(qū)域。
圖1 分排渦扇發(fā)動(dòng)機(jī)空氣和燃?xì)馊萸?/p>
容腔內(nèi)的質(zhì)量存儲(chǔ)效應(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和密度ρ來(lái)表征。假設(shè)空氣在標(biāo)稱(chēng)工作點(diǎn)附近為理想氣體,容腔內(nèi)壓力的微小變化量為:
式中,V是控制面包圍的體積,為常量。上述方程對(duì)時(shí)間微分,得到:
因而容腔內(nèi)壓力的變化率與溫度和質(zhì)量的變化率相關(guān)。然而,由于上述方程右邊第一項(xiàng)要遠(yuǎn)小于第二項(xiàng),可以舍去第一項(xiàng),這樣,上式中壓力項(xiàng)就近似正比于容腔內(nèi)質(zhì)量的變化率。
即:
分排渦扇發(fā)動(dòng)機(jī)工作過(guò)程中,各部件之間滿足以下共同工作條件:
壓氣機(jī)與高壓渦輪之間空氣流量連續(xù);
高壓渦輪與低壓渦輪之間燃?xì)饬髁窟B續(xù);
低壓渦輪與噴管之間空氣流量連續(xù);
風(fēng)扇進(jìn)口與發(fā)動(dòng)機(jī)出口總流量連續(xù);
高壓轉(zhuǎn)子功率平衡;
低壓轉(zhuǎn)子功率平衡。
根據(jù)容腔動(dòng)力學(xué)描述,將每個(gè)流量連續(xù)共同工作條件對(duì)應(yīng)一個(gè)容腔模塊,計(jì)算時(shí)根據(jù)容腔進(jìn)出口截面流量差更新模型中相對(duì)應(yīng)的壓力狀態(tài);另外,轉(zhuǎn)子的功率平衡約束是一組描述發(fā)動(dòng)機(jī)轉(zhuǎn)速動(dòng)態(tài)變化的微分方程,根據(jù)轉(zhuǎn)子功率平衡方程,計(jì)算過(guò)程中可采用歐拉積分的方法對(duì)高低壓轉(zhuǎn)子轉(zhuǎn)速進(jìn)行更新。
分排渦扇發(fā)動(dòng)機(jī)動(dòng)態(tài)模型可以表示成如下形式:
式(6)中,模型狀態(tài)x1,輸入變量u和輸出變量y分別為:
在給定的狀態(tài)、輸入條件下,通過(guò)發(fā)動(dòng)機(jī)部件計(jì)算可以求得對(duì)應(yīng)的發(fā)動(dòng)機(jī)整機(jī)及各個(gè)部件性能參數(shù)。部件計(jì)算后,根據(jù)部件容腔模塊更新發(fā)動(dòng)機(jī)動(dòng)態(tài)模型中的壓力狀態(tài)參數(shù),根據(jù)功率平衡方程更新高低壓轉(zhuǎn)子轉(zhuǎn)速。
分排渦扇發(fā)動(dòng)機(jī)動(dòng)態(tài)模型中將風(fēng)扇后作為1號(hào)容腔模塊,由1號(hào)容腔更新風(fēng)扇后總壓Pt2;在壓氣機(jī)與燃燒室之間加入2號(hào)容腔,由2號(hào)容腔更新壓氣機(jī)后總壓Pt3;在高低壓渦輪之間加入3號(hào)容腔,由3號(hào)容腔更新高壓渦輪后總壓Pt5;在低壓渦輪后與尾噴管之間加入4號(hào)容腔,由4號(hào)容腔更新低壓渦輪出口總壓Pt6。部件容腔的放置位置如圖3中所示。
圖3 分排渦扇發(fā)動(dòng)機(jī)模型部件容腔位置示意圖
容腔計(jì)算:
1)風(fēng)扇容腔
對(duì)風(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)子的動(dòng)力學(xué)方程。
根據(jù)動(dòng)能定理,高、低壓渦輪輸出功率與壓氣機(jī)、增壓級(jí)、風(fēng)扇部件的消耗功率之間滿足:
上式中,η為渦輪軸的機(jī)械效率,J為轉(zhuǎn)子轉(zhuǎn)動(dòng)慣量。
轉(zhuǎn)子轉(zhuǎn)速n與轉(zhuǎn)動(dòng)角速度ω之間的換算關(guān)系為:
將式(16)帶入式(14)、式(15)中,整理可得描述發(fā)動(dòng)機(jī)轉(zhuǎn)子功率平衡的微分方程:
為保證動(dòng)態(tài)模型程序的實(shí)時(shí)性,并使軟件在不同環(huán)境下具有良好的兼容性,采用標(biāo)準(zhǔn)C語(yǔ)言進(jìn)行代碼編寫(xiě)。根據(jù)發(fā)動(dòng)機(jī)動(dòng)態(tài)數(shù)學(xué)模型功能計(jì)算需求,按照“高內(nèi)聚、低耦合”的模塊劃分原則,采用面向?qū)ο笏枷氪罱òl(fā)動(dòng)機(jī)動(dòng)態(tài)數(shù)學(xué)模型。發(fā)動(dòng)機(jī)模型程序主要由16個(gè)子功能模塊組成,分別為氣動(dòng)函數(shù)計(jì)算模塊、矩陣運(yùn)算模塊、插值計(jì)算模塊、發(fā)動(dòng)機(jī)尺寸模塊、發(fā)動(dòng)機(jī)轉(zhuǎn)子模塊、發(fā)動(dòng)機(jī)輸入模塊、進(jìn)氣道計(jì)算模塊、風(fēng)扇計(jì)算模塊、壓氣機(jī)計(jì)算模塊、燃燒室計(jì)算模塊、高壓渦輪計(jì)算模塊、低壓渦輪計(jì)算模塊、外涵道計(jì)算模塊、噴管計(jì)算模塊、容腔計(jì)算模塊和發(fā)動(dòng)機(jī)上層組織模塊。
分排渦扇發(fā)動(dòng)機(jī)動(dòng)態(tài)模型程序模塊之間的組織關(guān)系如圖4所示。
圖4 分排渦扇發(fā)動(dòng)機(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ù)與模型輸出對(duì)比如圖7~圖9所示。
圖7 高壓轉(zhuǎn)子轉(zhuǎn)速對(duì)比曲線
圖8 低壓轉(zhuǎn)子轉(zhuǎn)速對(duì)比曲線
圖9 壓氣機(jī)出口總壓對(duì)比曲線
綜上,模型仿真數(shù)據(jù)與試驗(yàn)數(shù)據(jù)存在一定的差異(其中壓氣機(jī)出口總壓相差較大,約6%左右,高、低壓轉(zhuǎn)子轉(zhuǎn)速相差較小,約1%左右),其差異存在的主要原因是部件特性存在一定的誤差。
分排渦扇發(fā)動(dòng)機(jī)動(dòng)態(tài)數(shù)學(xué)模型是進(jìn)行發(fā)動(dòng)機(jī)性能仿真分析和進(jìn)行發(fā)動(dòng)機(jī)控制系統(tǒng)設(shè)計(jì)、分析和數(shù)值驗(yàn)證的重要基礎(chǔ)。本文詳細(xì)介紹了基于容腔動(dòng)力學(xué)約束的非線性動(dòng)態(tài)數(shù)學(xué)模型建模方法,該方法可以有效避免發(fā)動(dòng)機(jī)性能計(jì)算中的迭代過(guò)程,具有更好的計(jì)算穩(wěn)定性。文中通過(guò)仿真實(shí)例驗(yàn)證了模型的收斂性。該模型是發(fā)動(dòng)機(jī)實(shí)時(shí)動(dòng)態(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ā)動(dòng)機(jī)組態(tài)建模仿真技術(shù)研究[D].南京:南京航空航天大學(xué),2009.12.
[4]鄒先權(quán).渦軸發(fā)動(dòng)機(jī)自適應(yīng)模型的建立與魯棒控制[D].南京:南京航空航天大學(xué),2008.12.
[5]于龍江,樸英.一種簡(jiǎn)化算法的航空發(fā)動(dòng)機(jī)全狀態(tài)數(shù)學(xué)模型[J].航空動(dòng)力學(xué)報(bào),2008,23(3):510-515.
[6]楊剛,孫健國(guó),黃向華,等.一種不需要迭代的發(fā)動(dòng)機(jī)輔助變量建模方法[J].航空動(dòng)力學(xué)報(bào),2003,18(2):289-294.
[7]夏超,王繼強(qiáng),商國(guó)軍,等.基于Matlab/Simulink的航空發(fā)動(dòng)機(jī)部件級(jí)建模與分析[J].航空發(fā)動(dòng)機(jī),2012,38(4):31-33.
[8]張新國(guó),譯.飛機(jī)發(fā)動(dòng)機(jī)控制——設(shè)計(jì)、系統(tǒng)分析和健康監(jiān)視[M].北京:航空工業(yè)出版社,2011.