Modelling nitrogen oxide emissions in turbulent flames with air dilution: Application to LES of a non-premixed jet-flame

Combustion and Flame 161 (2014) 496–509

Contents lists available at ScienceDirect

Combustion and Flame j o u r n a l h o m e p a g e : w w w . e l s e v i e r . c o m / l o c a t e / c o m b u s t fl a m e

Modelling nitrogen oxide emissions in turbulent flames with air dilution: Application to LES of a non-premixed jet-flame François Pecquery a,b, Vincent Moureau a,⇑, Ghislain Lartigue a, Luc Vervisch a, Anthony Roux b a b

CORIA – CNRS, Normandie Université, INSA de Rouen, 76801 Saint-Etienne-du-Rouvray, France Turbomeca (Groupe Safran), 64510 Bordes, France

a r t i c l e

i n f o

Article history: Received 26 December 2012 Received in revised form 9 September 2013 Accepted 18 September 2013 Available online 21 October 2013 Keywords: Nitrogen oxide Chemistry tabulation Large-eddy simulation

a b s t r a c t The aeronautical industry faces the challenge of designing new engines with a better specific fuel consumption while decreasing pollutant emissions such as nitrogen oxide (NO). To this aim, many Computational Fluid Dynamics (CFD) models have been developed to predict NO emissions. Most of these models are based on a simplified description of the nitrogen chemistry, while few of them attempt to take into account the whole NO chemistry and its interactions with carbon oxidation. This article presents a new model, called NOMANI (Nitrogen Oxide emission model with one-dimensional MANIfold), which predicts NO emissions from all chemical processes present in detailed mechanisms and which is well suited to aeronautical burners. To have a good estimation of the NO mass fraction produced in the flame front and in the burnt gases, two progress variables are used. The first one is built from major carbon species, and once carbon chemistry is completed the NO mass fraction is used as a second progress variable with possibility of accounting for dilution by secondary air. This modeling is first assessed thanks to a phenomenological study based on the nitrogen chemistry under dilution by fresh air. The model is then applied and validated in highly resolved Large-Eddy Simulations (LES) of the TNF Sandia Flame D. Ó 2013 The Combustion Institute. Published by Elsevier Inc. All rights reserved.

1. Introduction Large-Eddy Simulations (LES) is a very promising tool in the design chain of aeronautical combustion chambers, because it provides access to unsteady flow features, essential to the description of scalar mixing and turbulent combustion. The major challenge in combustion LES is due to the large spectrum of turbulent and chemical scales, the stiffness of the chemical reactions and the high number of involved species, which imply that the reaction zones cannot be fully resolved and must be modeled. A related issue in the modeling of pollutant emissions arises from the fact that nitrogen oxide (NO) production involves several chemical processes. Most of them have to be taken into account to get accurate predictions of NO emissions in different regions of the flame. In this paper, a strategy is proposed to compute NO in turbulent flames, once a detailed chemical mechanism chosen; the mechanism is a given input of the model that is not the focus of the study. Here the methane–air GRI-3.0 [1] is used as the reference. The NO production processes have been intensively investigated in the past, and they can be grouped into two families that evolve with two different characteristic times. The first one corresponding to the prompt NO combines all processes that have the ⇑ Corresponding author. Fax: +33 2 32959780.

same characteristic time scale as carbon chemistry. The second major part of the NO production occurs once the exothermic reactions are completed, featuring much larger characteristic time scales. For fuels having no nitrogen atoms in their composition, the prompt NO production in the flame front area is known as the Fenimore mechanism [2]. The prompt NO results from the coupling between the nitrogen and carbon chemistry and has been extensively investigated. Recent works [3–5] highlighted the role of the NCN pathway in the prompt NO formation. Many authors [6,7] have emphasized the importance of evaluating the prompt NO for rich mixtures, even if the temperature remains moderate, i.e. T < 1100 K. The burnt gas zone combines all the other chemical processes, which produce 90% or more of NO mass fraction. These mechanisms include the Zeldovich mechanism, and those via N2O and NNH [8–12]. The prediction of NO formation in gas turbines is an active research field due to its environmental impact [13]. In the exhaust of combustion chambers, NO is always far from its equilibrium state, because the residence time of these devices is small compared to the time scale needed by NO to reach its equilibrium value [14]. Prompt NO may therefore be significant within the aeronautical context, moreover dilution of burnt gases by secondary air also impacts NO levels. Numerous RANS (Reynolds Averaged Navier Stokes) and LES models for NO emissions exist. Most of them, which have been

E-mail address: [email protected] (V. Moureau). 0010-2180/$ - see front matter Ó 2013 The Combustion Institute. Published by Elsevier Inc. All rights reserved. http://dx.doi.org/10.1016/j.combustflame.2013.09.018

F. Pecquery et al. / Combustion and Flame 161 (2014) 496–509

validated by comparison with experimental data, consist in the solving of a transport equation for the nitrogen oxide mass fraction. This equation needs to be closed since the nitrogen oxide production rate depends on pressure, enthalpy and all species mass fractions of the mixture. Several strategies have been investigated to close the transport equation of the nitrogen oxide mass fraction. Löffer et al. [15] developed a model that is coupled with a RANS approach. Applying quasi-steady-state (QSS) approximation for intermediate species leads to a linear equation system. It is then possible to close the NO balance equation, but this assumption fails in the flame front; hence, this model is limited to NO produced in the burnt gases. In the flame front, where prompt NO is produced, radicals are not at the equilibrium, and due to the QSS assumption, the model cannot track the prompt NO. As a result, Löffer closure is well suited to predict NO levels for configurations where the prompt NO stays negligible, mostly situations where the equivalence ratio stays below 0.75. A second approach to close the transport equation of the nitrogen oxide mass fraction is to resort to PDF and Conditional Moment Closure (CMC) methods to handle detailed chemistry in the simulation of turbulent flames [16–19]. An alternative to QSS, PDF and CMC models is the tabulated chemistry approach, which takes into account complex chemistry effects at low cost and complexity. The first tabulated model of nitrogen oxide emission in RANS calculations is presented by Bradley et al. in 1998 [20]. It is dedicated to the prediction of prompt NO formation in lean mixtures. The model consists in solving a transport equation for NO mass fraction, the source term being extracted from a laminar look-up table as a function of a progress variable, which is the normalized temperature. Bradley et al. also proposed to replace the normalized temperature by the NO mass fraction, but stressed that this progress variable is poorly coupled to the temperature field. A mathematical analysis of chemistry reduction is proposed in the Intrinsic Low-Dimensional Manifold (ILDM) concept, which was introduced in 1992 by Maas and Pope [21]. It is based on the analysis of the local time scales of the chemical system in a phase space. The evolution of the system is decomposed in characteristic times associated with eigenvalues of the Jacobian matrix of the species reaction rates. Through this analysis, Maas and Pope shed light on the existence of a manifold. The evolution along this manifold corresponds to the slowest processes, and the evolution outside the manifold is with fast relaxation times, bringing the mixture composition onto the manifold. Because of this property, it is possible to consider only the low-dimensional manifold instead of working in the high dimensional phase space. Maas and Pope [21] shown that a 1D or a 2D manifold provides accurate results for a system at high temperature, i.e. for a progress variable that is closed to its equilibrium value. Nafe and Maas [22] studied the nitrogen chemistry with the ILDM concept. They found that nitrogen chemistry adds one characteristic time, so that the dimension of the tabulated manifold has to be increased by one. To overcome this difficulty, they decomposed the kinetics into carbon and nitrogen chemistries. Nitrogen chemistry is considered as frozen during carbon oxidation. Under this assumption, the computation of nitrogen oxide is provided by the transport of an extra progress variable, which is NO itself. The drawback of this method comes from the decoupling of the flame reaction zone from the NO production. As a result, it is mostly dedicated to configurations where the prompt NO is negligible. Since 2007, new approaches based on tabulated flamelets appeared for NO. They are based on canonical combustion problems, as homogeneous reactors or one-dimensional premixed or diffusion flames [23–25]. These approaches are successful to predict the heat release and some species that evolve with the same

497

characteristic time as the progress variable. For species evolving differently from the usual progress variables, however, standard flamelet models may fail. Nitrogen chemistry and the progress variable Yc, usually based on major species Y CO ; Y CO2 and Y H2 O , feature different characteristic time scales. These differences are illustrated in Fig. 1 where strong variations in NO are observed in burnt gases for a very small variation in Yc. Due to this property, nitrogen oxide cannot be tabulated as a single function of Yc. To improve the capability of tabulated chemistry for predicting the slow evolution of NO, Godel et al. [26] and van Oijen et al. [27] suggested to include nitrogen species in the progress variable composition. This approach is very attractive but the filtered density function of this unique progress variable is very different from the classical progress variable based on CO and CO2. This is due to the fact that a unique progress variable has to describe the evolution of the chemical system over very different time scales. When using presumed-conditional moment closures, this implies to modify the modeling of the unmixedness factor of the progress variable [28]. Another solution is proposed by Ihme et al. [29], who also discussed a closure based on tabulated chemistry; the NO formation rate is then approximated by adopting a linear scaling of backward reaction rates tabulated versus temperature, pressure and enthalpy. This model also assumes that all species are in steady state, except NO. Three chemical pathways are then considered: prompt NO, N2O and thermal. Ketelheun et al. [30] compare the Bradley’s formalism with the Ihme’s formalism used in a premixed regime. Both of these models provided comparable results. The work of Zöller et al. [31] shows that the linear approximation is not suitable when the complexity of the chemical kinetics is increased. Due to some species, NO2 for instance, strong non-linear dependencies exist. To ensure an efficient closure of NO reaction rate, Zöller et al. integrate the full mechanisms, tabulating nitrogen chemistry with adiabatic and isobaric Perfectly Stirred Reactors (PSR). The NO mass fraction is used as a progress variable for the nitrogen chemistry. Some PSR are computed for the range of equivalence ratios (or mixture fractions) that are encountered during the CFD computation. Their compositions are constructed from the equilibrium state at the considered mixture fraction and temperature, but all nitrogen species are set to zero except the N2 mass fraction. During the PSR computation only the nitrogen species evolve. According to the authors this approach is accurate for flames where moderate NO concentrations are present. Introducing the tabulation of relaxation times as a function of perturbation amplitudes, Vervisch et al. [32] present a new model for the emission of thermal NO. This relaxation is tabulated as a function of equivalence ratio, pressure, temperature and dilution mass fraction. This method provides good results when NO is close

Fig. 1. Evolution of YNO in a 1D unstretched premixed methane/air flame as a function of the progress variable Y c ¼ Y CO þ Y CO2 þ Y H2 O for various mixture fractions Z.

498

F. Pecquery et al. / Combustion and Flame 161 (2014) 496–509

to its equilibrium state. This model was recently extended by adding alternative chemical routes and considering N2O and NO2 [33]. Biagioli and Güthe [34] showed the importance of the prompt NO pathway and underlined the effect of fuel–air unmixedness on NO emissions in a gas turbine burner. They developed a model suited to track the prompt NO and the NO produced in burnt gases. The model consists in solving two decoupled transport equations: one for the deviation of the prompt NO produced in the flame front from the reference value that would be obtained in a laminar premixed flame in the same operating conditions, and one transport equation for the NO produced in the post flame zone. To solve the NO in the post flame zone, the chemistry is assumed to be away from equilibrium, and produced at a reaction rate independent of the NO mass fraction. Within the large variety of these approaches, the purpose of the present work is to design an accurate model for nitrogen oxide emissions in aeronautical configurations, including air dilution. To this aim, two difficulties have to be overcome. The first is to track the NO mass fraction considering all the chemical pathways. From this non-exhaustive review of NO modeling in turbulent flames, it appears that the two different characteristic time scales of NO production may require different modeling approaches, but relying on mechanisms that are strongly coupled, making them difficult to compute independently. The second difficulty comes from aeronautical configurations. To increase the lifespan of aeronautical gas turbine combustors, the walls of the flame tube must be protected against the high surrounding temperature. One way to achieve this goal is to dilute the burnt gases with fresh air in order to decrease the temperature both at the walls and at the turbine inlet following the combustion chamber. Dilution strongly impacts NO formation due to the high sensibility of NO chemistry to temperature, O2 concentration and radical concentrations. It is therefore mandatory to develop a model able to track NO production in the presence of air dilution. In the second section of this paper, a study of the nitrogen chemistry under dilution is conducted. It consists in the solving of the NO formation in a constant pressure vessel with dilution by fresh air, in order to reproduce as closely as possible the conditions observed in the air-dilution zone of an aeronautical burner. The final objective of this section is to provide valuable information on the driving parameters of thermal NO formation in the burnt gases with dilution. In particular, two aspects have to be investigated to design a model:  Does the initial mixture fraction before dilution have an importance to evaluate the amount of NO in burnt gases?  Does the mixing dilution intensity play an important role on the NO emissions? In the third section, a new closure for nitrogen oxide emission is presented and validated in one-dimensional flames. The section four is dedicated to the extension of this NOMANI model to the LES formalism. Finally, it is applied and verified in a highly resolved LES of the TNF Sandia flame D.

2. Study of nitrogen oxide formation in a constant pressure vessel with air dilution

Fig. 2. Schematic representation of the canonical problem. Z is the mixture fraction. rdil and sdil are the dilution rate and time, respectively.

oxide formation, a constant pressure vessel, where diffusion effects are neglected, is simulated. In this Perfectly Stirred Reactor (PSR), auto-ignition of an initial methane/air mixture occurs until carbon chemistry is finished and dilution is then performed homogeneously by adding fresh air. Finally, when dilution is finished, the PSR continue to evolve to reach its equilibrium. These three steps are illustrated in Fig. 2 and defined as follows: 1. Mixture auto-ignites and evolves until carbon chemistry is complete. The evolution of the PSR is obtained by solving _ i , where x _ i is the chemical source term and V dqVY i =dt ¼ qV x the PSR volume which varies in time to keep the pressure constant. The total enthalpy H is also kept constant. 2. Mixture is diluted with fresh air by solving _ i þ fi where fi is the dilution mass flow rate of dqVY i =dt ¼ qV x species i. Dilution occurs until the mixture fraction of the PSR equals a value Zmin. The dilution is performed in such a manner to keep the dilution rate rdil constant, with rdil = dZ/dt. During this dilution step, reactions are not frozen and continue to happen. The total enthalpy H varies according to a mixing law: P dqVH=dt ¼ i¼1;nspecies Hi fi , where Hi is the total enthalpy of each species i at the fresh gases conditions. 3. Mixture evolves until reaching its equilibrium state at constant _ i as in the first mixture fraction Zmin by solving dqVY i =dt ¼ qV x step. In the cases presented in this section, Zmin = 0.035 which corresponds to an equivalence ratio of 0.623. The total enthalpy H is kept constant during this step. To follow the evolution of carbon chemistry, a progress variable is defined as Y c ¼ Y CO þ Y CO2 þ Y H2 O . This latter is monotonic in time and ensures a one-to-one relation between the temperature and the progress variable for the considered mixture fraction values. The progress variable can then be normalized with the equilibrium condition C ¼ Y c =Y eq c ðZÞ. The normalized progress variable is equal to zero in fresh gases and unity in burnt gases. The first step of the PSR ends when the normalized progress variable C reaches 99%, a value that will be investigated in a subsequent section. This threshold allows to distinguish carbon chemistry from other slowly varying processes in the burnt gases. Table 1 Fuel and oxidizer compositions, T0 is the initial temperature.

2.1. Description of the canonical problem Aeronautical burners are mainly a two-inlet problem, where fuel and oxidizer are brought separately. In these burners, fuel and oxidizer mix and burn rapidly in the primary zone and the burnt gases are then diluted by fresh air before exiting the combustor. To perform the analysis of fresh air dilution on the nitrogen

Y CH4 Y O2 Y N2 YAR T0[K] H[J/g]

Fuel

Oxidizer

1.0 0.0 0.0 0.0 850.0 4659.3

0.0 0.2321 0.7555 0.0124 850.0 7.1

499

F. Pecquery et al. / Combustion and Flame 161 (2014) 496–509

The solver used to compute this constant pressure vessel is Cantera [35], and the chemical kinetic scheme is the GRI-Mech 3.0 mechanism [1]. The pressure in the vessel is equal to 10 bars which is close to operating conditions in mid-size jet engines. The initial condition of the PSR is defined by a linear combination of oxidizer and fuel states in terms of mass versus mixture fraction with the inflow conditions given in Table 1.

12

8

4

2.2. Results The results obtained for the steps two and three, i.e. when C has reached the value 0.99, are presented in Figs. 3–5. In the present case, nitrogen chemistry subjected to dilution is influenced by two parameters: the mixing rate dilution intensity and the initial mixture fraction, the final mixture fraction being kept constant. The influence of each parameter is studied in this section. 2.2.1. Influence of the initial mixture fraction In this section, five cases with different initial compositions but with the same ending mixture fraction Zmin are investigated. The aim is to evaluate the influence of the initial conditions on the NO production in the burnt gases. Four PSR at different initial mixture fraction values have been computed and diluted with different amounts of fresh air, to end up at the same final composition featuring Zmin = 0.035, i.e. an equivalence ratio of 0.623. Their initial mixture fractions Z are 0.055; 0.05; 0.045 and 0.04. The value Z = 0.055 corresponds to a stoichiometric mixture. A fifth PSR is computed at an initial mixture fraction equal to Zmin, thus without _ Y NO of these five PSR is shown in the any dilution. The evolution of x phase space (Z, YNO) in Fig. 3. Figure 4 is a 2D projection of Fig. 3, which shows that once the final mixture fraction Zmin is reached, all the PSR evolutions in _ Y NO ðY NO Þ, collapse onto a single curve whatever the phase space, x initial conditions. 2.2.2. Impact of the dilution intensity To study the influence of the dilution intensity, NO mass fraction is computed for different mixtures at the same initial mixture fraction Z = 0.055 diluted up to the same Zmin with different dilution rates rdil. The evolution of NO mass fraction is then compared to the one obtained for auto-igniting constant pressure vessels at constant mixture fraction. Results are given in Figs. 5 and 6. In

11 10 9 8 7 6 5 4 3 2 1 0 0

0 0

0.001

0.002

0.003

0.004

0.005

0.006

0.007

Fig. 4. Projection of Fig. 3 in the (YNO, Z) plane.

11 10 9 8 7 6 5 4 3 2 1 0 0.002 0.004 0.006 0.008 0.01

0.055

0.05

0.045

0.04

0.035

Fig. 5. PSR responses in state phase YNO, Z, with the same initial and final mixture fractions but different dilution rates: from 1 s1 to 20 s1 in steps of 1 s1 and from 20 s1 to 50 s1 in steps of 5 s1. Lines: response under dilution. Dashed lines: flame response at constant mixture fraction. Dotted-dashed line: equilibrium.

Fig. 5, the mixture fraction and NO mass fraction for the diluted constant pressure vessels (continuous lines) are given only for the steps 2 and 3, i.e. when the carbon chemistry is complete. The surface defined by the diluted responses overlaps the one of non-diluted mixtures at the same (Z,YNO). It is even clearer in Fig. 6. So it is possible to compute the YNO evolution in the burnt gases just by tabulating flame responses without any dilution whatever the dilution rate of the mixture. 2.3. Characteristic time scales of NO production under fresh air dilution

0.002 0.004 0.006 0.008 0.055

0.05

0.045

0.04

0.035

Fig. 3. Evolution of NO source term for five different configurations in (YNO, Z) phase space. The dotted line corresponds to a mixture at Z = 0.055 diluted to Zmin. The dash-dot-dot line corresponds to a mixture at Z = 0.05 diluted to Zmin. The dash-dot line corresponds to a mixture at Z = 0.045 diluted to Zmin. The dashed line corresponds to a mixture at Z = 0.04 diluted to Zminsx. The continuous line corresponds to a mixture at Zmin. The four diluted mixture are performed with a dilution rate equal to 10 s1.

From the study of the canonical problem, it appears that NO production in the burnt gases with air dilution is piloted by a single time scale, which is large compared to those of the carbon chemistry. In a turbulent flame, a broad range of mixing time scales are present and it is important to know whether the NO production can still be simplified and modeled using a single time scale in the burnt gases. To this aim, the characteristic time scales of NO production are investigated by performing a perturbation analysis. The thermo-chemical state of an adiabatic constant–pressure system follows:

8 @h > ¼0 > < @t @p ¼0 @t > > : @Y i ¼ x _ Y i ; ði ¼ 1; . . . ; ns Þ @t

ð1Þ

500

F. Pecquery et al. / Combustion and Flame 161 (2014) 496–509 2.5 1 2 0.75 1.5 0.5 1 0.25 0.5

0

0 0

0.002

0.004

0.006

0.008

0

0.01

(a)

6

4

2

0 0.002

0.004

0.006

0.008

(b) Fig. 6. (a) Cut of Fig. 5 in the Z = 0.04 plane. (b) Cut of Fig. 5 in the Z = 0.05 plane. Dashed line: response of the PSR at constant mixture fraction. Dots: responses of diluted mixtures conditioned by the chosen mixture fraction.

where h is the total enthalpy, p the pressure and ns the number of species. This system evolves in a state space of ns + 2 dimensions. It can be written in vector form

@~ w wÞ; ¼ Fð~ @t

ð2Þ T

1

1.5

2

Fig. 7. Evolution of the normalized progress variable C versus C+ for a PSR in the conditions of Table 1 and Z = 0.055 at 10 bar.

8

0

0.5

w ¼ ðh; p; Y 1 ; . . . ; Y ns Þ . The system is highly non-linear and the with ~ chemical time scales may be extracted using the ILDM methodology [21], performing a first order approximation of the system. It follows that each eigenvalue of the Jacobian matrix of F, Jw, is associated with a time scale or mode of the locally linear solution to the chemical system, Eq. (1). After the removal of the four eigenvalues associated with the conservation of atoms C, H, O and N, the eigenvalues of Jw have a negative real part and may be ordered according to the increasing values of their real parts kr, such that kr1 6 kr2 6 . . . 6 krns . The real part represents the relaxation frequency of the system to a small perturbation. Relaxation times are linked to eigenvalues such that si ¼ 1=kri . A relaxation frequency cut-off parameter ktol is introduced. The dimensionality of the manifold, noted nm, is set equal to the number of eigenvalues, which have a real part smaller than  ktol, i.e. a relaxation time faster than 1/ktol. This dimensionality is strongly time dependent. To understand the weak dependence of the nitrogen oxide production to dilution once examined in the Z, YNO space, the computation of the manifold for a PSR in the conditions of Table 1 is performed. Three values of ktol are investigated: 105 s1, 106 s1 and 107 s1. To clarify the results, a þ eq new progress variable is introduced: C þ ¼ Y c =Y eq is c þ Y NO =Y NO ; C

equal to zero in fresh gases and is equal to two at equilibrium; C+ evolves monotonically both in the fresh and burnt gases. Figure 7 shows the evolution of the conventional progress variable Y c =Y eq c as a function of C+. Only the time scales of the kinetics are studied here: the eigenvalues of pressure and enthalpy are not computed. In practice, the first two columns and lines of Jw are thus ignored. Five eigenvalues are discarded because of the conservation of the five chemical elements: C, O, H, N and Ar. The number of remaining eigenvalues below ktol is shown in Fig. 8. It is seen that for the smallest ktol value, the manifold becomes one-dimensional when carbon chemistry is almost at equilibrium, i.e. when C+ P 0.8 and that at equilibrium, the number of dimensions is zero. For larger relaxation cut-off parameters, the number of dimensions increases. To know whether YNO can be considered as a progress variable for the nitrogen chemistry in the burnt gases, i.e. C+ P 1, one has to check that the direction of the eigenvector ~ v associated with the first eigenvalue is not orthogonal to the YNO axis in the phase space. If ~ Y NO is defined as the unitary vector of the YNO axis, one may compute the scalar product ~ v =j~ v j. This scalar product is Y NO  ~ always greater than 0.7 in the burned gases for various conditions. It validates the choice to tabulate the 1D manifold in burnt gases as a function of the NO mass fraction. It may be noticed that in the case of a 1D manifold, if the progress variable is chosen in the direction of the manifold, the Flamelet-Prolongation of ILDM (FPI) method [36,37] and the ILDM method lead to the same formalisms. Tabulating YNO as a function of Z is an easy task, be_ Y NO ðY NO ; ZÞ and its derivative are continuous. For a mixture cause x

Fig. 8. Evolution of the dimensionality nm as a function of C+ and the cut-off parameter ktol for a PSR in the conditions of Table 1 and Z = 0.055 at 10 bar.

501

F. Pecquery et al. / Combustion and Flame 161 (2014) 496–509

without dilution, the chemical composition evolves necessarily in the state space in a plane of constant mixture fraction. The 1D manifold evolution is driven by the smallest eigenvalue, i.e. the largest time scale. This slow evolution may be explicitly resolved with turbulent mixing and flow transport phenomena. Then, the relaxation of the chemical system to the 1D manifold as a response to any perturbation is driven by the other eigenvalues. The largest relaxation time corresponds to the second eigenvalue of the system: srelax ¼ 1=kr2 . If the dilution by fresh air is considered as a perturbation of the system, this time can be r linked to a maximum dilution rate rmax dil ¼ k2 due to the dilution of the burnt gases. If the burnt gases undergo smaller dilution rates, the chemical system relaxes quickly to the manifold. For typical conditions of an aeronautical mid-size gas turbine, i.e. p = 10 bar and T 0 ¼ 850 K;1=rmax dil ¼ 1=17; 000 s. An estimation of the mean dilution rate of the burnt gases in typical aeronautical engines at full load conditions can be performed. Assuming a residence time of 5ms in the burnt gases, a stoichiometric primary zone, and a global fuel/air ratio of the combustor around 30‰, one obtains a mean dilution rate hrdili of 5 s1. This mean dilution rate is several orders of magnitude smaller than the dilution rate corresponding to the second eigenvalue of the system. Obviously, the value of rmax is a function of the thermodynamic conditions: dil p, h and Yi, and the local dilution rate may be larger than the mean dilution rate in some regions close to the dilution holes and the cooling films but r max dil remains very large for the different investigated conditions, validating the 1D manifold assumption. 2.4. Conclusion The results of the canonical problem highlight two characteristics of the nitrogen chemistry under dilution. For the dilution rates observed in practical systems, the NO evolution may be determined from mixtures responses without dilution whatever the dilution intensity. Flame responses may be then tabulated as functions of YNO. This method takes into account all chemical processes that are considered in the chemical scheme used to compute the manifold, so that complex chemistry effects are included in the CFD computation. This property is only valid in the burnt gases. According to Schmidt et al. [38], an ILDM computed for a given range of element compositions, enthalpy and pressure can be used arbitrarily for laminar flames or turbulent flows, if the same thermodynamic properties are encountered. Therefore, results obtained with the PSR may be extrapolated to 1D premixed flames, and the 1D manifold may be tabulated from a collection of 1D premixed flames for a range of mixture fractions.

hand, extracting the NO source term as a function of the progress variable Yc in the flame front ensures a coupling (location, dynamic response) between the flame front chemistry and nitrogen chemistry. On the other hand, tabulating nitrogen chemistry of burnt gases as a function of YNO gives access to the nitrogen evolution in the _m dilution regime, as explained in Section 2. The source term x Y NO is expressed with a carbon chemistry table in the flame front _ FF (x Y NO ðC; ZÞ) and with a nitrogen chemistry table in the burnt gases _ BG (x Y NO ðY NO ; ZÞ). The parameter used to switch from the carbon chemistry table to the nitrogen chemistry table is the normalized progress variable C. The source term may be written: Prompt NO

m Y NO

x_

zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl}|fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{ _ FF _ BG ¼x Y NO ðC; ZÞHð1    CÞ þ xY NO ðY NO ; ZÞHðC  ð1  ÞÞ; |fflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflfflffl}

ð4Þ

NO produced in burnt gases

where H is the Heaviside function. This model makes no assumption on the chemistry itself, as all the mechanisms that are considered during the a priori calculation of the chemical table, are considered during the CFD computation. The parameter  is used to discriminate between the flame front and burnt gases. Furthermore, NO must be used in the burnt gases where dilution occurs. As a result,  has been chosen to be as small as possible and to have a NO mass fraction evolution as smooth as possible in the carbon chemistry table. To ease the choice of , the variation in YNO as a function of Yc and the variation in the NO source term as a function of YNO when C = (1  ) are plotted in Fig. 9. It is seen that a range of  values between 0.01 and 0.05 give dYNO/dYc small enough to allow for tabulation. Here,  is chosen as small as possible, i.e.  = 0.01 to make a clear distinction between the flame front and the burnt gases. As a preliminary test of the model, 1D premixed flames are simulated to verify that the model captures both prompt and thermal NO. The results with the NOMANI model are shown in Fig. 10 and compared to a detailed chemistry solver. The NO levels are well predicted in both the flame front and in the burnt gases for a wide range of mixture fraction values. These results validate the interpolations in the carbon and nitrogen look-up tables and also the transition between the two parametrizations. At this transition, where  = 0.01, source terms are retrieved at the boundaries of the lookup tables and interpolation errors may arise. 4. NOMANI in LES The LES formalism relies on the filtering of the Navier–Stokes equations. The filtering is denoted by an overbar and it is associated with the filter width D. In LES, the filter width is usually determined by the local mesh size. For the transported quantities, the

3. NOMANI model: nitrogen oxide emission model with one-dimensional MANIfold 3.1. Description of the NOMANI model The YNO balance equation reads:

  @ qY NO @ qui Y NO @ l @Y NO _ Y NO ; þ qx þ ¼ @xi Sc @xi @t @xi

ð3Þ

where q, ui, l, Sc and YNO represent density, velocity components, dynamic viscosity, Schmidt number (Sc = 0.75) and NO mass fraction. The novelty of this work is the source term closure, which is well suited to track NO both in the flame front and in the burnt _m gases. The model source term is named x Y NO and the fully detailed _ Y NO . The choice has been made to track the NO evochemistry one x lution with two progress variables, one for the carbon chemistry and another for the nitrogen chemistry. This method has two advantages in comparison with the previous models. On the one

Fig. 9. Chemical source term evolution in a 1D premixed flame for several values of _ NO =dY NO . the mixture fraction. Continuous lines: dYNO/dYc. Dotted lines: dx

502

F. Pecquery et al. / Combustion and Flame 161 (2014) 496–509

Fig. 10. Spatial evolution of YNO(Z, Yc) in 1D premixed flames for several mixture fraction values. Symbols: detailed mechanism. Lines: CFD solver with the NOMANI model.

e ¼ q/=q. This section focuses common Favre averaging is used: / e_ m . The filon the closure of Eq. (4), and particularly on the term x Y NO tered source term is approximated by the following expression:

_ BG qx_ mYNO ¼ x_ FF Y NO ðC; ZÞHð1    CÞ þ qxY NO ðY NO ; ZÞHðC  ð1  ÞÞ e_ FF ð C; e C v ; Z; e e Z v ÞHð1    CÞ x ’ q Y NO e_ BG ð Y e  ð1  ÞÞ; e NO ; ZÞHð e x þq C Y NO

ð5Þ

where Cv and Zv are the subgrid-scale variances of C and Z. The assumption of this closure is that the filtered Heaviside convoluted with a function is equal to the non-filtered Heaviside function of a e_ FF filtered parameter times the filtered function. The source term x Y NO is obtained from the PCM-FPI model [39], where the sub-grid filtered density function of the progress variable is assumed to be a beta function. Previous works [40] highlighted the efficiency of this e_ Y . It is important to stress that the formalism for the source term x c model relies on the retrieving of the NO source term and not the NO mass fraction directly. Therefore, if the progress variable variance is high, the beta-pdf function tends toward a two-delta function and more weight will be given to the extremal values C = 0 and C = 1. As shown in Fig. 11, the NO source term has a shape which is similar to the one of the progress variable, i.e. that the extremal values are negligible compared to the values reached in the flame front. In the burnt gases, mixture fraction and NO mass fraction variances are small compared to those in the flame front due to the strong increase in the diffusivities. As a result, a delta function is e_ BG . used for x Y NO To give more insight into these modeling assumptions, a stoichiometric 1D premixed unstrained methane/air flame at atmospheric pressure is spatially filtered with the following Gaussian filter:

e /ðxÞ ¼

R þ1

0

qðxÞ/ðxÞGD ðx0  xÞdx ; qðxÞGD ðx0  xÞdx0

1 R þ1 1

ð6Þ

where 2 6 2 Gðx; DÞ ¼ pffiffiffiffi exp6x =D :

p

ð7Þ

Fig. 11. Chemical source term evolution in a 1D flame for several mixture fraction _ Y c . Dotted lines: x _ Y NO . values. Continuous lines: x

Then, above closures may be compared with the exact filtered terms at different filter widths. An important parameter is the ratio between the filter width D and the flame thickness is obtained from the gradient of the progress variable, dY c ðZÞ ¼ Y eq c =jrY c jmax , which is equal to 386 lm. To validate the assumption on the Heaviside function, two ratios n1 and n2 are defined:

R þ1

e_ FF ð CÞHð1 e e Þdx C q x Y NO ; qx_ FF Y NO ðCÞHð1    CÞdx 1

n1 ðD; Þ ¼ R1 þ1

ð8Þ

and

R þ1

e_ BG ð CÞHð e e  ð1  ÞÞdx C q x Y NO : BG qx_ Y NO ðCÞHðC  ð1  ÞÞdx 1

n2 ðD; Þ ¼ R1 þ1

ð9Þ

The function n1(D,) is the ratio of the mass of NO that is evaluated by the modeled source term in the flame front and the exact mass of NO that is produced in a filtered flame. Results are shown in Fig. 12(a). The parameter that is the most influential is . If  is about 1%, the relative error is negligible even if the filter width is

503

F. Pecquery et al. / Combustion and Flame 161 (2014) 496–509

presented in Fig. 14(a) and in Fig. 14(b) for the same filter width, the second figure being a zoom of the first figure. Additionally, the NO mass fraction source term is given in Fig. 14(c) as a function of the normalized NO mass fraction. These figures highlight the fact that in the burnt gases, i.e. when Y þ NO > 1, all the filtered fields coincide with the unfiltered fields. This validates the assumption to use a delta pdf in the burnt gases.

1.2 1.18 1.16 1.14 1.12 1.1 1.08 1.06

5. Validation of the NOMANI model: LES of the Sandia flame D

1.04

5.1. Description of the configuration

1.02 1 0

2

4

6

8

10

12

14

16

(a) 1.002 1.0018 1.0016 1.0014 1.0012 1.001 1.0008 1.0006 1.0004 1.0002 1 0

2

4

6

8

10

12

14

16

(b)

LES results are compared with experimental data of the Sandia flame D studied experimentally by Barlow and Frank [41]. This amply documented flame, was previously used to validate nitrogen oxide emission models [31,29] and extinction effects [42,43]. It is a piloted partially premixed turbulent jet flame with a negligible degree of local extinctions. The Reynolds number at the jet nozzle is 22,400, the jet diameter D is 7.2 mm and the jet stream consists of a mixture of 25% methane and 75% air by volume. The flame is stabilized by a hot pilot stream surrounding the jet with a diameter of 2.62 D. The hot pilot flow is a lean burnt mixture of C2H2, air, CO2, H2 and N2 with a bulk velocity of 11.4 m/s. The coaxial burner is further surrounded by co-flowing air with a bulk velocity of 0.9 m/s. To describe the mixing of the main and pilot jets with the air coflow, a mixture fraction Z is defined, which is taken equal to zero in the air and one in the jet. In the pilot stream, the inlet conditions correspond to the equilibrium conditions of a methane/air perfectly stirred reactor at Z = 0.271. The thermodynamic properties used to compute the chemical table are referenced in Table 2. 5.2. Governing equations

Fig. 12. Representation of the two ratios: (a) n1(D, ) and (b) n2(D, ).

In the simulations, tabulated chemistry based on FPI and the Presumed-Conditional Moments (PCM) formalism [39] is used. Since the Mach number in the Sandia flames is quite small, the Navier–Stokes equations are not solved in the compressible form, but in the low-Mach number formalism with the massively parallel flow solver YALES2 [44]. Assuming that the unclosed terms can be written in terms of eddy viscosity and eddy diffusivity, equations solved for these simulations read:

1

0.8

0.6

 ej @q ei u ej u u @q @P @  þ ¼ þ2 l þ lt eS ij ; @xj @xi @t @xi

0.4

!   e c @q ec ei Y Y u @q @ l lt @ Ye c e_ ; x þ ¼ þ þq Yc @xi @t @xi Sc Sct @xi

0.2

0

ð10Þ

0

0.005

0.01

0.015

0.02

Fig. 13. Stoichiometric premixed flame filtered at 20 different filter widths by steps of 80 lm. Dashed line: unfiltered flame response. Continuous line: filtered flame response.

important. The function n2(D, ) is the ratio of the mass of NO that is evaluated by the modeled source term in the burnt gases and the exact mass of NO that is produced in the burnt gases of a filtered flame. Results are shown in Fig. 12(b). The relative error is negligible whatever D. Then, to study the model assumptions in the burnt gases, i.e. the e_ BG , various fields of filtered flames are use of a delta function for x Y NO presented in Figs. 13 and 14(c). The filtered progress variable in Fig. 13 allows to compare the filter width to the laminar flame thickness. The normalized NO mass fraction Y þ NO ¼ Y NO =Y NO ðC ¼ 0:99Þ is

e i Y cv  Y cv @ q u @q @ þ ¼ @xi @t @xi 



l

Sc

lt D2 Sct



þ

Y cv

lt @Y cv



ð11Þ

2  lt @ Ye c  þ2  2 

@xi Sct  @xi 

e_ Y  x ecx _g  Y þ 2q Yc Y c ; c

Sct

ð12Þ

  ! e @q e ei Z Z u @q @ l lt @ Ze ; þ ¼ þ @t @xi @xi Sc Sct @xi  Zv @q u e i Zv @q @ þ ¼ @xi @t @xi



l

Sc



þ

lt @Z v Sct

@xi



ð13Þ  2 lt @ Ze  l þ2    2 2 t Zv ; Sct @xi 

D Sct

ð14Þ ff f ee where Y cv ¼ Yg c Y c  Y c Y c and Z v ¼ ZZ  Z Z. Because of the lowMach formalism and the NOMANI closure, two additional equations

504

F. Pecquery et al. / Combustion and Flame 161 (2014) 496–509 25

2.5

20

2

15

1.5

10

1

5

0.5

0

0 0

0.2

0.4

0.6

0.8

1

0

0.01

0.02

(a)

0.03

0.04

(b) 0.012 0.01 0.008 0.006 0.004 0.002

4

8

12

16

20

24

28

(c) Fig. 14. Stoichiometric premixed flame filtered at 20 different filter widths by steps of 80 lm. Dashed line: Unfiltered flame response. Continuous lines: filtered flame response. Dotted lines: Y þ NO ¼ 1. (a and b) Normalized NO mass fraction in the physical phase space. (c) NO mass fraction source term as a function of the normalized NO mass fraction Y þ NO .

Table 2 Fuel and oxidizer compositions corresponding to the Sandia flame D.

Z Y CH4 Y O2 Y N2 YAR T[K] H[J/g]

Fuel

Oxidizer

1.0 0.1576 0.198 0.6428 0.0016 291.0 737.74

0.0 0.0 0.2321 0.7555 0.0124 294.0 7.13

are solved: one for the Poisson equation and one for the nitrogen oxide mass fraction. In Eq. (10), the resolved velocity gradient tensor reads:

  ei @ u ej 1 @u e ; þ S ij ¼ 2 @xj @xi

ð15Þ

and lt is closed with the Smagorinsky model [45], using the Germano dynamic procedure [46] to determine the constant. In Eq. (11) to Eq. (14), Sc is equal to 0.75 and Sct is equal to 1. The fuel injection in the main jet has a critical influence on the development of the flame. Therefore, a mean velocity profile corresponding to a pipe flow at a Reynolds number of 22,400 is imposed assuming that the turbulent flow in the inlet pipe is fully developed:

v v max

¼

1 1 y n ; b 0:5D

ð16Þ

where n = 6.6, b = 0.807, and v and vmax are the mean and maximum averaged velocities respectively. Homogeneous isotropic turbulence is added to this inlet condition with a level of fluctuations that

matches the experimental data and with approximately seven energetic scales in the jet diameter. The computational domain extends to one diameter upstream of the burner. A slip wall boundary condition is imposed at the walls of the burner to get a good agreement with the experimental data.

5.3. Computational procedure Most of the results shown in the literature are obtained with structured grids and do not exceed tenth of million cells. In order to assess the mesh resolution and turbulent closures dependence, two different unstructured grids with a high resolution have been designed. The first grid, M1 is composed of 43 million cells, with 40 cells in the inlet diameter and three cells at the lip of the burner. The second grid, M2 is a homogeneous refinement of M1. All the tetrahedra of M1 are split into eight parts, so that M2 is composed of 347 million cells. Some parameters of these two meshes are detailed in Table 3. Another important point in the tabulated chemistry approach is the choice of canonical flames used to parametrize the temperature and species as functions of the driving parameters, i.e. the progress variable Yc, the mixture fraction Z and their variances SC and SZ, respectively. Vreman et al. [25] and Fiorina et al. [47] have shown that one-dimensional premixed flames with a unity Lewis number are able to reproduce quite accurately the Sandia flame D. In the present work, one-dimensional premixed flames with unity and variable Lewis numbers are computed with the Cantera package and tabulated. Including the Lewis number effects in the flamelets allows to take into account differential diffusion effects in the progress variable gradient direction. Other differential diffusion effects are not considered in this methodology. However, the

505

F. Pecquery et al. / Combustion and Flame 161 (2014) 496–509 Table 3 Parameters of the two meshes.

2000 1500

Mesh

M1

M2

Cells, 106 Dmin, in the lip of the burner Cells in the diameter of the burner

43 82 lm 40

347 41 lm 80

1000 500 0 0

20

40

60

80

0

20

40

60

80

400

Table 4 Computational cases.

300

Computation

Case 1

Case 2

Case 3

Case 4

200

Mesh Transport model Accumulation time

M1 Lek – 1 50 ms

M1 Lek = 1 50 ms

M2 Lek – 1 80 ms

M2 Lek = 1 80 ms

100

differential diffusion effects in the progress variable gradient direction have a strong impact on the laminar flame speed and on the heat release rate. The laminar premixed flame computations are done with the detailed GRI-Mech 3.0 chemical mechanism [1]. Any local property of one-dimensional premixed flames /(C, Z, SC, SZ) is stored in a 4D-table (100  100  16  16). The progress variable direction is equally spaced unlike the other e_ BG is extracted from a 2D-table with directions. The source term x Y NO 300 points in the direction of C and 400 points in the direction of Z. Four different computations are performed with the two meshes M1 and M2 and with unity and variable Lewis numbers. The parameters of these four cases are detailed in Table 4.

0

Fig. 15. Axial mean and rms temperature at the centerline. Symbols: experiment, dashed lines: variable Lewis number, continuous lines: unity Lewis number, bold: mesh M2, thin: mesh M1.

1 0.75 0.5 0.25 0 0

20

40

60

80

0

20

40

60

80

0.2

5.4. Statistics 0.15

Before proceeding to the analysis of the flame topology, the statistics of the flame are first discussed. A common way to assess the mesh resolution in LES is to look at the ratio of turbulent to molecular viscosities lt/l. Computations performed with the mesh M2 exhibit a ratio, which is comprised between 0 and 2, which guarantees a good resolution of the flow field and the flame. Statistics of temperature, mixture fraction, velocity and NO mass fractions are presented on the jet axis in Figs. 15–19. The figures show that a good agreement is found between the experiment and the calculations especially for the finest mesh M2. Figure 16 shows that the quality of the mixing prediction fully depends on the mesh resolution. It may be explained by the aerodynamic instabilities that occur just after the lips of the burner, which play an important role in the mixing further streamwise and which are well resolved on the finest mesh M2. Compared to mixture fraction statistics, temperature statistics shown in Fig. 15 are more sensitive to the transport model. The mean mixture fraction and streamwise velocity are very similar for the cases 3 and 4 (Lek = 1) so that the flame length difference for these two cases is attributed to the source term of the progress variable. It is known that the source term of the progress variable is 25% smaller using the unity Lewis number assumption than using the complex transport model. The flame length is also sensitive to the mesh resolution as only the finest mesh M2 reproduces well the temperature statistics. Results are in agreement with Kemenov et al. [48] who have already shown the great influence of the grid resolution on this piloted flame. The NO emission results are more sensitive to the transport model, even if the mesh resolution determines the position of the YNO maximum. The NO level is smaller in all the cases with unity Lewis number assumption. To conclude, on the coarse mesh, the Lewis number choice has an important impact on the predictions and a unity Lewis number is preferable, when premixed flames are used for tabulation. For all

0.1 0.05 0

Fig. 16. Axial mean and rms mixture fraction at the centerline. Symbols: experiment, dashed lines: variable Lewis number, continuous lines: unity Lewis number, bold: mesh M2, thin: mesh M1.

8e-05 6e-05 4e-05 2e-05 0 0

20

40

60

80

0

20

40

60

80

2.5e-05 2e-05 1.5e-05 1e-05 5e-06 0

Fig. 17. Axial mean and rms NO mass fraction at the centerline. Symbols: experiment, dashed lines: variable Lewis number, continuous lines: unity Lewis number, bold: mesh M2, thin: mesh M1.

statistics the better results are obtained with the finest mesh M2 and the unity Lewis number assumption.

506

F. Pecquery et al. / Combustion and Flame 161 (2014) 496–509 70 6e-05

Y_NO Mean

60 50 40 30 20

4e-05 3e-05 2e-05 1e-05

10 0

5e-05

0 0

20

40

60

0

20

0

20

40

60

80

40

60

80

80

8

Y_NO RMS

2.5e-05

6 4

2e-05 1.5e-05 1e-05 5e-06

2

0

0 0

20

40

60

80

Fig. 18. Axial mean and rms streamwise velocity at the centerline. Symbols: experiment, dashed lines: variable Lewis number, continuous lines: unity Lewis number, bold: mesh M2, thin: mesh M1.

5e-05

5e-05

4e-05

4e-05

3e-05

3e-05

2e-05

2e-05

1e-05

1e-05

0

0

1

2

0

3

5e-05

5e-05

4e-05

4e-05

3e-05

3e-05

2e-05

2e-05

1e-05

1e-05

0

0

3

x/D

6

Fig. 21. Axial mean and rms NO mass fraction at the centerline for Case 2. Symbols: experiment, dashed line: GRI-2.11 mechanism, continuous line: GRI-3.0 mechanism.

9

Fig. 22. Instantaneous temperature fields (left) and instantaneous field of YNO (right), computed on the finest mesh M2.

0 0

1

2

3

0

3

6

9

5.5. Influence of the chemical mechanism Fig. 19. Radial mean NO mass fraction at various downstream positions. Symbols: experiment, dashed lines: variable Lewis number, continuous lines: unity Lewis number, bold: mesh M2, thin: mesh M1.

Mean Temperature

2000 1500 1000 500 0

0

20

40

60

80

RMS Temperature

400 300 200 100

In all the calculations presented in the previous sections, the methane-air GRI-3.0 [1] mechanism is used. It has been shown by Cao and Pope [49] that the GRI-3.0 mechanism may overpredict NO emissions in similar configurations while the GRI-2.11 mechanism gives better estimates. To confirm these results, the Case 2 with the mesh M1 and the unity Lewis number has been carried out with the GRI-2.11 mechanism. Both calculations are compared in Figs. 20 and 21. As highlighted by the temperature statistics in Fig. 20, all the fields are very similar. The only noticeable differences appear on the NO mass fraction as shown in Fig. 21. In accordance with similar studies, the overprediction of the mean NO mass fraction with the GRI-3.0 mechanism between 25D and 50D is corrected with the GRI-2.11 mechanism and a good agreement is obtained over the whole axis. An impact is also visible on the rms values with lower values below 40D and an overprediction beyond. The overprediction may be attributed to a lack of statistical convergence in this zone where the flow evolves very slowly. 5.6. Flame topology

0 0

20

40

60

80

x/D Fig. 20. Axial mean and rms temperature at the centerline for Case 2. Symbols: experiment, dashed line: GRI-2.11 mechanism, continuous line: GRI-3.0 mechanism.

The case 4 is in good agreement with experimental data, therefore this section is dedicated to the study of the instantaneous results of this case. The instantaneous flow field and flame structure are shown in Fig. 22. The flame stabilization provided by the pilot flames is clearly visible. The shear layers surrounding the jet core

F. Pecquery et al. / Combustion and Flame 161 (2014) 496–509

Fig. 23. Scatter plot of temperature versus mixture fraction Z at x = 7.5D. Black: experiment, red: simulation with the M2 mesh and unity Lewis number. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

507

Fig. 26. Conditional mean of temperature and rms versus mixture fraction Z at x = 30D. Black thick line: experiment, red thin line: simulation with the M2 mesh and unity Lewis number. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Fig. 24. Conditional mean of temperature and rms versus mixture fraction Z at x = 7.5D. Black thick line: experiment, red thin line: simulation with the M2 mesh and unity Lewis number. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.) e_ > 0:01, colored by C e (left), the Fig. 27. Instantaneous fields conditioned by x Yc e (right). An iso-contour of the mixture fraction equal to Takeno index (middle) and Z 0.6 is represented by a red line (middle). (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Fig. 25. Scatter plot of temperature versus mixture fraction Z at x = 30D. Black: experiment, red: simulation with the M2 mesh and unity Lewis number. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

region are dominated by quasi laminar structures, which can be attributed mainly to the increase in the molecular viscosity due to the burnt gases. Figures 23–26 are scatter plots and conditional means of temperature versus mixture fraction at two different heights. In the scatter plots, the same number of samples is plotted for experimental and numerical data. In the scatter plots, both experiments and simulations exhibit temperature values far from the equilibrium highlighting the fact that infinitely fast chemistry would not provide accurate results at this mesh resolution. In the conditional mean plots, Figs. 24 and 26, the experimental and numerical conditional means are in a good overall agreement with some under prediction of the mean at x = 7.5D. The conditional rms are overpredicted at x = 7.5D but well captured at x = 30D. The highly resolved LES of the case 4 allows the study of the flame stabilization. In the Sandia flames, the main jet is a rich premixed mixture and the flame is stabilized by burnt gases of a lean premixed mixture. It is therefore not obvious to extract the driving

508

F. Pecquery et al. / Combustion and Flame 161 (2014) 496–509

Fig. 28. Illustration of the flow asymmetry on temperature fields.

Fig. 29. Instantaneous normalized progress variable fields (left) and instantaneous field of YNO (right), computed on the finest mesh M2.

parameters of the flame stabilization. A first step in this direction is to analyze the normalized Takeno index [50], which can been defined as

rY f  rY o G¼ ; jrY f  rY o j

ð17Þ

where Yf is the fuel mass fraction and Yo the oxidizer mass fraction. This value is equal to 1 in the non-premixed regime and is equal to 1 in the premixed regime. Figure 27 shows the Takeno index for an instantaneous field. The particularity of the Sandia flame D is the long length core, where the mixture is premixed with mixture fraction values between 0.6 and 1. In this area, the combustion regime is exclusively a premixed regime, but when the major part of methane is burnt, the dominant regime is the non-premixed regime. The main jet has a great influence on the topology of the flame, which justifies in some sense the use of premixed look-up tables for the computation of this flame. Another feature of this Sandia flame D is the high amount of tangential asymmetry. This asymmetry is highlighted plotting

the mean temperature at different heights in Fig. 28. These lowfrequency modes, which are typical of jets and which could not be observed with RANS or coarse LES, seem to evolve very slowly during the computations and are insensitive to the turbulence injection parameters. The numerical set-up has been checked several times and the results are the same with different meshes. A particular attention is paid to these modes in the paper because they may cause some difficulties when comparing the experimental results to the computational ones. The experimental and numerical results may be averaged in the tangential direction to limit the potential discrepancies. As a consequence, all the radial statistics in this work are tangentially averaged. With the NOMANI model, it is possible to estimate the amount of the NO that is produced in the flame front and the amount of the NO that is produced in the burnt gases. The ratio of these two quantities is called c and is defined by the relation (18). Figure 29 shows instantaneous fields of C and YNO. Where C is greater than 0.99, i.e. in the red area, the mixture is exclusively composed of burnt gases. For the Sandia flame D and for a height below 80D, c is equal to 0.36. The important value of c means that it is mandatory to take into account the NO that is produced in the flame front.

R

e_ j C e < 0:99ÞdV x ðq NO e_ NO j C e > 0:99ÞdV x ðq V

c ¼ RV

ð18Þ

6. Conclusions In this paper, a novel model for NO formation in aeronautical burners is derived. The design of the model starts from the analysis of a canonical problem based on a constant pressure vessel with dilution by fresh air. This canonical problem is used to study the influence of the dilution on the NO emissions. The results exhibit a one-dimensional manifold in the burnt gases, so that the NO production can be known tabulating this one-dimensional manifold: x_ NO ðY NO ; ZÞ. The only assumption in this tabulation technique is that the chemical system evolves over a 1D-manifold, with a relaxation toward this manifold faster than any other process, as air

F. Pecquery et al. / Combustion and Flame 161 (2014) 496–509

dilution. This is in order for the chemical system to stay on the tabulated 1D manifold. Based on these results, a new model called NOMANI is designed. This model keeps track of the NO formation with two progress variables, the first one is Yc, the conventional progress variable for the carbon chemistry. And the second one is YNO, which parametrizes the 1D manifold in the burnt gases. This model is able to evaluate the fast production of NO in the flame front and the slower one in the burnt gases taking into account complex chemistry effects. The model has been extended to LES based on the PCM-FPI approach. Finally, it is applied and validated in highly resolved LES of the Sandia flame D experiment. Very good agreement is obtained for the mixture fraction, the temperature and the NO mass fraction when using the most refined mesh and a unity Lewis number assumption for the transport coefficients. During the analysis of the Sandia flame D calculations, it was observed that the flame features a strong asymmetry in the radial and tangential directions at certain heights. This asymmetry is attributed to low-frequency aerodynamic modes that are typical of jets. These modes may be an issue for the comparison of experimental and numerical data, especially when only a few planes orthogonal to the streamwise direction are extracted. In the present model, radiative heat losses, which may be important in some situations, are not taken into account. However, proposing an extension of the model with the enthalpy defect as an additional control parameter is possible in the present framework.

[8] [9] [10] [11] [12] [13] [14] [15]

Acknowledgments

[33] [34] [35] [36]

The authors acknowledge that the results in this paper have been obtained using the PRACE Research Infrastructure resource Curie based in France at the CEA-TGCC center under the allocation x2012026880. This work was also granted access to the HPC resources of IDRIS under the allocation 2012-6880 made by GENCI (Grand Equipement National de Calcul Intensif). References [1] M. Frenklach, H. Wang, M. Goldenberg, G.P. Smith, D.M. Golden, C.T. Bowman, R.K. Hanson, W.C. Gardiner, V. Lissianki, Tech. rep., Gas Research Institute, 1995. [2] C. Fenimore, H. Fraenkel, Proc. Combust. Inst. 18 (1981) 143–149. [3] A.E. Bakali, L. Pillier, P. Desgroux, B. Lefort, L. Gasnot, J. Pauwels, I. da Costa, Fuel 85 (7–8) (2006) 896–909. [4] A. Konnov, Combust. Flame 156 (11) (2009) 2093–2105. [5] N. Lamoureux, P. Desgroux, A.E. Bakali, J. Pauwels, Combust. Flame 157 (10) (2010) 1929–1941. [6] W. Bartok, V. Engleman, R. Goldstein, E. Del Valle, AIChE Symp. Ser. 68 (126) (1972) 30–38. [7] P. Glarborg, J. Miller, R. Kee, Combust. Flame 65 (2) (1986) 177–202.

[16] [17] [18] [19] [20] [21] [22] [23] [24] [25] [26] [27] [28]

[29] [30] [31] [32]

[37] [38] [39] [40] [41] [42] [43] [44] [45] [46] [47] [48] [49] [50]

509

Y. Zeldovich, Acta Physicochim., USSR 21 (1946) 577–628. G. De Soete, Symp. (Int.) Combust. 15 (1) (1975) 1093–1102. J. Wolfrum, Chem. Ing. Tech. 44 (10) (1972) 656–659. P. Malte, D. Pratt, Combust. Sci. Technol. 9 (5-6) (1974) 221–231. J. Bozzelli, A. Dean, Int. J. Chem. Kinet. 27 (11) (1995) 1097–1109. S. Correa, Combust. Sci. Technol. 87 (1–6) (1993) 329–362. A. Lefebvre, Gas Turbine Combustion, CRC, 1999. G. Loffer, R. Sieber, M. Harasek, H. Hofbauer, R. Hauss, J. Landauf, Ind. Eng. Chem. Res. 44 (17) (2005) 6622–6633. V. Raman, R.O. Fox, A.D. Harvey, Combust. Flame 136 (3) (2004) 327–350. H. Wang, Y. Chen, Chem. Eng. Sci. 59 (16) (2004) 3477–3490. S.H. Kim, K.Y. Huh, Combust. Flame 138 (4) (2004) 336–352. A. Kronenburg, M. Kostka, Combust. Flame 143 (4) (2005) 342–356. D. Bradley, P. Gaskell, X. Gu, M. Lawes, M. Scott, Combust. Flame 115 (4) (1998) 515–538. U. Maas, S. Pope, Symp. (Int.) Combust. 24 (1) (1992) 103–112. J. Nafe, U. Maas, Proc. Combust. Inst. 29 (1) (2002) 1379–1385. H. Pitsch, H. Steiner, Phys. Fluids 12 (10) (2000) 2541. A. Odedra, W. Malalasekera, Combust. Flame 151 (3) (2007) 512–531. A. Vreman, B. Albrecht, J. van Oijen, L. de Goey, R. Bastiaans, Combust. Flame 153 (3) (2008) 394–416. G. Godel, P. Domingo, L. Vervisch, Proc. Combust. Inst. 32 (1) (2009) 1555– 1561. J.A. van Oijen, P.H. de Goey, Proc. Eur. Combust. Meeting (2009) 810248. G. Godel, Turbulent combustion Sub-Grid Scale modeling. Towards pollutant emissions predictive tools for aeronautical chambers, Ph.D. thesis, INSA of Rouen, France, 2010. M. Ihme, H. Pitsch, Phys. Fluids 20 (2008) 055110. A. Ketelheun, C. Olbricht, F. Hahn, J. Janicka, Proc. Combust. Inst. 33 (2) (2011) 2975–2982. B. Zoller, J. Allegrini, U. Maas, P. Jenny, Combust. Flame 158 (8) (2011) 1591– 1601. P. Vervisch, O. Colin, J. Michel, N. Darabiha, Combust. Flame 158 (8) (2011) 1480–1490. V. Knopp, A. Nicolle, O. Colin, Proc. Combust. Inst. 34 (1) (2013) 667–675. F. Biagioli, F. Güthe, Combust. Flame 151 (1–2) (2007) 274–288. D. Goodwin, Chem. Vapor Deposition XVI EUROCVD 14 (2003) 2003–2008. O. Gicquel, N. Darabiha, D. Thévenin, Proc. Combust. Inst. 28 (2) (2000) 1901– 1908. J.A. van Oijen, F.A. Lammers, L.P.H. de Goey, Combust. Sci. Technol. 127 (2001) 2124–2134. D. Schmidt, T. Blasenbrey, U. Maas, Combust. Theory Modell. 2 (2) (1998) 135– 152. P. Domingo, L. Vervisch, D. Veynante, Combust. Flame 152 (3) (2008) 415–432. V. Moureau, P. Domingo, L. Vervisch, Combust. Flame 158 (7) (2011) 1340– 1357. J. Frank, R. Barlow, Symp. (Int.) Combust. 27 (1) (1998) 759–766 (Twentyseventh Symposium (International) on Combustion Volume One). M. Ihme, H. Pitsch, Combust. Flame 155 (1) (2008) 70–89. A. Garmory, E. Mastorakos, Proc. Combust. Inst. 33 (1) (2011) 1673–1680. V. Moureau, P. Domingo, L. Vervisch, C.R. Mec. 339 (2-3) (2011) 141–148. J. Smagorinsky, Mon. Weather Rev. 91 (1963) 99–164. M. Germano, U. Piomelli, P. Moin, W. Cabot, Phys. Fluids 3 (7) (1991) 1760– 1765. B. Fiorina, O. Gicquel, L. Vervisch, S. Carpentier, N. Darabiha, Combust. Flame 140 (3) (2005) 147–160. K. Kemenov, S. Pope, Combust. Flame 158 (2) (2011) 240–254. R.R. Cao, S.B. Pope, Combust. Flame 143 (4) (2005) 450–470. H. Yamashita, M. Shimada, T. Takeno, Symp. (Int.) Combust. 26 (1) (1996) 27– 34.