Jing CHEN ,Hua ZHOU ,Zhuyin REN

a Center for Combustion Energy,Tsinghua University,Beijing 100084,China

b School of Aerospace Engineering,Tsinghua University,Beijing 100084,China

KEYWORDS Bluff-body stabilized flames;Large Eddy Simulation(LES);Large-scale flow structures;Proper Orthogonal Decomposition(POD)

Abstract Large Eddy Simulations(LES)in conjunction with the Flamelet Progress Variable(FPV)approach have been performed to investigate the flame and large-scale flow structures in the bluff-body stabilized non-premixed flames,HM1 and HM3.The validity of the numerical methods is first verified by comparing the predicted velocity and composition fields with experimental measurements.Then the evolution of the flame and large-scale flow structures is analyzed when the flames approach blow-off.The analysis of instantaneous and statistical data indicates that there exists a shift of the control mechanism in the recirculation zone in the two flames.In the recirculation zone,HM1 flame is mainly controlled by the mixing effect and ignition mainly occurs in the outer shear layer.In HM3 flame,both the chemical reactions and mixing are important in the recirculation zone.The Proper Orthogonal Decomposition(POD)results show that the fluctuations in the outer shear layer are more intense in HM1,while the flow structures are more obvious in the outer vortex structure in HM3,due to the different control mechanism in the recirculation zone.It further shows that the flow structures in HM1 spread larger in the intense mixing zone due to higher temperature and less extinction.

1.Introduction

Large-scale flow structures that often appear as the vortex structures are naturally produced in many practical turbulent combustion applications.These flow structures can affect the flame stability and possibly lead to extinction/reignition1through enhancing the large-scale fuel/air mixing2and rolling up/stretching the flame surfaces3etc.In return,due to the baroclinic4and thermal expansion effects5,flame activities can also significantly alter flow structures.For example,the flow structures in reactive flows can be more broken and the vorticity magnitude is usually smaller6,7.The investigation of the subtle flow-flame interaction in turbulent reactive flows has attracted much attention8-12, especially for near-limit flames exhibiting extinction/reignition that have significance for practical engine applications13.

For the non-premixed flames,on one hand,the fundamental flame-vortex interaction has been widely investigated under low Reynolds number flows14,15.For example,the extinction and reignition dynamics in the interactions of flames with counter-rotating vortex pairs have been analyzed in the Ref.15,with the unsteadiness and curvature effects on extinction/reignition being summarized in Reynolds-Damko¨hler regime diagrams.Thevenin et al.14investigated the effects of Damko¨hler number on the interaction of flame and flow structures,and illustrated eight different interaction types in a gas turbine combustor.It is shown15,16that local flame extinction occurred frequently when the flames are subject to excessive strain rates,which may lead to flame instability.On the other hand,for highly turbulent flames,large-scale flow structures were found to have substantial impact on flame stabilization.Substantial efforts have been made to visualize and analyze these energetic large-scale structures, e.g. recirculation zone and coherent structure.For example,the flame-vortex interaction after the triangular flame holder and the blow off scaling have been reported in the Ref.13;the blow-off behavior in the bluffbody stabilized turbulent flames has been visualized in the work of17;the flame dynamics and the blow-off in the asymmetrical bluff-body stabilized flames have been captured in the work of18;the flame-vortex interaction in the turbulent swirl flames has been studied in the Ref.19with the flow structure being identified.

The POD method20,a powerful tool to extract the coherent structures in turbulent flow field,decomposes fluctuations into the linear sum of the orthogonal eigen-functions.It has been widely employed to visualize and analyze the large-scale flow structures in bluff-body flames21,swirling flames22,and flame holder23etc.Blanchard et al.23compared the POD modes from experimental measurements and numerical results downstream of the flame holder,demonstrating the overall good agreement between the shapes and magnitudes of the energetic modes.For bluff-body stabilized flames,analysis showed that obvious changes could be found for the large-scale flow structures,which play a key role of entraining mixtures into the recirculation zone,when the flames approach blow-off8.

In this study,Large Eddy Simulations(LES)of the Sydney bluff-body stabilized flames24HM1 and HM3 have been performed to investigate the evolution of the flame and large-scale flow structures approaching blow-off.These two flames,having relevance to many practical systems,are target flames for model development and validation25.In these flames,the fuel with high velocity is injected from the central jet surrounded by a bluff-body, and the coflow air enters from the outside of the bluff-body with lower velocity. The HM1 and HM3 CH4/H2(1:1 in volume) fueled flames have Reynolds numbers of 15800 and 28700 respectively,and HM3 features certain levels of local extinction.In these flames,the large scale recirculation zone is of crucial importance to flame stabilization since it entrains the downstream combustion products into the upstream. The neck region further downstream may have a high level of local extinction and subsequent reignition.Most of the previous studies focused on the validation of turbulence and combustion models in the prediction of flow and composition statistics.For example,RANS/LES simulations in conjunction with various combustion models have been reported e.g., RANS/transport PDF model26, LES/transport PDF model27,LES/conditional moment closure model28and the LES/flamelet approach29.It is observed that large eddy simulations are preferred for capturing the complex vortex shedding and recirculation,and combustion models accounting for finite rate kinetics are needed to capture local extinction/reignition in the near blow-off HM3 flame.

Fewer studies however have been reported on the analysis of flame and large-scale flow interaction26,27,30. Kim and Pitsch30studied the complex flow structures interaction in the intense mixing zone of the flame HM1 via the instantaneous plots and the analysis of the conditional scalar dissipation rates.However,the detailed analysis of flow structure through POD and the evolution of flow and flame structure with enhanced turbulence-chemistry interaction, especially when approaching blow-off,have not been reported.This study aims to investigate the flow and flame structures in HM1 and HM3 with POD analysis to reveal the subtle difference in the flow structures when the flames approach blow-off.The rest of the paper is organized as follows.In Section 2,the LES/Flamelet Progress Variable (FPV) approach and the POD method are first described,followed by a description of the computational configuration and settings.In Section 3,the predicted flow and composition statistics are first compared to experimental measurements,followed by a detailed analysis of flow and flame structures. Conclusions are in Section 4.

2.Methodology

2.1.LES-FPV approach

In the LES-FPV approach,with the assumption of incompressibility,the mass and momentum conservation equations after applying the filter operator are given by

where the bar and tilde denote the spatial filtering and densityweighted spatial filtering,respectively.andare the filtered density, pressure, velocity and kinematic viscosity,respectively.I and ?are the unit vector and the del operator.The sub-grid stressis modeled with the gradient transport assumption as

where the subgrid turbulent viscosity νtis modeled with the dynamic Smagorinsky model31.With the FPV approach,the following transport equations of the mean mixture fractionits varianceand the mean progress variableare transported and solved

in which the sub-grid scalar flux terms,i.e.andare also modeled based on the gradient transport assumption with the diffusion coefficient being modeled using the dynamic Smagorinsky model.The sub-grid scalar dissipation rateis modeled using a linear relaxation model as in the work of Ref.32,i.e.where τtrepresents subgrid turbulent time scale and the model parameter Cχis taken to be the typical value of 2.033.The progress variable Ycis defined based on the sum of the species mass fraction of CO,CO2,H2,and H2O,i.e.Yc=yCO+yCO2+yH2+yH2O,withbeing its source term.

In the FPV approach32,the density,viscosity and diffusivity,the composition(e.g.,species mass fraction and temperature)and the source term ~˙ωYcare evaluated from the pretabulated flamelets,which is parameterized in terms of mixture fraction and progress variable.Specifically,the steady flamelet equations are solved,

where ρ,T and cpare density,temperature,specific heat at constant pressure,respectively;are the mass fraction,net creation rate,and specific enthalpy of species i,respectively;χ is the scalar dissipation rate of mixture fraction dissipation rate.Then,species mass fraction,temperature and chemical source terms are parameterized by mixture fraction and scalar dissipation rate,i.e.,

where χstis the scalar dissipation rate at stoichiometry.Since,this representation may result in multiple solutions at certain scalar dissipation rates due to the S-shaped solution curve34,the progress variable Yc,a quantity independent of mixture fraction,is introduced in the FPV approach to parameterize the flamelet solution as

In the LES-FPV approach, the sub-grid turbulencechemistry interaction is modeled with the presumed probability density function(PDF)approach.Specifically,the filteredare obtained through

2.2.POD analysis

In POD,the coherent structures can be regarded as a weighted linear sum of the eigen-functions of the two-point correlation matrix, based on the assumption that different types of large-scale coherent motions in flow field lead to different POD eigen-functions.The mode with the largest eigenvalue represents the most energetic flow structure.In this study,the snapshots POD proposed by Sirovich20is adopted to solve the eigenvalue problem.Specifically,for the simulated bluffbody stabilized turbulent flames,two-dimensional flow fields u(x,t )in the axial-radial plane are saved as snapshots with a total number of N series at different time instants.The fluctuating component is computed as

where x represents a two-dimensional flow field in the axialradial plane,and 〈u〉represents time-averaged velocity.Then the auto-correlation matrix Cmncan be formed as

The eigenvalues represent the strength of each mode,i.e.,kinetics energy for the velocity field,with each POD mode being computed as

The corresponding POD mode coefficientis computed as

and the fluctuating velocity can be reconstructed from the POD modes via

In the bluff-body stabilized flames,the axial velocity component contributes most to the turbulent kinetic energy.Therefore,only the axial velocity component is applied in the POD analysis as in the work of Ref.35.

2.3.Computational configuration and settings

The schematic structure of the burner and the computational domain is shown in Fig.1.In the burner,the inner and outer diameters of the bluff-body are Dj=3.6 mm and Db=50 mm(radius Rb=25 mm),respectively.The inlet boundaries are listed in Table 1.Note that as the central jet is composed of the methane-hydrogen mixture(1:1 in volume)at 298 K,in which the speed of sound is around 610 m/s,and therefore,the low-Mach assumption still holds for the high speed HM3. The domain is 0.3 m×0.12 m×2π (6.0Db×2.4 Db×2π)in the axial(x),radial(y)and circumferential(z)directions.The mesh adopted is 320×216×64.The mesh is uniform in the circumferential direction,but non-uniform in the axial and radial directions.It is refined axially near the end of the recirculation zone and the inlet zone and refined radially near the inner and outer shear layers.Note that the large scale mixing was reported to be reasonably resolved using a grid of 256×152×64 (in axial, radial and circumferential directions respectively) for a domain of 5.5Db×1.5Db×2π in the work of Ref.30.In this work,a finer grid was applied to resolve a slightly larger domain,and therefore,the current resolution is sufficient to predict the flow field with adequate accuracy for the POD analysis.For the flamelet table in the FPV approach,100 nodes are uniformly employed in each of theanddirection.

The filtered governing equations,e.g.,Eqs.(1)and(2)and Eqs.(4)-(6)are solved using finite difference schemes with the SIMPLE algorithm for pressure-velocity coupling. The momentum equations are discretized by a second-order energy conserving scheme,while the scalar equations are discretized using the 3rd Weighted Essentially Non-Oscillatory(WENO)scheme.The temporal discretization is via the second-order semi-implicit Crank-Nicolson scheme.In the simulation,the convective boundary condition is adopted for the axial outlet while the zero gradient boundary condition for the other outlets.The Dirichlet boundary condition is applied for the inlet.A 1/7 exponential distribution law is given to the mean inlet velocity and the while noise is added to model the Root Mean Square(RMS)value of the inlet velocity.The simulation was carried out for more than two flow-through times to reach statistically stationary,followed by another six flow-through times for statistics.The time step is set by limiting the Courant-Friedrichs-Lewy (CFL) number to be 0.5. Four hundred three-dimensional instantaneous data(N=400)were sampled with a time interval of 5×10-4s for post-processing.The sampling process lasts more than 100 turn-over time of the largest eddy in the region for POD analysis,and a smaller value of N=300 hardly made any difference to the first three POD modes.These demonstrate that N=400 satisfies the weak correlation condition for POD analysis.

Fig.1 Schematic structure of the burner and the corresponding computational domain.

3.Results and discussion

3.1.Statistics of flow and composition

The time-averaged streamlines of HM1 is shown in Fig.2 to illustrate the locations of the recirculation zone and the shear layers.In Fig.2,"OVS"and"IVS"represent the outer/inner vortex structures in the recirculation zone,"OSL"and"ISL"the outer/inner shear layers,"Ep"the end point of the recirculation zone.The left box in the figure marks the"ignition zone"which covers parts of the outer shear layer and the outer vortex structure.The right box marks the intense mixing zone where the mixing of the flow structures from the vortex shedding and the central jet is intense.The mixing zone locates at the tail of the recirculation zone and is also pointed out via the scalar dissipation rate analysis in the work of Kim30.In fact,extinction in this zone is also intense in HM3 flame and this can also be marked as‘‘extinction zone”.The flow structures in these two zones are studied in the following sections.

Figs.6 and 7 show the radial profiles of the time-averaged mean and the RMS of the mixture fraction,temperature,CO and OH profiles at three axial locations x/Db=0.26,0.90,1.80 from HM3 flame.The length of the recirculation zone in HM3 is 1.5Db,which is close to the 1.6Dbin HM1.In this sense,the selected axial locations(x/Db=0.26,0.90,1.80)are comparable to HM1.In Fig.6(a),the mixture fractions at x/Db=0.26,x/Db=0.90 are close to the stoichiometric combustion state (Zst=0.05), which indicates that the flame may fill the front of the entire recirculation zone in HM3.In Fig.6(b),the predicted temperature inside the recirculation zone agrees well with the experimental data,only downstream of the recirculation zone(x/Db=1.80)the predicted values are slightly larger.According to the experimental analysis by Dally et al.36,HM3 is extinguished downstream of the neck region of the recirculation zone.This demonstrates that the present LES-FPV approach is able to capture the local extinction to some extent.

The time-averaged temperature in HM3 upstream of the recirculation zone(shown in Fig.6(b))is about 2100 K,whichis much higher than 1600 K in HM1(Fig.3(c)),in which the jet flow velocity is much lower and the fuel residence time is longer.For HM1,the mixture is in fuel-rich state inside the recirculation zone and the flame occurs at the outer shear layer and downstream of the recirculation zone.The high temperature in the recirculation zone is mainly caused by the hightemperature products transported from the outer shear layer.While in HM3,little fuel enters the recirculation zone and the mixture is stoichiometric upstream of the recirculation zone where stable combustion takes place.The high temperature in the recirculation zone is produced by both mixing and chemical reactions,resulting in much higher temperature than HM1.Meanwhile,it is found in Fig.7 that downstream of the recirculation zone,the temperature of HM3 is lower than that of HM1 due to the local extinction.

Table 1 Inlet parameters of HM1and HM3 flames.

Fig.2 Time-averaged streamlines in longitudinal section from HM1 flame.

The comparisons of the time-averaged mean and RMS of components CO,CO2,H2,and H2O show that the distributions are all in good agreement with the experiment data.In HM1,the OH peak distributes along the outer edge of the recirculation zone while its mass fraction is almost zero at other positions(see Fig.5(c)).However,in HM3,the OH peak distributes in the middle of the recirculation zone and almost fills the entire recirculation zone(see Fig.6(d)).All these indicate the different mechanism in the recirculation zone in HM1 and HM3 flames.Meanwhile,the distribution of OH in the intense mixing zone in HM3 is lower than that of the HM1 flame,indicating certain level of local extinction.In general,the comparison of the time-averaged profiles with the experimental data demonstrates that the LES-FPV approach predicts large scale mixing and local extinction with adequate accuracy for the POD analysis in the subsequent section.

3.2.Large-scale flow structures in HM1 and HM3

Fig.3 Radial profiles of time-averaged mean axial velocity,mixture fraction and temperature at three locations x/Db=0.26(upper row),0.90(middle row),1.80(lower row)from HM1 flame compared with experimental data.

Fig.4 Radial profiles of RMS of the axial velocity,mixture fraction and temperature at three locations x/Db=0.26(upper row),0.90(middle row),1.80(lower row)from HM1 flame compared with the experimental data.

Fig.5 Radial profiles of the time-averaged mean mass fractions of the species CO2,CO and OH at three locations x/Db=0.26(upper row),0.90(middle row),1.80(lower row)from HM1 flame compared with the experimental data.

The large-scale flow structures in HM1 and HM3 flames are studied based on the instantaneous and the statistical data in this section.Fig.8(a)shows that the recirculation zone is in fuel-rich state and the stoichiometric mixture fraction iso-line locates at the outer shear layer of the recirculation zone in HM1.The peak temperature in Fig.8(b)distributes along the iso-line of the stoichiometric mixture fraction,indicating that the HM1 flame features the typical combustion characteristics of a diffusion flame.The similar feature can be also found in HM3 flame.In the recirculation zone,the flame locates at the outer shear layer and moves toward the center jet as it goes downstream.This indicates that the ignition mainly appears along the outer shear layer. Meanwhile,upstream of the recirculation zone,the distribution of mixture fraction is rather uniform,with low scalar dissipation rate.The uniformly distributed mixture fraction and the low scalar dissipation rate produced by the recirculation zone helps to stabilize the flame.Fig.8(c)shows the distribution of instantaneous scalar dissipation rates.It can be seen that there exist many strips,which gradually increases downstream along the iso-line of the stoichiometric mixture and reaches a maximum at the tail of the recirculation zone,which again illustrates the low and uniform scalar dissipation rate in the recirculation zone,as in the work of Ref.30.

Fig.6 Radial profiles of time-averaged mean mixture fraction,temperature,CO and OH at three axial locations x/Db=0.26(upper row),0.90(middle row),1.80(lower row)from HM3 flame compared with experimental data.

Fig.7 Radial profiles of the RMS of mixture fraction,temperature,CO and OH at three axial locations x/Db=0.26(upper row),0.90(middle row),1.80(lower row)from HM3 flame.

The time-averaged streamlines of the HM1 and HM3 flames are shown in Fig.9.The dash-dot lines are the isolines of the axial velocity U=0.With the increase of the central jet velocity,the length of the recirculation zone becomes slightly shorter,from 1.6Dbin HM1 to 1.5Dbin HM3.The change is not obvious and the double vortex structures do not change significantly.However,the distribution of the mean mixture fraction in the recirculation zone has changed drastically as shown in Fig.10.The HM1 flame has good mixing inside the recirculation zone and the mixture fraction inside the recirculation zone is higher than 0.15.On the contrary,the HM3 flame has little mixing and the time-averaged mixture fraction upstream of the recirculation zone is between 0.05 and 0.10.This results in the time-averaged temperature distributions in Fig.11,where the maximum temperature in HM1 locates near the outer shear layer of the recirculation zone,while the maximum value appears upstream of the recirculation zone in HM3.This indicates that the mixing effect leads to the uniform and high temperature in the recirculation zone in HM1,while the higher temperature in the recirculation zone in HM3 is produced by both mixing and chemical reactions.Fig.12 shows the time averaged mass fraction of OH where the difference in the flame position near the recirculation zone between HM1 and HM3 flames can also be observed.

Fig.8 Contour plots of mixture fraction,temperature and scalar dissipation rate in HM1 flame.

Fig.9 Time-averaged flow fields in HM1 and HM3 flames.Red line denotes U=0 iso-line.

Fig.10 Time-averaged mixture fraction in HM1 and HM3 flames.

Fig.11 Contour plots of time-averaged temperature in HM1 and HM3 flames.

Fig.12 Contour plots of the time-averaged OH mass fraction in HM1 and HM3 flames.

Fig.13 Comparison of the first three POD modes in the ignition zone in HM1 and HM3 flames.

3.3.Large-scale flow structures from POD analysis

The large-scale flow structures from POD analysis in the‘‘ignition zone”and‘‘extinction zone”(intense mixing zone)are analyzed in this section.The first three eigenmodes in the‘‘ignition zone”from HM1 and HM3 flames are identified in Fig.13,which shows the projected axial kinetic energy in the corresponding mode in two-dimensional space.The first mode(Fig.13(a)and(d))shows many similarities while the second mode(Fig.13(b)and(e))and third mode(Fig.13(c)and(f))exhibit large difference.The large-scale flow structures in the outer vortex structure in HM3 are more obvious,which is due to more intensive chemical reactions in this region.Meanwhile,the fluctuations in the outer shear layer are more intense in HM1.Fig.14 shows the first three eigenmodes in the‘‘extinction zone”.For the first mode,HM1 and HM3 are very similar.For the second mode,the paired vortical regions appearing spread much larger in HM1.For the third mode,the peaks are much nearer to the centerline in HM3.This indicates that the structures spread larger in HM1 than in HM3.The larger structures in HM1 is caused by the high temperature.While in HM3,the vortex shedding from the outer shear layer and the intensive mixing causes local extinction,resulting in smaller structure.

Fig.14 Comparison of the first three POD modes in the extinction zone in HM1 and HM3 flames.

The POD analysis illustrates that the large-scale flow structures in the ignition and extinction zones are quite different.Within the ignition zone,the fluctuations in the outer shear layer are more intense in HM1 while the large-scale flow structures in the outer vortex structure are more obvious in HM3.This is caused by the different control mechanism in the recirculation zone.Within the extinction zone,the flow structures in HM1 spread larger in the intense mixing zone due to higher temperature and less extinction compared to HM3.

4.Conclusions

Large eddy simulations with FPV approach have been performed to investigate the flame and large-scale flow structures in the bluff-body stabilized non-premixed flames,especially when approaching blow-off.

(1)Good agreement between the experimental data and the simulation has been obtained for the flow and composition statistics in HM1 and HM3 flames,demonstrating the capability of the LES/FPV approach to reasonably predict large scale mixing and local extinction.

(2)The instantaneous data and the time-averaged data of the two flames show different flow and mixing characteristics,and there exists a shift of the control mechanism in the recirculation zone:HM1 is mainly controlled by the mixing,while HM3 is governed by both the chemical reactions and mixing.

(3)The POD analysis shows that the fluctuations in the outer shear layer are more intense in HM1,while the flow structures are more obvious in the outer vortex structure in HM3.This is closely related to the different control mechanism in the recirculation zone.The POD analysis also illustrates that the flow structures in HM1 spread larger in the intense mixing zone due to higher temperature and less extinction compared to HM3.

Acknowledgements

The work was supported by the National Natural Science Foundation of China(Nos.91441202 and 51476087).Simulations are done with the computational resources from the Tsinghua National Laboratory for Information Science and Technology.