Large-eddy simulations of turbulent methane jet flames with filtered mass density function

International Journal of Heat and Mass Transfer 53 (2010) 2551–2562

Contents lists available at ScienceDirect

International Journal of Heat and Mass Transfer journal homepage: www.elsevier.com/locate/ijhmt

Large-eddy simulations of turbulent methane jet flames with filtered mass density function M. Yaldizli, K. Mehravaran, F.A. Jaberi * Department of Mechanical Engineering, Michigan State University, East Lansing, MI 48824-1226, USA

a r t i c l e

i n f o

Article history: Received 1 June 2009 Received in revised form 13 December 2009 Accepted 13 December 2009 Available online 20 February 2010 Keywords: LES Filtered mass density function PDF methods Monte-Carlo simulations Methane jet flames

a b s t r a c t The filtered mass density function (FMDF) model (Jaberi et al. 1999 [1]) is employed for large eddy simulations (LES) of ‘‘high speed” partially-premixed methane jet flames with the ‘‘flamelet” and ‘‘finite-rate” kinetics models. The FMDF is the joint probability density function (PDF) of the scalars and is determined via the solution of a set of stochastic differential equations. The LES/FMDF is implemented using a highly scalable, parallel hybrid Eulerian–Lagrangian numerical scheme. The LES/FMDF results are shown to compare well with the experimental data for all flow conditions when ‘‘appropriate” reaction and mixing models are employed. Ó 2009 Elsevier Ltd. All rights reserved.

1. Introduction The experimental data generated by the Sandia National Laboratory researchers and their collaborators for turbulent jet flames [2,3] have been widely employed for the assessment of ReynoldsAveraged Simulation (RAS) and Large-Eddy Simulation (LES) [4– 6] models. For example, Jones and Kakhi [7] used the RAS based probability density function (PDF) method with linear mean square estimation (LMSE) and Curl’s coalescence-dispersion (CCD) mixing closures to simulate Sandia’s non-piloted non-premixed jet flames at various jet speeds. Their results indicate that at lower jet speed of 41 m=s (flame L), there is no significant extinction, but at higher speed of 48 m=s (flame B), there are some. In the extinction–reignition region ð15 6 x=D 6 35Þ of the flame B, the predicted temperatures by the LMSE model are shown to be lower than the experimental values, while the CCD model overpredicts the experimental data at the same locations. Muradoglu et al. [8] have conducted similar jet simulations with flamelet chemistry and LMSE models and have found good agreement with the low-speed flame L data. The assumed PDF method was used by Landenfeld et al. [9] to simulate the Sandia’s near-equilibrium piloted jet flame (flame D). The reaction was implemented with the intrinsic low-dimensional manifolds (ILDM) algorithm. Additionally, transport equa-

* Corresponding author. Tel.: +1 517 432 4678; fax: +1 517 353 1750. E-mail addresses: [email protected] (M. Yaldizli), [email protected] (K. Mehravaran), [email protected] (F.A. Jaberi). 0017-9310/$ - see front matter Ó 2009 Elsevier Ltd. All rights reserved. doi:10.1016/j.ijheatmasstransfer.2009.12.061

tions for the mean and variance of mixture fraction and species mass fractions are solved, together with a Beta PDF for the mixture fraction and a Reynolds stress closure for the turbulence. The numerical results were found to compare well with the flame D data. Xu and Pope [4] also simulated the piloted jet flames at various jet speeds with a joint velocity-composition-turbulent frequency RAS/PDF model. The 16 species augmented methane mechanism have been implemented with the insitu adaptive tabulation (ISAT) [10] algorithm for the chemistry. The Euclidian minimum spanning tree (EMST) model [11], which treats the mixing locally in the composition space, is used for the mixing model. Overall, the results were found to be in good agreement with the experimental data. The LES method has also been used for predictions of turbulent flames with various subgrid-scale (SGS) models [5,12,14,15]. One of these models is the filtered density function (FDF) [16]. The FDF is the counterpart of the PDF method in RAS [16,17]. The fundamental advantage of the LES/FDF is that it accounts for the effects of chemical reaction exactly while calculating the resolved turbulence field, thus allowing a more reliable prediction of turbulent flames. The literature on SGS closures via FDF has been growing at a relatively fast pace since its first introduction [16]. The scalar FDF is considered in Refs. [18–21], the scalar filtered mass density function (FMDF) in Refs. [22–24], the velocity FDF in Ref. [25] and the velocity–scalar FDF in Ref. [26]. Some recent data on experimental validation of FDF are also available [27]. In the present work, the scalar FMDF methodology is employed as a SGS closure for LES of high-speed turbulent methane jet

2552

M. Yaldizli et al. / International Journal of Heat and Mass Transfer 53 (2010) 2551–2562

Nomenclature

Symbols filtered value h i‘ Favre filtered value h iL conditional Favre filtered value h j iL ^ quantities which depend only on the scalar composition, i.e. b Sð/; x; tÞ  Sð/ðx; tÞÞ 1 ambient Subscript k dummy index Superscripts + properties of the stochastic Monte Carlo particles (n) index of the Monte Carlo particles

x z y

Methane + Air

flames. The main objectives are (1) to further assess the applicability and the extent of validity of the FMDF method for modeling of complex turbulent flames experiencing local flame extinction, (2) to investigate the accuracy of mixing models, (3) to study the effects of chemical kinetics on LES predictions of near-equilibrium and non-equilibrium flames. For (3), the reaction is implemented following two different approaches. In the first approach, the chemical species are obtained by direct solution of the reaction equations. In the second approach, only the FMDF of the mixture fraction is considered; the temperature and species concentrations are obtained from a set of flamelet tables. The tables are generated by solving the steady-state one-dimensional opposed jet equations with a detailed (GRI) mechanism. Considering the limitations of steady-state flamelet assumption in the second approach, it is expected that the local extinction not to be predicted by the FMDF/ flamelet method. This study is focused on Sandia’s ‘‘high speed” piloted turbulent methane jet flames D and F [2]. Flame D involves limited regions of local extinction, while flame F tends towards total extinction. The existence of different levels of local extinction in these flames provides a good means of assessing the capabilities of the models to predict realistic combustion systems. The geometrical configurations in these two flames are the same (see Fig. 1), but the jet inlet velocity in flame F is twice of that in flame D. The formulation of LES/FMDF methodology is presented in Section 2, followed by Section 3 which describes the consistency issues, the parallelization algorithm, the implementation of chemistry, and the computational parameters. Results obtained by the LES/FMDF for flames

Greek symbols r gradient operator D grid spacing in LES Kronecker delta dij DE ensemble domain width DH grid level filter width DH0 secondary level filter width l molecular viscosity, l ¼ qm Xm SGS mixing frequency / scalar field the compositional values of scalar a /a /þ the compositional values of stochastic scalar a a q density r number of scalars, r ¼ Ns þ 1 sij molecular stresses s/ scalar mixing timescale

Pilot

specific heat of the mixture at constant pressure molecular diffusion coefficient subgrid diffusion coefficient joint scalars filtered mass density function enthalpy enthalpy of species a enthalpy of formation of species a ith component of the diffusive flux of scalar a thermal conductivity of the mixture molecular Lewis number ith component of the subgrid scalar flux of species a number of species molecular Prandtl number universal gas constant Reynolds number SGS Schmidt number reaction rate of species a strain rate tensor temperature reference temperature time ith component of the velocity vector molecular weight of species a position vector ith component of the position vector Lagrangian position of the particles streamwise coordinate mass fractions of species a coordinates defining the plane normal to x mixture fraction

Coflow

cp D Dt FL h ha 0 ha a Ji k Le Mai Ns Pr R0 Re Sct Sa Sij T T0 t ui Wa x xi Xi x Ya y; z Z

Fig. 1. A schematic view of Sandia’s piloted methane jet flame. Main and pilot jet diameters are 7:2 mm and 18:2 mm, respectively.

D and F are discussed in Section 4, and the main conclusions are stated in 5.

2553

M. Yaldizli et al. / International Journal of Heat and Mass Transfer 53 (2010) 2551–2562

@h/a iL ; @xi

2. Formulation

Mai ¼ hqi‘ Dt

In the hybrid Eulerian–Lagrangian LES/FMDF methodology, the filtered velocity equations are solved together with the joint scalar FMDF equation. These equations are presented in two different sections. In LES, the resolved variables are obtained by application of the spatial filtering operator,

It must be indicated here that this model is not used directly in the FMDF approach but modeled FMDF transport equation is constructed to be consistent with it as discussed below.

hf ðx; tÞi‘ ¼

Z

a ¼ 1; 2; . . . ; r:

ð9Þ

2.2. The filtered mass density function (FMDF)

þ1

f ðx0 ; tÞHðx0 ; xÞ dx0 ;

ð1Þ

1

where H is a positive, spatially and temporally invariant, and localized filter function with the properties, HðxÞ ¼ HðxÞ, and R1 HðxÞ dx ¼ 1. The filtered value of the transport variable f ðx; tÞ 1 is represented by hf ðx; tÞi‘ , where x ¼ ðx; y; zÞ denotes the spatial coordinates, with x axis parallel to the jet and y; z perpendicular to x. In variable density flows it is more convenient to consider the Favre filtered quantity, hf ðx; tÞiL ¼ hqf i‘ =hqi‘ : 2.1. Filtered equations With the application of the filter function to the transport equations one can derive the following filtered equations:

@hqi‘ @hqi‘ hui iL ¼ 0; þ @t @xi @hqi‘ huj iL @hqi‘ hui iL huj iL @hpi‘ @hsij i‘ @Tij ¼ þ  ; þ @t @xi @xj @xi @xi @hqi‘ h/a iL @hqi‘ hui iL h/a iL @hJ a i @Mai ¼ i ‘ þ hqSa i‘ ; þ @t @xi @xi @xi a ¼ 1; 2; . . . ; r:

ð2Þ ð3Þ

Ns X h/a iL Wa 1

  @hui iL @huj iL 2 @huk iL ; hliL ¼ Prhk=cp iL ; hsij i‘  hliL þ  dij 3 @xj @xi @xk   @h/a iL 1 k hJ ai i‘  hqi‘ hDiL ; hDiL ¼ : hqi‘ Le cp L @xi

ð5Þ

Z

þ1

qðx0 ; tÞf½w; /ðx0 ; tÞHðx0  xÞ dx0

ð10Þ

1

f½w; /ðx; tÞ ¼ d½w  /ðx; tÞ 

r Y

d½wa  /a ðx; tÞ;

ð11Þ

a¼1

where d denotes the delta function and w denotes the composition domain of the scalar array /ðx; tÞ. The term f½/; wðx; tÞ is the ‘‘finegrained” density, and Eq. (10) implies that the FMDF is the mass weighted spatially filtered value of the fine-grained density. To ensure that the FMDF has all the properties of the PDF, we will only consider ‘‘positive” filter functions for which all the moments R1 m x HðxÞ dx exist for m P 0. The following exact deterministic 1 equation describes the variations of FMDF in space and time:

  @F L @hui iL F L @ @f ¼ qD þ @xi @xi ‘ @t @xi    2 @ @/ @/  qD a b jw F L =q^ @wa @wb @xi @xi ‘ 

@½hui jwi‘  hui iL F L @½b S a ðwÞF L   : @wa @xi

ð12Þ

In Eq. (12), hAjBi‘ denotes the filtered value of the variable A, ‘‘conditioned” on B, and the hat is used for the quantities which are dependent only on the scalar field. In Eq. (12), the last term on the right-hand side represents the effects of chemical reaction and is in a closed form. The second and the third terms on the right-hand side are unclosed. They represent the effects of SGS mixing and convection, respectively. The convective flux is modeled here as

@ðF L =hqi‘ Þ ; @xi

ð6Þ

½hui jwi‘  hui iL F L ¼ hqi‘ Dt

ð7Þ

where mt is the SGS viscosity and is to be determined by hydrodynamic LES and Sct is the SGS Schmidt number. The first ‘‘moment” of Eq. (13) recovers the model given in Eq. (9). The closure for the SGS mixing can be via any of the ones used in RAS/PDF methods [29]. The simplest one is the LMSE model,

In Eqs. (2)–(4), the hydrodynamic SGS closure problem is associated with Tij ¼ hqi‘ ðhui uj iL  hui iL huj iL Þ and Mai ¼ hqi‘ ðhui /a iL  hui iL h/a iL Þ, denoting the SGS stresses and the SGS scalar fluxes, respectively. In reacting flows, an additional model is required for the filtered reaction term hqSa i‘ ¼ hqi‘ hSa iL . It will be shown in the next section that both hqi‘ and hSa iL are determined exactly with the knowledge of the FMDF. For Tij , the model used by Jaberi et al. [1] is considered:

  1 2 Tij ¼ 2hqi‘ mt hSij iL  hSkk iL dij þ C I hqi‘ Edij ; 3 3

F L ðw; x; tÞ 

ð4Þ

In these equations, the scalar (composition and energy) field is P s ha /a , in which denoted by /a  Y a ; a ¼ 1; 2; . . . ; N s ; /r  h ¼ Na¼1 RT 0 0 0 ha ¼ ha þ T 0 cpa ðT ÞdT . Eqs. (2)–(4) are closed by the following constitutive relations:

hpi‘  hqi‘ R0 hTiL

The scalar filtered mass density function (FMDF) is defined as [19,1,24,20]

ð8Þ

where hSij iL is the resolved strain rate and E ¼ jhu0i iL hu0i iL  pffiffiffi hhu0i iL i‘0 hhu0i iL i‘0 j, mt ¼ C R DH E; u0i ¼ ui  Ui (Ui is a reference velocity in the xi direction and the subscript ‘0 denotes filtering at the secondary level with a characteristic size larger than that at grid level). We have found the performance of this model to be better than the Smagorinsky type closures, and we did not even have to implement the ‘‘dynamic” procedure [28] for evaluation of the model coefficients. For the closure of the subgrid mass flux Mai , the model most often used in LES of non-reacting flows is used

@2 @wa @wb



qD

@/a @/b jw @xi @xi

Dt ¼ mt =Sct ;

ð13Þ



 @ ^ ¼ F L =q ½Xm ðwa  h/a iL ÞF L ; @wa ‘ ð14Þ

where Xm ðx; tÞ is the frequency of mixing within the subgrid and is modeled as Xm ¼ C X ðhDiL þ Dt Þ=ðD2H Þ. The mixing model coefficient C X represents the velocity-to-scalar timescale ratio [4]. It needs to be either empirically specified [4,30] or directly obtained via a dynamic model [31,32]. In the present work, a range of constant empirical values and two different dynamic models are employed to calculate C X : 1. In this approach, the LMSE model is used in its standard form, while the mixing model coefficient C X is assumed to be constant [30]. Here, C X ¼ 5 is used as the base value for the reference case since it generates the closest predictions to the experimental data. However, other values of C X are also considered.

2554

M. Yaldizli et al. / International Journal of Heat and Mass Transfer 53 (2010) 2551–2562

2. In this approach, the mixing coefficient C X is computed as a function of time and space by a dynamic method. The dynamic method is based on the generalized subgrid variance, r2a ¼ h/2ðaÞ iL  h/ðaÞ i2L and the total subgrid scalar dissipation rate, a . The latter is defined as



a ¼ 2hqi‘ ðhDiL þ Dt Þ

@h/ðaÞ iL @h/ðaÞ iL @xi @xi

 ð15Þ

¼ 2Xm hqi‘ ðh/2ðaÞ iL  h/ðaÞ i2L Þ; where the subscripts in parentheses are excluded from the summation convention. The model for a in Eq. (15) is obtained by assuming that the dissipation term balances the production term in the generalized subgrid variance equation [19,31]. By using Xm ; r2a and a , the model coefficient can be computed locally as

CX ¼

D2H

rh/ðaÞ iL  rh/ðaÞ iL

r2a

;

a ¼ 1; 2; . . . ; Ns :

ð16Þ

Evidently, this formulation allows different coefficients for different species, meaning ‘‘Differential diffusion” effects may be included. However, the computational time needed for the calculation of C X could become excessive when realistic chemistry models with numerous species are employed. Therefore, to reduce the computational overhead, a modified version of Eq. (16) is considered in which the mixture fraction variable Z is used,

CX ¼

D2H rhZiL  rhZiL : hZ 2 iL  hZi2L

ð17Þ

In the numerical implementation of the model, the gradient of the Favre-filtered mixture fraction rhZiL is calculated from the fixed finite difference grid points with a purpose of decreasing the numerical noise, while the subgrid scalar variance is determined from the Monte Carlo particles. An additional averaging over the homogeneous (azimuthal) direction is also performed to prevent numerical instabilities. In the discussion below, this dynamic model is referred to as DM1. 3. In this approach, a dynamic method proposed by Raman and Pitsch [31] is employed for calculating the coefficient C X . This model is referred to as DM2. In the DM2 model, the subgrid scale mixing frequency Xm is again evaluated from the subgrid scalar variance and the dissipation rate. However, in contrast to DM1, the subgrid scalar variance in DM2 is obtained by a model. The equations [31] describing the conserved scalar varg iance Z 002 and the dissipation rate v are

g e  r Z; e Z 002 ¼ C Z D2H r Z e v ¼ 2ð De þ Dt Þr Ze  r Z;

ð18Þ ð19Þ

where the operator ð.g . .Þ denotes the Favre-filtering operator h. . . iL . Based on these equations, the scalar mixing timescale s/ (which is the inverse of subgrid mixing frequency ð1=Xm Þ) is computed as

s/

C Z D2H g ; ¼ Z 002 =v ¼ e þ DT Þ 2ð D

ð20Þ

where the model coefficient C Z ¼ 2=C X is obtained by a dynamic procedure [33]. The IEM model has some limitations [29] which are not discussed here, but per results obtained in previous studies [29,34], it can be safely indicated that while the LMSE model is not quite satisfactory in RAS/PDF, it functions reasonably well in LES/FMDF, and its numerical implementation is computationally convenient.

With the closures for the SGS convection and mixing, the modeled FMDF transport equation may be written as

  @F L @½hui iL F L  @ @ðF L =hqi‘ Þ ¼ hqi‘ ðhDiL þ Dt Þ þ @xi @xi @xi @t @ @½b SaFL þ ½Xm ðwa  h/a iL ÞF L   : @wa @wa

ð21Þ

The transport equations for the SGS moments can be obtained by direct integration of the above FMDF equation over the composition domain. The equations for the first subgrid Favre moment, h/a iL is identical with that in Eq. (4). 3. Numerical solution The numerical method used for LES/FMDF has two major components. For the solution of hydrodynamic field, a high-order accurate finite difference (FD) method which has proven effective for LES [17] is employed. The FMDF is obtained by a Lagrangian Monte Carlo (MC) method. The FD discretization procedure is based on the ‘‘compact parameter” scheme [35] which yields up to 6th order spatial accuracy. All finite difference operations for obtaining the filtered values are performed on fixed and uniform grid points. The most convenient means of solving the FMDF transport equation is via the ‘‘Lagrangian Monte Carlo” procedure [38]. The basis of this procedure is the same as that described in recent papers [19,1,24,25]. Therefore, here only some basic features of the procedure are described. With the Lagrangian procedure, the FMDF is represented by an ensemble of computational ‘‘stochastic elements” (or ‘‘particles”) which are transported in the ‘‘physical space” by the combined actions of large scale convection and diffusion (molecular and subgrid). In addition, transport in the ‘‘composition space” occurs due to chemical reaction and SGS mixing. In doing so, the notional particles evolve via a ‘‘stochastic process,” described by the following stochastic differential equations (SDEs)

dX i ðtÞ ¼ Di ðXðtÞ; tÞ dt þ EðXðtÞ; tÞ dWi ðtÞ;

d/þa ðtÞ ¼ Ra ð/þ ; tÞ dt; ð22Þ

where X i is the Lagrangian position of the particles, D and E are known as the ‘‘drift” and ‘‘diffusion” coefficients, and Wi denotes the Wiener–Lévy process [39]. /þ a denotes the scalar value of the particle with the Lagrangian position vector X i . Eq. (22) defines what is known as the general ‘‘diffusion” process; thus the PDFs of the stochastic processes (X i ðtÞ; /þ a ðtÞ) are governed by the Fokker–Planck equation. A comparison between the standard Fokker– Planck equation corresponding to Eq. (22) with the FMDF equation (Eq. (21)) under consideration identifies the parameters of Eq. (22),

qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi E  2ðhDiL þ Dt Þ; 1 @½hqi‘ ðhDiL þ Dt Þ ; hqi‘ @xi Ra   Xm ð/þa  h/a iL Þ þ b S a ð/þ Þ: Di hui iL þ

ð23Þ

With this analogy, the FMDF is represented by an ensemble of Monte Carlo particles, each with a set of scalars /ðnÞ a ðtÞ ¼ /a ðXðnÞ ðtÞ; tÞ and Lagrangian position vector XðnÞ . A splitting operation then can be employed in which the transports in physical and compositional domains are treated separately. The simplest means of simulating the spatial transport in Eq. (22) is via the Euler–Maruyamma approximation [40]. The transfer of information from the fixed finite difference points to the location of the Monte Carlo particles are conducted via (fourth and second order) interpolation. The filtered scalar quantities are calculated by weighted averaging of the particle values over space with average volumetric size of DE . Ideally, DE has to be very small.

2555

M. Yaldizli et al. / International Journal of Heat and Mass Transfer 53 (2010) 2551–2562

The SGS empirical ‘‘constants” in LES/FMDF include C I ; C R ; Sct . Based on our experience, in which a variety of different flows (2D and 3D, constant and variable density, different chemistry schemes, etc.) are considered, we have determined C I  0:01; C R  0:02; Sct  0:4  0:7. Interestingly, the range of some of these values is the same as that typically used in equivalent models in RAS [36]. The magnitudes of the molecular parameters are close to those typically used for hydrocarbon–air flames [37]. For methane hk=cp iL  2:58  105 ðhTiL =298Þ0:7 ; Pr  0:75 and cp a is specified through polynomial fits as functions of the temperature.

3.1. Consistency As stated before, the equation governing the first subgrid Favre-moment of the scalar /a obtained from the FMDF equation is identical with the filtered scalar equation (4), indicating that the filtered temperature, scalar and density may be alternatively obtained from the Eulerian, finite difference (FD) or Lagrangian, Monte Carlo (MC) solutions. This implies a mathematical consistency between FD and MC parts of the hybrid scheme. However, due to finite grid size in FD and limited number of MC particles, consistency may not be achieved in practice. In the following, the possible sources of inconsistency are identified and conditions leading to a consistent solution are discussed. Four preliminary simulations with constant values of C / ¼ 2C X ¼ 8; 16 and DE ¼ D (D is grid spacing) and DE ¼ 2D have been performed to study the consistency of LES/FMDF for conditions of Sandia’s piloted jet flames. The reaction is turned off in these simulations but the variable density/temperature effects are still important due to pilot flame at inflow. The instantaneous temperature profiles for cases with DE ¼ 2D and DE ¼ D (not shown) indicates that the difference between the MC and FD solutions becomes negligible when DE ¼ D. Moreover, the mixing model coefficient, which does not appear explicitly in the Favre-filtered equation, does not seem to have a significant effect on the filtered temperature when DE ¼ 2D, even though the temperature profiles seem to be slightly more diffused for higher C / . For the reacting simulations considered in the next section, DE was chosen to be equal to 2D. In all of these simulations we have found the difference between the FD and MC solutions to be less than 3% for instantaneous temperatures and even lower for time averaged values.

1. In the first approach, all processors are solving the same FD equations for the whole domain. In the MC part, particles are equally divided among the processors. Since the particles are not transferred between the processors, there is no communication load for the MC calculations and the load distribution is exactly uniform. Inter-processor communications are only required in the averaging operations where a local particle average value is calculated on each FD cell. 2. The second approach is considerably more complex than the first one. However, it is proven to be scalable and much more efficient. In this approach, the FD domain and the MC particles are distributed among the processors. A limited number of communications are required when a subgroup of particles are located in the boundaries of the FD subdomain defined for each processor. With this, the communication required for the averaging process described in the first parallelization method is eliminated. Limited communications are also required for calculating the ensemble averages on the boundary points and for interpolating a FD quantity to a particle location.

3.3. Chemistry In quantitative comparison with the laboratory data, the role that the chemical kinetics model plays may become important. In this work, the chemistry model is based on (i) non-equilibrium (finite-rate) and (ii) near-equilibrium models. In (i), the finite rate kinetics effects are modeled with a one-step global mechanism, or a 12-step reduced mechanism. In (ii), the FMDF transport equation for the mixture fraction is solved together with a set of flamelet tables, generated by laminar opposed jet flame simulations and detailed kinetics models. All other thermo-chemical variables are constructed from the flamelet data. Additionally, in (i), the transport equation for the sensible enthalpy RT 0 ðhs ¼ T 0 cp ðT 0 Þ dT Þ is solved

@hqi‘ hhs iL @hqi‘ hui iL hhs iL @hJ r i @Mri 0 ¼ i ‘ þ hqSa ha i‘ ; þ @t @xi @xi @xi

ð24Þ

where

hJ ri i‘  

1 k @hhs iL h i ; Le cp L @xi

Mri ¼ hqi‘ Dt

@hhs iL @xi

ð25Þ

0

3.2. Parallelization As mentioned before, in the hybrid LES/FMDF methodology the filtered continuity, momentum and energy equations (Eqs. (2)–(4)) are solved by a finite difference (FD) method over a fixed Eulerian grid system. On the other hand, all scalars are obtained from the FMDF by solving its transport equation with a Lagrangian stochastic Monte Carlo (MC) method. The employed MC method involves grid-free particles interacting with the fixed background Eulerian grid. However, there is no inter-particle interactions, suggesting that MC calculations are potentially efficient in a parallel environment. To obtain a good statistical representation, a significant number of particles are required within each FD cell. Typically, about ten particles are required for each FD grid for DE ¼ 2D. A typical hybrid (FD–MC) simulation with ten particles per cell is about 5 times more expensive than its corresponding FD simulation with no MC particles. Hence, the usage of ‘‘complex” multi-step kinetics mechanisms in these simulations would only be possible by proper application of parallel processing techniques. Two different parallelization schemes may be used:

and the term hqSa ha i‘ is obtained from the FMDF and Monte Carlo particles. In (ii), the transport equation for the hEiL ¼ hRTiL is solved (in the flamelet table RT is a function of mixture fraction). This equation, as derived by multiplying the modeled FMDF equation by E and integrating over the mixture fraction space, is

@hqi‘ hEiL @hqi‘ hui iL hEiL @hJr i @Mri ¼ i ‘ þ @t @xi @xi @xi þ hqi‘ Xm ½hZGiL  hZiL hGiL ;

ð26Þ

where

G¼

dEðZÞ ; dZ

hJ ri i‘  

  1 k @hEiL Le cp L @xi

Mri ¼ hqi‘ Dt

@hEiL : @xi

ð27Þ

The reduced mechanism considered in (i) is the 12-step mechanism of Sung et al. [41]. This mechanism is developed from the GRI 1.2 detailed mechanism and involves 16 species ðH2 ; H; O2 ; OH; H2 O; HO2 ; H2 O2 ; CH3 ; CH4 ; CO; CO2 ; CH2 O; C2 H2 ; C2 H4 ; C2 H6 ; N2 Þ. It contains more unsteady intermediates than the conventional 4 and 5-steps mechanisms, and it has proven effective in a range of applications including auto-ignition, laminar flame propagation, and variety of premixed, non-premixed and partially-premixed laminar opposed jet flames [42]. In all these

2556

M. Yaldizli et al. / International Journal of Heat and Mass Transfer 53 (2010) 2551–2562

Table 1 Important parameters of Sandia’s piloted turbulent methane jet flames D and F.

Rejet Main jet diameter Pilot jet diameter Main jet bulk velocity Main jet peak velocity Main jet temperature Pilot bulk velocity Pilot temperature Pilot mixture fraction Co-flow bulk velocity

Flame D

Flame F

22,400 7:2 mm 18:2 mm 49:6 m=s 63:1 m=s 300 K 11:4 m=s 1880 K 0.27 0:9 m=s

44,800 7:2 mm 18:2 mm 99:2 m=s 126:2 m=s 300 K 22:8 m=s 1880 K 0.27 0:9 m=s

applications, the 12-step results were found to be almost indistinguishable from the GRI results. The detailed mechanism in (ii) is based on GRI-3. 3.4. Jet parameters and boundary conditions Fig. 1 shows a schematic view of the Sandia’s piloted methane jet flame, and the coordinate system used in our simulations. For these simulations, a FD mesh with 160  161  161 grid points was considered for a domain of 16  12  12 jet diameters in the x, y, and z directions, respectively. The approximate number of MC elements per each FD cell is 8 and DE ¼ 2D. This provides 64 MC particles per ensemble domain for averaging. The main jet composition is 25% CH4 and 75% air and the Reynolds number is 22,400 for the Flame D. Flame F has the same parameters, except the jet speed or Reynolds number, which is doubled. Detailed specifications of the flames and measurement methods are available elsewhere [43,44], and are not discussed here. However, for convenience the main parameters of flames D and F are listed in Table 1. Non-reflecting boundary conditions [45] are considered for the outlet boundary and zero-gradient conditions are used for the lat-

eral boundaries. For the inflow, non-reflecting, turbulent boundary conditions based on experiment are used. It should be noted that the temperature is not measured at the nozzle exit, and x=d ¼ 1 is the closest distance to the nozzle where such measurements are conducted. The sensitivity of the model predictions to the uncertainty in the pilot boundary conditions is important in LES [46], specially with regards to flame F, which is at a state very close to global extinction. Here, we use values close to those suggested by Barlow [43]. 4. Results and discussion Isosurfaces and isocontours of the instantaneous vorticity magnitude, and temperature for flames D and F simulations are shown in Fig. 2. Evidently, there is a transition to turbulence at x=D  5 in both flames. Nevertheless, the flow field in flame F appears to be more turbulent due to increased jet speed, and lesser overall heat release effects particularly when the 12-step mechanism is employed. The results (not shown) for flame D indicate that the flow field is not very different when flamelet and 1-step models are used but the turbulence is stronger when 12-step chemistry is employed. This is due to damping effect of the reaction on turbulence, which is more significant when 1-step and flamelet models are employed. The temperature isosurfaces in flame D, as obtained by the 1-step and flamelet models (not shown) seem to be continuous without a noticeable sign of extinction. This is expected from 1-step reaction and flamelet models which do not allow significant or any local extinction. However, the temperature isosurface as obtained by the 12-step mechanism in flame D (Fig. 2c) exhibits discontinuities at x=D  5 and x=D  15; suggesting limited flame extinction at some regions of the flow. For flame F, the local extinction is much more significant and clearly visible in Fig. 2d, where the temperature isosurfaces are shown to be broken after x=D  5 when the 12-step model is used. This is consistent with Fig. 2b that shows a ‘‘stronger” turbulence in flame F and a lesser

Fig. 2. Three-dimensional and two-dimensional isosurfaces and isocontours of vorticity magnitude X (a and b), and filtered temperature hTiL (c and d), obtained by LES/FMDF with the 12-step chemistry model; (a), (c) flame D and (b), (d) flame F.

2557

M. Yaldizli et al. / International Journal of Heat and Mass Transfer 53 (2010) 2551–2562

effect of the combustion on the turbulence due to enhanced local flame extinction.

40

⎮ ⎮

0

1

2

80 40

20 0

(b)

0

3

0

1

20

(c) RMS(〈u〉L )

RMS(〈u〉L )

12 9 6 3 0

0

1

2

r/D

2

3

r/D

r/D

3

(d)

15 10 5 0

0

1

2

3

r/D

Fig. 3. The time-averaged (mean) values and the RMS values of the filtered axial velocity field at x=D ¼ 15; (a) mean velocity for flame D, (b) mean velocity for flame F, (c) velocity RMS for flame D and (d) velocity RMS for flame F. Symbols, experimental data; dashed-dot line, 1-step model; dashed line, flamelet model; solid line, 12-step model.

〈T〉L

〈T〉L

〈u〉L

〈u〉L

⎮ ⎮

120

(a)

(b)

1500 ⎮ ⎮

1000

1000

500

500 0

1

2

3

0

1

r/D 600

RMS(〈T〉L )

(c)

400 200 0

0

1

2

r/D

2

3

r/D

600

RMS(〈T〉L )

Fig. 3 shows the radial variations of the mean axial velocity and the RMS of axial velocity for flames D and F at x=D ¼ 15. In this figure, the LES results obtained with the flamelet, 1-step and 12-step chemistry models are compared to the Favre-averaged experimental data. The measurement error is estimated to be below 5% for the mean velocity and about 10% for the RMS [44]. Overall, the agreement between the calculated and the measured mean values is good for both flames for all chemistry models, suggesting that the chemistry effects on the mean axial velocity is not significant. For both flames, the predicted RMS values are lower and higher than the experimental values at r=D  1 and r=D  0, respectively. The agreement is better at r=D > 1. For r=D < 1, the 12-step model predicts lower RMS values in comparison to other reaction models. This is understandable since in the jet core region the effect of combustion on the turbulence is more significant. In this region, the turbulence damping is less when 12-step model is used. The 1-step and flamelet models predict higher heat release and mean temperatures (Fig. 4a). A similar trend is observed in flame F (Fig. 3d). For this flame, the predicted values by the flamelet model are the lowest due to the highest heat release effects. The LES/ FMDF results in Fig. 3 are consistent with those shown in Figs. 2 and 4. The mean values of the temperature for flames D and F are shown in Fig. 4a and 4b, respectively. The measured temperatures in flame F are much lower than those in flame D; at x=D ¼ 15, the peak temperature is about 1100 K in flame F compared to 1750 K in flame D. Nevertheless, the radial variation of the mean temperature is observed to be reasonably well predicted by the LES/FMDF when the 12-step model is employed. The predicted mean temperatures are higher than the corresponding experimental values for both flames when 1-step or flamelet models are used. The 1-step model is based on an irreversible and fast reaction that generally overpredicts the flame temperature, particularly in the rich side of the flame. However, the mean temperatures as calculated by the flamelet model are higher than those obtained by the 1-step in flame F, since the flamelet model does not allow extinction.

2000

(a)

1500 ⎮ ⎮

4.1. Turbulence and flame statistics

60

2000

3

(d)

400 200 0

0

1

2

3

r/D

Fig. 4. The time-averaged (mean) temperature and RMS of temperature at x=D ¼ 15; (a) mean temperature for flame D, (b) mean temperature for flame F, (c) RMS of temperature for flame D and (d) RMS of temperature for flame F. Symbols, experimental data; dashed-dot line, 1-step model; dashed line, flamelet model; solid line, 12-step model.

These results show the importance of chemical kinetics model particularly when flame F is simulated. It is possible to improve the 1-step model predictions by adjusting the reaction parameters. However, we have decided to use the same original values for these parameters. Fig. 4c and d shows the predicted and measured values of the RMS of temperature. Overall, the computed results for flames D and F are not very different than the experimental data in the jet core region when the 12-step model is used. The difference is much more significant away from the core at r=D > 1:5. This can be attributed to the differences in the physical structures of the turbulence. The results shown below indicate a better comparison between the experiment and LES when flow statistics are plotted in the mixture fraction domain. The RMS of temperature is affected by the combined effects of turbulence and heat conduction on one hand, and the chemical heat release on the other. The mean values of the temperature, mixture fraction and axial velocity along the jet centerline (not shown) are in good agreement with the experimental data for both flames D and F when the 12step model is used. The 1-step and flamelet model predictions are less accurate, as expected. Also, the temperature and mixture fraction profiles as obtained from the finite difference part of the LES/ FMDF model are close with those calculated from the Monte-Carlo part, indicating the accuracy of the numerical methods. The radial profiles of CO2 mass fraction for flames D and F as obtained by the 12-step model (Fig. 5a) again show good agreement between the experiments and LES/FMDF in the jet core region. The agreement is less away from the centerline. Again, this can be primarily attributed to the differences in the turbulence structure. The radial profiles of the intermediate species H2 , CO and OH in Fig. 5b–d also indicate a good agreement. The results for other species (e.g. CH4 and O2 ) are similar to those shown in Fig. 5. The LES/ FMDF predictions with 1-step and flamelet models are much less accurate, specially for flame F. Radial profiles of the mean heat release in the energy equation as calculated by the flamelet and 12-step models (not shown) indicate that the LES/FMDF results for these two models are relatively close to each other, which is consistent with the temperature plots. However, for flame F, the predicted heat release with flamelet model is substantially higher than that with the 12-step model. This is also consistent with the mean and isocontours of the filtered temperature in Fig. 4b and c.

2

0

0

0.5

(d)

0

1

r/D

1

Cond.MeanYO 2

0

2

3

4.2. Compositional flame structure The conditional mean and RMS values of various flame variables, conditioned on the mixture fraction, in Figs. 6 and 7 show the structure of the piloted jet flames D and F [43]. The numerical results are obtained by the LES/FMDF with the 12-step chemistry model and C / ¼ 10. As commonly done in experimental and computational studies [4,30], the conditional averages are computed from the data in a plane perpendicular to the jet axis at several different times. Fig. 6a and c shows that for flame D the peak values of conditional mean temperature are slightly overpredicted by the LES/FMDF at both locations from the nozzle, which is consistent with the underprediction of O2 mass fraction. The slightly higher values of mean temperature in these figures are consistent with less scatter in the temperature data in Fig. 8c and d. This could be due to SGS stress, SGS scalar flux, or mixing models. The conditional RMS of temperature, as shown in Figs. 6a and c, are also slightly underpredicted by the LES/FMDF.

0.2

Temp.

1000

0.1

O2 rms Temp.

0

0.5

1

0

2500

0.1 O2 rms Temp.

0

0.5

1

CO

0.04

0

0.04 H2

0 0

Mixture Fraction

0.5 Mixture Fraction

x/D = 15

Temp.

0.2

1000

0.1 O2 rms Temp.

0

0

0.5 Mixture Fraction

1

0

Cond. Mean YH2 ×10

0.08

Cond. Mean YO2

Cond. Mean Temp.

x/D = 15 2000 (c)

0.08

(d) CO

0.04

0.04 H2

0 0

0.5

0

x/D = 15 0.08

(d)

0.04 CO

0

H2

0

0.5

1

0

MixtureFraction

2500

2000 1500 1000 500

(b)

2000 1500 1000 500

0

0.5

1

0

Mixture Fraction 2500

0.5

1

Mixture Fraction 2500

(c)

2000 1500 1000

(d)

2000 1500 1000 500

0.5

1

0

0.5

1

Mixture Fraction

Fig. 8. Scatter plots of temperature vs. mixture fraction at r=D ¼ 0:83 for flame D; (a) measurements at x=D ¼ 7:5 (975 sample points); (b) measurements at x=D ¼ 15 (780 sample points); (c) calculations with the LES/FMDF and 12-step mechanism at x=D ¼ 7:5 (973 sample points); (d) calculations with the LES/FMDF and 12-step mechanism at x=D ¼ 15 (762 sample points).

0

1

1

Cond. Mean YCO

0

0.08

(b)

1

0.04

(a)

Mixture Fraction Cond. Mean YCO

1000

Cond. Mean YH2 ×10

0.2

Cond. Mean YO2

Cond. Mean Temp.

0.08

Temp.

0.5

Fig. 7. Conditional mean temperature and O2 ; CO; H2 mass fractions for flame F as predicted by the LES/FMDF with 12-step kinetics mechanism and with C / ¼ 10; (a) and (b) show the results, at x=D ¼ 7:5; (c) and (d) show the results, at x=D ¼ 15 (lines, calculations; symbols, measurements).

0 2000 (a)

0

Mixture Fraction

500

x/D = 7.5

H2

0

0.08

2000 (c)

0

0.04 CO

Mixture Fraction

x/D = 15

Fig. 5. Profiles of various species mass fractions as predicted by the LES/FMDF with the 12-step model at x=D ¼ 15; (a) CO2 (b) H2 , (c) CO, (d) OH. Symbols, measurement; dashed line, flame D; solid line, flame F.

x/D = 7.5

0.04

Mixture Fraction

0.01

3

rms Temp.

0

0.08

Cond. Mean YCO

3

0.1

x/D = 7.5

(b)

Cond. Mean YCO

0.02

r/D

r/D

2

O2

Cond. Mean YH2 ×10

⎮ ⎮ ⎮ ⎮

1

1

0.02

(c)

0

0

1000

Temperature

〈YCO〉L

0.04

0

0

3

0.2

Temp.

0.08

Temperature

r/D

2

2000 (a)

Cond. Mean YO2

1

0.01

Cond. Mean Temp.

0

0.02

Cond. Mean Temp.

⎮ ⎮ ⎮ ⎮

0.05

(b)

Temperature

0.1

0

⎮ ⎮ ⎮ ⎮

〈YH2〉L ×10

⎮ ⎮ ⎮ ⎮

x/D = 7.5

0.03

(a)

〈YOH〉L ×10

〈YCO2〉L

0.15

Cond. Mean YH2 ×10

M. Yaldizli et al. / International Journal of Heat and Mass Transfer 53 (2010) 2551–2562

Temperature

2558

0

Mixture Fraction

Fig. 6. Conditional mean temperature and O2 ; CO; H2 mass fractions for flame D as predicted by the LES/FMDF with 12-step kinetics mechanism and with C / ¼ 10; (a) and (b) show the results, at x=D ¼ 7:5; (c) and (d) show the results, at x=D ¼ 15 (lines, calculations; symbols, measurements).

Fig. 6b and d shows the conditional mean mass fractions of CO and H2 species for flame D. Evidently, the H2 concentration is very well predicted at both lean and rich sides of the flame away from the stoichiometric location, where it is over predicted. The predicted CO mass fractions also agree well with the experimental data on the lean side of the flame, but they tend to be higher than the measured data on the rich side. For flame D, the conditional averages of the temperature, O2 ; CO and H2 mass fractions exhibit trends similar to those reported in the literature [5,47]. The conditional mean and RMS of the temperature and the O2 mass fraction as obtained by the LES/FMDF with the 12-step kinetics model for flame F are shown in Figs. 7a and c. At x=D ¼ 7:5, the LES/FMDF results are in good agreement with the experimental data. However, the average temperature is slightly overpredicted on the rich side of the flame at x=D ¼ 15, which is again consistent

M. Yaldizli et al. / International Journal of Heat and Mass Transfer 53 (2010) 2551–2562

1000

500

500

0

0.5

1

0

0.5

Mixture Fraction 2500

(c) Temperature

Temperature

2500

1

Mixture Fraction

2000 1500 1000

(d)

2000 1500 1000

500

500

0

0.5

1

0

Mixture Fraction

0.5

1

Mixture Fraction

Fig. 9. Scatter plots of temperature vs. mixture fraction at r=D ¼ 0:83 for flame F; (a) measurements at x=D ¼ 7:5 (1550 sample points); (b) measurements at x=D ¼ 15 (976 sample points); (c) calculations with the LES/FMDF and 12-step mechanism at x=D ¼ 7:5 (1347 sample points); (d) calculations with the LES/FMDF and 12-step mechanism at x=D ¼ 15 (1097 sample points).

with the underestimation of the O2 mass fraction. Interestingly, the RMS values of the temperature are well predicted at x=D ¼ 15. The increase in the magnitude of the temperature RMS in flame F is due to stronger strain field in this flame in comparison to that for flame D. However, the higher level of local extinction in flame F leads to lower mean temperature in this flame in comparison to that in flame D. The conditional CO and H2 mass fraction for flame F are shown in Fig. 7b and d. Again, the experimental results are well predicted by the LES/FMDF with the 12-step model for all mixture fraction values at x=D ¼ 7:5, while there is a slight discrepancy between the experimental and computational data for Z P 0:45 at x=D ¼ 15. The error in the conditional mean CO mass fraction and temperature on the rich side of the flame are consistent with those reported in the literature [4,30].

2000

0.2

0

0.1⎮

O2

⎮

rms Temp. 0

0.5

1

⎮ ⎮ 1000

0

0

rms Temp. 0

0.04

0.5

0.05 ⎮

⎮ ⎮

H2 0.5

⎮ ⎮

0

〈Z〉L

1

0

1

0

⎮ ⎮ ⎮

0.1

(d)

〈CO〉L

0.1

CO

0.02

0

⎮

〈Z〉 L

〈YH2〉L×10

〈CO〉L

⎮ ⎮ ⎮

(c)

⎮

0.1 ⎮

Temp.

〈Z〉 L

0.04

0.2

O2

〈T〉L

⎮

1000

(b)

Temp.

〈YO2〉L

⎮ ⎮

(a)

⎮ ⎮

〈T〉L

2000

〈YO2〉L

1000

1500

〈YH2〉L×10

1500

2000

⎮ ⎮

2000

The scatter plots of temperature in the mixture fraction space as obtained by the LES/FMDF with the 12-step mechanism for both flames D and F are compared with the experimental data in Figs. 8 and 9. The steady-state results obtained from low-strain laminar opposed jet flame simulations are also shown for comparison. For flame D, as Fig. 8a and b indicates, there is a limited scatter in the experimental data suggesting that the local extinction is insignificant [2]. For flame F, there are considerable local extinction and significant scatter in the data at x=D ¼ 7:5, and 15 (Fig. 9a and b). Also at both locations, the LES/FMDF results for flame F are considerably lower than the laminar results, while they are close to the laminar results for flame D. Overall, the computed temperatures show a reasonably good agreement with the experimental data at different locations. For flame D, there are some finite-rate effects in the experimental data that is not fully captured by the LES/FMDF. This could be due to SGS turbulence and mixing models. At x=D ¼ 7:5, a limited local extinction is predicted by the LES/FMDF for flame F when the 1-step model is employed, which is consistent with the isosurface contours of the temperature. Fig. 10 shows the predicted (Favre time-averaged) mean and RMS values of the temperature, and species mass fractions in flames D and F as a function of the mean mixture fraction in comparison to the corresponding experimental data. These values are obtained by cross-referencing the radial mean and rms values of these variables with the radial mixture fraction. The computed data are obtained by the 12-step chemistry model and the mixing model constant of C / ¼ 10. For flame D, both the peak and the shape of scalar and temperature profiles are well predicted by the LES/FMDF (Fig. 10a). At axial location of x=D ¼ 15, the RMS of temperature appears to be slightly underpredicted for hZiL 6 0:125. Considering the highly sensitive and oscillatory behavior of flame F, the mean and RMS of temperature are well predicted by the LES/FMDF (Fig. 10b). The RMS of temperature for flame F exhibit similar trend to that for flame D. It is also shown in Fig. 10c and d that the mean mass fractions of CO and H2 are well predicted for both flames at x=D ¼ 15. However, the CO profiles for flame F at x=D ¼ 7:5 are slightly different on the rich side (not shown), mainly due to finite rate effects. This is consistent with the underpredictions of the CO mass fraction in the rich side of the flame that has been reported in the literature [47]. Additional error might be resulted from the constant mixing coefficient in the LMSE mixing model. However,

(b)

0.05 ⎮

CO

0.02

⎮ ⎮

H2 0

0

0.5

1

0

⎮ ⎮

2500

(a) Temperature

Temperature

2500

2559

〈Z〉L

Fig. 10. Comparison between LES/FMDF and experimental data in the mixture fraction space for the 12-step model with C / ¼ 10 at x=D ¼ 15; (a) mean temperature and O2 mass fraction for flame D; (b) mean temperature and O2 mass fraction for flame F; (c) mean CO and H2 mass fractions for flame D; (d) mean CO and H2 mass fractions for flame F (lines, calculations; symbols, measurements).

Temp.

1500 1000 500 0

Cond. Mean YH ×10 2

0.08

0.06 CO

0.04

0.03 H2

0

0

1000 500 0

0.5 1 Mixture Fraction

(b)

0.5 1 Mixture Fraction

0

Temp.

1500

rms Temp. 0

0.08

Cond. Mean YH ×10 2

0

rms Temp.

(c)

2000

0.5 1 Mixture Fraction

(d)

CO

0.04

0.06

0.03 H2

0

0

0.5 1 Mixture Fraction

0

Cond. Mean YCO

Cond. Mean Temp.

(a)

2000

Cond. Mean YCO

Cond. Mean Temp.

Fig. 11. Conditional mean and RMS of various flame variables as obtained by the LES/FMDF with the 12-step model and different C / values for flame D; (a) mean and RMS of temperature at x=D ¼ 7:5; (b) mean and RMS of temperature at x=D ¼ 15; (c) mean of CO and H2 mass fractions at x=D ¼ 7:5; (d) mean of CO and H2 mass fractions at x=D ¼ 15 (symbols, measurements; long-dashed line, C / ¼ 6; dasheddot line, C / ¼ 8; solid line, C / ¼ 10; dashed-dot-dot line, C / ¼ DM1; dashed line, C / ¼ DM2).

1000 rms Temp.

500 0

0

0.08

0.5 1 Mixture Fraction 0.06

(b) CO

0.04

0.03

H2 0

0

0.5 1 Mixture Fraction

0

(c)

2000

Temp.

1500 1000 rms Temp.

500 0

0

0.08

0.5 1 Mixture Fraction 0.06

(d) CO

0.03

0.04 H2 0

0

0.5 1 Mixture Fraction

0

Cond. Mean YCO

Temp.

1500

Cond. Mean YH ×10 2

In this section we focus on the sensitivity of the LES/FMDF calculations to the mixing model coefficient C / which is either empirically prescribed or dynamically evaluated in our simulations. All of the following simulations are conducted with the 12-step reaction model. The radial profiles of the mean temperature for flame F as obtained by the LES/FMDF with various C / (not shown) indicates that by increasing C / the mixing/combustion enhances and the average temperature increases, albeit not linearly, which is expected and is consistent with the previous studies [4,30]. However, the influences of C / on the temperature do not seem to be very significant in flame D. This is because of the less sensitivity of flame D, which is close to equilibrium, to the molecular mass and thermal diffusivity coefficients. The mean temperatures predicted by the dynamic models DM1 and DM2 are also close to those obtained via constant C / for flame D. However, comparison of temperatures obtained with the DM1 and DM2 models with those computed with constant C / ¼ 6; 8; 10 for flame F indicates that the flame structure is affected more by the mixing model in this flame. The higher mean temperatures predicted with the DM2 for flame F are due to higher level of mixing. With the DM1, the mixing seems to be relatively lower than that via constant coefficient of C / ¼ 10. The dynamic model DM2 predicts the highest RMS temperature. The results for flame F indicate that the mean heat release is significantly lower for lower C / . Combined with the stronger turbulence in flame F, these lower heat release values yield substantially lower mean temperatures. Also, consistent with the mean temperature profiles shown in the preceding figures, substantially higher mean heat release values are generated when the DM2 model is employed, particularly for flame F. The effects of mixing model on the compositional structure of flames D and F are shown in Figs. 11 and 12, where the variations of various flame variables in the mixture fraction domain for different mixing model coefficients are considered. In Fig. 11a and b, the

Cond. Mean YH ×10 2

4.3. Subgrid-scale mixing

(a)

2000

Cond. Mean YCO

similar trends have been observed in the calculated mean profiles of Y CO by Pitsch et al. [31,47], who used a dynamic method to compute the mixing model coefficient. This issue is discussed further in the following section.

Cond. Mean Temp.

M. Yaldizli et al. / International Journal of Heat and Mass Transfer 53 (2010) 2551–2562

Cond. Mean Temp.

2560

Fig. 12. Conditional mean and RMS of various flame variables as obtained by the LES/FMDF with the 12-step model and different C / values for flame F; (a) mean and RMS of temperature at x=D ¼ 7:5; (b) mean and RMS of temperature at x=D ¼ 15; (c) mean of CO and H2 mass fractions at x=D ¼ 7:5; (d) mean of CO and H2 mass fractions at x=D ¼ 15 (symbols, measurements; long-dashed line, C / ¼ 6; dasheddot line, C / ¼ 8; solid line, C / ¼ 10; dashed-dot-dot line, C / ¼ DM1; dashed line, C / ¼ DM2).

conditional mean and RMS of temperature for flame D as obtained with different constant C / and dynamic models at axial locations of x=D ¼ 7:5 and x=D ¼ 15 are compared with the experimental data. Overall, the mean temperature is well predicted by the LES/FMDF, despite some differences at x=D ¼ 7:5. The conditional RMS of temperature is shown to be also well predicted with the LES/FMDF when C / ¼ 6; 8 or DM1 is used, while it is underestimated with the DM2 model at all mixture fraction values. This could be due to subgrid scalar variance model used in the DM2 model rather than the local equilibrium assumption in the dynamic model as Raman and Pitch [31] suggested, since the dynamic model DM1 does not show such an underestimation. Consistent with our earlier results and with those obtained by others [4], the RMS of temperature decreases and the mean temperature increases as C / increases. The results at x=D ¼ 15 are somewhat similar to those at x=D ¼ 7:5, but more visibly the conditional mean temperature are overestimated by the LES/FMDF in the stoichiometric and the fuel-rich regions. Also, the deviations between the experimental and LES/FMDF values of the temperature RMS are higher at x=D ¼ 15, mainly due to stronger turbulence at this location. Fig. 11c and d shows the conditional averages of CO and H2 mass fractions for flame D at two different downstream locations. At x=D ¼ 7:5, the numerical values are generally in agreement with the reported experimental data at all mixture fraction values for all C / values. There seems to be a slight overprediction of the experimental data around the stoichiometric point when C / ¼ 8; 10 or DM2 is used, which is mainly due to differences in heat release. This effect is more noticeable at further downstream locations (Fig. 11d). The mean profiles of H2 mass fraction are also overpredicted by the LES/FMDF in the stoichiometric region but are well predicted at other mixture fractions. Fig. 12a and b shows the computed and measured conditional mean and RMS of temperature for flame F at two downstream locations. At x=D ¼ 7:5, the conditional mean temperature and the extent of local extinction present in flame F are well predicted for C / ¼ 6; 8 and DM1 model. However, the DM2 model overpredicts the peak temperature, even though the predictions on the lean and rich sides of the flame are still in fairly good agreement with the experiment. Also at this location, the RMS fluctuations of the conditional mean temperature are well predicted for all C / s. At downstream location of x=D ¼ 15, where the effects of turbulence is more significant the conditional mean temperature is

2561

M. Yaldizli et al. / International Journal of Heat and Mass Transfer 53 (2010) 2551–2562

2500

(a) Temperature

Temperature

2500 2000 1500 1000 500

(b)

2000 1500 1000 500

0

0.5

1

2500

Temperature

Temperature

0.5

2000 1500 1000 500

(d)

2000 1500 1000 500

0

0.5

Mixture Fraction

1

Temperature

Temperature

2000 1500 1000 500

(b)

2000 1500 1000 500

0

0.5

1

0

Mixture Fraction 2500

0.5

1

Mixture Fraction 2500

Temperature

(c)

2000 1500 1000

(d)

2000 1500 1000 500

500 0

0.5

Mixture Fraction

1

0

0.5

1

Mixture Fraction

Fig. 14. Scatter plots of temperature vs. mixture fraction as obtained by the LES/ FMDF with the 12-step model and different C / values for flame F x=D ¼ 7:5, and r=D ¼ 0:83; (a) C / ¼ 6 (1418 sample points); (b) C / ¼ 8 (1350 sample points); (c) C / ¼ DM1 (1376 sample points); (d) C / ¼ DM2 (1134 sample points).

ment, by collecting the data from the cells that are at radial distance of x=d ¼ 0:83. Only the results at x=D ¼ 7:5 are shown. The trends are similar at different axial locations. The solid line in Figs. 13 or 14 shows the temperature for a low-strain ð10=sÞ laminar opposed jet flame. Figs. 13 shows the scatter plots of temperature for flame D for various C / s. Evidently, the lowest scatter in the data is for the DM2 model, which is consistent with the overestimation of the peak conditional mean temperature in Fig. 11a by this model. Simulations with the dynamic model DM1 and with constant coefficients of C / ¼ 6; 8; 10 exhibit more significant finite rate chemistry effects and better prediction of the conditional mean temperature. The results at x=D ¼ 15 (not shown) exhibit similar trends. For flame F, it can be safely stated that the challenging task of capturing the increasing level of local extinction is fairly accomplished by the LES/FMDF when the 12-step model is employed. As shown in Fig. 14, temperatures predicted with different C / coefficients are not very different. Simulations with the dynamic models DM1 and DM2 show less localized extinction and consequently higher average temperatures. Nevertheless, they seem to be able to capture the trend in the experimental data (Fig. 9).

1

Mixture Fraction

(c)

2500

(a)

5. Summary and conclusions 0

Mixture Fraction 2500

2500

Temperature

underpredicted at stoichiometric and lean side of the flame by the LES/FMDF when mixing coefficient is C / ¼ 6; 8. The simulations conducted with the dynamic models overpredict the measured mean temperatures for the mixture fraction values of Z P 0:25. The poorest agreement is for the DM2 model, mainly due to higher heat release. At x=D ¼ 15, the computed conditional RMS of temperatures are also in good agreement with the experimental data for C / ¼ 6; 8 and DM1. With the DM2 model, there is an underprediction of data at lean side of the flame, possibly due to sensitivity of the flame F to the oscillations in the dynamic model coefficient. Similar to what have been shown for flame D, the CO and H2 mass fractions profiles agree well with the experimental data for flame F at x=D ¼ 7:5 as observed in Fig. 12c. Only the peak values of these species are overestimated when DM2 is used, which is consistent with the expository temperature profiles in Fig. 12a. Similar results have been observed at lower jet speed for the (piloted jet) flame E by Raman and Pistch [31]. Close examination of the mass fraction profiles at x=D ¼ 15 (Fig. 12d) indicates that the LES/FMDF results with constant mixing coefficients of C / ¼ 6; 8 are better than those with DM1 and DM2 models. Both dynamic models tend to underpredict the rate of change of the mean CO mass fraction profile with the Reynolds number, which is consistent with the results obtained by others for jet flames [4]. Nevertheless, the strong influences of the turbulence on the flame and the level of local flame extinction as seen in the conditional variables are captured by the LES/FMDF. For flame D, the computed mean temperatures with different mixing model coefficients (not shown) are found to be close to each other, and generally in agreement with the experimental data. Similarly, the predicted values of mean CO and H2 mass fractions agree well with the experimental data. For flame F, the measured mean temperature profiles and the level of local flame extinction are well predicted by the simulations when C / ¼ 10. There are some minor differences between the LES/ FMDF and the experiment when the DM2 model is employed. This is consistent with the overprediction of CO mass fraction by this model. The scatter plots of temperature in the mixture fraction domain for various mixing models and for both flames D and F are shown in Figs. 13 and 14. These results are obtained, similar to experi-

0

0.5

1

Mixture Fraction

Fig. 13. Scatter plots of temperature vs. mixture fraction as obtained by the LES/ FMDF with the 12-step model and different C / values for flame D at x=D ¼ 7:5, and r=D ¼ 0:83; (a) C / ¼ 6 (1102 sample points); (b) C / ¼ 8 (1022 sample points); (c) C / ¼ DM1 (1042 sample points); (d) C / ¼ DM2 (973 sample points).

This paper is concerned with the application and assessment of the filtered mass density function (FMDF) for large-eddy simulation (LES) of Sandia’s ‘‘high speed” piloted turbulent methane jet flames [2]. The LES/FMDF calculations are conducted with both finite-rate and flamelet-type reaction models. The flamelet model employs detailed (GRI-3) kinetics mechanism. The finite-rate model employs 1-step and 12-step mechanisms. Various subgrid-scale mixing models are also used and tested. A scalable algorithm for parallelization of the hybrid (Eulerian–Lagrangian) LES/FMDF methodology is presented and the consistency of its Monte Carlo – finite difference solutions are discussed. The parallel algorithm is shown to yield superlinear speedups with respect to single node calculations. The Favre-averaged temperatures for the two jet speeds (referred to as flames D and F) are shown to be in fairly good agreement with the experimental data. For the lower speed

2562

M. Yaldizli et al. / International Journal of Heat and Mass Transfer 53 (2010) 2551–2562

jet (flame D), the numerical results are consistent with the experiment exhibiting near-equilibrium flame structure with limited local extinction. The higher degrees of local extinction observed in the higher speed jet (flame F) is successfully captured and reproduced with the LES/FMDF when the 12-step chemistry model is employed. As expected, the flamelet and 12-step models generate results close to each other and comparable to experiment for flame D. For flame F, the results are more sensitive to the reaction and mixing models. For this flame, the higher heat release predicted by the 1-step and flamelet models suppress the effects of turbulence on the flame. The sensitivity of the calculations to the subgrid mixing model is also investigated by performing simulations with different model coefficients and with two dynamic models. Even though, the results with different mixing coefficients and models are not that much different for flame D, the calculations conducted with constant coefficient of C / ¼ 10 yield the most accurate overall agreements for both flames. Further improvements in the predictions might be possible with better mixing and subgrid stress/scalar flux closures. Acknowledgments This work was sponsored by the National Science Foundation under grant CTS0092665. Additional support was provided by the NASA under grant NNX07AC50A. References [1] F. Jaberi, P. Colucci, S. James, P. Givi, S. Pope, Filtered mass density function for large-eddy simulation of turbulent reacting flows, J. Fluid Mech. 401 (1999) 85–121. [2] R. Barlow, J. Frank, Effects of turbulence on species mass fractions in methane/ air jet flames, in: 27th International Symposium on Combustion, 1998, pp. 1087–1095. [3] A. Masri, R. Bilger, R. Dibble, Turbulent nonpremixed flames of methane near extinction – mean structure from Raman measurements, Combust. Flame 71 (1988) 245–266. [4] J. Xu, S. Pope, PDF calculations of turbulent nonpremixed flames with local extinction, Combust. Flame 123 (2000) 281–307. [5] H. Pitsch, H. Steiner, Large-Eddy simulation of a turbulent piloted methane/air diffusion flame (Sandia Flame D), Phys. Fluids 12 (10) (2000) 2541–2554. [6] M. Sheikhi, T. Drozda, P. Givi, F. Jaberi, S. Pope, Large eddy simulation of a turbulent nonpremixed piloted methane jet flame (Sandia Flame D), Proc. Combust. Inst. 30 (2005) 549–556. [7] W. Jones, M. Kakhi, PDF modeling of finite-rate chemistry effects in turbulent nonpremixed jet flames, Combust. Flame 115 (1998) 210–229. [8] M. Muradoglu, S. Pope, D. Caughey, The hybrid method for the PDF equations of turbulent reactive flows: consistency conditions and correction algorithms, J. Comput. Phys. 172 (2001) 841–878. [9] T. Landenfeld, A. Sadiki, J. Janicka, A turbulence–chemistry interaction model based on a multivariate presumed beta-PDF method for turbulent flames, Flow, Turbulence Combust. 68 (2002) 111–135. [10] S. Pope, Computationally efficient implementation of combustion chemistry using in situ adaptive tabulation, Combust. Theory Model. 1 (1997) 41–63. [11] S. Subramaniam, S. Pope, A mixing model for turbulent reactive flows based on Euclidean minimum spanning trees, Combust. Flame 115 (1998) 487–514. [12] M. Yaga, H. Endo, T. Yamamoto, H. Aoki, T. Miura, Modeling of Eddy characteristic time in LES for calculating turbulent diffusion flame, Int. J. Heat Mass Transfer 45 (11) (2002) 2343–2349. [14] L.Y. Hu, L.X. Zhou, J. Zhang, Large-eddy simulation of the sydney swirling nonpremixed flame and validation of several subgrid-scale models, Numer. Heat Transfer 53 (1) (2008) 39–58. [15] H.W. Wu, S.W. Perng, LES analysis of turbulent flow and heat transfer in motored engines with various SGS models, Int. J. Heat Mass Transfer 45 (11) (2002) 2315–2328. [16] S.B. Pope, Computations of turbulent combustion: progress and challenges, in: Proceedings of 23rd Symp. (Int.) on Combustion, The Combustion Institute, Pittsburgh, PA, 1990, pp. 591–612. [17] T. Poinsot, D. Veynante (Eds.), Theoretical and Numerical Combustion, R.T. Edwards Inc., Philadelphia, PA, 2001.

[18] F. Gao, E.E. O’Brien, A large Eddy scheme for turbulent reacting flows, Phys. Fluids A 5 (6) (1993) 1282–1284. [19] P.J. Colucci, F.A. Jaberi, P. Givi, S.B. Pope, Filtered density function for large Eddy simulation of turbulent reacting flows, Phys. Fluids 10 (2) (1998) 499– 515. [20] S.C. Garrick, F.A. Jaberi, P. Givi, Large Eddy simulation of scalar transport in a turbulent jet flow, in: D. Knight, L. Sakell (Eds.), Recent Advances in DNS and LES, Fluid Mechanics and its Applications, vol. 54, Kluwer Academic Publishers., The Netherlands, 1999, pp. 155–166. [21] F.A. Jaberi, P.J. Colucci, Large Eddy simulation of heat and mass transport in turbulent flows. Part 2: scalar field, Int. J. Heat Mass Transfer 46 (10) (2003) 1827–1840. [22] F.A. Jaberi, P.J. Colucci, S. James, P. Givi, S.B. Pope, Filtered mass density function for large Eddy simulation of turbulent reacting flows, J. Fluid Mech. 401 (1999) 85–121. [23] F.A. Jaberi, Large Eddy Simulation of Turbulent Premixed Flames via Filtered Mass Density Function, AIAA Paper 99-0199, AIAA, 1999. [24] S. James, F. Jaberi, Large scale simulations of two dimensional non-premixed methane jet flames, Combust. Flame 123 (4) (2000) 465–487. [25] L.Y.M. Gicquel, P. Givi, F.A. Jaberi, S.B. Pope, Velocity filtered density function for large Eddy simulation of turbulent flows, Phys. Fluids 14 (2002) 1196– 1213. [26] M.R.H. Sheikhi, T.G. Drozda, P. Givi, S.B. Pope, Velocity–scalar filtered density function for large Eddy simulation of turbulent flows, Phys. Fluids 15 (8) (2003) 2321–2337. [27] C. Tong, Measurements of conserved scalar filtered density function in a turbulent jet, Phys. Fluids 13 (10) (2001) 2923–2937. [28] M. Germano, U. Piomelli, P. Moin, W.H. Cabot, A dynamic subgrid-scale eddy viscosity model, Phys. Fluids A 3 (7) (1991) 1760–1765. [29] R.S. Miller, S.H. Frankel, C.K. Madnia, P. Givi, Johnson–Edgeworth translation for probability modeling of binarymixing in turbulent flows, Combust. Sci. Technol. 91 (1–3) (1993) 21–52. [30] R. Cao, S. Pope, The influence of chemical mechanisms on PDF calculations of nonpremixed piloted jet flames, Combust. Flame 143 (2005) 450–470. [31] V. Raman, H. Pitsch, A consistent LES/filtered-density function formulation for the simulation of turbulent flames with detailed chemistry, Proc. Combust. Inst. 31 (2007) 1711–1719. [32] C.M. Cha, P. Trouillet, A subgrid-scale mixing model for large-eddy simulations of turbulent reacting flows using the filtered density function, Phys. Fluids 15 (6) (2003) 1496–1504. [33] C. Pierce, P. Moin, A dynamic model for subgrid-scale variance and dissipation rate of a conserved scalar, Phys. Fluids 10 (12) (1998) 3041–3044. [34] F.A. Jaberi, P. Givi, Inter-layer diffusion model of scalar mixing in homogeneous turbulence, Combust. Sci. Technol. 104 (4–6) (1995) 249. [35] M.H. Carpenter, A high-order compact numerical algorithm for supersonic flows, in: K.W. Mortron (Ed.), Proc. 12th International Conference on Numerical Methods in Fluid Dynamics, Lecture Notes in Physics, vol. 371, Springer-Verlag, 1990, pp. 254–258. [36] S.B. Pope, Turbulent Flows, Cambridge University Press, Cambridge, UK, 2000. [37] R.J. Kee, F.M. Rupley, J.A. Miller, The Chemkin Thermodynamic Database, Tech. Rep. SAND87-8215B, Sandia National Laboratory, Livermore, CA, 1994. [38] S.B. Pope, Lagrangian PDF methods for turbulent flows, Ann. Rev. Fluid Mech. 26 (1994) 23–63. [39] S. Karlin, H.M. Taylor, A Second Course in Stochastic Processes, Academic Press, New York, NY, 1981. [40] P.E. Kloeden, E. Platen, Numerical Solution of Stochastic Differential Equations, Applications of Mathematics, Stochastic Modelling and Applied Probability, vol. 23, Springer-Verlag, New York, NY, 1995. [41] C.J. Sung, C.K. Law, J.-Y. Chen, An augmented reduced mechanism for methane oxidation with comprehensive global parametric validation, in: Proceedings of 27th Symp. (Int.) on Combustion, The Combustion Institute, Pittsburgh, PA, 1998, pp. 295–303. [42] H. Mohammad, Numerical Simulation of Partial Premixed Methane Combustion in Laminar and Turbulent Flows, Masters Thesis, MSU. [43] R. Barlow, Sandia/TUD Piloted CH4/Air Jet Flames, Available from: . [44] C. Schneider, A. Dreizler, J. Janicka, Flow field measurements of stable and locally extinguishing hydrocarbon-fuelled jet flames, Combust. Flame 135 (2003) 185–190. [45] T. Poinsot, S. Lele, Boundary conditions for direct simulations of compressible viscous flows, J. Comput. Phys. 101 (1992) 104–129. [46] Q. Tang, J. Xu, S. Pope, Probability density function calculations of local extinction and no production in piloted-jet turbulent methane/air flames, in: 28th Symposium (Int.) on Combustion, The Combustion Institute, Pittsburgh, PA, 2000, pp. 133–139. [47] H. Pitsch, Improved pollutant predictions in large-eddy simulations of turbulent non-premixed combustion by considering scalar dissipation rate fluctuations, Proc. Combust. Inst. 29 (2002) 1971–1978.