L Ma

F Nmira

J L Consalvi

Exploring subgrid-scale variance models in LES of lab-scale methane fire plumes

Keywords: Lab-scale fire plume, subgrid-scale variance, flamelet model, large eddy simulation

Large Eddy Simulation (LES) of the lab-scale methane fire plumes investigated experimentally by McCaffrey are performed using the steady laminar flamelet/presumed beta filtered density function model on grids of different resolution ranging from the Taylor length scale to about six times the Kolmogorov length scale. This work focuses on investigating existing subgrid (SGS) mixing models for mixture fraction variance prediction. Three different models based on the local equilibrium assumption, the variance transport equation (VTE) and the second moment transport equation (STE) are assessed. In the non-equilibrium modeling (VTE and STE), the scalar dissipation rate is modeled with an algebraic expression involving a SGS mixing time-scale. The comparison of the solutions is based on the convergence properties of LES statistics for mixture fraction, temperature and axial velocity with respect to the filter width. The simulations show that the equilibrium algebraic model is not suitable for purely buoyant flows. On the other hand, simulations performed with the transport models show that grids coarser than 1 cm cannot resolved adequately the natural laminar instability near the edge of the plume that governs the formation of large scale vortex and, therefore, underestimate the mixing process, especially in the lower part of the continuous flame. For grid resolutions finer than 1 cm, the STE model is less sensitive to grid refinement than the VTE formulation and differences between the two models are reduced with grid refinement. The STE model predicts also a stronger mixing, resulting in a slightly larger lateral expansion of the fire plume. Predicted solutions by the two models are in quantitative agreement with the experimental data in terms of axial temperature, velocity and temperature fluctuations.

Introduction

The purely buoyant fire plume is an important canonical problem in fire science since its dynamics governs problems related to fire detection, fire heating of building structures, smoke-filling rates, and fire venting. This has motivated a significant amount of experimental studies over the years for different scales and fuels [START_REF] Mccaffrey | Purely buoyant diffusion flames some experimental results[END_REF][START_REF] Cox | A study of the deterministic properties of unbounded fire plumes[END_REF][START_REF] Cox | Some stochastic properties of fire plumes[END_REF][START_REF] Cox | Some source-dependent effects of unbounded fires[END_REF][START_REF] Crauford | Experimental and numerical simulation of a buoyant fire[END_REF][START_REF] Fischer | The structure and radiation of an ethanol pool fire[END_REF][START_REF] Weckman | Experimental investigation of the turbulence structure of medium-scale methanol pool fires[END_REF][START_REF] Koseki | Large scale pool fires: results of recent experiments[END_REF][START_REF] Tieszen | Experimental study of the effect of fuel mass flux on a 1-m-diameter methane fire and comparison with a hydrogen fire[END_REF][START_REF] Blanchat | Well-characterized open pool experiment data and analysis for model validation and development[END_REF][START_REF] Blanchat | Hydrocarbon characterization experiments in fully turbulent fires -results and data analysis[END_REF]. These experimental data were widely used to validate numerical studies.

Large Eddy Simulations (LES) has become the standard for numerical modeling of turbulent fire over the last twenty years [START_REF] Branley | Large eddy simulation of a turbulent non-premixed flame[END_REF][START_REF] Ma | Numerical simulation of axi-symmetric fire plumes: accuracy and limitations[END_REF][START_REF] Desjardin | Modeling of conditional dissipation rate for flamelet models with application to large eddy simulation of fire plumes[END_REF][START_REF] Xi | Fire dynamics simulation of a turbulent buoyant flame using a mixture-fraction-based combustion model[END_REF][START_REF] Cheung | A fully-coupled simulation of vortical structures in a largescale buoyant pool fire[END_REF][START_REF] Wang | Large eddy simulation of fire plumes[END_REF][START_REF] Chatterjee | Application of a subgrid soot radiation model in the numerical simulation of a heptane pool fire[END_REF][START_REF] Maragkos | Advances in modelling in cfd simulations of turbulent gaseous pool fires[END_REF][START_REF] Maragkos | Towards predictive simulations of gaseous pool fires[END_REF][START_REF] Sikic | Non-grey radiative heat transfer modelling in les-cfd simulated methanol pool fires[END_REF] where the large scales of the flow are resolved while the small scales are modeled. Since chemical reactions in non-premixed combustion occur only by molecular mixing of fuel and oxidizer, which in practical applications occurs only on the dissipative turbulent scales, the combustion process occurs essentially at the smallest scales of the sub-filter level, and has to be modeled entirely. A significant part of these studies used the eddy dissipation concept (EDC) proposed by [START_REF] Ertesvåg | The eddy dissipation turbulence energy cascade model[END_REF], to model the turbulence-combustion interaction in different fire scenarios [START_REF] Chatterjee | Application of a subgrid soot radiation model in the numerical simulation of a heptane pool fire[END_REF][START_REF] Maragkos | Advances in modelling in cfd simulations of turbulent gaseous pool fires[END_REF][START_REF] Maragkos | Towards predictive simulations of gaseous pool fires[END_REF][START_REF] Sikic | Non-grey radiative heat transfer modelling in les-cfd simulated methanol pool fires[END_REF]. In EDC, the reactions are assumed to occur only in the smallest eddies along the turbulence energy cascade, which are called fine structures. To account for detailed chemical kinetics, these fine structures are typically treated as a Perfectly Stirred Reactors (PSRs) in which chemical reactions depend on the molecular mixing between the reactants. In fire research, the detailed kinetics is generally not retained and an infinitely fast chemistry single-step irreversible reaction is assumed. Due to nonlinearity of the combustion process, the filtered reaction rate is modeled by introducing model constants. However, no universal values of these constants exist [START_REF] Panjwani | Turbulence combustion closure model based on the eddy dissipation concept for large eddy simulation[END_REF][START_REF] Bösenhofer | The eddy dissipation concept-analysis of different fine structure treatments for classical combustion[END_REF]. On the other hand, only few LES of fire plumes have reported the use of mixture fraction based combustion models. Noticeable examples are the works of Desjardin and Cheung and Yeoh [START_REF] Desjardin | Modeling of conditional dissipation rate for flamelet models with application to large eddy simulation of fire plumes[END_REF][START_REF] Cheung | A fully-coupled simulation of vortical structures in a largescale buoyant pool fire[END_REF] where a one-step reaction based steady laminar flamelet model (SLFM) was considered and that of Wang et al. [START_REF] Wang | Large eddy simulation of fire plumes[END_REF] who used Burke-Shumann state relationships. These studies considered a classical presumed Filtered Density Function (FDF), parametrized by the filtered mixture fraction and SGS mixture fraction variance, to characterize the mixture fraction statistics at the subgrid scale.

Mixture fraction based combustion models, including flamelet models and conditional moment closure [START_REF] Echekki | Turbulent Combustion Modeling. Advances, New Trends and Perspectives[END_REF], are state-of-the art turbulent combustion models and their development for LES of non-premixed turbulent diffusion flames has retained a considerable attention over the last twenty years. In these approaches, the effects of small-scale mixing processes are typically characterized by subgrid mixture fraction variance and its dissipation rates. The subgrid scalar variance quantifies the level of unmixedness at the subgrid scales and the filtered scalar dissipation rate describes the rate at which scalars relax toward a fully mixed state. The modeling of both quantities has received significant attention over the last years and several models using different assumptions can be found in the literature. Cook and Riley [START_REF] Cook | A subgrid model for equilibrium chemistry in turbulent flows[END_REF] proposed a model based on the scale similarity assumption for the subfilter variance. Jiménez et al. [START_REF] Jimenez | A priori testing of subgrid models for chemically reacting non-premixed turbulent shear flows[END_REF] showed that the model constant depends on the exponent of the scalar spectrum function in the large Reynolds number limit. Cook [START_REF] Cook | Determination of the constant coefficient in scale similarity models of turbulence[END_REF] derived a method for calculating the coefficient and showed that it varies with the test filter scale and the Reynolds number. Another approach to model the subfilter variance was developed by Pierce and Moin [START_REF] Pierce | Dynamic model for subgrid-scale variance and dissipation rate of a conserved scalar[END_REF] who assumes a local equilibrium between the production and dissipation of subfilter variance. The model coefficient is calculated by a dynamic procedure, thus eliminating the need to specify it a priori. The use of this equilibrium model leads to a simple model for the subfilter scalar dissipation rate. Balarac et al. [START_REF] Balarac | Development of a dynamic model for the subfilter scalar variance using the concept of optimal estimators[END_REF] modified the dynamic procedure of Pierce and Moin to include certain leading-order terms in the Taylor series expansion of the left hand side of the dynamic closure. Although being computationally efficient, it was reported that the dynamic local equilibrium-based models produce inaccurate estimation of scalar mixing and hence scalar dissipation rate [START_REF] Kaul | Numerical errors in the computation of subfilter scalar variance in large eddy simulations[END_REF][START_REF] Kaul | A posteriori analysis of numerical errors in subfilter scalar variance modeling for large eddy simulation[END_REF][START_REF] Knudsen | Modeling scalar dissipation and scalar variance in large eddy simulation: Algebraic and transport equation closures[END_REF][START_REF] Kaul | Large eddy simulation of a lifted ethylene flame using a dynamic nonequilibrium model for subfilter scalar variance and dissipation rate[END_REF].

An alternative to equilibrium models is to consider a transport equation for the subgrid scalar variance. This approach, initiated by Jiménez et al. [START_REF] Jiménez | Subgrid scale variance and dissipation of a scalar field in large eddy simulations[END_REF], will be referred to as the VTE model hereafter. Kaul et al. [31,[START_REF] Kaul | A posteriori analysis of numerical errors in subfilter scalar variance modeling for large eddy simulation[END_REF][START_REF] Kaul | Large eddy simulation of a lifted ethylene flame using a dynamic nonequilibrium model for subfilter scalar variance and dissipation rate[END_REF] proposed to solve the transport equation for the filtered second moment of mixture fraction (STE), referred to as the STE model hereafter, instead of directly solving the VTE. After the STE is solved, the subgrid mixture fraction variance is then computed by subtracting the square of filtered mixture fraction from the filtered second moment of mixture fraction. The two models are identical at the level of their continuous equations, but give different results when discretized. In Ref. [START_REF] Kaul | Numerical errors in the computation of subfilter scalar variance in large eddy simulations[END_REF], Kaul et al. claimed that the calculation of production term of VTE introduces discrepancies due to the chain-rule of product operation at the discrete level, especially in high Reynolds number flows. The STE model does not consider this term and thus is not affected by this discrete chain-rule caused error. They also mentioned that another advantage of STE is to recover the maximum subgrid variance in the case where dissipation is negligible. Based on these reasons, they argued that the use of the STE model is more accurate than the VTE model. Kemenov et al. [START_REF] Kemenov | Modelling effects of subgrid-scale mixture fraction variance in les of a piloted diffusion flame[END_REF] and Jain et al. [START_REF] Jain | On the non-equilibrium models for subfilter scalar variance in large eddy simulation of turbulent mixing and combustion[END_REF] also addressed the differences between both models at the discrete level. Kemenov et al. [START_REF] Kemenov | Modelling effects of subgrid-scale mixture fraction variance in les of a piloted diffusion flame[END_REF] performed a comparative grid refinement study for the two models on the Sandia Flame D. The VTE model was found to provide a better convergence behavior with respect to grid refinement than the STE model. The discrepancy of the STE model was attributed to the conservation issue of the square of the filtered mixture fraction at the discrete level, leading to the appearance of an extra source term, depending on grid spacing, in the subgrid scalar variance evolution for the STE model. The VTE model is inherently consistent and capable to enforce the conservation law for the square of filtered mixture fraction. Jain et al. [START_REF] Jain | On the non-equilibrium models for subfilter scalar variance in large eddy simulation of turbulent mixing and combustion[END_REF] demonstrated in a planar jet flame study, that the VTE model predictions agree better with the DNS data, and the STE model may cause sharp oscillations of subgrid scalar variance.

Both two models require modeling for the scalar dissipation rate. Closure may be achieved through algebraic expressions using a turbulent time-scale [START_REF] Knudsen | Modeling scalar dissipation and scalar variance in large eddy simulation: Algebraic and transport equation closures[END_REF][START_REF] Jiménez | Subgrid scale variance and dissipation of a scalar field in large eddy simulations[END_REF]. Alternatively, closure can also be performed by solving a transport equation for the filtered square gradient of mixture fraction [START_REF] Knudsen | Modeling scalar dissipation and scalar variance in large eddy simulation: Algebraic and transport equation closures[END_REF]. However, studies of these models for subgrid scalar variance and scalar dissipation rate were limited to momentum-driven turbulent jet flames and investigations of these models on purely buoyant fire plumes were not reported to date.

In this context, the article aims to study different models for the subgrid variance and scalar dissipation for LES of methane fire plumes. The non-adiabatic steady laminar flamelet model will be used as combustion model and local equilibrium-based algebraic models and non-equilibrium models involving either the VTE or STE formulations will be compared. Following Kemenov et al. [START_REF] Kemenov | Modelling effects of subgrid-scale mixture fraction variance in les of a piloted diffusion flame[END_REF] the comparison of LES solutions will be based on the convergence properties of their statistics with respect to the filter size, ∆. As discussed by Pope and co-workers [START_REF] Kemenov | Modelling effects of subgrid-scale mixture fraction variance in les of a piloted diffusion flame[END_REF][START_REF] Pope | Ten questions concerning the large-eddy simulations of turbulent flows[END_REF], characterizing the dependence of predicted statistics on ∆ must be part of the overall LES methodology and gain importance in combustion problems since the rate-controlling processes such as reactant mixing and chemical reactions are modeled at the subgrid-scale level. Such convergence studies are rare in the LES modeling of fire plumes and were limited to mean temperature and/or axial velocity [START_REF] Sikic | Non-grey radiative heat transfer modelling in les-cfd simulated methanol pool fires[END_REF], [START_REF] Koo | Model sensitivities in les predictions of buoyant methane fire plumes[END_REF]. The present study reports a more exhaustive convergence study involving mean and rms of mixture fraction, temperature and velocity. The methane fire plumes investigated experimentally by McCaffrey are considered for this exercise. These purely buoyant diffusion flames belong to the list of target flames identified by the workshop MacFP [START_REF] Brown | Proceedings of the first workshop organized by the iafss working group on measurement and computation of fire phenomena (macfp)[END_REF] recently initiated by members of the fire community. The goals of this workshop are to develop and validate predictive models of fire plume, to identify well-defined target flames that are suitable for modeling, and to archive detailed data sets for these target flames.

This article is organized as follows. The second section presents the numerical methods and LES fire modeling. Section 3 describes the experimental case and computational setup. The results are discussed in section 4. Finally, the concluding remarks are drawn.

Numerical and physical models

Governing equations

LES is based on a separation of scales. This separation is commonly introduced by a filtering operation which decomposes the velocity and scalar fields into a resolved (filtered) part and an unresolved (subgrid-scale) part. In addition, for variable density flow, the resolved quantities are density-weighted (Favre filtered) as defined by ρ φ = ρφ. The Favre-filtered Navier-Stokes equations in low Mach number formulation supplemented with transport equations for the filtered enthalpy, h and mixture fraction, Z, are:

∂ ρ ∂t + ∂ ρũ i ∂x i = 0 (1)
∂ ρũ j ∂t + ∂ ρũ i ũj ∂x i = - ∂ p ∂x j + ∂ ∂x i (μ + µ t ) Sij + (ρ -ρ ∞ )g j (2)
∂ ρ Z ∂t + ∂ ρ ũi Z ∂x i = ∂ ∂x i ρ( D + D t ) ∂ Z ∂x i (3)
∂ ρh ∂t + ∂ ρ ũi h ∂x i = ∂ ∂x i ρ( D + D t ) ∂ h ∂x i -∇ • q R (4)
The subgrid (SGS) contribution to the momentum stress and scalar flux is computed using a dynamic Smagorinsky model and a dynamic eddy diffusivity model, respectively [START_REF] Moin | A dynamic subgrid-scale model for compressible turbulence and scalar transport[END_REF]:

τ sgs uiuj = -ρ( u i u j -ũi ũj ) ≈ 2µ t Sij (5)
τ sgs uiZ = -ρ( u i Z -ũi Z) ≈ ρD t ∂ Z ∂x i (6)
where µ t = C s ρ∆ 2 | S| is the turbulent eddy viscosity and ρD t = C z ρ∆ 2 | S| is the SGS diffusivity, with | S| = 2 Sij Sij being the norm of the resolved strain rate tensor, Sij , and ∆ the filter width. Here, we take the filter width ∆ equal to the grid-spacing. The coefficients C s and C z are calculated using dynamic procedure according to Refs. [START_REF] Moin | A dynamic subgrid-scale model for compressible turbulence and scalar transport[END_REF][START_REF] Lilly | A proposed modification of the germano subgrid-scale closure method[END_REF]. (•) and (•) represent filtered and density-weighted filtered quantities, respectively.

The thermochemical variables, such as density, ρ, molecular diffusivity, D, and viscosity, µ, are provided by the combustion model described in the next section.

Combustion model

The combustion model is based on the non-adiabatic steady laminar flamelet model [START_REF] Peters | nlaminar diffusion flamelet models in non-premixed turbulent combustio[END_REF]. State relationships for the different reactive scalars can then be generated, such as the temperature or the species mass fraction, as a function of a reduced number of parameters describing the local flow conditions (mixture fraction, strain rate or scalar dissipation rate, radiative loss). The flamelet library was generated by solving the governing equations of counterflow diffusion flames in physical space at a series of specified strain rates. Starting with a low strain rate of 0.5 s -1 , flamelet was calculated by gradually increasing the strain rate until the local quenching occurs at 260 s -1 .

The strain rate values encountered during the CFD calculations were found to range from 1.5 to 40 s -1 . The link between the physical space and the mixture fraction coordinate is given by a properly defined mixture fraction. Once the solutions of these counterflow diffusion flames are obtained in the physical space using CHEMKIN code [START_REF] Lutz | Oppdif: a fortran program for computing opposed-flow diffusion flames[END_REF], the mixture fraction distributions are obtained by solving transport equations for the mixture fraction, as suggested by Pitsch and Peters [START_REF] Pitsch | A consistent flamelet formulation for non-premixed combustion considering differential diffusion effects[END_REF]. Radiative loss was incorporated in the flamelet library by using the methodology described by Carbonell et al. [START_REF] Carbonell | Flamelet mathematical models for non-premixed laminar combustion[END_REF]. The idea consists in generating flamelet profiles for each strain rate with different degrees of heat losses. In practice, a volumetric radiative heat sink term based on the optically-thin approximation was introduced in the flamelet energy equation and was multiplied by a constant δ to allow variation in the degrees of radiative loss.

For each strain rate, the flamelet temperature and species equations were solved for a set of prescribed δ factors ranging from 0 (adiabatic) to the maximum value which was adjusted to be near quenching conditions. Figure 1 illustrates flamelet profiles for the temperature for different δ and two strain rates. The full chemical kinetic scheme developed by Qin et al. [START_REF] Qin | Combustion chemistry of propane: A case study of detailed reaction mechanism optimization[END_REF] was used.

The thermochemical states of the flamelet are then parametrized by a mixture fraction, Z, scalar disspation rate, χ, and enthalpy defect, X R = h -h ad where h and h ad are the enthalpy and adiabatic enthalpy, respectively, as:

φ = φ f l (Z, χ, X R ) (7)
In LES of non-premixed combustion, many important mixing and reaction processes occur at scales that are unresolved on the computational grid. Combustion-related variables such as the mixture fraction fluctuating quantities and their statistical distribution need to be considered. Subgrid fluctuations of the combustion variables can have an important impact on the filtered properties due to the sensitivity and nonlinearity of combustion processes. To account for subgrid fluctuations in the combustion variables, filtered combustion variables are obtained by integrating Eq. ( 7) over the joint FDF of Z, χ and X R :

φ = φ f l (Z, χ, X R ) P (Z, χ, X R )dZdχdX R (8)
Z, χ and X R are assumed to be statistically independent and the marginal FDFs are modeled by a beta distribution for Z and delta distributions for χ and X R :

P (Z, χ, X R ) = β(Z; Z, V Z )δ(χ -χ)δ(X R -XR ) (9)
Once the flamelet library is computed and assumed FDF integrals are evaluated, a look-up table can be generated to provide the filtered thermochemical variables as functions of the quantities readily available from LES ( Z, V Z , χ, and XR ) [START_REF] Pierce | Progress-variable approach for large-eddy simulation of nonpremixed turbulent combustion[END_REF]:

φ = φ( Z, V Z , χ, XR ) (10)
2.3. Subgrid scalar variance and dissipation modeling

Non-equilibrium models

The subgrid scalar variance V Z , required to retrieve solutions from the filtered flamelet look-up table, is defined in terms of a FDF as [START_REF] Jiménez | Subgrid scale variance and dissipation of a scalar field in large eddy simulations[END_REF]:

V Z (x, t) = Z 2 (x, t) -Z2 (x, t) (11)
where the appearance of the second moment Z 2 makes the variance unclosed. The subgrid scalar variance can be computed by using its transport equation (VTE) [START_REF] Jiménez | Subgrid scale variance and dissipation of a scalar field in large eddy simulations[END_REF]:

∂ ρV Z ∂t + ∂ ρũ i V Z ∂x i = ∂ ∂x i ρ D + D t ∂V Z ∂x i + P -ρ χ (12)
where the production term is defined as:

P = 2ρ D + D t ∂ Z ∂x i ∂ Z ∂x i (13)
Another way to compute the subgrid-scale variance consists in solving the transport equation for the second moment of the mixture fraction Z 2 , (STE) [START_REF] Kaul | Numerical errors in the computation of subfilter scalar variance in large eddy simulations[END_REF][START_REF] Kaul | A posteriori analysis of numerical errors in subfilter scalar variance modeling for large eddy simulation[END_REF]. Once the second moment Z 2 and Z are resolved, the subgrid scalar variance V Z thus can be computed from its definition, Eq. [START_REF] Blanchat | Hydrocarbon characterization experiments in fully turbulent fires -results and data analysis[END_REF].

∂ ρ Z 2 ∂t + ∂ ρ ũi Z 2 ∂x i = ∂ ∂x i ρ( D + D t ) ∂ Z 2 ∂x i -ρ χ (14)
In both VTE and STE models, the filtered scalar dissipation rate term, χ, is unclosed and requires modeling. Jiménez et al. [START_REF] Jiménez | Subgrid scale variance and dissipation of a scalar field in large eddy simulations[END_REF] proposed to model the filtered scalar dissipation rate by relating the subgrid-scale scalar mixing time scale to the subgrid-scale turbulent time scale. The characteristic mixing time is assumed proportional to the turbulent characteristic time. In LES, the subgrid-scale scalar mixing time can be defined as:

1 τZ = χ V Z (15)
An equivalent subgrid-scale turbulent characteristic time τ is introduced as the ratio between subgrid-scale kinetic energy, κ sgs = 1/2( u i u i -ũi ũi ), and the filtered kinetic energy dissipation rate, ε = ν ∂ui ∂xj ∂ui ∂xj . Given a proportionality between both time scales, the model for χ is derived as:

χ V Z = 1 τZ ≈ C τ = C ε κ sgs (16)
The filtered kinetic energy dissipation rate and the subgrid-scale kinetic energy are modeled as:

ε = 2 (ν + ν T ) Sij Sij ( 17
)
κ sgs = 2C I ∆ 2 Sij Sij ( 18
)
As proposed by Jiménez et al. [START_REF] Jiménez | Subgrid scale variance and dissipation of a scalar field in large eddy simulations[END_REF], C = 1/Sc, and C I can be computed by dynamic approach. Therefore, the filtered scalar dissipation rate χ can be written as [START_REF] Jiménez | Subgrid scale variance and dissipation of a scalar field in large eddy simulations[END_REF]:

χ = D + D T C I ∆ 2 V Z (19)

Algebraic equilibrium model

If production and destruction of SGS the mixture fraction variance is assumed to be in equilibrium [START_REF] Pierce | Dynamic model for subgrid-scale variance and dissipation rate of a conserved scalar[END_REF], the scalar dissipation rate can be computed as:

χ = 2 D + D t ∂ Z ∂x i ∂ Z ∂x i (20)
With the local equilibrium assumption (LEA), the SGS mixture fraction variance is computed through a scale similarity model [START_REF] Pierce | Dynamic model for subgrid-scale variance and dissipation rate of a conserved scalar[END_REF], denoted hereafter as the algebraic model:

V Z = C v ∆ 2 ∂ Z ∂x i ∂ Z ∂x i (21)
where C v is calculated dynamically.

Radiation modeling

The spectral coverage range in terms of wavenumber, η, is 0-25000 cm -1 and H 2 O and CO 2 are considered as the only radiating species since the contributions of soot and CO can be neglected. The Rank-Correlated Full-Spectrum k-distribution (RC FSK) method is used as gas radiative property model [START_REF] Solovjov | The rank correlated fsk model for prediction of gas radiation in non-uniform media, and its relationship to the rank correlated slw model[END_REF]. As in the classical FSK [START_REF] Modest | Assembly full spectrum k-distribution from a narrow band database: effects of mixing gases, gases and non-gray absorpting particles and non-gray scatters in non-gray enclosures[END_REF], the FS cumulative k-g distribution function, is defined as

g(k, φ, T p ) = ∞ 0 H[k - κ η (φ)]I bη (T p )dη/I b (T p )
, where H is the Heaviside function, κ η is the spectral absorption coefficient, φ = {x CO2 , x H2O , T } is an array of thermodynamic variables affecting κ η .

x CO2 and x H2O represent the mole fractions of CO 2 and H 2 O, respectively. I bη and I b are the spectral and total blackbody intensities at the blackbody temperature, T p , respectively. The main advantage of the RCFSK is that it does not require any specification of a reference state [START_REF] Modest | Assembly full spectrum k-distribution from a narrow band database: effects of mixing gases, gases and non-gray absorpting particles and non-gray scatters in non-gray enclosures[END_REF]. Mixed FS k-g distributions (for mixtures of H 2 O and CO 2 ) are constructed from HITEMP 2010 [START_REF] Rothman | Hitemp: the high-temperature molecular spectroscopic database[END_REF] by using the procedure proposed by Modest and Riazzi [START_REF] Modest | Assembly full spectrum k-distribution from a narrow band database: effects of mixing gases, gases and non-gray absorpting particles and non-gray scatters in non-gray enclosures[END_REF]. The FSK radiative transfer equation (RTE) can be written as:

dI g0 ds = -k * (g 0 )I g0 + k * (g 0 )a(g 0 )I b (T p ) ( 22
)
where g 0 corresponds in the present study to a quadrature-point of a 10-point Gauss-Legendre quadrature scheme and I g0 is the radiative intensity at this quadrature point [START_REF] Solovjov | The rank correlated fsk model for prediction of gas radiation in non-uniform media, and its relationship to the rank correlated slw model[END_REF]. The RCFSK scheme determines the absorption coefficient by solving g(k * , φ, T p ) = g 0 whereas the stretching function is computed as a = ∂g[k(g 0 , φ, T p ), φ, T ]/∂g 0 [START_REF] Solovjov | The rank correlated fsk model for prediction of gas radiation in non-uniform media, and its relationship to the rank correlated slw model[END_REF]. The total radiative intensity, I, and the total incident radiation, G, are computed as I = 1 0 I g0 dg 0 and G = 4π IdΩ, respectively. The divergence of the radiative flux is then calculated from the following equation:

∇ • q R = 1 0 4πk * aI b dg 0 - 0 1 k * Gdg 0 (23)
Predictions were found insensitive to the choice of T p . In the present simulations, T p was set equal to 1500K. The filtered RTE and divergence of the radiative flux are obtained by applying the filtering operation to Eqs. [START_REF] Ertesvåg | The eddy dissipation turbulence energy cascade model[END_REF][START_REF] Panjwani | Turbulence combustion closure model based on the eddy dissipation concept for large eddy simulation[END_REF]:

dI g0 ds = -k * I g0 + k * aI b (24)
∇ • q R = 1 0 4πk * aI b dg 0 - 0 1 k * Gdg 0 (25)
The filtered absorption terms, k * I g0 or k * G, are closed by neglecting the subgridscale absorption, leading to k * I g0 ≈ k * I g0 and k * G ≈ k * G. The filtered absorption coefficient and emission terms are closed by using the presumed filtered PDF approach:

k * = ρ (k * ) f l (Z, χ, χR ) ρ f l (Z, χ, χR ) β(Z; Z, V Z )dZ (26) k * aI b = ρ (k * aI b ) f l (Z, χ, χR ) ρ f l (Z, χ, χR ) β(Z; Z, V Z )dZ ( 27
)
where the superscript f l refers to the flamelet library. The Filtered RTE is solved by using the Discrete Ordinates Method with a S8 quadrature scheme [START_REF] Modest | Radiative Heat Transfer[END_REF].

Numerical solution

The Favre-filtered transport equations are solved by using the second-order iterative variable-density solver developed by Ma et al. [START_REF] Ma | Verification and validation of a variable-density solver for fire safety applications[END_REF] and implemented in the finite volume code Code Saturne v5.0.9 [START_REF]Code Saturne[END_REF], employing cell-centered collocated grids.

The second-order Crank-Nicolson scheme is used for time advancement. For scalar transport, numerical oscillations must be prevented through non-oscillatory schemes. We employ a TVD (Total Variation Diminuishing) superbee scheme that blends a second order central scheme with a first order upwind scheme in a way that combines good accuracy with limited numerical oscillation. Avoiding numerical oscillation is particularly important for reacting scalars (here the mixture fraction), as the density depends on them in a strongly non-linear manner. A second-order central difference scheme (CDS) is used for diffusion. For the momentum equation both convective and diffusive terms are also discretized by using a second-order CDS. Numerical details and a method of manufactured solutions (MMS) verification of the predictor-corrector approach can be found in Ref. [START_REF] Ma | Verification and validation of a variable-density solver for fire safety applications[END_REF].

Computational domain and boundary conditions

Five methane fire plumes in open conditions, investigated experimentally by Mc-Caffrey [START_REF] Mccaffrey | Purely buoyant diffusion flames some experimental results[END_REF], were simulated. These flames were generated from a burner of 0.3 m in diameter with heat release rates (HRR, Q) of 14.4, 21.7, 33, 45 and 57.5 kW, respectively. Simulations were performed in the rectangular domain of 3 × 3 × 3 m 3 in a Cartesian coordinate. Consistently with the experiments, the burner is lifted by 0.3 m above the floor to prevent the influence of the ground.

In order to examine the effect of the grid resolution on the LES results, five grids G1, G2, G3, G4 and G5 have been used with progressively increasing the resolution in the three directions from 0.909 to 6.32 million cells. These grids are uniform in the region 0.3 × 0.3 × 1.0 m 3 which results in a minimal grid spacing equal to 1.5, 1.25, 1.0, 0.75 and 0.5 cm, respectively (see Table 1). Outside this region, all grids are progressively stretched in all the direction. Throughout all simulations, time step is set to 5 × 10 -4 s which corresponds to an averaged maximum CFL of 0.6. Simulations were run for 25 s and the time-averaged mean and root mean square values were collected over the last 19 s. The first 6 s of simulation were used to establish a statistically stationary flow.

The convergence of statistics will be investigated for the 57.5 kW fire plume at heights above the burner of 0.1, 0.2, 0.4 and 0.6 m, respectively. The corresponding normalized heights are z * = z/ Q2/5 = 0.019, 0.039, 0.079, and 0.119 m • kW -2/5 , showing that they are located close to the burner, in the continuous flame region (CF, z * ≤ 0.08 m • kW -2/5 ), at the transition between the CF and the intermittent flame regions (IF) and in the IF region (0.08 < z * ≤ 0.2 m • kW -2/5 ), respectively [START_REF] Mccaffrey | Purely buoyant diffusion flames some experimental results[END_REF]. These locations cover all the reactive parts of the fire plume where SGS mixing and resulting combustion processes are expected to be the most important. In addition, it should be pointed out that the positions z = 0.1 and 0.2 m are located in the CF where the mixing is controlled by the formation and growth of laminar instabilities [START_REF] Tieszen | Chapter 2. Transport phenomena that affect heat transfer in fully turbulent fires[END_REF]. The Kolmogorov length scale, η k , the Taylor length scale, λ, and the diffusive layer thickness, l d , at the CF tip can be estimated from the data of McCaffrey [START_REF] Mccaffrey | Purely buoyant diffusion flames some experimental results[END_REF] and Cox and Chitty [START_REF] Cox | Some stochastic properties of fire plumes[END_REF] following the analysis proposed in Refs. [START_REF] Williamson | Developing extinction criteria for fires[END_REF][START_REF] Wu | Detailed modeling of a small-scale turbulent pool fire[END_REF]. The integral length scale, L t , is assumed to be one-half of the equivalent burner diameter [START_REF] Williamson | Developing extinction criteria for fires[END_REF][START_REF] Wu | Detailed modeling of a small-scale turbulent pool fire[END_REF], i.e. L t = 0.17 m. Based on the turbulent intensity reported by Cox and Chitty [START_REF] Cox | Some stochastic properties of fire plumes[END_REF], the fluctuating velocity, w , is assumed to be 30 % of the mean axial velocity at the CF tip, w CF . This latter is estimated from the correlation of McCaffrey [START_REF] Mccaffrey | Purely buoyant diffusion flames some experimental results[END_REF], leading to w CF = 1.93 Q1/5 = 4.33 m/s and w = 1.3 m/s, respectively. The kinematic viscosity is computed from ν = ν ∞ (1 + ∆T CF /T ∞ ) 1.7 = 187 × 10 -6 m 2 /s [START_REF] Williamson | Developing extinction criteria for fires[END_REF] where ν ∞ = 15×10 -6 m 2 /s is the kinematic viscosity at the ambient temperature, T ∞ = 293 K, and ∆T CF = 1000 K is the temperature rise at the continuous flame tip [START_REF] Cox | Some stochastic properties of fire plumes[END_REF]. The turbulent Reynolds number is then computed as Re t = w L t /ν = 1175, leading to a Kolmogorov length scale, η k = L t Re t -3/4 , of about 0.82 mm. The Taylor length scale can then be deduced from λ = √ 10L t Re t -1/2 ≈ 15 mm [START_REF] Pope | Turbulent flows[END_REF]. In addition, the diffusive layer thickness is estimated as l d = D st /χ st ≈ 14 mm, where D st and χ st are the molecular diffusivity and the scalar dissipation rate at stoichiometry, respectively. The scalar dissipation rate is determined from

χ st = a t ϕexp[-2(erf c -1 (2Z st )) 2 ] with ϕ = 3[(T st /T ∞ ) 1/2 + 1] 2 /(4[2(T st /T ∞ ) 1/2 + 1] 2 )
. T st = 2000 K and Z st = 0.055 are the temperature and the mixture fraction at the stoichiometry, respectively. The strain rate is calculated from a t = 0.28τ k [START_REF] Williamson | Developing extinction criteria for fires[END_REF] where τ k = η k w k and w k = w Re -1/4 t are the Kolmogorov time scale and velocity scales, respectively. Table 1 shows that the present LES are resolved beyond the Taylor microscale and the diffusive layer thickness for all the grids. In addition, the finest filter size, G5, is about 6 times the Kolmogorov length scale.

Concerning the boundary conditions, an inlet velocity was imposed at the burner to ensure the specified HRR for each methane flame. Fixed values of fuel mass flow rate and enthalpy flow rate are maintained according to the specified HRR. Both convective and diffusive mass and enthalpy fluxes are accounted for at the inlet. Typical entrainment boundary conditions are used for lateral sides. At the domain exit, a convective condition was used:

∂ φ ∂t + u c ∂ φ ∂n = 0 ( 28
)
where φ = {u j , Z, Z 2 , V Z , h}. Following [START_REF] Craske | Robust and accurate open boundary conditions for incompressible turbulent jets and plumes[END_REF], the convective velocity u c is given by a Gaussian profile, u c = u out 0 exp(-r 2 /b 2 ). The mean axial velocity on the plume centreline, u out 0 , and the plume radius, b, at the exit height were obtained from the Heskestad correlations [START_REF] Heskestad | Engineering relations for fire plumes[END_REF]. Homogeneous Neumann condition is used for the pressure. In the rest domain, the classical wall boundary condition is imposed.

Results and discussions

Impact of grid spacing

The performance of algebraic equilibrium, VTE and STE models are examined based on the convergence of the LES statistics with respect to the filter width, ∆, in the 57.5 kW fire plume. For a consistent LES, time-averaged values must converge, root-meansquare (RMS) resolved values must increase and the SGS variance diminishes when ∆ decreases [START_REF] Boudier | Effects of mesh resolution on large eddy simulation of reacting flows in complex geometry combustors[END_REF]. These three models differ by the calculation of the SGS mixture fraction variance. Consequently, the convergence properties of mean mixture fraction and mixture fraction variance are first analyzed [START_REF] Kemenov | Modelling effects of subgrid-scale mixture fraction variance in les of a piloted diffusion flame[END_REF]:

Z = Z + Z = Z ( 29
)
Z 2 = Z 2 -Z 2 = ( Z2 -Z 2 ) ZV,ReS + ( Z 2 -Z2 ) VZ (30)
For sufficiently fine LES, it is reasonable to expect that the residual mixture fraction, Z , is equal to 0 [START_REF] Kemenov | Modelling effects of subgrid-scale mixture fraction variance in les of a piloted diffusion flame[END_REF]. Note that this assumption has been used to express Eq. (30). The SGS variance reaches on the whole maximum values along the plume axis at all the heights (see Fig. 3). As a consequence, the dependence of mean mixture fraction, SGS mixture fraction variance and total mixture fraction variance on ∆ along the plume centerline is further investigated in Figs. 6 in order to provide more quantitative information.

Figure 2(g) and (h) and the diagram (b 1 ) in Fig. 6 show clearly that the STE model exhibits convincing convergence for mean mixture fraction on the different grids at z = 0.4 m and 0.6 m. In the CF region (z = 0.1 m and 0.2 m), the STE solutions are more grid-dependent, especially along the plume wings (see Fig. 2(e) and (f)). Nevertheless, a reasonable convergence state is reached along the plume axis for these heights on the grids G3-G5 (see Fig. 6(b 1 )). On the other hand, the convergence of VTE is more complicated whatever the height at the vicinity of the fire plume axis (see Fig. 2 (a)-(d) and 6(a 1 )). In addition, similar grid effects as for STE are observed along the plume wings in the CF region (see Figs. 2(a) and (b)). Figure 2 shows also that the VTE model predicts higher mixture fraction close to the fire plume axis (see also Fig. 6) whereas the STE model predicts slightly wider radial distributions.

The radial profiles of SGS mixture fraction variance are plotted in Fig. 3 for all the grids. In accordance with the criterion of consistent LES, V Z decreases continuously as ∆ is reduced for STE. This behavior is further highlighted in the diagram (b 2 ) of Fig. 6. The rate of decrease for VTE is in overall slower than that of STE up to G4 and accelerates between G4 and G5 as illustrated in the diagram (a 2 ) of Fig. 6. Figure 3 shows also that, for a given grid up to G4, the VTE model predicts lower V Z than the STE model close to the centreline at z = 0.1 m and z = 0.2 m (see also the diagrams (a 2 ) and (b 2 ) of Fig. 6) and in the profile wings at all heights. This behavior is consistent with the observations made in other flame configurations [START_REF] Kaul | Numerical errors in the computation of subfilter scalar variance in large eddy simulations[END_REF][START_REF] Kaul | A posteriori analysis of numerical errors in subfilter scalar variance modeling for large eddy simulation[END_REF][START_REF] Kemenov | Modelling effects of subgrid-scale mixture fraction variance in les of a piloted diffusion flame[END_REF][START_REF] Jain | On the non-equilibrium models for subfilter scalar variance in large eddy simulation of turbulent mixing and combustion[END_REF] and is attributed to numerical errors stemming from an underresolved discrete representation of the squared-gradient term in the production rate of the SGS mixture fraction variance (see Eq. 13). This production term is explicitly involved in the VTE model (see Eq. 12) but not in the STE model that does not suffer from these approximations. For the finest grid G5, a detailed examination of the results show that both models predict the similar radial profiles of V Z at all the heights. This latter result suggests that the aforementioned numerical errors affecting the VTE model are reduced when the grid is sufficiently fine and that both models tend toward similar solutions.

As expected, Fig. 4 shows that Z V,ReS increases as the grid becomes finer and finer up to G4 for both models. In addition, it can be observed that the corresponding rate of increase is higher close to the source (z = 0.1 m and 0.2 m). For G5, Z V,ReS keeps on increasing for both models at z = 0.1 m and 0.2 m whereas the shape of the radial profiles changes at z = 0.4 m and 0.6 m, resulting in lower values close to the centerline and wider radial extensions. STE reaches an approximate convergent state on G3-G5 grids for Z V,tot as shown in Fig. 5(e)-(h). This is confirmed by the diagram (b 3 ) of Fig. 6. On the other hand, Fig. 5(a)-(d) and the diagram (a 3 ) of Fig. 6 show that the convergence is less convincing for VTE whatever the height, especially close to the fire plume axis. A careful examination of the results shows that the differences between VTE and STE are substantially less important on G5 than on G4 and G3, which reinforces the feeling that both models converge toward the same solution when the grid becomes sufficiently fine. These differences between VTE and STE are the following. First at z = 0.1 m, the VTE model exhibits a more pronounced double-peak radial profile than the STE model. This behavior is also observed for V Z and Z V,ReS (see Figs. 3 and4). Second, close to the plume axis, the VTE model predicts slightly higher Z V,tot than the STE model (see the diagrams (a 3 ) and (b 3 ) of Fig. 6). Third, it appears also that the STE model predicts a slightly higher lateral spreading.

In order to complete the previous discussion, time-averaged temperature rise above the ambient, ∆T , rms of total temperature fluctuations, σ T,tot , defined in accordance with Eq. 30, time-averaged axial velocity, w , and rms of axial velocity fluctuations, σ w , are plotted for the two models on the different grids in Figs. 7 to 10, respectively.

These figures indicate that the grids G1 and G2, that have a resolution larger than 1 cm, predict a significantly narrower radial spreading of the fire plume than the grids G3-G5, especially in the CF region (z = 0.1 m and 0.2 m). These mesh effects are particularly remarkable for ∆T , σ T,tot and σ w whereas they are less pronounced for w . In addition, the radial extensions of w and σ w are insensitive to further grid refinement from G3 to G5 whereas those of ∆T and σ T,tot keep on being enhanced in the CF region. This latter behavior is less pronounced for STE than for VTE. The mixing in the lower part of the CF region is driven by the formation and growth of initially laminar non-dissipative instabilities in the near source region of the flow [START_REF] Tieszen | Chapter 2. Transport phenomena that affect heat transfer in fully turbulent fires[END_REF]. These instabilities growth by vortex dynamic to become nonlinear and energy bearing. This process is illustrated in Fig. 11 that shows instantaneous snapshots of the density field in the centerline x-z plane over a sequence of 0.4 s. The formation and the growth of initially laminar bubble structures is clearly observed for the grids G3 to G5 whereas G1 and G2 are not sufficiently fine to capture these details. In addition, this figure shows that these instabilities become better resolved as the grid is refined from G3 to G5. As discussed in Ref. [START_REF] Tieszen | Toward the development of a LES-SGS closure model for buoyant plumes[END_REF], the dynamic SGS turbulence model considered in the present study is dissipative in nature and is not capturing the mixing resulting from these structures. As a consequence, the present results evidence that grid resolution has a profound impact on the dynamic of the flow in the CF region and suggest that grid resolution coarser than 1 cm are not able to describe adequately the formation of the initial base instability near the edge of the fire plume. These observations are in accordance with those of Tieszen et al. [START_REF] Tieszen | Toward the development of a LES-SGS closure model for buoyant plumes[END_REF] for 1m diameter helium plume, who also showed that too coarse grids cannot capture the formation of the instabilities and their subsequent growth.

In addition, Figs. 7 to 10 emphasize the main observations made previously. i) STE exhibits a better convergence behavior than VTE on G3 to G5 for the temperature statistics (see Figs. 7 and8) . On the other hand, both models exhibits a reasonable convergence on G3-G5 for the axial velocity statistics (see Figs. 9 and 10). ii) A careful examination of Figs. 7 to 10 shows that the differences between VTE and STE are significantly more pronounced on G3 and G4 than on G5 where both STE and VTE predictions tend to become similar. This behavior is further highlighted in Figs. S1 to S4 in the Supplementary Materials that compare VTE and STE solutions on G4 and G5. The observed differences are consistent with those described for mean mixture fraction and mixture fraction variances and are evidenced by the comparison of VTE and STE solutions on G4 reported in the Figs. S1 to S4 of the Supplementary Materials. STE exhibits in overall a wider fire plume spreading at all heights. In addition, as for mixture fraction variances, the double peak radial profiles, observed for ∆T , σ T,tot and σ w at z = 0.1 m and resulting from the fact that close to the burner the high-temperature combusting sheets occur predominately close to the plume edge, are more pronounced for VTE than for STE. Finally, it can be also observed that σ T,tot is higher for STE than for VTE (see Fig. S2).

The differences between STE and VTE originate from different formulations at the discrete level of the transport equation of the SGS mixture fraction variance. Numerical errors stemming from finite differences lead to an underestimation of the SGS variance production term, which appears explicitly in the VTE formulation as discussed in [START_REF] Kaul | Numerical errors in the computation of subfilter scalar variance in large eddy simulations[END_REF][START_REF] Kaul | A posteriori analysis of numerical errors in subfilter scalar variance modeling for large eddy simulation[END_REF]. For the present fire plume, as discussed previously, this results in SGS variance lower for VTE than for STE in the CF region (z = 0.1 m and z = 0.2 m), which, in turn, affects the combustion process that occurs at the SGS level. In particular, in the CF region (z = 0.1 m and z = 0.2 m), the temperature and, as a result, the temperature gradient that drives the formation and the growth of the laminar instabilities in the near source region of the flow, are higher for STE than for VTE. This is evidenced in Fig. S1 of the Supplementary Materials that compares the temperature rise above the ambient for STE and VTE on G4. Therefore, the resolved-scale mixing is stronger for STE than for VTE. This is illustrated in the Fig. S5 in the Supplementary Materials that compares the resolved temperature fluctuation predicted by STE and VTE on G4. This explains that STE predicts in overall wider radial profiles and lower centerline mixture fraction. This explains also why the double peak radial profiles observed for V Z , Z V,ReS , Z V,tot , ∆T , σ T,tot and σ w at z = 0.1 m are less pronounced for STE than for VTE. In addition, these interactions between SGS mixing and resolved-scale mixing in the CF region through the combustion process may explain that the differences between VTE and STE are higher in the present configuration than for the SANDIA flame D where, as stated by Kemenov et al. [START_REF] Kemenov | Modelling effects of subgrid-scale mixture fraction variance in les of a piloted diffusion flame[END_REF], the flow is relatively insensitive to the incurred errors on the SGS variance.

Table 2 shows the puffing frequency, f , predicted by both VTE and STE for the different grids. These predictions can be compared with the correlation of Cetegen and Ahmed [START_REF] Cetegen | Experiments on the periodic instability of buoyant plumes and pool fires[END_REF], i.e. f = 1.5 1/D eq = 2.58 Hz where D eq = 0.34 m is the equivalent burner diameter. It can be observed that the puffing frequencies computed with STE and VTE converge toward 2.5 Hz on G3-G5.

The dependence of the radiant fraction and the optical thickness, defined as the part of the radiative emission reabsorbed within the flame [START_REF] Nmira | Assessment of subfilter-scale turbulence-radiation interactions in non-luminous pool fires[END_REF], on ∆ are reported in Table 3. The radiant fraction computed with STE is nearly constant on G3-G5 with a value around 0.27 whereas that computed with VTE is more sensitive to the grid. It can be also observed that the differences between the two models are reduced as the grid is refined from G3. The flame optical thickness is less sensitive to the grid resolution than the radiant fraction and is around 0.47-0.48 on G3-G5 for the two models. These results show clearly that this fire plume cannot be assumed as optically-thin.

Impact of grid spacing on the algebraic model

The convergence behavior of the algebraic model is investigated in Fig. 6(c). Figure 6c 1 shows that mean mixture fraction exhibits a good convergence for all the grids at z = 0.4 m and z = 0.6 m. On the other hand, the convergence is not established in the CF region (z = 0.1 m and z = 0.2 m) where the mean mixture fraction increases continuously as the turbulence resolution scale is reduced (see Fig. 6(c 1 )). The SGS mixture fraction variance does not behave as expected for consistent LES along the plume axis at z = 0.1 m, z = 0.2 m and z = 0.6 m where it increases as the grid is refined up to G3 before reaching a nearly constant value as ∆ is further reduced (see Fig. 6(c 2 )). The mean total variance increases continuously as the grid size is refined and, therefore, does not convergence whatever the grid resolution more particularly in the CF region (see Fig. 6(c 3 )). In addition, Table 2 shows that the convergence is not established for the puffing frequency whatever the grid.

It can be also observed by comparing the diagram (c 2 ) and the diagrams (a 2 ) and (b 2 ) of Fig. 6 that, whatever the grid, the algebraic model predicts substantially lower centerline SGS scalar variance at z = 0.1 m and z = 0.2 m than both STE and VTE models. This underprediction is in line with observations made in other simpler flames [START_REF] Kaul | Numerical errors in the computation of subfilter scalar variance in large eddy simulations[END_REF][START_REF] Knudsen | Modeling scalar dissipation and scalar variance in large eddy simulation: Algebraic and transport equation closures[END_REF][START_REF] Jiménez | Subgrid scale variance and dissipation of a scalar field in large eddy simulations[END_REF][START_REF] Jain | On the non-equilibrium models for subfilter scalar variance in large eddy simulation of turbulent mixing and combustion[END_REF] and, for the present fire plume configuration, can be related to the highly non-equilibrium nature of buoyancy-driven flows as described by Bakosi and Ristorcelli [START_REF] Bakosi | Exploring the beta-pdf in variable-density turbulent mixing[END_REF]. Figure 12 compares the time-averaged production and dissipation rates of the SGS variance. It indicates clearly that the local equilibrium assumption is not satisfied in the CF region although it seems to become more valid as the height increases in the IF region and in the plume. These results evidence the limitation of algebraic equilibrium based models for the simulation of fire related flows.

Comparisons with experimental data

Figures 13 shows the axial profiles of mean axial velocity (first line, index 1), mean temperature rise (second line, index 2), resolved temperature fluctuation (third line, index 3) and total temperature fluctuation (fourth line, index 4) as a function of the normalized height for the five plumes investigated experimentally by McCaffrey [START_REF] Mccaffrey | Purely buoyant diffusion flames some experimental results[END_REF]. Model predictions obtained on G3 with VTE (first column, diagrams a), STE (second column, diagrams b) and the algebraic model (third column, diagrams c) are compared with experimental data.

McCaffrey [START_REF] Mccaffrey | Purely buoyant diffusion flames some experimental results[END_REF] found that scaling the height and axial velocity by Q2/5 and Q1/5 , respectively, allows to collapse the experimental data for mean temperature rise above the ambient and mean axial velocity on a single curve independent of the HRR. Let us start to discuss the solutions obtained with VTE and STE (Figs. 13(a 1 ), (a 2 ), (b 1 ) and (b 2 )). Both predictions are similar and follow the same scaling as the experiments. For velocity, the deviation of these models from the correlation near the burner surface shows a similar trend and magnitude as the deviation of the measurement points to the correlation. The temperatures reported by McCaffrey were not corrected for thermocouple radiation which explains the over-prediction by both VTE and STE models, especially in the CF and IF regions. A much better agreement is observed with the correlation of Cox and Chitty, obtained also from 30 cm diameter methane fire plumes of different HRR [START_REF] Cox | A study of the deterministic properties of unbounded fire plumes[END_REF] and corrected for thermocouple radiation [START_REF] Cox | Some source-dependent effects of unbounded fires[END_REF], and the data reported by Crauford et al. [START_REF] Crauford | Experimental and numerical simulation of a buoyant fire[END_REF] for a 28 kW methane fire plume generated by circular burner with a diameter of 25 cm. The algebraic model, meanwhile, provides a less satisfactory agreement with the experimental correlations and data than the STE and VTE models. First, it introduces more scatter for the data corresponding to the different HRR than observed experimentally and predicted by both VTE and STE. Second, it overestimates the temperature and the axial velocity in IF region.

The axial resolved and total temperature fluctuations are compared with the data of Cox and Chitty [START_REF] Cox | Some stochastic properties of fire plumes[END_REF] and Crauford et al. [START_REF] Crauford | Experimental and numerical simulation of a buoyant fire[END_REF] (see the last two lines of Fig. 13). Cox and Chitty [START_REF] Cox | Some stochastic properties of fire plumes[END_REF] measured temperature fluctuations by using uncompensated 13 µm thermocouples. An analysis performed by these authors showed that these probes can adequately follow signals with frequency up to 10 Hz. Crauford et al. [START_REF] Crauford | Experimental and numerical simulation of a buoyant fire[END_REF] also measured temperature fluctuations by using 50 µm thermocouple with a compensation technique based on a thermocouple time constant of 30 ms thorough the flame zone. They suggested that their data may significantly underestimate the temperature fluctuations. The data of Crauford et al. [START_REF] Crauford | Experimental and numerical simulation of a buoyant fire[END_REF] follow quantitatively the same trends as those of Cox and Chitty [START_REF] Cox | Some stochastic properties of fire plumes[END_REF]. Both VTE and STE predictions reproduce correctly similar scaling and trends as the experiments (see Figs. 13(a 3 ), (b 3 ), (a 4 ) and (b 4 )). The total rms of temperature fluctuation overestimate the experimental data which can be attributed, as discussed above, to the measurement techniques (see Fig. 13(a 4 ) and (b 4 )). Nevertheless, it is interesting to note that the resolved temperature fluctuations predicted by both VTE and STE models are in better agreement with the experimental data (see Fig. 13(a 3 ) and (b 3 )). As for mean temperature and axial velocity, the algebraic model introduces more scatter for the different HRR than observed experimentally and predicted by both VTE and STE, especially in the CF region (see Figs. 13(c 3 ) and (c 4 )). In addition, as expected from the discussion in section 3.1.2, the resolved and total fluctuating temperatures are significantly underestimated in the CF region, especially for HRRs higher than 33 kW (see Figs. 13(c 3 ) and (c 4 )).

Figure 14 compares VTE and STE radial profiles of mean temperature on G3 and G4 with data extracted from the experimental results of McCaffrey for the 45 kW fire plume at different heights. The correlations of Cox and Chitty [START_REF] Cox | A study of the deterministic properties of unbounded fire plumes[END_REF] for a 47 kW fire plume are also reported. The two sets of experimental data are in good accordance at heights of 0.1, 0.5 and 0.6 m but present large discrepancies for heights of 0.2 to 0.4 m. On G3, the largest discrepancies between the model and the experimental data are observed at 0.2 and 0.3 m that correspond to regions where the intermittency is the largest. At these heights, the model underestimates significantly the measurements for r > 0.05 m and the maximum deviations between the calculations and the data of McCaffrey are of about 200K at z = 0.2 m and r ≈ 0.05 m. As the grid is refined to G4, a much better agreement is observed at these heights. As expected the STE predicts on the whole slightly wider radial profiles than the VTE.

Conclusions

LES of lab-scale methane fire plumes were performed using the SLF model coupled to a presumed beta FDF approach. Three modeling approaches of the SGS variance, namely an algebraic dynamic model and two non-equilibrium transport models based on either VTE or STE formulation, were analyzed. These models were assessed based on the convergence properties of their statistics for mixture fraction, temperature and axial velocity with respect to the filter width ranging from the Taylor length scale to about six times the Kolmogorov length scale. The following conclusions can be drawn from the present study:

(1) The simulations confirm the limitations of equilibrium algebraic model for purely buoyant flows. (2) The simulations with non-equilibrium transport models show that computational grids with resolutions coarser than about 1 cm cannot adequately resolve the initial base instability near the edge of the fire plume that plays an important role to describe accurately the mixing process, especially in the continuous flame. (3) For resolutions finer than about 1 cm, the STE model is less sensitive to grid refinement than the VTE formulation and predicts a higher mixing that results in a slightly larger lateral expansion of the fire plume. The differences between the two models tend to be reduced with the turbulent resolution scale and comparable predictions are obtained on the finest grid. (4) For grids finer than 1 cm, both VTE and STE are in quantitative agreement with available experimental data in terms of axial velocity, temperature and rms of temperature fluctuations.

List of figure captions

Figure 1. Temperature profiles as a function of the mixture fraction for different radiative losses and two strain rates of a=2 s -1 and a=38 s -1 . -

Figure 2 .

 2 Figure 2. Radial distribution of time-averaged mixture fraction with VTE model (a)-(d) and STE model (e)-(h) at different heights for the 57.5 kW fire plume.

Figure 3 .

 3 Figure 3. Radial distribution of subgrid mixture fraction variance with VTE model (a)-(d) and STE model (e)-(f) at different heights for the 57.5 kW fire plume.

Figure 4 .

 4 Figure 4. Radial distribution of resolved mixture fraction variance with VTE model (a)-(d) and STE model (e)-(h) at different heights for the 57.5 kW fire plume.

Figure 5 .

 5 Figure 5. Radial distribution of total mixture fraction variance with VTE model (a)-(d) and STE model (e)-(h) at different heights for the 57.5 kW fire plume.

Figure 6 .

 6 Figure 6. Convergence of mean mixture fraction (index 1), subgrid mixture fraction variance (index 2) and total mixture fraction variance (index 3) with grid spacing along the plume axis for: (a) VTE, (b) STE and (c) the algebraic model.

Figure 7 .

 7 Figure 7. Radial distribution of time-averaged temperature with VTE model (a)-(d) and STE model (e)-(h) at different heights for the 57.5 kW fire plume.

Figure 8 .

 8 Figure 8. Radial distribution of total rms temperature with VTE model (a)-(d) and STE model (e)-(h) at different heights for the 57.5 kW fire plume.

Figure 9 .

 9 Figure 9. Radial distribution of time-averaged axial velocity with VTE model (a)-(d) and STE model (e)-(h) at different heights for the 57.5 kW fire plume.

Figure 10 .

 10 Figure 10. Radial distribution of rms resolved axial velocity with VTE model (a)-(d) and STE model (e)-(h) at different heights for the 57.5 kW fire plume.

Figure 11 .

 11 Figure 11. Instantaneous snapshots of the density field, computed with the STE model on the different grids, in the center-line (x-z) plane over a sequence of 0.4 s for the 57.5 kW fire plume.

Figure 12 .

 12 Figure 12. Axial evolution of the time-averaged production and dissipation rates of the subgrid mixture fraction variance computed by VTE on the grid G5 for the 57.5 kW fire plume.

Figure 13 .

 13 Figure 13. Axial distribution of axial velocity (index 1), mean temperature (index 2), resolved temperature fluctuation (index 3) and total temperature fluctuation (index 4) for the five fire plumes as a function of the normalized height with: (a) VTE, (b) STE and (c) the algebraic model.

Figure 14 .Figure 1 .Figure 2 .

 1412 Figure 14. Radial distribution of mean temperature at different heights for the 45 kW fire plume.

Figure 4 .

 4 Figure 4. Radial distribution of resolved mixture fraction variance with VTE model (a)-(d) and STE model at different heights for the 57.5 kW fire plume.

Figure 5 .Figure 6 .Figure 9 . 32 -

 56932 Figure 5. Radial distribution of total mixture fraction variance with VTE model (a)-(d) and STE model (e)-(h) at different heights for the 57.5 kW fire plume.

Figure 10 .

 10 Figure 10. Radial distribution of rms resolved axial velocity with VTE model (a)-(d) and STE model (e)-(h) different heights for the 57.5 kW fire plume.

Figure 11 .

 11 Figure 11. Instantaneous snapshots of the density field, computed with the STE model on the different grids, in the center-line (x-z) plane over a sequence of 0.4 s for the 57.5 kW fire plume.

Figure 12 .

 12 Figure 12. Axial evolution of the time-averaged production and dissipation rates of the subgrid mixture fraction variance computed by VTE on the grid G5 for the 57.5 kW fire plume.

Figure 13 .

 13 Figure 13. Axial distribution of axial velocity (index 1), mean temperature (index 2), resolved temperature fluctuation (index 3) and total temperature fluctuation (index 4) for the five fire plumes as a function of the normalized height with: (a) VTE, (b) STE and (c) the algebraic model.

  3.1.1. Impact of grid spacing on VTE and STE Let us start by investigating the effects of the turbulent resolution scale ∆ on VTE and STE. Radial profiles of time-averaged mixture fraction, Z , SGS variance, V Z , resolved-scale variance, Z V,ReS , and total variance, Z V,tot , at different heights of 0.1, 0.2, 0.4 and 0.6 m, computed with the VTE and STE models on the different grids, are presented in Figs. 2 to 5, respectively.

Table 1 .

 1 Resolution parameters for the LES calculations.

	Grids	G1	G2	G3	G4	G5
	Number of cells
	in	20 × 20 × 66 24 × 24 × 80 30 × 30 × 100 40 × 40 × 133 60 × 60 × 200
	0.3 × 0.3 × 1.0 m 3
	∆ (cm)	1.5	1.25	1.00	0.75	0.50
	λ/∆	1.04	1.25	1.56	2.08	3.12
	l d /∆	0.95	1.14	1.42	1.90	2.85
	∆/η k	17.8	14.8	11.9	8.9	5.9

Table 2 .

 2 Puffing frequency computed by VTE, STE and algebraic closure models for the 57.5 kW fire plume on the different grids.

	Grids	G1	G2 G3	G4	G5
	VTE	3.3	3.5	2.3	2.5	2.5
	STE	2.5	3.0	2.5	2.5	2.5
	algebraic model	3.4	3.3	3.0	1.7	2.4

Table 3 .

 3 Radiant fraction and optical thickness computed by VTE and STE for the 57.5 kW fire plume on the different grids.

	Grids	Model	G1	G2	G3	G4	G5
	Radiant fraction	VTE	0.288 0.293 0.272 0.269 0.271
		STE	0.276 0.261 0.255 0.261 0.269
	Optical thickness	VTE	0.470 0.473 0.477 0.477 0.482
		STE	0.470 0.468 0.474 0.473 0.480

Acknowledgements

JLC wishes to express his gratitude to Electricité de France (EDF) for financial supports.