2. 浙江大学 航空航天学院, 杭州 310027
2. School of Aeronautics and Astronautics, Zhejiang University, Hangzhou 310027, China
航空发动机与燃气轮机被誉为装备制造业“皇冠上的明珠”,是保证国防与能源安全,保持工业竞争力的核心战略产业。我国已启动的“两机专项”被列为国家“十三五”期间百个大项目之一,是工程研究的重中之重。目前我国航空发动机和燃气轮机的发展仍有较大不足,而燃烧室作为航空发动机和燃气轮机实现能量转化的核心部件之一,通常工作在高温、高压等极端环境之下,其设计水平直接影响到整机的性能[1], 对其开展深入研究具有极大的现实意义。
近年来,燃烧室内旋流火焰结构和机理的研究越来越引起相关科研工作者的重视。在实验分析方面,Watanabe[2]和Zhou[3]等在常温常压下对旋流预混火焰进行了实验研究,分析了单旋流燃烧器中的火焰结构,阐释了组分扩散与热扩散间的作用关系,并对湍流流场进行了测量和分析;Li[4]和Wang[5]等对航空发动机贫燃预混预蒸发模型燃烧室进行了实验,研究了高温高压下航空发动机内的污染物生成机制;Worth[6]和Singh[7]等则着重分析了环形燃烧室内火焰和流场的相互作用关系。在数值模拟方面,苗文博[8]等采用二维化学非平衡Navier-Stokes方程分析了空间发展的超声速燃烧反应剪切层的放热效应;Butz[9]采用大涡模拟(Large Eddy Simulation, LES)方法对湍流预混旋流火焰进行了模拟,发现了内部回流区对燃尽产物的阻滞作用;Tanaka[10]和Aoki[11]等采用直接数值模拟(Direct Numerical Simulation, DNS)方法对常温常压下预混甲烷/空气旋流火焰进行了计算,从机理层面分析了旋流火焰和流场的结构以及湍动能的来源。
然而由于受高压强湍流条件下的实验诊断技术和DNS计算量的限制,大部分研究均未考虑高压强条件下湍流特征量与实际航空发动机燃烧室特征量的内在关系。与常压燃烧最大的区别在于,高压预混燃烧的火焰面极薄。在湍流预混燃烧理论中,常使用湍流均方根速度u′、湍流积分长度尺度lt、湍流Kolmogorov尺度η、层流火焰速度SL和层流火焰厚度δL等特征量来定义无量纲特征量如u′/SL、lt/δL、Karlovitz数(Ka)等。对于高压火焰,极薄的火焰厚度直接影响着上述参数的量级。对于航空发动机燃烧室内的湍流预混燃烧,Ka大致为100~102,在这一参数范围内,预混燃烧模式位于薄反应区和破碎火焰区中[12]。
近期也有学者针对不同Ka下火焰开展了研究工作。Carlsson[13]通过对各向同性湍流燃烧进行直接数值模拟,研究了高Ka下的破碎火焰结构特征;MacArt[14]同样采用DNS方法模拟分析了不同Ka下平板预混燃烧射流中标量通量和梯度分布;Wang[15]则对圆柱预混射流进行了直接数值模拟,分析了较高Ka下的火焰拉伸和火焰增厚现象。然而基于不同Ka的旋流火焰的研究还有所不足,对燃烧室内相应旋流火焰特征的理解仍需加强。
基于此,为了进一步理解航空发动机燃烧室内旋流预混燃烧的火焰结构特征和机理,本文设计构造了高温高压环形预混旋流模型燃烧室,并采用直接数值模拟方法研究了不同Ka的火焰和湍流结构特征,分析了旋流火焰中标量通量和梯度的对齐关系,讨论了对相关模型的影响,为发展适用于航空发动机燃烧的高精度模型和方法提供基础。
1 直接数值模拟 1.1 计算模型实际航空发动机燃烧室多为环形结构,且入口燃料温度高,室内压力和湍流强度都极高。考虑到直接数值模拟计算量的限制,对实际结构等比例建模较为困难,所以基于典型工况(p、T、U0)及其对应的典型湍流预混火焰参数(u′/SL、lt/δL、Ka)设计了简化的计算模型。本研究中将航空发动机内的环形燃烧室简化为两个相邻的预混旋流。图 1所示为计算域示意图,在z方向上采用周期性边界条件,而在y方向上采用壁面边界条件以模拟实际环形燃烧器。单个环状旋流入口的内径Din设置为0.1 mm,外径Dout设置为0.3 mm,两相邻旋流中心间距为2Dout,计算域整体尺寸为5.5Dout×2Dout×4Dout(x×y×z)。
![]() |
图 1 计算域示意图(涡量等值面采用温度着色,左侧切片图为热释率) Fig.1 Schematic illustration of the computational domain (the isosurface colored by temperature presents vorticity and the slice on the left presents heat release rate) |
为了研究不同Ka下的湍流火焰特征,设置了两个工况,其入口速度U0分别为40 m/s和15 m/s,分别对应于航空发动机燃烧室内的高/低负荷工况,其入口雷诺数Re0分别约为4000和1500,如表 1所示。其余参数将在后文详细介绍。
表 1 计算工况 Table 1 Calculation conditions |
![]() |
为了获得具有实际意义的边界条件,首先单独计算了环形管道湍流,待其充分发展稳定后切片保存生成了入口文件,然后在正式计算时对入口文件进行读取从而得到流场入口信息。计算中,旋流数为1.0,燃料为航空煤油,当量比为0.7,入口预混气体温度为860 K,计算域压力为3 MPa。基于以上条件,通过计算一维层流火焰可以得到其层流火焰速度SL=0.738 m/s, 层流火焰厚度δF=18.2 μm。
1.2 数值方法本研究采用了变密度、低马赫数的燃烧求解器,求解了连续性方程、动量守恒方程、能量守恒方程和组分守恒方程[16]。在气体预混燃烧直接数值模拟中,影响计算网格的因素主要有两个方面:从火焰角度,为了准确解析火焰面处组分扩散和燃烧放热效应,在火焰面处需要保持至少8~10个网格[17];从湍流角度,网格尺度需要能解析湍流最小特征尺度(η)的涡结构,通常而言,需要满足η>Δxmin/2[18]。本文模拟中最小涡尺度在入口附近产生,约为1.5 μm。综合考虑火焰和湍流的精度要求,本研究中网格数量(Nx×Ny×Nz)设置为768×256×512≈9.3×107。其中在x方向的前半段采用2 μm的均匀网格,在后半段采用拉伸率为1.05的拉伸网格,在y、z方向均采用大小为2.3 μm的均匀网格。为了降低计算量,采用高温高压条件下两步航空煤油/空气化学反应机理[19]。
2 旋流预混火焰分析 2.1 旋流预混燃烧状态为了确定本文湍流预混燃烧火焰的参数特征,需在火焰面位置统计各湍流燃烧参数。本文中,火焰面根据过程变量C来定义:
$ C=\left(Y_{\mathrm{O}_{2}, u}-Y_{\mathrm{O}_{2}}\right) /\left(Y_{\mathrm{O}_{2}, u}-Y_{\mathrm{O}_{2}, b}\right) $ | (1) |
式中YO2, u为未燃气体中氧气初始质量分数,YO2, b为燃尽物中氧气质量分数,YO2为当地氧气质量分数,C=0表示未燃物,C=1表示燃尽物。
图 2为工况1的y=0 mm截面上时间平均后的过程变量
![]() |
图 2 工况1中y=0 mm截面上的过程变量分布 Fig.2 Distribution of progress variable at y=0 mm of case 1 |
如前文所述,在实际航空发动机燃烧器中,湍流预混燃烧特征参数处于一定的分布区间,其中具有代表性意义的无量纲特征量有u′/SL、lt/δL和Ka等。其中Ka从物理意义上表示层流火焰特征时间尺度与最小涡特征时间尺度的比值,定义为:
$ K a=\tau_{L} / \tau_{\eta}=\left(\delta_{F} / S_{L}\right)(\varepsilon / \nu)^{1 / 2} $ | (2) |
式中τL为层流火焰特征时间尺度,τη为最小涡特征时间尺度,ε为湍流耗散率,ν为运动黏度。
图 3为高/低负荷工况中不同流向位置上基于过程变量
![]() |
图 3 不同x位置截面Ka基于过程变量的分布 Fig.3 Variation of Ka based on progress value at different x location |
在这两个工况中,Ka的值均在100~102范围内,与航空发动机燃烧室工况的参数一致。对于其他参数,图 4显示了湍流预混燃烧模式分区图[20],其中横纵坐标分别为湍流预混燃烧无量纲特征值lt/δL和u′/SL。可以看到,本文设计的两个工况参数均分布于“Aircraft engines”工况范围内。
![]() |
图 4 湍流预混燃烧模式分区图(三角形:工况1;矩形:工况2) Fig.4 Regime diagram for premixed turbulent combustion(triangle: case 1; rectangle: case 2) |
图 5、图 6分别显示了某个时刻两个工况在y=0 mm处的瞬时速度和瞬时涡量云图。旋流预混气体从喷口进入燃烧室后,由于旋转离心导致产生张角,流场呈双“V”字形。在燃烧室内,由于旋转射流内外剪切层处存在较为强烈的剪切作用,从而产生速度小于0的回流区。回流区主要分布在旋流中心区域和两旋流间的角落处,如图 5所示。相比于单旋流火焰,双旋流火焰在火焰下游处由于存在相邻旋流的相互作用,导致下游湍流强度显著增强,另外由于环形燃烧室内相邻旋流之间有限的空间,导致双旋流火焰下游的扩张效应明显弱于单旋流火焰[10],即其内部回流区小于单旋流火焰,且外部回流区相比于单旋流而言也更小。相比于低Ka工况,高Ka工况下速度场褶皱更强烈,这与图 6中的涡量场相对应。在工况1中,湍流涡的分布更加细小,对火焰的影响更大。
![]() |
图 5 y=0 mm截面处的瞬时速度云图(单位: m/s) Fig.5 Instantaneous velocity contours at y=0 mm (unit: m/s) |
![]() |
图 6 y=0 mm截面处的瞬时涡量云图 Fig.6 Instantaneous vorticity contours at y=0 mm |
图 7所示为瞬时热释率云图。火焰面在湍流作用下弯曲褶皱,表明在当前工况下湍流预混燃烧模式处于薄反应区中,如图 4的燃烧模式分区图所示。与常压旋流火焰相比,高压下的火焰面厚度极薄,因此存在于更小的湍流尺度范围内,火焰面受小尺度湍流的影响更显著,对应地,在预混湍流火焰燃烧模式分区图中,高压火焰的横坐标lt/δL的值较常压火焰更大,火焰褶皱也就更剧烈。对比图 5可以发现,在靠近入口的上游燃烧区,火焰面沿着旋流内外剪切层双侧分布,而在入口的下游燃烧区,火焰面相互作用,变得更加褶皱甚至破碎,这在x=0.3 mm截面上更加突出。特别是工况1中,在入口速度更高的高负荷情形下,由于其较高的Ka,下游火焰的破碎更加明显。
![]() |
图 7 瞬时热释率云图 Fig.7 Instantaneous contours of heat release rate |
如图 5、图 7所示,预混旋流喷射入燃烧室后,在上游稳定燃烧,而在下游发生两旋流间的相互作用。为了理解旋流湍流预混燃烧中Ka的影响,分别对x=0.15 mm(上游)和x=0.30 mm(下游)两处的截面进行了速度分布的统计比较。
图 8显示了轴向、径向和周向的平均速度在过程变量空间上的分布。在上游x=0.15 mm且
![]() |
图 8 平均速度统计 Fig.8 Statistics of mean velocity |
从平均轴向速度来看,在x=0.15 mm、
图 9所示为轴向、径向和周向的脉动速度统计。在上游处(x=0.15 mm),由于工况1的雷诺数更大,湍流扰动效应也更强,所以脉动速度均大于工况2。而在下游情况则有所不同,对轴向和径向脉动而言,虽然在
![]() |
图 9 脉动速度统计 Fig.9 Statistics of fluctuating velocity |
在雷诺平均(Reynolds Averaged Navier-Stokes, RANS)方法和大涡模拟(Large Eddy Simulation, LES)方法中,质量加权的标量通量Fj, k是守恒方程的重要项,其定义为:
$ F_{j, k}=\widetilde{u_{j} Y_{k}}-\widetilde{u_{j}} \widetilde{Y_{k}}=\widetilde{u_{j}^{\prime \prime} Y_{k}^{\prime \prime}} $ | (3) |
式中,Yk为组分k的质量分数。
然而,在RANS或LES中,根据式(3)所得到的标量通量项不能封闭,通常采用梯度扩散假设来建模,如:
$ F_{j, k} \approx-\frac{\nu_{t}}{S c_{t}} \frac{\partial \widetilde{Y_{k}}}{\partial x_{j}} $ | (4) |
式中:νt为涡黏度,Sct为湍流施密特数。
从式(4)中可以看出,梯度扩散假设是假设组分通量
为研究航空发动机燃烧室内典型梯度输运假设的有效性,本节对组分扩散和组分梯度的对齐性进行了统计。
$ \cos \theta_{k}=(\widetilde{u^{\prime \prime}_{j} Y_{k}^{\prime \prime}} \frac{\partial \widetilde{Y_{k}}}{\partial x_{j}}) /\left(|\widetilde{u_{j}^{\prime \prime} Y_{k}^{\prime \prime}}|\left|\frac{\partial \widetilde{Y_{k}}}{\partial x_{j}}\right|\right) $ | (5) |
cosθk为组分k通量与组分梯度间夹角的余弦值。当cosθk→-1时,通量与梯度为反对齐关系,即梯度扩散占主导地位;当cosθk→1时,通量与梯度为对齐关系,即逆梯度扩散占主导地位。
图 10为基于主要反应物、生成物和中间组分的标量通量和梯度的对齐性统计。可以发现,反应物(O2)和生成物(H2O)的对齐性在各位置的过程变量空间中基本保持一致,夹角的余弦值均小于-0.4,且均在火焰面附近取得最大值,表明采用梯度扩散假设来对主要反应物和生成物通量进行建模是基本合理的。且由于工况1的Ka更大,其梯度扩散效应更明显。
![]() |
图 10 基于反应物、生成物和中间组分的标量通量和梯度的对齐性统计 Fig.10 Alignment statistics of scalar flux and gradient for fuel, product and radial mass fraction |
然而对于中间产物,虽然同样当Ka较大时,夹角的余弦值更小,但是在火焰面附近,cosθCO的最大值约等于0.2,表明组分通量更偏向于与组分梯度呈对齐关系,即逆梯度扩散效应更为显著,传统梯度扩散模型失效,如图 10(c)所示。这表明传统上采用的基于梯度扩散的模型对于发动机燃烧室中间产物的模拟存在问题,在发展适用于航空发动机燃烧模拟的高精度模型和方法时,需要对相关通量项的封闭进行进一步的改进和优化。
4 结论基于航空发动机燃烧室湍流预混燃烧的特征参数,构造了环形旋流预混模型燃烧室,并研究了不同Ka条件下湍流火焰的区别及其内在机制,加强了对航空发动机燃烧室内预混湍流火焰和流场结构的理解。结果表明,发动机燃烧室内湍流预混燃烧模式位于薄反应区,且Ka越大时火焰的弯曲褶皱越强烈。火焰面沿着剪切层内外呈双侧分布,与单旋流燃烧不同的是,双旋流火焰在下游由于相邻旋流的碰撞而产生较为破碎的火焰结构,这在高Ka的工况中尤为明显。由于旋流剪切层的作用,导致在下游产生内部回流区,且Re和Ka越大,回流效应越显著。在剪切层附近,脉动速度最大,同样Re和Ka越大时,湍流扰动效应也越强。
在典型的RANS或LES方法中,常采用梯度扩散假设对组分通量进行建模,为了验证这一方法在发动机燃烧室内的适用性,对梯度输运效应进行了统计分析。分析表明,在航空发动机燃烧室工况下,反应物和生成物基本满足梯度输运效应,且Ka越大,梯度输运效应越显著,这与此前基于各向同性湍流和平板射流的研究结果相符。然而对于中间产物,如CO,梯度扩散假设在火焰面附近严重失效。因此在发展适用于航空发动机燃烧室的高精度数值模拟方法时,需要对标量通量项的封闭进行进一步的改进。
[1] |
GICQUEL L Y M, STAFFELBACH G, POINSOT T. Large eddy simulations of gaseous flames in gas turbine combustion chambers[J]. Progress in Energy and Combustion Science, 2012, 38(6): 782-817. DOI:10.1016/j.pecs.2012.04.004 |
[2] |
WATANABE H, SHANBHOGUE S J, TAAMALLAH S, et al. The structure of swirl-stabilized turbulent premixed CH4/air and CH4/O2/CO2 flames and mechanisms of intense burning of oxy-flames[J]. Combustion and Flame, 2016, 174: 111-119. DOI:10.1016/j.combustflame.2016.09.015 |
[3] |
ZHOU R G, BALUSAMY S, SWEENEY M S, et al. Flow field measurements of a series of turbulent premixed and stratified methane/air flames[J]. Combustion and Flame, 2013, 160(10): 2017-2028. DOI:10.1016/j.combustflame.2013.04.007 |
[4] |
LI L, LIN Y Z, FU Z B, et al. Emission characteristics of a model combustor for aero gas turbine application[J]. Experimental Thermal and Fluid Science, 2016, 72: 235-248. DOI:10.1016/j.expthermflusci.2015.11.012 |
[5] |
WANG Z C, LIN Y Z, WANG J C, et al. Experimental study on NO emission correlation of fuel staged combustion in a LPP combustor at high pressure based on NO-chemiluminescence[J]. Chinese Journal of Aeronautics, 2020, 33(2): 550-560. DOI:10.1016/j.cja.2019.09.004 |
[6] |
WORTH N A, DAWSON J R. Cinematographic OH-PLIF measurements of two interacting turbulent premixed flames with and without acousticforcing[J]. Combustion and Flame, 2012, 159(3): 1109-1126. DOI:10.1016/j.combustflame.2011.09.006 |
[7] |
SINGH P, VELAMATI R K, PRATHAP C, et al. Study of flow patterns and impingement heat transfer for an annular array of eight co-rotating dual-swirlingflames[J]. International Journal of Heat and Mass Transfer, 2019, 144: 118657. DOI:10.1016/j.ijheatmasstransfer.2019.118657 |
[8] |
苗文博, 程晓丽, 王强. 超声速空间发展燃烧反应剪切层放热效应分析[J]. 空气动力学学报, 2008, 26(3): 339-343. MIAO W B, CHENG X L, WANG Q. Analysis of the energy release effects in supersonic reactingflows[J]. Acta Aerodynamica Sinica, 2008, 26(3): 339-343. DOI:10.3969/j.issn.0258-1825.2008.03.012 (in Chinese) |
[9] |
BUTZ D, GAO Y, KEMPF A M, et al. Large Eddy Simulations of a turbulent premixed swirl flame usingan algebraic scalar dissipation rate closure[J]. Combustion and Flame, 2015, 162(9): 3180-3196. DOI:10.1016/j.combustflame.2015.05.003 |
[10] |
TANAKA S, SHIMURA M, FUKUSHIMA N, et al. DNS of turbulent swirling premixed flame in a micro gas turbine combustor[J]. Proceedings of the Combustion Institute, 2011, 33(2): 3293-3300. DOI:10.1016/j.proci.2010.07.034 |
[11] |
AOKI K, SHIMURA M, NAKA Y, et al. Disturbance energy budget of turbulent swirling premixed flame in a cuboid combustor[J]. Proceedings of the Combustion Institute, 2017, 36(3): 3809-3816. DOI:10.1016/j.proci.2016.08.033 |
[12] |
杨越, 陈正, 周豪, 等. 高雷诺数预混湍流火焰的数值模拟与结构表征[J]. 中国科学:物理学力学天文学, 2017, 47(7): 070005. YANG Y, CHEN Z, ZHOU H, et al. Numerical simulations and structure characterizations in high-Reynolds-number premixed turbulent flames[J]. SCIENTIA SINICA Physica, Mechanica & Astronomica, 2017, 47(7): 070005. (in Chinese) |
[13] |
CARLSSONH, YU R X, BAI X S. Flame structure analysis for categorization of lean premixed CH4/air and H2/air flames at high Karlovitz numbers:Direct numerical simulation studies[J]. Proceedings of the Combustion Institute, 2015, 35(2): 1425-1432. DOI:10.1016/j.proci.2014.09.002 |
[14] |
MACART J F, GRENGA T, MUELLER M E. Effects of combustion heat release on velocity and scalar statistics in turbulent premixed jet flames at low and highKarlovitz numbers[J]. Combustion and Flame, 2018, 191: 468-485. DOI:10.1016/j.combustflame.2018.01.022 |
[15] |
WANG H O, HAWKES E R, CHEN J H, et al. Direct numerical simulations of a high Karlovitz number laboratory premixed jet flame-an analysis of flame stretch and flame thickening[J]. Journal of Fluid Mechanics, 2017, 815: 511-536. DOI:10.1017/jfm.2017.53 |
[16] |
DESJARDINS O, BLANQUART G, BALARAC G, et al. High order conservative finite difference scheme for variable density low Mach number turbulentflows[J]. Journal of Computational Physics, 2008, 227(15): 7125-7159. DOI:10.1016/j.jcp.2008.03.027 |
[17] |
HAWKES E R, CHATAKONDA O, KOLLA H, et al. Apetascale direct numerical simulation study of the modelling of flame wrinkling for large-eddy simulations in intense turbulence[J]. Combustion and Flame, 2012, 159(8): 2690-2703. DOI:10.1016/j.combustflame.2011.11.020 |
[18] |
POPE S B. Turbulent flows[M]. Cambridge: Cambridge University Press, 2000.
|
[19] |
FRANZELLI B, RIBER E, SANJOSÉ M, et al. A two-step chemical scheme for kerosene-air premixed flames[J]. Combustion and Flame, 2010, 157(7): 1364-1373. DOI:10.1016/j.combustflame.2010.03.014 |
[20] |
POINSOT T, VEYNANTE D, Theoretical and numerical combustion[M]. Philadelphia: R.T. Edwards, Inc., 2005. ISBN-13: 978-1930217102 ISBN-10: 1930217102
|
[21] |
MOSS J B. Simultaneous measurements of concentration and velocity in an open premixed turbulentflame[J]. Combustion Science and Technology, 1980, 22(3/4): 119-129. DOI:10.1080/00102208008952377 |
[22] |
VEYNANTE D, TROUVÉ A, BRAY K N C, et al. Gradient and counter-gradient scalar transport in turbulent premixed flames[J]. Journal of Fluid Mechanics, 1997, 332: 263-293. DOI:10.1017/s0022112096004065 |