基于動態(tài)增厚模型的甲烷/空氣預混旋流燃燒回火大渦模擬

夏 昊，王金華，衛(wèi)旭濤，韓 旺，吳鳳榮，黃佐華

(1.西安交通大學動力工程多相流國家重點實驗室，西安 710049；2.愛丁堡大學工程學院，愛丁堡EH8 3JL，蘇格蘭，英國；3.中國石油天然氣管道工程有限公司，廊坊131000)

回火是一種存在于燃氣輪機中常見的燃燒不穩(wěn)定問題.火焰向上游傳播會損傷燃燒器頭部，降低燃燒器壽命.主要的回火機理主要有自點火、邊界層回火、燃燒誘導渦破碎導致的回火、燃燒不穩(wěn)定性引起的回火等.燃氣輪機中燃燒器結構以鈍體旋流燃燒器為主.鈍體旋流燃燒器回火機理以邊界層回火和燃燒誘導渦破碎機理為主.邊界層回火主要是由于在邊界層內(nèi)，流場速度小于在壁面熄火距離以上范圍內(nèi)的火焰?zhèn)鞑ニ俣?，從而導致回?這個模型被稱為臨界速度梯度模型.但是，在湍流條件下，臨界速度梯度無法準確預測邊界層回火.Putnam[1]在臨界速度梯度模型的基礎上提出了無量綱的 Peclet 數(shù)模型.Peclet 數(shù)模型考慮燃料組分以及燃燒室內(nèi)壓力等其他因素對回火的影響.Kurdyumov 等[2]以及Eichler 等[3]提出了Damk?hler 數(shù)模型，在臨界速度梯度模型的基礎上，考慮了火焰速度和火焰厚度變化對回火預測的影響.Gruber 等[4]通過直接數(shù)值模擬(DNS)對低雷諾數(shù)湍流槽道流中的邊界層回火進行了模擬，結果表明在湍流邊界層中，火焰誘發(fā)了逆流區(qū)的產(chǎn)生，而逆流區(qū)又促進火焰的上游傳播，當逆流區(qū)消失時，火焰的上游傳播速度減小.邊界層回火是流場與火焰相互作用的結果，而Peclet 數(shù)模型和Damk?hler 數(shù)模型都忽略了火焰與流場的相互作用，不能夠真實地描述邊界層回火中發(fā)生的物理過程.燃燒誘導渦破碎是旋流火焰發(fā)生回火的主要原因.燃燒誘導渦破碎最先由Fritz 等[5]在研究旋流燃燒器預混管中的回火現(xiàn)象時提出.在實驗中，由于火焰始終沿著中心線向上游傳播，而穩(wěn)態(tài)時中心流的軸向速度比火焰?zhèn)鞑ニ俣雀叩枚?，同時燃燒室內(nèi)并未出現(xiàn)燃燒不穩(wěn)定現(xiàn)象.因此，他們提出一種新的機理，即燃燒誘導渦破碎促使旋流火焰發(fā)生回火.反應區(qū)域的斜壓轉矩在回流泡中產(chǎn)生了負渦量，進而產(chǎn)生負的軸向速度，回流泡到達了一個更上游的位置，同時整個流場與反應達到一個穩(wěn)定狀態(tài)；隨著當量比的增加，反應區(qū)會變得更加緊湊，火焰跟隨著回流區(qū)運動.此時，由斜壓轉矩產(chǎn)生的負的周向渦量也會增加，當量比超過臨界值時，就可能引起回火.但是對于鈍體旋流火焰，由于中心鈍體的存在，中心區(qū)域流場受鈍體壁面影響，僅僅用燃燒誘導渦破碎無法完全描述回火過程.

目前，國內(nèi)外已開展一些對回火的數(shù)值模擬研究.Lietz 等[6]通過大渦模擬對純甲烷和甲烷/氫氣鈍體旋流火焰回火過程進行預測.雖然在宏觀結構上能夠較好預測，但是在回火速度等方面仍存在差異.雖然鈍體旋流火焰回火機理方面已有一些研究，但對回火過程中邊界層回火和燃燒誘導渦破碎機理的作用機理認識仍不夠清楚.因此需要進一步實驗和模擬研究鈍體旋流火焰回火過程中火焰-流場-壁面相互作用機理.筆者通過大渦模擬結合增厚火焰模型對得克薩斯大學旋流燃燒器回火過程進行模擬.國內(nèi)關于鈍體旋流火焰回火機理的研究相對較少.Xu 等[7]采用計算和實驗相結合的方法研究甲烷與富氫及純氫合成氣貧預混燃燒的回火極限，并分析了不同燃料、不同穩(wěn)定方式的回火特性.數(shù)值和實驗結果均表明，回火極限可以關聯(lián)為Peclet 數(shù)模型.Tian 等[8]選用雷諾應力模型耦合有限速率/渦黏模型對二維氫氣軸對稱旋流預混火焰進行模擬，研究了預混區(qū)出口周向渦量對燃燒誘導渦破碎回火的影響.綜合來看目前國內(nèi)針對富氫及純氫燃料在燃氣輪機背景鈍體旋流燃燒回火方面系統(tǒng)性實驗研究和高精度模擬研究非常缺乏.對于回火的非定常性，大渦模擬由于具有較低的計算成本和較高的湍流解析度而被廣泛使用.本文利用增厚火焰模型大渦模擬耦合甲烷兩步機理對甲烷預混火焰在得克薩斯大學旋流燃燒器中的回火過程進行模擬.分別采用絕熱邊界條件及熱損失邊界條件對回火過程進行模擬并對比熱損失對火焰結構及回火速度的影響，并對預混鈍體旋流火焰回火過程進行分析.

1 數(shù)值模擬方法與計算設置

1.1 數(shù)值模擬方法

本文采用大渦模擬的方法對甲烷預混旋流火焰回火進行模擬.湍流模型采用Dynamic Smagorinsky 模型[9].燃燒模型采用動態(tài)增厚模型.增厚F 倍的火焰面可以通過求解方程(1)獲得：

式中：xk為組分k 的摩爾分數(shù)；Dkj為組分j 和 k 的二元擴散系數(shù)，通過Chapman-Enskog 表達式計算.對式(1)濾波后得到式(3)：

式中：n 為火焰面增厚之后包含的網(wǎng)格數(shù)，本文n＝5；Δmesh為網(wǎng)格尺寸；δL為層流火焰厚度；f(Ω) ＝tanh(αΩ)，α控制過渡區(qū)域厚度[12].最終，組分方程重整為

1.2 計算設置

本文模擬的對象為德克薩斯大學旋流燃燒器，其主要結構如圖1 所示.混合段長度為150 mm，出口直徑為52 mm.混合段中心為直徑25.4 mm 的圓柱形中心體，中心體底部安裝軸向旋流器.燃燒室結構為圓柱形，直徑為100 mm，高度為150 mm.軸向旋流器由8 個與軸向呈60°夾角的葉片構成，旋流數(shù)約為0.9[13].

圖1 德克薩斯大學旋流燃燒器[13]Fig.1 UT swirl combustor[13]

計算域設置如圖2 所示，計算域入口為燃燒器預混段入口，在出口上方設置延伸段以防止出口回流.入口處速度邊界為預先計算的湍流流場.初始時刻預混段內(nèi)為預混氣，在出口上方設置高溫區(qū)域.壁面為無滑移邊界.本文計算的回火工況為當量比0.8，平均速度為2.5 m/s，化學反應機理采用包含6 組分、兩步反應的2S_CH4_BFER[14].在計算回火工況之前，首先計算當量比為0.6 的穩(wěn)定工況.之后將上游當量比改至回火當量比.

圖2 計算域及邊界設置Fig.2 Computational domain and boundary settings

本文采用的求解器基于開源 CFD 平臺OpenFOAM 編寫.壓力-速度耦合通過PIMPLE 算法計算.采用隱式歐拉格式對時間推進進行求解.整個計算域網(wǎng)格數(shù)約為1 千萬.對混合段內(nèi)尤其是鈍體壁面附近的網(wǎng)格進行加密，壁面附近網(wǎng)格尺度約為0.4 mm.通過網(wǎng)格敏感性評估，發(fā)現(xiàn)當前網(wǎng)格可以對流場統(tǒng)計預測.

2 結果與分析

2.1 冷態(tài)流場對比

由于回火主要發(fā)生在預混段內(nèi)，本文主要對比預混段內(nèi)不同高度的冷態(tài)流場.分別提取預混段出口上游5 mm、29 mm、54 mm、73 mm 處的軸向平均速度與周向平均速度與實驗值進行對比，如圖3 所示.圖中虛線為模擬值，散點為實驗值.由圖3 可知，周向速度峰值出現(xiàn)在中心體附近，而軸向速度峰值出現(xiàn)在外壁面附近.除了在預混段出口上游5 mm 處靠近中心體處計算得到的平均軸向速度和周向速度略微低于實驗值以及周向速度在上游73 mm 處整體略低于實驗值，整體冷態(tài)流場可以較好地反映實驗流場.由于下文火焰結構分析區(qū)域主要集中在遠離出口處以及旋流器的中間位置，因此可以認為這些偏差對后續(xù)結果分析影響不大.

圖3 冷態(tài)流場對比Fig.3 Cold flow comparison

2.2 溫度邊界對回火速度的影響

中心體壁面分別采用絕熱以及考慮熱損失兩種熱邊界條件對回火過程進行計算.考慮熱損失通過在溫度邊界設置壁面導熱系數(shù)，外部環(huán)境換熱系數(shù)及外界環(huán)境溫度實現(xiàn).本文中考慮熱損失工況環(huán)境溫度設置為350 K，環(huán)境換熱系數(shù)為10 W/(m2·K)，壁面導熱系數(shù)為10 W/(m·K).圖4 為兩種不同熱邊界條件下計算得到的回火過程中瞬時熱釋放.在絕熱邊界條件下，火焰在中心體壁面附近熱釋放較強，而在考慮熱損失情況下，火焰在壁面處熱釋放減弱.這主要是由于額外的熱損失使壁面溫度低于絕熱邊界條件下的壁面溫度，從而導致反應減弱.而在絕熱邊界條件下，由于邊界層內(nèi)速度低，未燃氣駐留時間更長，且沒有熱損失，因此熱釋放率更高.此外，在考慮熱損失后，火焰整體熱釋放率低于絕熱邊界條件下得到的火焰熱釋放率.

圖4 不同熱邊界條件下火焰熱釋放分布Fig.4 Flame heat release rate under different thermal boundary conditions

本文提取溫度1 400 K 等勢面以表征進展變量為0.7 的火焰面，如圖5 所示.在預混火焰回火過程中，火焰面出現(xiàn)火舌結構.由于存在熱損失，火舌頂點徑向位置與中心體表面有一定距離，而在絕熱邊界條件下火舌完全附著在中心體表面.根據(jù)一定時間內(nèi)火舌頂點傳播的軸向距離，可以計算得到平均回火速度.不同邊界條件下計算得到的平均回火速度及實驗得到的平均回火速度如表1 所示，在絕熱邊界條件下，平均回火速度約為 0.71 m/s，遠高于實驗值0.58 m/s[15].在考慮熱損失后，平均回火速度為0.57 m/s，接近于實驗值.考慮熱損失后，回火速度準確度相比于絕熱邊界條件有明顯改善.

圖5 瞬時溫度1 400 K 等勢面Fig.5 Instantaneous isosurface of 1 400 K

表1 不同熱邊界條件下回火速度Tab.1 Flashback speed under different thermal boundary conditions

2.3 回火過程分析

基于上一部分的分析，本部分根據(jù)考慮熱損失得到的計算結果對回火過程進行進一步的分析.圖6為回火過程中火焰面瞬時結構的變化過程.在回火過程后期，由于燃燒器出口下游為高溫已燃氣，溫度較高，仍然能得到1 400 K 等勢面，但已沒有熱釋放區(qū)域.由圖可知，在回火過程中，火焰出現(xiàn)火舌結構，火舌隨預混段內(nèi)部流場作周向運動，并沿中心體向下傳播.圖7 為回火過程中某一瞬時t 時刻截面軸向速度及二維渦量分布.圖中黑線表示進展變量0.7 的火焰面.圖7(b)中渦量場分布為當?shù)販u量在截面法向上的分量.在回火過程中，上游受火焰面的作用形成的渦量場促進了邊界層內(nèi)的回流，使火焰更容易向上游傳播.而當火焰面?zhèn)鞑サ缴嫌螀^(qū)域，當?shù)亓鲌鍪艿交鹧嬗绊懚a(chǎn)生回流，進一步促進了火焰向上游傳播.因此，火焰不斷向上游傳播.

圖6 回火過程Fig.6 Flashback process

圖7 軸向速度場與二維渦量場Fig.7 Axial flow velocity and two-dimensional vorticity

提取每一時刻的火焰面頂點，可以得到圖8 所示的火舌頂點軌跡.由圖8(a)可得，火舌周向運動方向與預混段內(nèi)部流場方向相同，沿中心體壁面向下傳播.在火舌傳播一定距離后，在某一位置會出現(xiàn)新的火舌，因此火舌頂點軌跡會出現(xiàn)斷點.火舌軸向運動軌跡如圖8(b)所示.在回火過程中，火舌軸向運動并不是勻速的，從出口至上游150 mm 處，火舌軸向傳播速度逐漸減小.藍線為對軌跡線性擬合得到的平均速度.當火焰?zhèn)鞑サ叫髌鞲浇鼤r，火舌頂點幾乎不再向下傳播，僅沿著鈍體周圍作圓周運動.這主要由于當火焰?zhèn)鞑サ叫髌鞲浇鼤r，截面積減小，局部流速高于火焰?zhèn)鞑ニ俣?鈍體旋流火焰主要的回火機理為邊界層回火與旋流燃燒誘導渦破碎回火.火舌作圓周運動表明此時旋流燃燒誘導渦破碎占主導.由圖8(a)可得，火舌的軌跡并非完全連續(xù).當火舌傳播一定距離后，在某一位置形成新的火舌繼續(xù)傳播.而在火舌傳播過程中，有些時刻火舌頂點圓周運動并不明顯，表明此時火舌頂點相對來流的運動較為明顯，受邊界層回火機理影響較大.

圖8 回火火舌頂點三維軌跡及軸向軌跡Fig.8 Three dimensional trajectory and axial trajectory of flame tongue

3 結論

本文通過動態(tài)增厚火焰大渦模擬結合兩步機理對德克薩斯大學旋流燃燒器中甲烷預混火焰回火進行模擬，分別計算了絕熱條件下及考慮熱損失下的回火過程.主要結論如下：

(1) 動態(tài)增厚模型耦合兩步反應機理能夠對甲烷預混旋流火焰邊界層回火宏觀結構及回火速度進行有效預測.熱損失對回火速度預測有很大影響.相比于絕熱邊界條件，考慮熱損失后回火速度的預測更接近實驗值.而采用絕熱邊界會使回火速度提高約22%.

(2) 預混鈍體旋流火焰回火過程中會形成火舌結構.火舌順著預混段內(nèi)流場方向作周向運動，同時沿著中心體壁面向上游傳播.

(3) 對于預混鈍體旋流火焰，火舌結構的傳播并非完全連續(xù)，而是會出現(xiàn)新的火舌.當邊界層回火機理占主導時，火舌頂點相對來流傳播變明顯，周向運動會減弱.