龐 杰, 張 健
(清華大學 工程力學系,北京 100084)
許多條件下燃燒火焰具有較低的流向速度或Froude數(Fr=u2/g l),如各種層流火焰、火災中可燃物表面附近的湍流火焰、旋風爐和四角切圓鍋爐爐膛中的湍流火焰等.在這類火焰中,氣體流動和燃燒過程以及火焰結構與輸運特性等都會受到浮力的作用.對這類火焰的理論描述和數值模擬需要考慮浮力的影響.
已有學者對受浮力作用的湍流火焰進行了理論與數值模擬研究.Jeng等[1]采用k-ε-g模型對靜止空氣環境中受浮力作用的甲烷湍流擴散火焰進行了數值模擬,得到的氣體溫度場和組分體積分數場基本與試驗數據相符.Liu等[2]分別應用代數Reynolds應力/熱流通量模型和低Re數k-ε模型對受浮力作用的天然氣湍流擴散火焰進行了數值模擬,前者得到的氣體速度場和溫度場與試驗數據相符.Xin等[3]采用快速反應的混合物分數湍流燃燒模型對靜止空氣中受浮力作用的甲烷湍流火焰進行了大渦模擬,得到的氣體速度、溫度和混合物分數分布基本與試驗數據相符.
旋流燃燒在各種工程燃燒反應裝置(如四角切圓鍋爐和旋風爐)中有許多應用[4].然而,在現有對受浮力作用的湍流火焰的數值模擬研究中,較少涉及浮力對湍流旋流火焰的影響.為合理地描述并預測湍流旋流火焰的特性,筆者對低Fr條件下的湍流旋流火焰進行了數值模擬,研究了浮力對具有不同初始切向動量或旋流數的湍流旋流火焰的影響,并將模擬結果與試驗數據進行了對比.
對受浮力作用的湍流旋流燃燒的模擬,已有的湍流模型包括浮力修正的k-ε模型、Reynolds應力輸運方程模型和代數Reynolds應力模型等[5],但后兩種模型應用尚不多.筆者采用浮力修正的k-ε模型[6]對受浮力作用的湍流旋流火焰中湍流輸運進行了模擬,除在氣體平均動量方程中引入浮力作用項外,還在湍流參數如湍流動能k及其耗散率ε方程中引入了浮力修正項.
浮力修正的k-ε模型中k方程和ε方程的具體形式分別為:


對湍流燃燒采用渦團耗散(EDC)模型,對氣體輻射傳熱采用四熱流通量模型.在軸對稱坐標系中建立受浮力作用的湍流旋流火焰的平均控制方程組,略去方程中各一階矩上的平均號“ˉ”,各控制方程的通用形式可表示為:

當 αφ=ρ,φ=1、u、v 、w 、k、ε、h 和Y s時 ,式(4)分別對應連續、軸向動量、徑向動量、切向動量、湍能、湍能耗散率、焓和組分質量分數方程;當αφ=0,φ=Rx和Ry時,式(4)分別對應軸向與徑向輻射熱流方程.
為求解上述受浮力作用的湍流旋流火焰平均控制方程組,采用有限容積法[5]對各控制方程進行離散化.在求解區域上布置交錯網格系,各控制方程中的對流-擴散項采用混合格式離散,對源項進行負斜率線性化處理,k和ε方程浮力產生項中的密度梯度采用中心格式離散.對各變量離散化方程組的求解采用沿流動方向從上游到下游的逐線 TDMA低松弛迭代算法,對壓力與速度的耦合求解采用Simp le算法.
對旋流燃燒室內受浮力作用的湍流火焰進行了數值模擬.圖1所示為圓柱形旋流燃燒室,該燃燒室豎直向上布置.在燃燒室進口處,空氣分為旋流一次風和直流二次風,分別通過同軸的內、外環形通道進入燃燒室內,氣體燃料則通過同軸的中心管噴入燃燒室.中心燃料管、旋流一次風管、直流二次風管和燃燒室的內直徑分別為 8mm、40 mm、58 mm和160mm,外直徑分別為 14 mm、44 mm、68 mm 和236 mm,燃燒室長度為1 000 mm.

圖1 圓柱形旋流燃燒室示意圖Fig.1 Schematic diagram of the cylindrical swirl combustion chamber
數值模擬采用自行編制的湍流多相流動與燃燒計算程序CCVC進行,選取與文獻[7]和文獻[8]試驗條件相同的三組工況.表1給出了三組工況下燃燒室各進口的參數,包括燃料、一次風和二次風進口處的流量和軸向速度,一次風進口處的旋流數和切向速度.燃燒室各進口處的氣體溫度均按實測的環境溫度取為282 K,燃燒室外壁面溫度按試驗測得的平均值取348 K,燃燒室底壁內側表面溫度也取348 K.從表1可以看到,工況1和工況2的差異僅在于一次風進口處的切向速度和旋流數不同.對比這2組工況的模擬結果可以得出浮力對具有不同旋流數的湍流火焰的影響.

表1 工況參數Tab.1 Parameters for calculation
由于所研究的對象具有軸對稱性,僅對燃燒室的半個縱剖面進行計算.計算時沿軸向和徑向共選取了62×44個網格節點,且均采用非均勻網格布置.在燃料與空氣進口、軸線、壁面及燃燒室出口附近網格布置相對較密.
圖2~圖6給出了三組工況下受浮力作用的湍流旋流火焰的模擬結果.為對比浮力的影響,圖中還同時給出了平均動量方程中不含重力項,且采用不含浮力修正的k-ε湍流模型的模擬結果.將工況1和工況2的氣體溫度及O2和CO2體積分數分布的計算結果與文獻[7]的試驗數據進行了比較,將工況3的氣體軸向速度和軸向脈動速度均方根值分布計算結果與文獻[8]的試驗數據進行了對比.
圖2(a)給出了工況1下氣體溫度分布的模擬結果,并與文獻[7]的試驗數據進行了對比.在燃燒室前部,浮力修正的k-ε模型與k-ε模型的計算結果差異較大.在靠近壁面的區域內,采用浮力修正的k-ε模型預測得到的溫度與試驗數據相符,而采用k-ε模型計算出的溫度則高于試驗值.在中心區域,采用浮力修正的k-ε模型比k-ε模型計算得出的溫度峰值更靠近軸線,高溫區的徑向范圍更窄,與試驗結果更接近.在燃燒室后部,兩種模型計算得到的溫度分布趨于均勻,且二者逐漸接近,并均與試驗結果相符.圖2(b)給出了工況1下O2體積分數分布的模擬結果與試驗結果的比較.在燃燒室前部,采用浮力修正的k-ε模型計算出的O2體積分數在近壁區域與試驗值相符,在軸線附近較試驗值偏低,而采用k-ε模型預測得出的O2體積分數比試驗值和浮力修正的k-ε模型的計算結果均偏低較多.在燃燒室后部,兩種模型的計算結果逐漸趨于一致,均與試驗結果相符,O2體積分數分布趨于均勻.

圖2 工況1下氣體溫度和O2體積分數模擬結果與試驗數據的比較Fig.2 Comparison between calculated and experim ental results for gas temperature and O2 volume fraction in case 1
圖3(a)給出了工況1下CO2體積分數分布的模擬結果,并與文獻[7]中的試驗數據進行了對比.在燃燒室前部除中心軸線附近以外的區域,采用浮力修正的k-ε模型計算出的CO2體積分數與試驗值相符,而采用k-ε模型得到的CO2體積分數則比試驗值偏高,只在近壁區與試驗值相符.在中心軸線附近的區域,兩種模型的計算結果較接近,在燃燒室進口處均比試驗值偏低,在遠離進口處則與試驗值相符.在燃燒室后部,浮力修正的 k-ε模型與k-ε模型的計算結果相差較小,均與試驗值相符,CO2體積分數分布趨于均勻.圖3(b)給出了工況1下CH4體積分數分布的模擬結果.CH4體積分數在燃燒室中心區域較高,近壁區域較低,峰值位于中心軸線處.采用浮力修正的k-ε模型與k-ε模型計算得到的結果有一定差異,前者預測得到的CH 4體積分數峰值比后者高一些,且體積分數較高區域的范圍窄一些.從圖3(b)中還可以看出,伴隨燃燒的進行,軸線處CH4體積分數峰值沿軸線逐漸下降,到燃燒室后部,燃料已消耗完畢.

圖3 工況1下CO2和CH4體積分數模擬結果及其與試驗數據的比較Fig.3 Comparison between calculated and experim ental results for CO2 and CH4 volume fraction in case 1
圖4(a)給出了工況2下氣體溫度分布的模擬結果及其與試驗數據的比較.在燃燒室前部,氣體溫度在中心區域較高,近壁區域較低.在近壁區域,采用浮力修正的k-ε模型和k-ε模型計算得到的氣體溫度分布均與試驗結果相符.在中心區域,兩種模型計算得到的氣體溫度均比試驗值偏高,且采用浮力修正的k-ε模型計算得到的溫度分布與試驗結果更接近,高溫區范圍更窄.在燃燒室后部,兩種模型得到的溫度分布仍比試驗值偏高,但趨勢上與試驗結果相符,溫度分布逐漸趨于均勻,且浮力修正的k-ε模型比k-ε模型的計算結果更接近試驗值.圖 4(b)給出了工況2下O2體積分數分布模擬結果與試驗結果的比較.在燃燒室前部的中心區域,浮力修正的k-ε模型和k-ε模型計算出的O2體積分數均低于試驗值,但前者的計算結果與試驗值更接近.在燃燒室前部的近壁區域,兩種模型的計算結果一致,均與試驗數據相符.在燃燒室后部,兩種模型預測得出的O2體積分數均比試驗值偏低,其中采用浮力修正的k-ε模型計算出的O2體積分數分布較平緩,趨勢上與試驗值更接近.燃燒室內O2體積分數分布的特點是中心區域較低,壁面附近較高,在燃燒室后部趨于均勻.

圖4 工況2下氣體溫度和O2體積分數模擬結果與試驗數據的比較Fig.4 Comparison between calculated and experimental results for gas temperature and O2 volume fraction in case 2
圖5(a)給出了工況2下CO2體積分數模擬結果及其與試驗數據的比較.在燃燒室前部靠近壁面的較大區域及中心軸線附近,采用浮力修正的k-ε模型和k-ε模型計算出的CO2體積分數分布很接近,在近壁區域與試驗數據相符,在中心軸線附近則比試驗值偏低.在偏離軸線的區域,兩種模型計算得到的CO2體積分數均比試驗值偏高,其中浮力修正的k-ε模型的計算結果與試驗值更接近,CO2體積分數較高的區域也略窄一些.在燃燒室后部,CO2體積分數分布趨于均勻,兩種模型預測得到的CO2體積分數均與試驗值相符.圖5(b)給出了工況2下CH4體積分數分布的模擬結果.采用浮力修正的k-ε模型與k-ε模型預測得出的CH4體積分數的分布趨勢一致.由于燃料的消耗,CH4體積分數沿軸向和徑向均較快地降低到接近零,體積分數峰值位于中心軸線上.采用浮力修正的 k-ε模型計算出的CH4體積分數比k-ε模型的結果低一些,燃料沿徑向的消耗要快一些.

圖5 工況2下CO2和CH4體積分數模擬結果及其與試驗數據的比較Fig.5 Comparison between calculated and experimental results for CO2 and CH4 volume fraction in case 2
圖6(a)給出了工況3下氣體軸向速度分布的模擬結果,并與文獻[8]的試驗數據進行了比較.在靠近壁面的區域,浮力修正的k-ε模型和k-ε模型的計算結果較接近,均與試驗值相符,其中前者與試驗值符合得更好.在燃燒室中心區域,浮力修正的k-ε模型比k-ε模型預測得出的軸向速度偏高,而后者與試驗值更接近.圖6(b)給出了工況3下氣體軸向脈動速度均方根值分布的模擬結果,并與試驗數據進行了對比.在靠近燃燒室進口的2個截面上,浮力修正的k-ε模型和k-ε模型的計算結果相差不大,均與試驗值相符.在x/R=1.81的截面上,浮力修正的k-ε模型預測得出的軸向脈動速度均方根值比k-ε模型預測得出的大,且與試驗值更接近.而在中心軸線處,k-ε模型的預測結果與試驗值更接近.
(1)采用浮力修正的k-ε模型計算得到的不同工況下的O2和CO2體積分數分布均與試驗結果相符,得到的氣體溫度、軸向速度和軸向脈動速度均方根值分布基本與試驗結果相符.
(2)采用浮力修正的k-ε模型計算得到的氣體組分體積分數場、溫度場和軸向脈動速度均方根值分布比k-ε模型的計算結果有較明顯的改進.在燃燒室前部,各物理量變化較劇烈,浮力修正的k-ε模型和k-ε模型的計算結果差別較大.在燃燒室后部,各物理量變化趨于平緩,浮力修正的k-ε模型與k-ε模型的模擬結果相差不大.
(3)采用浮力修正的k-ε模型計算得到的氣體溫度峰值更靠近軸線,高溫區的徑向范圍更窄,氣體組分體積分數變化較大的區域也更窄且更集中在軸線附近.

圖6 工況3下氣體軸向速度和軸向脈動速度均方根值模擬結果與試驗數據的比較Fig.6 Comparison between calculated and experimental results for gasaxial velocity and root mean square of axial fluctuating velocity in case 3
(4)在燃料和空氣的初始軸向動量相同的條件下,對于初始切向速度或旋流數相對較高的湍流火焰,采用浮力修正的k-ε模型與k-ε模型得到的模擬結果差異更明顯,說明浮力對切向速度或旋流數較高的湍流火焰影響更大.
[1] JENG SM,FAETH G M.Species concentrations and turbulence properties in buoyant methane diffusion flames[J].ASME JHeat Transfer,1984,106(4):721-727.
[2] LIU F,W EN JX.The effect of turbulence modeling on the CFD simulation of buoyant diffusion flames[J].Fire Safety Journal,2002,37(2):125-150.
[3] XIN Y,GORE JP,MCGRATTAN K B,et al.Fire dynamics simulation o f a turbulent buoyant flame using a mixture-fraction-based combustion model[J].Combustion and Flame,2005,141(4):329-335.
[4] 楊震,莊恩如,張建文,等.大型電站鍋爐采用切向燃燒方式燃用無煙煤的研究[J].動力工程,2006,26(6):766-772.YANG Zhen,ZHUANG Enru,ZHANG Jianwen,et a l.Research on tangential firing of anthracite in large capacity station boilers[J].Journal of Power Engineering,2006,26(6):766-772.
[5] 陶文銓.數值傳熱學[M].2版.西安:西安交通大學出版社,2001.
[6] RODIW.環境問題的紊流模型[M].賀益英,譯.北京:水利電力出版社,1987.
[7] 普勇,張健,周力行.一次風旋流數對燃燒室內湍流燃燒與NOx生成的影響[J].工程熱物理學報,2005,26(增刊):265-268.PU Yong,ZHANG Jian,ZHOU Lixing.Effects of swirl number for the primary air on the turbulent combustion and NOxformation in a swirl combustor[J].Journal of Engineering Thermophysics,2005,26(s):265-268.
[8] 普勇,張健,周力行.旋流燃燒室內湍流燃燒速度場的實驗研究[J].力學學報,2003,35(3):341-347.PU Yong,ZHANG Jian,ZHOU Lixing.Measurements of the velocity fields for turbulent combustion in a swirl combustor[J].Acta Mechanica Sinica,2003,35(3):341-347.