徐寶鵬，曾佑杰，馬宏宇，趙凱嵐，金 戈

（1.大連理工大學(xué)能源與動力學(xué)院，遼寧大連116024；2.中航工業(yè)沈陽發(fā)動機設(shè)計研究所，沈陽110015）

燃氣輪機模型燃燒室的大渦模擬

徐寶鵬1，曾佑杰1，馬宏宇2，趙凱嵐2，金 戈2

（1.大連理工大學(xué)能源與動力學(xué)院，遼寧大連116024；2.中航工業(yè)沈陽發(fā)動機設(shè)計研究所，沈陽110015）

燃燒室內(nèi)的燃油霧化、蒸發(fā)以及和空氣進行混合過程對燃燒過程有重要影響。提出1種基于大渦模擬的數(shù)學(xué)模型來模擬燃燒室內(nèi)燃料噴射、蒸發(fā)和混合過程。被空間濾波掉的亞網(wǎng)格尺度渦對大尺度渦的影響由求單方程SGS湍流模型進行模擬。采用拉格朗日法和蒙特卡洛技術(shù)對流場中的噴霧粒子進行采樣跟蹤，采樣噴霧粒子在流場中作為點源項與氣相進行質(zhì)量、動量和能量的雙向耦合。提出1個基于SGS湍流動能的雙向耦合模型來模擬SGS脈動速度對噴霧粒子運動的影響以及噴霧相對SGS湍流動能的影響。通過對1個同軸模型燃燒室中的噴霧蒸發(fā)及混合過程的大渦模擬，將預(yù)測結(jié)果和試驗值進行了比較，預(yù)測值和試驗值吻合良好，驗證了模型的可靠性。

燃燒室；燃氣輪機；大渦模擬；雙向耦合；燃油霧化

0 引言

試驗研究表明，燃油在燃氣輪機燃燒室內(nèi)的霧化、蒸發(fā)及與空氣進行混合的過程對直噴或預(yù)混燃燒室中的燃燒過程起至關(guān)重要的作用[1-3]。Fric[3]通過試驗研究發(fā)現(xiàn)，燃料混合在空間上的不均勻性和時間上的脈動性對氮氧化物的排放量有顯著影響；此外,燃燒室內(nèi)的流動狀態(tài)以及燃油噴霧粒子的大小和速度分布對燃燒過程同樣有重要影響。

光學(xué)測試和數(shù)值模擬是目前研究燃燒室內(nèi)燃油混合及燃燒過程的2種主要方法。光學(xué)測試已被成功應(yīng)用于燃油噴射及混合過程的研究中[4-6]，但使用數(shù)值模擬方法對燃燒前噴霧的分布情況進行理論研究的文獻較少，且現(xiàn)有的數(shù)值模擬工作大都采用基于雷諾時間平均的RANS方法。以往的數(shù)值研究表明，RANS方法無法準(zhǔn)確預(yù)測燃燒室內(nèi)的回流流動，且其穩(wěn)態(tài)特性也不能準(zhǔn)確預(yù)測燃料的混合及燃燒過程[7]。大渦模擬（LES）技術(shù)介于直接數(shù)值模擬（DNS）和RANS之間，尤其適合捕捉復(fù)雜的非穩(wěn)態(tài)流動。對簡單形狀的燃燒室中湍流混合及燃燒過程進行的數(shù)值模擬研究表明LES明顯優(yōu)于RANS方法[7]。

本文提出1種基于LES模擬燃燒室內(nèi)燃料噴射、蒸發(fā)和混合過程的數(shù)學(xué)模型。該模型采用單方程的亞網(wǎng)格尺度（SGS）湍流模型，考慮了氣液兩相的雙向耦合作用。將模型植入到開源CFD程序OpenFoam中，對1個同軸模型燃燒室中的噴霧蒸發(fā)及混合過程進行大渦數(shù)值模擬，并將預(yù)測結(jié)果與試驗值進行比較。

1 數(shù)學(xué)模型

與基于時間平均的RANS方法不同， LES采用空間濾波方法將小于網(wǎng)格尺度的湍流渦從控制方程中濾掉，而與問題相關(guān)的大尺度湍流渦通過經(jīng)空間濾波的控制方程進行直接求解。亞網(wǎng)格湍流渦對控制方程求解出大尺度湍流渦的影響由亞網(wǎng)格湍流模型進行模擬。本文所求解經(jīng)過空間濾波的可壓縮、多組分、具有噴霧源項的控制方程組如下。

1.1 控制方程組

式中：量符號上的“—”表示空間濾波變量，（fxi，t）=∫ΩG（xi-z）i（fzi，t）dzi；量符號上的“～”表示Favre空間濾波變量，~f=ρf ρ；ρ為密度；ρm為組分m的分密度；Ym為組分m的質(zhì)量分數(shù)；u為流體速度矢量；D為擴散系數(shù)；為噴霧質(zhì)量源項；p為流體壓力；σ為黏性應(yīng)力張量；Fs為噴霧動量源項；I為比內(nèi)能；J為熱傳導(dǎo)和焓擴散項；Q˙s為噴霧能量源項；?sgs、σsgs、Hsgs、ψsgs為由空間濾波引起的SGS源項，由SGS湍流模型確定。

1.2 SGS湍流模型

本文采用單方程SGS湍流模型[8]。此模型放棄了湍流平衡假設(shè)，因此允許采用相對稀疏的網(wǎng)格來加快計算速度；根據(jù)此模型求解出的SGS湍流動可以計算出SGS脈動速度，為氣液兩相雙向耦合模擬提供1個重要參數(shù)。

根據(jù)求解SGS湍流動能，控制方程中的SGS源項為

式中：Schmidt數(shù) Sct=1；Prandtl數(shù) Prt=0.9。

1.3 噴霧兩相流模型

本文采用歐拉-拉格朗日法和蒙特卡洛（Monte Carlo）法對噴霧粒子在流場中的采樣進行跟蹤[9]。采樣噴霧團在流場中作為點源項與氣相進行質(zhì)量、動量和能量的雙向耦合。本文提出1個基于SGS湍流動能的雙向耦合模型來模擬SGS脈動速度對噴霧粒子運動的影響以及噴霧相對SGS湍流動能的影響。

噴霧粒子加速度F為

式中：CD為阻力系數(shù)；g為重力加速度；u'為SGS湍流脈動速度，由下述方法確定。

考慮到SGS速度的各向同性，可以假定u'服從Gaussian分布

u'在SGS相關(guān)時間tsub內(nèi)僅采樣1次，SGS相關(guān)時間取為SGS湍流渦的破碎時間和噴霧粒子穿越SGS尺度渦所需時間的最小值

式中：cps為經(jīng)驗系數(shù)，取值為0.16。

u'確定后，控制方程中的噴霧源項為

式中：F'=F-g；f為噴霧粒子概率密度分布函數(shù)。

式（13）～（16）由噴霧方程進行求解[10]。粒子直徑變化率R由Frossling相關(guān)公式確定[11]；粒子變形y及其變化率y˙由粒子破碎ETAB模型[12]確定。

2 數(shù)值方法

為保證離散控制方程的守恒特性，采用有限容積法離散所求解的控制方程組。離散的控制方程滿足2階時間和空間精度：對時間項的離散采用2階隱式格式；在空間上對擴散項、壓力相關(guān)項以及動量方程中的對流項采用2階中心差分格式；為了保證數(shù)值穩(wěn)定性，對能量方程和組分方程中的對流項采用2階的TVD格式；為提高在低馬赫數(shù)下的求解效率，采用PISO算法對壓力進行迭代求解。以上數(shù)值方法被植入到開源CFD計算程序OpenFoam中。

3 問題描述

Sommerfeld等[13]對1個同軸燃氣輪機模型燃燒室（如圖1所示）中噴霧的蒸發(fā)及混合進行了詳細試驗研究，試驗結(jié)果被廣泛用于燃氣輪機燃燒室大渦模擬的計算程序考核中[14-15]。試驗采用多普勒測速儀測量了氣液兩相的速度及噴霧粒子的直徑。加熱空氣由筒形燃燒室頭部環(huán)形入口進入試驗段，具有高揮發(fā)性異丙醇通過位于燃燒室頭部中心位置的空錐形壓力霧化噴嘴進入燃燒室。具體試驗條件請參見文獻[13]。

圖1 同軸模型燃燒室

采用結(jié)構(gòu)化六面體網(wǎng)格進行計算，將計算區(qū)域劃分為2×107個網(wǎng)格單元。在環(huán)形入口平面上設(shè)定非穩(wěn)態(tài)入口邊界條件，其入口速度值取為平均速度型線上疊加5%的湍流隨機脈動值；在出口平面設(shè)定連續(xù)出口邊界條件；所有固體壁面設(shè)定為無滑移絕熱壁。噴霧相的入口條件取自試驗數(shù)據(jù)，噴霧粒子直徑通過試驗離散的概率密度函數(shù)進行采樣；噴霧的空錐角分別為 60°（外角）和 15°（內(nèi)角）；為了保證計算網(wǎng)格內(nèi)有足夠數(shù)目的噴霧團用于統(tǒng)計平均數(shù)，噴霧團的入口噴射率取為1.5×108個/s。

4 計算結(jié)果及分析

大渦模擬不同于RANS方法，是1種非穩(wěn)態(tài)數(shù)值模擬方法，可以直接求解出大于網(wǎng)格尺度的渦運動。圖2（a）中顯示壓力脈動預(yù)測值的等值面，數(shù)值模擬可以成功捕捉到這種模型燃燒室固有的渦結(jié)構(gòu)：由Kelvin-helmholtz剪切層引起的環(huán)形渦及其下游流線方向的蟲形渦。圖2（b）中顯示噴霧粒子在中截面上的瞬時分布，噴霧粒子進入燃燒室后，受環(huán)形入口氣流的限制，集中在軸線附近；當(dāng)噴霧粒子隨氣流移到下游時，在湍流彌散的作用下，其分布更加均勻，同時受蒸發(fā)影響，噴霧粒子逐漸變小。在圖2（c）中，截面中平均軸向速度云圖上有2個回流區(qū)（如藍色區(qū)域所示）：1個為圓環(huán)射流內(nèi)的小回流區(qū)，另1個為壁面附近相對較大的回流區(qū)。圓環(huán)入口射流在z=10 cm處匯合，并開始向壁面擴展。

圖2 總體預(yù)測結(jié)果

噴霧粒子軸向速度平均值及其脈動值與試驗數(shù)據(jù)的對比如圖3所示。預(yù)測的軸向速度脈動值和試驗值吻合良好；預(yù)測的平均軸向速度與試驗值基本吻合。在z=10、20、30 cm 3個截面上，平均軸向速度的預(yù)測值略大于試驗值；但在z=40 cm截面上，受蒸發(fā)影響，粒子直徑變小，更容易跟隨氣相流動，所以預(yù)測值和試驗值吻合得更好。

圖3 噴霧粒子軸向速度平均值及其脈動值與試驗數(shù)據(jù)的對比

噴霧粒子索特平均直徑與試驗值的對比如圖4所示。除z=30 cm截面外，其他5個截面上預(yù)測值和試驗值都吻合良好。預(yù)測的型線反映出典型的空錐噴嘴的特性，即噴霧進入到燃燒室后，小霧滴集中在中心區(qū)，大霧滴分布在噴霧邊緣。隨著霧滴遠離入口，霧滴的平均直徑受蒸發(fā)和湍流的混合共同影響，在橫截面上的分布趨于更加均勻。

圖4 噴霧粒子索特平均直徑與試驗值的對比

噴霧軸向粒子質(zhì)量通量與試驗值的對比如圖5所示。除了z=5 cm截面外，其他5個截面上的預(yù)測值和試驗值吻合良好。預(yù)測的型線在z=2.5和z=5 cm截面上表現(xiàn)出2個和空錐噴嘴特性相關(guān)的峰值。由于存在中心小回流區(qū)，在z=2.5 cm截面上出現(xiàn)負的質(zhì)量通量。在z=10 cm及其下游截面上質(zhì)量通量開始沿軸線分布，且受蒸發(fā)和湍流的彌散共同作用，其峰值迅速變小且分布趨于更加均勻，在z=40 cm截面上幾乎所有的霧滴蒸發(fā)成氣相。

圖5 噴霧軸向粒子質(zhì)量通量與試驗值的對比

5 結(jié)論

本文提出了1種基于LES模擬燃燒室內(nèi)燃油噴射蒸發(fā)和混合過程的數(shù)學(xué)模型，模型中首次提出基于SGS湍流動能的氣液兩相亞網(wǎng)格尺度雙向耦合作用。通過對Sommerfeld同軸模型燃燒室中的噴霧蒸發(fā)及混合過程的大渦數(shù)值模擬，對提出的數(shù)學(xué)模型進行了驗證。結(jié)果表明：本文提出的大渦模擬兩相流模型成功得捕捉了混合過程的非穩(wěn)態(tài)特性以及湍流渦結(jié)構(gòu)；預(yù)測的噴霧粒子平均速度及脈動速度、噴霧粒子大小以及質(zhì)量通量與試驗值吻合較好；驗證了模型的正確性和可靠性；噴霧兩相流的大渦模擬數(shù)學(xué)模型具有一定的普遍性，可以應(yīng)用到真實燃氣輪機和內(nèi)燃機燃燒室內(nèi)的燃油噴射、蒸發(fā)和混合過程。

[1]Greenhalgh D A.Laser imaging of fuel injection systems and combustors[J].Proceedings of the Institution of Mechanical Engineers，1991，214（4）：367-376.

[2]Tait N P,Greenhalgh D A.2D Laser induced fluorescence imaging o parent fuel fraction in nonpremixed combustion[J].Proceedings of Combustion Institute,1992，24：1621-1628.

[3]Fric T F.Effects of fuel-air unmixedness on NOx emissions[J].Journal of Propulsion Power,1993，9（5）：708-713.

[4]Harding SC,Greenhalgh D A.Fuel-airmixing and combustion in an optical,lean,premixed,prevaporised gas turbine combustor[R].ASME 98-GT-553.

[5]Mongia R,Tomita E,Hsu F K，et al.Use of an optical probe for time-resolved in-situmeasurementof local air-to-fuel ratio and extent of fuelmixing with applications to low NOx emissions in premixed gas turbines[J].Proceedings of Combustion Institute,1996，26：2749-2755.

[6]Farrugia N,Lockett R D,Harding SC,etal.LIF imaging of fuel distribution in gas turbine combustors[C]//Proceedings of ILASS-Europe,Florence,Italy,1997：299-305.

[7]Mahesh K,Constantinescu G,Moin P.A numericalmethod for large eddy simulations in complex geometries[J].Journal of Computational Physics，2004，197：215-240.

[8]Menson S，Yeung P K，Kim W W.Effect of sub-grid models on the computed interscale energy transfer in isotropic turbulence[J].Computers and Fluids，1996，25（2）：165-180.

[9]Dukowicz J K.A particle fluid numerical method for liquid sprays [J].Journal of Computational Physics，1980，35：229-253.

[10]Willians F A.Spray combustion and atomization[J].Physics of Fluids，1958，1（6）：541-546.

[11]Faeth GM.Current status of dropletand liquid combustion[J].Progress in Energy and Combustion Science，1977，3（4）：191-224.

[12]Tanner F X.Liquid jet atomization and droplet break-up modelling of non-evaporating diesel fuel sprays[J].SAE transactions：Journal of Engines，1998，106：127-140.

[13]Sommerfeld M,Qiu H H.Experimental studies of spray evaporation in turbulent flow[J].International Journal of Heat and Fluids Flows，1998，19（1）：10-22.

[14]Moin P，Apte SV.Large eddy simulation of realistic gas turbine combustors[J].AIAA Journal，2006，44（4）：698-708.

[15]James S，Zhu J，Anand M S.Large eddy simulations as a design tool for gas turbine combustion systems[J].AIAA Journal，2006，44（4）：674-686.

俄羅斯PD-14發(fā)動機研制概況

PD-14發(fā)動機于2006年由彼爾姆航空發(fā)動機公司牽頭，與俄羅斯聯(lián)合發(fā)動機制造集團（ODK）的多家單位以及烏克蘭伊夫琴科-進步設(shè)計局共同研制，其研制基礎(chǔ)為彼爾姆航空發(fā)動機公司于1999年研制的推力范圍為68.6～196 kN的PD-12發(fā)動機。

PD-14發(fā)動機基準(zhǔn)型推力為137.2 kN，風(fēng)扇直徑為1900mm，發(fā)動機干質(zhì)量為2870 kg，涵道比為8.5，總增壓比為37.2，爬升狀態(tài)總壓比為41.0，渦輪前燃氣溫度起飛狀態(tài)為1730K，巡航狀態(tài)（H=11 km，Ma=0.8）為1450 K，配裝MS-21-300飛機。改型包括：PD-14A，推力為122.5 kN，采用節(jié)流方案，用于配裝MS-21-200飛機；PD-14M，推力為152.9 kN，采用增推方案，配裝MS-21-400飛機；PD-18，推力為176.4～196 kN。

PD-14發(fā)動機的研制采用了大量新技術(shù)，廣泛使用復(fù)合材料，與目前世界上正在使用的民用發(fā)動機 V2527-A5、CFM56-5B4、CFM56-7B27 相 比 ，PD-14發(fā)動機耗油率降低約10%～17%，全壽命周期費用降低約11%～24%。

PD-14發(fā)動機及其改型機可用于配裝SSJ-100、安 -148、圖 -204、圖 -214、圖 -334、別-200、伊爾-76、伊爾-96等飛機及其改型機。（中航工業(yè)動力所 劉 靜）

Large Eddy Simulation of a Gas Turbine ModelCombustor

XU Bao-peng1,ZENG You-jie1,MA Hong-yu2,ZHAO Kai-lan2,JIN Ge2
（1.Schoolof Energy and Power Engineering,Dalian University of Technology,Liaoning Dalian 116024,China;2.AVIC Shenyang Engine Design and Research Institute,Shenyang 110015,China）

Fuelatomization,evaporation andmixingwith air in gas turbine combustors are vital to the subsequent combustion process.Numerical formulation based on large eddy simulation is proposed to model fuel injection,evaporation and mixing in a gas turbine combustor.The proposedmodeladopts a one-equation subgrid scale turbulentmodel to handle the effectof the filtered subgrid scale eddies on the solved large scale eddies.Spray droplets are tracked using both Lagrangian method and Mento Carlo technique,and the sampled spray particles are regarded as point sources to conduct two-way couplings ofmass,momentum and energy.The two-way couplingmodel based on SGS turbulent kinetic energy is used tomodel themutual influences between SGS fluctuating velocity and themovement of spray droplets.The proposed models are validated againsta large simulation of a co-axialmodel combustor and the predictions are compared to the experimental data.Good agreementsare obtained,which demonstrate the reliability of the proposedmodels.

combustor；gas turbine;large eddy simulation;two-way coupling;atomization

V211.3

A

10.13477/j.cnki.aeroengine.2014.03.003

2013-08-14 基金項目：燃氣輪機重大項目聯(lián)合培育基金（2011LH006）資助

徐寶鵬（1969）,男，博士，教授，主要研究方向為計算流體力學(xué)、兩相流和燃燒學(xué)；E-mail:xbp624@gmail.com。

徐寶鵬，曾佑杰，馬宏宇，等.燃氣輪機模型燃燒室的大渦模擬[J].航空發(fā)動機，2014，40（3）：14-18.XUBaopeng，ZENGYoujie，MAHongyu，etal.Large eddy simulation ofagas turbinemodelcombustor[J].Aeroengine，2014，40（3）：14-18.