基于MATLAB的微型燃?xì)廨啓C(jī)建模與仿真研究

鐘亞飛， 延偉東

(1. 北京航空航天大學(xué) 能源與動力工程學(xué)院， 北京 100191；2. 北京惠斯特卓越科技發(fā)展有限公司， 北京 100081)

隨著人類社會對能源的需求越來越高，傳統(tǒng)能源逐漸枯竭，環(huán)境污染問題日益嚴(yán)重，社會經(jīng)濟(jì)可持續(xù)發(fā)展目標(biāo)受到了嚴(yán)峻挑戰(zhàn)。為了解決當(dāng)前的能源危機(jī)，須要不斷發(fā)展新型、高效、清潔的能源技術(shù)，智能電網(wǎng)、分布式能源的概念應(yīng)運(yùn)而生[1-4]。微型燃?xì)廨啓C(jī)作為新能源技術(shù)的核心設(shè)備，其技術(shù)的成熟對推動智能電網(wǎng)、分布式能源的發(fā)展起著至關(guān)重要的作用。除此之外，由于微型燃?xì)廨啓C(jī)具備高效、安全、靈活、環(huán)保等優(yōu)點(diǎn)，在交通以及國防建設(shè)等方面均具有很好的發(fā)展前景和應(yīng)用價值。

微型燃?xì)廨啓C(jī)的工作條件具有地域性差異，建立環(huán)境溫度、海拔高度可變的總體仿真模型，對研制新型實(shí)用微型燃?xì)廨啓C(jī)具有重要意義。陳強(qiáng)[5]、李瑩瑩[6]在研究中均建立了環(huán)境溫度與微型燃?xì)廨啓C(jī)發(fā)電量和效率之間關(guān)系的仿真模型，由于在參數(shù)選擇過程中缺乏準(zhǔn)確的部件性能參數(shù)數(shù)據(jù)，致使計算結(jié)果的精度較低。同時，國內(nèi)關(guān)于介紹海拔高度變化對微型燃?xì)廨啓C(jī)發(fā)電量影響的文獻(xiàn)較少，為了彌補(bǔ)現(xiàn)有研究的不足，滿足微型燃?xì)廨啓C(jī)研發(fā)的指導(dǎo)要求，亟須建立更準(zhǔn)確、更全面的變工況參數(shù)仿真模型。

為了驗(yàn)證燃?xì)廨啓C(jī)性能是否達(dá)標(biāo)，試車試驗(yàn)是必不可少的環(huán)節(jié)，但地面試車往往是在不同大氣條件下進(jìn)行的，微型燃?xì)廨啓C(jī)的發(fā)電量、扭矩等參數(shù)會有很大差異。為了與標(biāo)準(zhǔn)微型燃?xì)廨啓C(jī)性能做比較，需要把不同大氣條件下的試車數(shù)據(jù)換算成標(biāo)準(zhǔn)大氣條件下的數(shù)據(jù)。目前，國內(nèi)關(guān)于介紹微型燃?xì)廨啓C(jī)性能換算模型的文獻(xiàn)較少，亟須進(jìn)一步開展微型燃?xì)廨啓C(jī)參數(shù)換算技術(shù)的研究，建立微型燃?xì)廨啓C(jī)參數(shù)換算模型。

筆者借鑒以往燃?xì)廨啓C(jī)分模塊建模經(jīng)驗(yàn)[7-8]，基于MATLAB軟件中GUI程序設(shè)計平臺，在氣體變比熱容條件下建立了包括微型燃?xì)廨啓C(jī)進(jìn)氣管、壓氣機(jī)、換熱器、燃燒室、透平、發(fā)電機(jī)等部件和性能參數(shù)換算模塊的整體模型，計算得到了不同環(huán)境溫度和海拔高度條件下發(fā)電效率和發(fā)電量，并且繪制了轉(zhuǎn)速特性曲線。將仿真結(jié)果與C30微型燃?xì)廨啓C(jī)的標(biāo)準(zhǔn)數(shù)據(jù)進(jìn)行了對比，驗(yàn)證了模型的準(zhǔn)確性，并對試驗(yàn)燃?xì)廨啓C(jī)的性能進(jìn)行了評估。

1 微型燃?xì)廨啓C(jī)數(shù)學(xué)模型

基于模塊化建模理念，建立了包括微型燃?xì)廨啓C(jī)進(jìn)氣管、壓氣機(jī)、換熱器、燃燒室、透平、發(fā)電機(jī)等部件和性能參數(shù)換算模塊的整體模型，其流程圖見圖1。

圖1 微型燃?xì)廨啓C(jī)模型流程圖

1.1 氣體熱力學(xué)性質(zhì)

在實(shí)際過程中，氣體經(jīng)過燃?xì)廨啓C(jī)每個部件，其溫度會有較大變化，致使氣體的比熱容也隨之變化。為了保證模型的準(zhǔn)確性，運(yùn)用MATLAB軟件擬合空氣和煙氣的比定壓熱容隨溫度變化的關(guān)系式。擬合過程中假設(shè)：(1)氣體的比定壓熱容只與溫度有關(guān)；(2)混合氣體的總壓是組成氣體分壓之和[9]。

根據(jù)空氣的物性參數(shù)表[10]，空氣比定壓熱容擬合式為：

cp,air=1.08-6.11×10-4T+1.65×10-6T2-

1.49×10-9T3+6.09×10-13T4-

9.23×10-17T5

(1)

式中：T為空氣熱力學(xué)溫度，K；cp,air為空氣的比定壓熱容，kJ/(kg·K)。

設(shè)燃料的化學(xué)式為CmHn，則燃料燃燒的方程式為：

(2)

根據(jù)方程式，燃燒產(chǎn)物可以認(rèn)為是由CO2、H2O、N2和空氣4部分組成，比定壓熱容計算公式[6]為：

cp,gas=∑(cp,i×yi)

(3)

式中：cp,gas為煙氣比定壓熱容，kJ/(kg·K)；yi為煙氣中i組分的摩爾分?jǐn)?shù)；cp,i為煙氣中i組分的比定壓熱容，kJ/(kg·K)。

1.2 壓強(qiáng)與海拔

隨著海拔的升高，環(huán)境壓力、空氣密度有所下降。根據(jù)質(zhì)量流量的計算公式qm=ρυA，假設(shè)燃?xì)廨啓C(jī)進(jìn)口速度v和截面面積A不變，則進(jìn)口空氣質(zhì)量流量qm,air和密度ρ成正比。表1為海拔高度和密度的參數(shù)表，根據(jù)表中數(shù)據(jù)通過MATLAB軟件擬合多項(xiàng)式，擬合出密度比α(對應(yīng)海拔高度的空氣密度與海平面的空氣密度的比值)隨海拔高度h的變化關(guān)系式見式(4)，擬合數(shù)據(jù)與實(shí)際數(shù)據(jù)的均方差為1.72×10-5，擬合效果良好。

α=-1.91×10-20h5+3.27×10-16h4-

1.96×10-12h3+7.79×10-9h2-

9.91×10-5h+1

(4)

表1 不同海拔高度空氣壓力、密度與密度比數(shù)據(jù)表

1.3 進(jìn)氣管

空氣流經(jīng)進(jìn)氣管吸收燃?xì)廨啓C(jī)余熱，起到對整機(jī)冷卻的作用，維持發(fā)電機(jī)在發(fā)電過程中溫度的穩(wěn)定，從而改善燃?xì)廨啓C(jī)的發(fā)電性能。由于進(jìn)氣管存在一定的流阻，空氣壓力會有所降低，進(jìn)氣管出口氣流溫度、壓強(qiáng)以及余熱表達(dá)式[5]如下：

(5)

p1=p0×(1-ξ0)

(6)

Q0=WT-WG-WC

(7)

式中：t1為進(jìn)氣管出口溫度，℃；Q0為進(jìn)氣管中空氣吸收的余熱，kW；qm,0air為流過進(jìn)氣管空氣的質(zhì)量流量，kg/s；cp,0air為進(jìn)氣管中空氣的比定壓熱容，kJ/(kg·K)；t0為進(jìn)氣管入口溫度，℃；p0為環(huán)境壓力，Pa；p1為進(jìn)氣管出口(即壓氣機(jī)進(jìn)口)壓力，Pa；ξ0為進(jìn)氣管壓力損失系數(shù)；WT為渦輪輸出功率，kW；WG為發(fā)電量，kW；WC為壓氣機(jī)所耗功率，kW。

1.4 壓氣機(jī)

由于存在摩擦損失以及二次流損失，氣體在壓氣機(jī)內(nèi)部的流動極其復(fù)雜，可以認(rèn)為壓氣機(jī)內(nèi)部的流動為熵增絕熱過程。引入壓氣機(jī)絕熱效率ηC，相關(guān)計算公式[5-6]如下：

p2=p1×πC

(8)

(9)

式中：πC為壓氣機(jī)增壓比；p2為壓氣機(jī)出口壓力，Pa；t21為在理想絕熱狀態(tài)下壓氣機(jī)出口處空氣的溫度，℃；t2為壓氣機(jī)出口空氣的實(shí)際溫度，℃。

(10)

(11)

式中：T1為壓氣機(jī)進(jìn)口熱力學(xué)溫度，K；k為空氣比熱比；cp,acom為壓氣機(jī)內(nèi)空氣的比定壓熱容，kJ/(kg·K)。

1.5 回?zé)崞?/h3>

回?zé)崞髯鳛橐粋€熱交換裝置，有效利用了煙氣的熱量，升高了進(jìn)入燃燒室的空氣總溫，提高了燃燒效率?；?zé)岫萢是衡量回?zé)崞餍阅艿闹匾獏?shù)，其擬合式[6]見式(12)。由于內(nèi)部壓阻的原因，流過回?zé)崞鞯臍饬鲏毫τ兴档?。相關(guān)計算公式[5-6]如下：

a=0.846 5-0.002 9t0+0.000 05t02

(12)

t4=a×(t6-t3)+t3

(13)

p4=p3×(1-ξr)

(14)

式中：t3為回?zé)崞鬟M(jìn)口空氣溫度，℃；t4為回?zé)崞鞒隹诳諝鉁囟?，℃；t6為渦輪出口煙氣溫度，℃；p3為回?zé)崞鬟M(jìn)口空氣壓力，Pa；p4為回?zé)崞鞒隹诳諝鈮毫?，Pa；ξr為回?zé)崞鲏毫p失系數(shù)。

1.6 燃燒室

燃料與高溫高壓空氣在燃燒室內(nèi)混合燃燒，釋放出的化學(xué)能轉(zhuǎn)換成煙氣的內(nèi)能，根據(jù)能量守恒定律，忽略進(jìn)口燃料的焓，可得燃燒室能量平衡關(guān)系式如下：

cp,5gas×(qm,4air+qm,fuel)×T5=cp,4air×qm,4air×

T4+ηb×qm,fuel×Hu

(15)

式中：qm,4air為進(jìn)入燃燒室的空氣質(zhì)量流量，kg/s；qm,fuel為燃料質(zhì)量流量，kg/s；Hu為燃料的低熱值，kJ/kg；T4為燃燒室進(jìn)口空氣熱力學(xué)溫度，K；T5為燃燒室出口煙氣熱力學(xué)溫度，K；cp,4air為燃燒室進(jìn)口空氣的比定壓熱容，kJ/(kg·K)；cp,5gas為燃燒室出口煙氣的比定壓熱容，kJ/(kg·K)；ηb為燃燒效率。

根據(jù)油氣比定義[10]可得理想氣體油氣比計算公式如下：

(16)

p5=p4×(1-ξn)

(17)

式中：p5為燃燒室出口壓力，Pa；p4為燃燒室入口壓力，Pa；ξn為燃燒室壓力損失系數(shù)。

1.7 渦輪

微型燃?xì)廨啓C(jī)渦輪一般為徑向結(jié)構(gòu)，與壓氣機(jī)結(jié)構(gòu)相似。作為做功部件，煙氣在渦輪內(nèi)部膨脹，向外部輸出機(jī)械功，其計算公式見式(18)。

(18)

式中：qm,5gas為燃燒室出口煙氣質(zhì)量流量，kg/s；kg為煙氣比熱比；Rg為煙氣的氣體常數(shù)。

與壓氣機(jī)相似，渦輪效率的近似理論公式為：

(19)

(20)

(21)

式中：πT為渦輪落壓比；ηT為渦輪效率；p6為渦輪出口壓力，Pa；t5為渦輪入口(即燃燒室出口)煙氣的溫度，℃；t61為理想絕熱狀態(tài)下渦輪出口煙氣的溫度，℃；cp,tur為渦輪內(nèi)部煙氣的比定壓熱容，kJ/(kg·K)。

1.8 發(fā)電機(jī)

渦輪產(chǎn)生的軸功一部分用于壓氣機(jī)的壓縮功，一部分轉(zhuǎn)變成電力輸出，最后一部分以熱能的形式耗散。引入發(fā)電機(jī)效率ηe，由轉(zhuǎn)軸上的能量平衡可得：

Wout=WT-WC

(22)

WG=Wout×ηe

(23)

(24)

式中：Wout為發(fā)電機(jī)理論輸出功率，kW；η為發(fā)電效率。

1.9 性能參數(shù)換算

在微型燃?xì)廨啓C(jī)試車過程中，實(shí)際條件下燃?xì)廨啓C(jī)工作的轉(zhuǎn)速特性與標(biāo)準(zhǔn)條件下燃?xì)廨啓C(jī)工作的轉(zhuǎn)速特性之間的差異是衡量發(fā)電機(jī)與發(fā)動機(jī)之間是否很好配合的重要信息，是確定微型燃?xì)廨啓C(jī)性能是否符合要求的標(biāo)準(zhǔn)之一。根據(jù)相似理論，把已知環(huán)境溫度為T0、環(huán)境壓力為p0條件下轉(zhuǎn)速、發(fā)電量、扭矩的試驗(yàn)數(shù)據(jù)換算成標(biāo)準(zhǔn)條件下對應(yīng)的數(shù)值。

根據(jù)轉(zhuǎn)速、溫度及壓力相似換算原理，可得渦輪功換算式如下[11]：

(25)

式中：WT,cor為渦輪膨脹功換算值，kW；cp,turcor為渦輪內(nèi)部煙氣的比定壓熱容換算值，kJ/(kg·K)。

由于比定壓熱容變化很小，所以得渦輪功的換算值與實(shí)際值之間的關(guān)系式為：

(26)

同理，壓氣機(jī)在標(biāo)準(zhǔn)條件下的耗功為：

(27)

根據(jù)發(fā)電量的計算公式，在標(biāo)準(zhǔn)條件下工作，燃?xì)廨啓C(jī)的發(fā)電量為:

(28)

式中：WC,cor為壓氣機(jī)耗功換算值，kW；WG,cor為發(fā)電量換算值，kW。

輸出扭矩Mcor換算式為：

(29)

式中：ncor為轉(zhuǎn)速的換算值，r/min。

2 模型功能驗(yàn)證

以C30微型燃?xì)廨啓C(jī)為驗(yàn)證對象，將模型計算結(jié)果對比廠家提供的標(biāo)準(zhǔn)數(shù)據(jù)來驗(yàn)證模型的準(zhǔn)確性并評估實(shí)際燃?xì)廨啓C(jī)性能的優(yōu)劣。

2.1 參數(shù)選擇

筆者使用NUMECA軟件完成了壓氣機(jī)、渦輪的三維計算，在設(shè)計轉(zhuǎn)速96 000 r/min的條件下，求得質(zhì)量流量、壓比以及效率等參數(shù)(見表2)。根據(jù)C30微型燃?xì)廨啓C(jī)廠家提供的尾氣質(zhì)量流量為0.308 kg/s,并且考慮引氣情況，選取設(shè)計點(diǎn)：(1)壓氣機(jī)質(zhì)量流量為0.298 kg/s、增壓比為3.24、效率為71.3%；(2)渦輪質(zhì)量流量為0.31 kg/s、落壓比為2.94、效率為82.8%。廠家提供的燃料熱值、發(fā)電機(jī)效率、燃燒效率以及各部件總壓恢復(fù)系數(shù)等參數(shù)見表3。

表2 壓氣機(jī)、渦輪性能參數(shù)表

表3 相關(guān)參數(shù)表

2.2 結(jié)果驗(yàn)證

將壓比、效率與壓損系數(shù)等填入軟件對應(yīng)的輸入框，改變環(huán)境溫度、海拔高度及轉(zhuǎn)速計算得到對應(yīng)的仿真結(jié)果，并繪制成曲線。

(1) 應(yīng)用該模型計算了環(huán)境溫度在-40～60 ℃變化時的發(fā)電效率與發(fā)電量，并繪制成關(guān)系曲線(見圖2)，為了驗(yàn)證計算的準(zhǔn)確性，繪制了計算結(jié)果相對于標(biāo)準(zhǔn)數(shù)據(jù)的相對誤差隨溫度變化的關(guān)系曲線(見圖3)。當(dāng)環(huán)境溫度低于15 ℃時，微型燃?xì)廨啓C(jī)接近滿負(fù)荷運(yùn)行，發(fā)電效率比較高，由于系統(tǒng)自身控制作用發(fā)電量保持不變；當(dāng)環(huán)境溫度高于15 ℃時，發(fā)電效率隨溫度升高逐漸降低，在影響發(fā)電水平的因素中起主導(dǎo)作用，致使發(fā)電量隨環(huán)境溫度升高而下降。

圖2 發(fā)電量、發(fā)電效率與環(huán)境溫度的關(guān)系曲線

圖3 發(fā)電量、發(fā)電效率相對誤差

(2) 應(yīng)用該模型計算了海拔高度在0～4 500 m變化時微型燃?xì)廨啓C(jī)的發(fā)電量比(某一海拔高度微型燃?xì)廨啓C(jī)的發(fā)電量與海平面上微型燃?xì)廨啓C(jī)發(fā)電量的比值)，并繪制成關(guān)系曲線(見圖4)；為了驗(yàn)證計算的準(zhǔn)確性，繪制了計算結(jié)果與標(biāo)準(zhǔn)數(shù)據(jù)的相對誤差隨海拔高度變化的關(guān)系曲線(見圖5)。

圖4 發(fā)電量比隨海拔高度變化的關(guān)系曲線

圖5 發(fā)電量比的相對誤差

隨著海拔高度的增加，空氣密度會降低，進(jìn)入微型燃?xì)廨啓C(jī)的氣體流量會減少，這是發(fā)電機(jī)發(fā)電量降低的重要因素。

(3) 地面試車環(huán)境溫度為298.15 K，壓強(qiáng)為101 325 Pa,得到轉(zhuǎn)速與發(fā)電量的測試數(shù)據(jù)，利用參數(shù)轉(zhuǎn)換模塊得到轉(zhuǎn)速、發(fā)電量及扭矩在標(biāo)準(zhǔn)大氣條件下的換算數(shù)據(jù)(見表4)。

表4 測試數(shù)據(jù)與換算數(shù)據(jù)

為了得到微型燃?xì)廨啓C(jī)性能的優(yōu)劣，將換算后的發(fā)動機(jī)性能參數(shù)與標(biāo)準(zhǔn)大氣條件下發(fā)動機(jī)轉(zhuǎn)速特性曲線相比較(見圖6)。根據(jù)式(29)可知，發(fā)電量、扭矩的換算值與標(biāo)準(zhǔn)值的相對誤差相同，繪制相對誤差隨轉(zhuǎn)速變化的關(guān)系曲線(見圖7)。

圖6 發(fā)動機(jī)轉(zhuǎn)速特性曲線

圖7 性能參數(shù)的相對誤差

結(jié)果表明：仿真結(jié)果與廠家提供的試驗(yàn)數(shù)據(jù)基本一致，驗(yàn)證了該仿真模型的準(zhǔn)確性。根據(jù)轉(zhuǎn)速特性曲線可知在高轉(zhuǎn)速情況下，微型燃?xì)廨啓C(jī)性能參數(shù)與標(biāo)準(zhǔn)值相似，但在低轉(zhuǎn)速條件下，參數(shù)值差別較大，部件結(jié)構(gòu)需要進(jìn)一步優(yōu)化。

3 結(jié)語

通過將仿真結(jié)果與C30微型燃?xì)廨啓C(jī)標(biāo)準(zhǔn)性能參數(shù)進(jìn)行對比得到：(1)在不同環(huán)境溫度下，該模型計算的發(fā)電量、效率與標(biāo)準(zhǔn)值的相對誤差低于0.7%；(2)不同海拔高度下，該模型計算的發(fā)電量比與標(biāo)準(zhǔn)值的相對誤差低于3%，驗(yàn)證了模型的準(zhǔn)確性。從換算值與標(biāo)準(zhǔn)值繪制成的轉(zhuǎn)速特性曲線中可以得到，在高轉(zhuǎn)速條件下?lián)Q算值與標(biāo)準(zhǔn)值的相對誤差低于8%，在低轉(zhuǎn)速條件下?lián)Q算值與標(biāo)準(zhǔn)值的相對誤差大于60%，說明試驗(yàn)的微型燃?xì)廨啓C(jī)還需要進(jìn)一步優(yōu)化。

在不同環(huán)境溫度與海拔高度下，使用該仿真模型計算得到的發(fā)電量、發(fā)電效率等參數(shù)，可以用來分析外界因素對微型燃?xì)廨啓C(jī)性能影響的程度。通過換算模塊計算得到的發(fā)電量、扭矩等與標(biāo)準(zhǔn)參數(shù)之間的比較，可以反映燃?xì)廨啓C(jī)性能的優(yōu)劣，也可以為部件優(yōu)化指明方向。綜上所述，該模型的建立為新型燃?xì)廨啓C(jī)的性能分析和部件優(yōu)化設(shè)計提供了參考。

致謝：感謝北京惠斯特卓越科技發(fā)展有限公司提供的C30微型燃?xì)廨啓C(jī)相關(guān)數(shù)據(jù)。