基于火焰?zhèn)鬟f函數(shù)的旋流預(yù)混火焰穩(wěn)定性分析

楊 旸, 余志健

(1. 中國(guó)科學(xué)院工程熱物理研究所 南京未來(lái)能源系統(tǒng)研究院，南京 210000；2. 中國(guó)科學(xué)院工程熱物理研究所 先進(jìn)燃?xì)廨啓C(jī)實(shí)驗(yàn)室，北京 100190；3. 中國(guó)科學(xué)院大學(xué)，北京 100190)

貧預(yù)混燃燒是現(xiàn)代燃?xì)廨啓C(jī)燃燒室實(shí)現(xiàn)低NOx排放的主要途徑之一。但該燃燒室燃燒處于貧預(yù)混狀態(tài)及火焰筒取消了各種氣膜孔，容易產(chǎn)生熱聲不穩(wěn)定問(wèn)題[1]。根據(jù)瑞利準(zhǔn)則，當(dāng)壓力波動(dòng)相位與燃燒室內(nèi)體積積分熱釋放率脈動(dòng)相位差小于90°，且系統(tǒng)的阻尼不大于注入系統(tǒng)的能量時(shí)，則會(huì)產(chǎn)生熱聲不穩(wěn)定[2]。根據(jù)先前的研究[3]，相位差可以轉(zhuǎn)換為從燃料噴射位置開始到達(dá)火焰鋒面處的燃料飛行時(shí)間。因此，火焰形狀和混合區(qū)長(zhǎng)度是產(chǎn)生熱聲不穩(wěn)定的兩個(gè)關(guān)鍵因素。為了控制火焰形狀，燃燒室頭部出口形狀是除旋流數(shù)之外的重要幾何控制參數(shù)。這是因?yàn)槿紵翌^部出口形狀會(huì)影響火焰根部的發(fā)散角，進(jìn)而導(dǎo)致火焰形狀變?yōu)閂形、M形或圓錐形等。但國(guó)內(nèi)外文獻(xiàn)對(duì)該結(jié)構(gòu)的影響特性報(bào)道較少。

同時(shí)，在燃?xì)廨啓C(jī)燃燒室設(shè)計(jì)階段如何采用數(shù)值方法分析不同燃燒室結(jié)構(gòu)下火焰動(dòng)態(tài)特性及預(yù)判火焰是否穩(wěn)定，這對(duì)燃燒室結(jié)構(gòu)篩選和確定具有一定作用。雷諾平均Naiver-Stokes(RANS)方法[4]和大渦模擬(LES)[5]是燃燒室設(shè)計(jì)階段中廣泛使用的兩種數(shù)值方法。有了準(zhǔn)確的燃燒模型，可以使用RANS方法在設(shè)計(jì)階段預(yù)測(cè)燃燒室總體性能，而LES方法可以更準(zhǔn)確地解析流動(dòng)-火焰的瞬態(tài)相互作用，例如吹熄，火焰不穩(wěn)定性等。但對(duì)于LES來(lái)說(shuō)，由于計(jì)算的物理時(shí)間較短和缺乏相應(yīng)阻抗邊界條件，在燃燒室的設(shè)計(jì)階段通常采用不可壓縮LES計(jì)算[6]。

不可壓LES計(jì)算忽略燃燒室出口的反向傳播波，且無(wú)法獲得潛在的自激熱聲特性。盡管如此，該方法仍可以結(jié)合系統(tǒng)辨識(shí)技術(shù)(system identification)，獲得火焰相應(yīng)的頻率響應(yīng)信息。通過(guò)在入口處施加適當(dāng)?shù)耐獠考?lì)，可以基于放熱速率的響應(yīng)時(shí)間序列，通過(guò)辨識(shí)技術(shù)提取火焰?zhèn)鬟f函數(shù)(FTF)或火焰描述函數(shù)(FDF)[7-8]。該方法將火焰視為一個(gè)線性時(shí)不變系統(tǒng)，該系統(tǒng)在弱激勵(lì)下可由辨識(shí)獲得火焰動(dòng)態(tài)響應(yīng)參數(shù)。該火焰?zhèn)鬟f函數(shù)可添加至低階熱聲網(wǎng)絡(luò)模型中或三維亥姆霍茲聲學(xué)模型中預(yù)估燃燒系統(tǒng)的穩(wěn)定性[9-12]，或單獨(dú)分析預(yù)判火焰發(fā)生不穩(wěn)定的風(fēng)險(xiǎn)。分析火焰?zhèn)鬟f函數(shù)的方法已成功應(yīng)用于西門子[13]及阿爾斯通[14]等燃?xì)廨啓C(jī)燃燒室設(shè)計(jì)中，并在國(guó)際學(xué)術(shù)界得到廣泛的報(bào)道。

本文提取和分析了不同燃燒室頭部結(jié)構(gòu)下的火焰形狀及火焰?zhèn)鬟f函數(shù)。采用試驗(yàn)設(shè)計(jì)DoE方法對(duì)四個(gè)頭部關(guān)鍵幾何控制參數(shù)進(jìn)行了分析。采用不可壓大渦模擬(LES)獲得精確的旋流火焰結(jié)構(gòu)。通過(guò)在入口處施加弱激勵(lì)，獲得燃燒室熱釋放率動(dòng)態(tài)響應(yīng)，進(jìn)而辨識(shí)出相應(yīng)火焰?zhèn)鬟f函數(shù)，并分析了各頭部結(jié)構(gòu)下火焰的響應(yīng)頻率帶及是否有不穩(wěn)定風(fēng)險(xiǎn)。最終分析了燃燒室頭部幾何結(jié)構(gòu)與火焰形狀及火焰?zhèn)鬟f函數(shù)的對(duì)應(yīng)關(guān)系，為設(shè)計(jì)燃燒不穩(wěn)定發(fā)生風(fēng)險(xiǎn)低的低污染燃燒室提供理論指導(dǎo)。

1 數(shù)值方法

1.1 CFD計(jì)算

1.1.1 幾何及邊界條件

本文對(duì)一個(gè)單頭部模型燃燒室進(jìn)行大渦模擬及火焰?zhèn)鬟f函數(shù)提取。采用的旋流器為自行設(shè)計(jì)的、可采用SLM(選擇性激光熔融)一體化打印成型的且適用于E級(jí)F級(jí)燃?xì)廨啓C(jī)燃燒室的軸向旋流器。該軸向旋流器葉片形狀如圖1所示，有8個(gè)長(zhǎng)弦葉片，葉片具有一定厚度可作為燃料通道，葉片壓力面開用一個(gè)燃料孔，吸力面開有兩個(gè)燃料孔，燃料孔直徑選為1.2 mm。

圖1 軸向旋流器

模型燃燒室完整域如圖2所示，實(shí)際計(jì)算的計(jì)算域?yàn)?5°扇形，采用周期性邊界條件降低數(shù)值計(jì)算成本。氣流方向?yàn)閺淖蟮接摇交靺^(qū)長(zhǎng)度為l，中心體末端半徑為R1_TE，擴(kuò)張角為α。旋流氣流隨擴(kuò)張角向燃燒室擴(kuò)散時(shí)，當(dāng)氣流離心力與徑向壓力梯度平衡時(shí)，根據(jù)渦破碎理論形成回流區(qū)。回流區(qū)的形成通常由標(biāo)準(zhǔn)旋流數(shù)控制，該旋流器速度基旋流數(shù)為0.4。邊界條件如表1所示。當(dāng)前幾何條件下，當(dāng)量比為0.514。

圖2 計(jì)算域和幾何結(jié)構(gòu)參數(shù)

表1 邊界條件

1.1.2 DoE設(shè)計(jì)表

選擇了控制燃燒室頭部結(jié)構(gòu)的四個(gè)關(guān)鍵幾何設(shè)計(jì)參數(shù)，如圖2所示，表述如下：

(1) 中心體末端到燃燒室前墻的距離(d)；

(2) 擴(kuò)張角(α)，控制氣流進(jìn)入燃燒室的方向；

(3) 擴(kuò)張比，燃燒室截面和頭部截面面積比，選擇燃燒室截面半徑R4為控制參數(shù)；

(4) 中心體末端半徑(R1_TE)。

采用中心體半徑R1作為歸一化尺寸，進(jìn)行二水平全因子實(shí)驗(yàn)設(shè)計(jì)(DoE)，共需計(jì)算16個(gè)工況，如表2所示。

表2 二水平全因子試驗(yàn)設(shè)計(jì)表

1.1.3 網(wǎng)格劃分和網(wǎng)格分辨率

網(wǎng)格采用Star ccm+生成的混合網(wǎng)格，邊界層采用三層棱柱層網(wǎng)格，其它區(qū)域采用多面體網(wǎng)格。燃料噴嘴、旋流器區(qū)、摻混區(qū)及燃燒室上游火焰位置區(qū)網(wǎng)格分辨率分別為0.2 mm、0.35 mm、0.5 mm和0.5 mm，以捕捉流動(dòng)火焰細(xì)節(jié)。燃燒室下游網(wǎng)格分辨率為1.1 mm。全局及局部網(wǎng)格劃分如圖3所示，總網(wǎng)格數(shù)約為470萬(wàn)，網(wǎng)格質(zhì)量大于85%。

圖3 全局及局部網(wǎng)格劃分

本文采用基于能量指數(shù)法的后驗(yàn)指數(shù)[15]來(lái)判斷劃分的網(wǎng)格是否足夠用于LES計(jì)算。當(dāng)求解的湍流動(dòng)能達(dá)到總湍流動(dòng)能的85%以上時(shí)，認(rèn)為網(wǎng)格求解精度達(dá)到[15]。由于采用邊界中心差分格式對(duì)動(dòng)量方程進(jìn)行離散，所以由數(shù)值格式引起的湍流動(dòng)能耗散忽略不計(jì)。截面后檢驗(yàn)指數(shù)M分布云圖如圖4所示，可以看出大部分區(qū)域后檢驗(yàn)指數(shù)均大于85%，僅燃燒室下游中心局部區(qū)域解析不夠，由于本研究集中在燃燒室頭部出口的火焰區(qū)，所以可以認(rèn)為該網(wǎng)格劃分可滿足本文研究要求。

圖4 截面后檢驗(yàn)指數(shù)M分布云圖

1.1.4 數(shù)值模型

不可壓大渦模擬的質(zhì)量、動(dòng)量和物質(zhì)傳輸Favre濾波守恒方程如下：

(1)

(2)

先進(jìn)行熱態(tài)RANS計(jì)算以獲得收斂場(chǎng)，而后激活LES計(jì)算。LES進(jìn)行0.1 s物理時(shí)間后(10倍以上的流通時(shí)間)，在空氣入口處施加激勵(lì)，激勵(lì)施加持續(xù)時(shí)間為0.25 s，期間同時(shí)采集數(shù)據(jù)。0.25 s的采集時(shí)間進(jìn)行頻域分析可獲得頻率分辨率為4 Hz。通過(guò)湍流火焰模型(TFC)計(jì)算湍流火焰速度[16-17]，該模型基于火焰團(tuán)內(nèi)小尺度湍流平衡假定，建立了湍流火焰速度與大尺度湍流參數(shù)關(guān)聯(lián)，適用于工業(yè)預(yù)混系統(tǒng)中的褶皺火焰。甲烷燃燒采用GRI 3.0反應(yīng)機(jī)理[18]。

1.2 火焰?zhèn)鬟f函數(shù)

火焰?zhèn)鬟f函數(shù)計(jì)算基于部分預(yù)混火焰是線性系統(tǒng)的假設(shè)。在線性時(shí)不變系統(tǒng)中，如果系統(tǒng)由弱激勵(lì)信號(hào)激勵(lì)，則可以識(shí)別系統(tǒng)的頻率響應(yīng)[19]。為了克服較短的LES仿真時(shí)間，采用離散隨機(jī)二進(jìn)制信號(hào)(DRBS)進(jìn)行速度入口激勵(lì)[5]。該辨識(shí)過(guò)程也被認(rèn)為是Winer-Hopf方法，其優(yōu)點(diǎn)是與掃頻激勵(lì)方法相比，單個(gè)LES仿真可以提供整個(gè)頻帶范圍，激勵(lì)幅值選為5%。圖5為激勵(lì)所用的離散隨機(jī)二進(jìn)制信號(hào)時(shí)域圖。

圖5 離散隨機(jī)二進(jìn)制信號(hào)(DRBS)

方程(3)為火焰?zhèn)鬟f函數(shù)定義式:

(3)

圖6 火焰?zhèn)鬟f函數(shù)計(jì)算步驟

2 結(jié)果與分析

2.1 燃燒場(chǎng)分析

上述16個(gè)DoE工況計(jì)算后的產(chǎn)物生成速率分布云圖如圖7所示，其中黑線表示軸向速度為零的等值線。由圖7中軸向速度零等值線可知，在所有工況下燃燒室內(nèi)部均可有效形成中心回流區(qū)和角回流區(qū)。形成的火焰結(jié)構(gòu)可分為三部分：外層火焰、內(nèi)層火焰和火焰前沿。外層火焰的根部顯示出較高的放熱速率，而內(nèi)層火焰的根部為相反的結(jié)果。同時(shí)若火焰不一直沿徑向傳播，則中心回流區(qū)可能會(huì)變形。圖7中，如工況6，12，13和15的產(chǎn)物生成速率輪廓所示，其火焰內(nèi)層未錨定在中心體末端附近。所有工況中，火焰外層均維持在頭部與燃燒室之間的連接點(diǎn)上。根據(jù)產(chǎn)物生成速率所展示的火焰形狀輪廓，以上16種工況可劃分為3種火焰類型，如表3所示。其中，A型火焰末端沖擊到燃燒室外殼上，如工況1，4，7和8所示。C型火焰與中心回流區(qū)有強(qiáng)烈的相互作用，如工況6,12,13,15和16所示。B型火焰同時(shí)具有A型和C型火焰的復(fù)合特性，如工況2，3，5,9,10,11和14所示。

圖7 各工況產(chǎn)物生成速率分布云圖(黑線表示軸向速度為零的等值線)

表3 火焰類型劃分

2.2 DoE分析

火焰形狀與特定幾何控制參數(shù)之間的關(guān)系可通過(guò)統(tǒng)計(jì)學(xué)中的方差分析(ANOVA)方法確定。表4列出了不同類型火焰的全因子回歸模型所對(duì)應(yīng)的p值，包括單個(gè)幾何控制參數(shù)的p值和參數(shù)間相互作用關(guān)系p值。該表中還列出了各工況下燃燒室的總壓損失，以綜合考慮燃燒室氣動(dòng)壓損特性。當(dāng)p值小于0.01時(shí)，表明此幾何參數(shù)對(duì)該類火焰形狀有非常大的影響。由表4可知，對(duì)于A型火焰來(lái)說(shuō)，沒(méi)有明顯的幾何控制參數(shù)。對(duì)于C型火焰，最重要的幾何控制參數(shù)為擴(kuò)張角α，且較小的擴(kuò)張角α將有利于C型火焰生成。B型火焰主要受擴(kuò)張角α與燃燒室半徑R4之間的相互作用影響，且中心體末端半徑R1_TE和燃燒室半徑R4也單獨(dú)對(duì)該類火焰有作用關(guān)系。燃燒室總壓損失大小主要受中心體末端半徑R1_TE影響。總壓損失隨R1_TE的增加而增加，表明燃燒室總氣動(dòng)阻力增加。R1_TE的增加使得混合區(qū)整體流通截面積減小，導(dǎo)致混合區(qū)域流速增加，且壓降與流速成正相關(guān)，因而流速增加后，壓降得到升高。

表4 不同類型火焰全因子回歸模型的p-value

總而言之，B型和C型火焰形成主要受擴(kuò)張角α及其與燃燒室半徑R4的相互作用控制。總壓損失主要受中心體末端尺寸R1_TE控制。

2.3 火焰?zhèn)鬟f函數(shù)(FTF)分析

可以從LES結(jié)果采用上述的激勵(lì)響應(yīng)方法辨識(shí)反應(yīng)火焰響應(yīng)特性的火焰?zhèn)鬟f函數(shù)。該火焰?zhèn)鬟f函數(shù)可添加至低階熱聲網(wǎng)絡(luò)模型中或三維亥姆霍茲聲學(xué)模型中預(yù)估燃燒系統(tǒng)的穩(wěn)定性。火焰?zhèn)鬟f函數(shù)還可以擴(kuò)展為弱非線性火焰描述函數(shù)，進(jìn)而預(yù)測(cè)燃燒室熱聲極限環(huán)特性。燃燒室設(shè)計(jì)階段也可單獨(dú)分析火焰?zhèn)鬟f函數(shù)，通過(guò)FTF的增益來(lái)判斷火焰是否有不穩(wěn)定的風(fēng)險(xiǎn)，F(xiàn)TF的相位信息反應(yīng)火焰的時(shí)滯。當(dāng)增益高于1時(shí)，一般可認(rèn)為火焰存在發(fā)生燃燒不穩(wěn)定的高風(fēng)險(xiǎn)。

圖8為16個(gè)工況下的火焰?zhèn)鬟f函數(shù)增益和相位信息。當(dāng)頻率高于400 Hz時(shí)，火焰?zhèn)鬟f函數(shù)增益小于1，火焰起低通濾波器的作用，火焰僅對(duì)低頻有響應(yīng)，而對(duì)高頻響應(yīng)不明顯。某些工況下，在零頻率附近具有較高的增益，而其它情況則顯示較低的幅值。表5列出了增益大于1的頻率帶，如41 Hz、82 Hz、122 Hz和204 Hz。特別工況火焰?zhèn)鬟f函數(shù)增益存在明顯的特征峰，如工況9的204 Hz、工況15的363 Hz及工況16的45 Hz。

表5 不同工況下特征頻率帶信息

圖8 各工況火焰?zhèn)鬟f函數(shù)增益和相位

比較工況9和工況11，可以看出中心體后緣的直徑(R1_TE/R1)和中心體末端到燃燒室前墻的距離(d)可能不會(huì)強(qiáng)烈影響參數(shù)該頻段的火焰?zhèn)鬟f函數(shù)。比較工況3與工況11，表明燃燒室的半徑R4對(duì)火焰響應(yīng)頻率特性影響不可忽略。比較圖8火焰?zhèn)鬟f函數(shù)結(jié)果與表3中所示的火焰類別劃分結(jié)果，可以看出兩種結(jié)果之間沒(méi)有太大的相關(guān)性。

總體來(lái)看，工況5，7和8只有一個(gè)響應(yīng)頻率帶，且位于100 Hz以下；工況9，14和15也均只有一個(gè)頻率響應(yīng)帶，位于200 Hz以上；其它工況響應(yīng)頻率帶均有2個(gè)以上。在進(jìn)行燃燒室和阻尼器的結(jié)構(gòu)設(shè)計(jì)時(shí)，可結(jié)合火焰?zhèn)鬟f函數(shù)結(jié)果，優(yōu)化設(shè)計(jì)結(jié)構(gòu)與方案。

3 結(jié)論

(1) 對(duì)控制燃燒室頭部形狀的四個(gè)關(guān)鍵幾何控制參數(shù)進(jìn)行DoE大渦模擬數(shù)值分析，可得生成的火焰形狀分為三種類型。A型火焰：火焰沖擊燃燒室外殼；C型火焰與中心回流區(qū)相互摻混；B型火焰：同時(shí)包含A和C種火焰特性。

(2) 對(duì)DoE結(jié)果方差分析可得：對(duì)A型火焰來(lái)說(shuō)，無(wú)明顯的幾何控制參數(shù)；對(duì)C型火焰來(lái)說(shuō)，幾何控制參數(shù)主要為擴(kuò)張角α；B型火焰主要受α與燃燒室半徑R4之間的相互作用的影響；燃燒室總壓損失主要受中心體末端半徑R1_TE影響。

(3) 火焰?zhèn)鬟f函數(shù)結(jié)果表明，火焰僅對(duì)低頻有響應(yīng)，而對(duì)高頻響應(yīng)不明顯。工況5，7和8只有一個(gè)響應(yīng)頻率帶，且位于100 Hz以下；工況9，14和15也均只有一個(gè)頻率響應(yīng)帶，位于200 Hz以上；其它工況響應(yīng)頻率帶均有2個(gè)以上。