{'fraction_non_alphanumeric': 0.03024331389503318, 'fraction_numerical': 0.012507540719887392, 'mean_word_length': 4.729493087557604, 'pattern_counts': {'":': 0, '<': 0, '<?xml version=': 0, '>': 0, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
A Study on turbulence-combustion interaction and Sub-grid Scale model in the simulation of Methane pool fire using LES

In this paper, the effect of the combustion and turbulence sub-grid scale (SGS) model on the simulation of pool fire turbulence field has been studied in open source CFD software, OpenFOAM. Two combustion models of Eddy Dissipation Model (EDM) and infinite fast chemistry, with the one-equation and Smagorinsky SGS model, is evaluated for a large-scale pool fire. In general, fast kinetic-based combustion models predict excessive heat release rate. The mean squared of the velocity fluctuations is over-predicted. In this simulation, the turbulence models have no significant effect on the results. In fact, the effect of the combustion model is dominant. The EDM combustion model is more compatible when used with the one-equation SGS model and improves the results compared to other cases. In addition, the infinite fast chemistry combustion model is not a suitable model for fire simulation.

In the early modeling of the fire, the Reynolds Averaged Navier Stokes (RANS) methods were used in order to model the turbulent flow in the fire scenarios [7,8]. In 1994, McGrattan [9] was developed a Large Eddy Simulation (LES) code for the fire simulation in open and closed (compartment) area. In 2000, Fire Dynamic Simulator (FDS) software was released based on the LES method [10,11].
Combustion models determine the reaction rate in the reacting flow and are essential in turbulent combustion regime [12,13]. Xue et al. [14] considered the results of three combustion models of Volumetric Heat Source (VHS), Eddy Brake-Up (EBU), and prePDF for three fire scenarios in large and small rooms, as well as a fire in a tunnel. The results of these three methods were compared by using the k-ԑ model. They preferred the prePDF combustion model with respect to others.. Huang et al. [15] also considered three combustion models of VHS, EBU, and prePDF with the k-ԑ model for a fire scenario in a room. As reported, the VHS model has the least computational time but is not a suitable model for estimating species because chemical reactions do not consider in this model. Some other studies have tested the effectiveness of various combustion models for the LES method [16]. Yeoh et al. [17] compared the results of two combustion models of EBU and Laminar Flamelet using the LES method for combustion in one room, two and more rooms. They reported that the results of these two models are almost similar but the computational cost of the Laminar Flamelet method is lower.
Yang et al. [18] examined the three combustion models of Eddy Dissipation Concept (EDC), infinite fast chemistry, and combustion model based on mixture fraction using FDS software for the fire scenario in a single room with an opening. In this study, the combustion model based on the mixture fraction with two-step kinetics was defined. The results showed that the prediction of the heat release rate in two-step kinetics was more accurate than one-step. Maragkos and Merci [19] compared the efficiency of the EDM combustion model in OpenFOAM and FDS software.
They predicted the mean flame temperature by the OpenFOAM of around 100 Kelvin higher than the experimental results package.
Pasdarshahri et al. [20] examined the simulation of a large pool fire with Smagorinsky and oneequation SGS using the infinite fast chemistry combustion model. It was observed that the model of one-equation and Smagorinsky, estimate the mean vertical velocity with 7% and 12% difference than the experimental results, respectively, but the computational time of one-equation SGS is about 16% higher than of Smagorinsky one.
Yuen et al. [21] studied the three sub-grid models of Smagorinsky, Wale, and Verman using the laminar flamelet model, and concluded that the WALE sub-grid scale model was more compatible with the laminar flamelet method. They also presented [22] the appropriate value for three constant coefficients of Smagorinsky, turbulence Prandtl and Schmidt numbers, which are more compatible with the laminar flamelet combustion model. Margkos et al. [23] corrected EDM model coefficients by comparing EDM and EDC results in a small pool fire. In this study, the EDM combustion coefficient changed from 1 to 8 and concluded that coefficient 1 predicts a more accurate result.
According to this review, the combustion model should be selected according to the fire scenario and the turbulence model [24].
In previous studies, no comparison has been made between different combustion model and subgrid scale methods in turbulence field of pool fire simulation. The contributions of the present study are mainly referred to: Investigation of the turbulence field and comparison of combustion model and sub-grid scale methods in pool fire modeling. So, in the present study, large scale pool fire was considered, to examine the effect of the combustion and the SGS models. The importance of using a suitable combustion model in modeling large scale pool fire turbulence field has been investigated. Considering the simple scenario, the efficiency of the EDM and the infinite fast chemistry combustion models in two different models of the Smagorinsky and one-equation SGS are examined, to simulate the turbulence field in the large pool fire.

Governing equation
The equations of non-premixed reacting flow such as continuity, momentum, energy, species and turbulence equations can be filtered by using the low-Mach-number Favre-filtering method; in this case, the required equations for the combustion modeling are equations 1-4 [25].
In these equations,  is the mixture density, ui, P and T are the velocity, pressure, and temperature respectively.  can be each scalar quantities (such as mass fraction) and T  is the rate of heat generation by the combustion. The terms rad S and S  are the radiation heat transfer and production term in the equations of energy and scalar quantity. In equations (2) to (4), viscous stress tensor, species diffusion vector, and heat flux vector are modeled using the Newton, Fourier and Fick rules, respectively [26].
In this simulation, Smagorinsky [27] and one-equation models [28] are used as SGS. Magnussen and Hjertager presented the EDM combustion model in 1976, with the idea proposed by Spalding [19]. In this combustion model, the fuel consumption rate is calculated by equation (5) In equation (6), sgs sgs k ε , are mean turbulent kinetic energy and turbulent dissipation rate which calculated from SGS parameter. α is thermal diffusivity and Δ is the filter width in LES method.
The infinite fast chemistry model is based on this assumption that the reactions occur very quickly; In other words, the reaction takes place by the mixing of fuel and oxidizing. This means that the time scale of chemical reactions are much faster than the time scale of diffusion and flow transfer.
In this case, the fuel reaction rate which appears in the form of the fuel mass transport source term, can be modeled with equation (7). Radiation heat transfer models with the discrete ordinate method (DOM). The interaction between the radiation and the SGS is ignored. In most studies of fire modeling, it is noted that for the use of the radiation model, it is better to use S4 (8×3 space angles) and to use the upwind method to discretize the radiation equation [29].

2-3. Numerical method:
For this simulation, the fireFoam solver of OpenFOAM (which is the free, open-source CFD software developed primarily by OpenCFD Ltd since 2004) has been used. FireFOAM is a C++, object-oriented LES code based on OpenFOAM, which developed by FM Global. This solver uses the finite volume method [30], and the PIMPLE algorithm is used for taking to account the coupling of velocity and pressure [31]. Also, Rhie-Chow interpolation is used to avoid odd-even decoupling, and the pressure equation is solved by a linear generalized geometric-algebraic multigrid (GAMG) solver.
The maximum local courant number is considered as 0.6, which the time step is about 0.0005 second in this case for 2000000 mesh. The second-order approximation is used for all convective terms in the equations of momentum, sub-grid turbulence kinetic energy, and energy equation. For the temporal term in all equations, the first order Eulerian method is chosen. The first order derivative of the radiation intensity in different directions is discretized by first-order upwind approximation.
Two combustion models of EDM and infinite fast chemistry with two SGS models of Smagorinsky and one-equation are evaluated, in order to investigate the effect of the combustion model as well as the compatibility of the combustion model with SGS model, according to Table 1.
A computer system with a shared memory system and 16 Intel-Core i5-2400 with frequencies of 3.1 GHz, which have four real kernel cores, and 4 GB of temporary memory is used. This machine could perform parallel processing with the MPI library of OpenFOAM code. In most of the computations, 12 cores are used for each study, at the same time. As an estimate, the computational time for the mentioned models is approximately 168 hours for 3 million grids for the modeling of pool fire behavior in 35 seconds (to achieve quasi-steady state) by 12 cores. Table 2 shows the computational cost of the implemented simulation for different grids (1000000, 2000000, 3000000 grids).

2-4. Description of the experiment
The modeling is based on the test done by Tieszen et al. [32] for a methane pool fire with an inlet diameter of one meter. In [32] instantaneous velocity was measured by the PIV method on the central fire plane, and the turbulence parameters were reported. Vertical velocity was reported instantly for a point at 0.505 m above the fuel inlet. Finally, the average velocity and turbulence parameters in the central plate were presented. The accuracies of 2 u , 2 v and uv values are about ± 35%.
As indicated in Figure 1, a computational domain of 3×3×7 m 3 was considered to simulate this test of a one-meter pool fire. The fuel is supplied from a circular inlet of one meter in diameter with The results should be averaged to compare with experimental data. Experimental results have been averaged over 11 cycles (equivalent to about 6 seconds). In this numerical study, the average period is about 20 seconds after the quasi-steady state behavior of the fire.
The proper length scale for fire plume and buoyancy flow is defined by equation (10).
This length scale is an estimate of the effective area around the pool fire which is affected by the combustion. This parameter is used as an estimate of the required computing grid. In general, large-scale values properly resolved when at least ten computational grid are considered for this length scale [19]. In this modeling, the length scale is about 1.3 meters which calculated from equation (8).
The flame plume resolution index (PRI) is defined for achieving the required computational grid in a parametric fire model. According to this parameter, one can talk about the grid used in LES fire modeling. This parameter is defined as equation (11).
In equation (11), x  is the size of the computational grid. For general modeling of the fire flow, it is usually recommended that the PRI value be between 5 and 15 [33]. However, for a detailed study of the oscillatory behavior, the PRI value should be greater than 40 [19]. In the present study, the PRI value is almost 55, so it presents the modeling details with high precision.

3-1. Grid sensitivity
One of the concerns raised in the use of LES method is the number of the required computational grid. In the previous section, a brief account of the computational grid number was stated. In this section, to investigate the effect of the computational grid on fire simulation results, the results of vertical velocity and vertical velocity fluctuation are evaluated.
In Figure 2 The PRI index is one of the fire-specific indexes. The other general index for verifying the computational grid is the LES_iq. This index expresses the ratio of the amount of modeled turbulence kinetic energy to the total turbulence kinetic energy of the flow. Pope [34] proposed that in the simulation of the LES method, at least 80% of the turbulent kinetic energy must be solved directly. LES_iq is defined as equation (12). indicates the optimal accuracy of these grids. In Figure 3, LES_iq is plotted at the z=0 for the two cases O-E and O-I which the mesh number is 3000000. As observed, the index value of LES_iq is higher than 80 % in most points the quality of the mesh is required for this case.

3-2. Puffing phenomena
One of the most important phenomena that occur in the fire is called puff phenomenon [35,36].
In fact, this phenomenon describes the appearance of vortices on the flame surface [37,38]. Due to the difference in density of the horizontal direction of the flame range, and on the other hand, the pressure difference, which is created in the vertical direction due to gravity, the baroclinic and gravity produce the Puffing vortex. This process can be seen in This phenomenon causes oscillating behavior and instability at the flame surface. Moreover, it gradually creates a fluctuation in velocity, temperature, and turbulence parameter. In Figure 4, by using the temperature contour, as well as the velocity vector, the cycle of formation and dissipation of the puffing vortex (after reaching the quasi-steady state) is depicted. As shown in this figure (Part A), a small vortex was initially created near the inlet due to the difference in flow density in the horizontal direction and the pressure difference in the vertical direction [39].
By decreasing the density of the mixture relative to the environment, the buoyancy force is activated which is corresponding to the difference of the density of the mixture and the surrounding environment. Thus the hot gases generated by the combustion move upward due to buoyancy forces and getting bigger along their path (part B). Vortices while moving upward, merged (Part C) and forming a larger vortex (part D). So, the combustion gradually disappears and the energy has been transferred to the vortices is dissipated, and gradually vortices disappear (the end of the fire).

3-3. Puffing frequency
In the modeling of pool fire with LES method, the results are changed periodically, because of the puffing phenomenon. The experimental results [32] have been reported the instantaneous velocity analyses to evaluate the puffing frequency. Figure 5 shows the frequency analysis of vertical velocity data by different methods (for the point at the height of 0.505 m). According to the figure, the dominant frequency which representing the puffing frequency, are given in Table 3.
The infinite fast chemistry combustion model has a lower accuracy in Smagorinsky and one-Equation SGS for prediction of puffing frequency. In contrast, the EDM combustion model has given puffing frequency more precisely. When EDM combustion model is used with a oneequation SGS, the puffing frequency has 11% relative error from experimental results; although the accuracy of this combustion model is lower in the case of using the Smagorinsky SGS.
In general, numerical models report higher frequencies than experimental results. Because these combustion models are based on fast kinetics and predict the amount of heat release more than reality; therefore, by releasing too much energy, the gas velocity increases, and the puffing vortex move up faster. So, the value of the puffing frequency increases accordingly.
The LES model is based on the assumption that the large vortices are solved directly and the effects of small vortices are modeled with SGS models [40]. The energy of large scale vortices gradually transfers into smaller vortices, until they eventually turn into vortices known as Kolmogorov and then, dissipate their energy. This vortex energy transformation from large scale vortices to Kolmogorov's ones called energy cascade. Then the PSD value reaches its maximum value at the puffing frequency. Therefore, the value of PSD and velocity amplitude are maximum at the puffing frequency.
Taylor-scale vortices is in the area where energy is transmitted from large vortices to small-scale vortices [41]. Within this area, the energy spectrum has a power functional of -5/3 of the wave number that is proportional to the inverse of the vortex length or the frequency. In Figure 6, the -5/3 slop is shown; as can be seen, the slopes of the PSD lines in two combustion models are almost close to the -5/3 slop, although somewhere the slope is more than -5/3. In Figure 7, the frequency analysis of these four points are shown. According to Figure 7(a), the dominant frequency of the point A (fA) is the frequency of vortex that forms near the fuel bed.
Then the vortex goes up to reach point B, while its frequency does not change in this path (i.e., fB = fA), but according to Figure 7(c), the amount of kinetic energy associated with this frequency increases, due to the energy received from the combustion process. As seen in Figure 7 According to Figure 7(c), from point A to C, the level of energy that received has increased. From point C to point D, the amount of intake energy decrease. At the height between point C and D, the amount of energy, which released by combustion, reduce and also the vortices begin to dissipate. In Figure 4, also this behavior could be observed schematically.

3-4. Comparison of the turbulence field
One of the issues in the investigation of the turbulence field by LES modeling is the calculation of parameters such as fluctuations velocity, mean squared of velocity fluctuations, turbulence kinetic energy. The main problem in computing these parameters is that in the LES modeling the filtered parameters are calculated and these parameters differ from instantaneous parameters such as instantaneous velocity, which results from the SGS.
In [32], the mean square curve of the horizontal and vertical velocity fluctuations and turbulence kinetic energy on the middle plate (Z=0), was reported for a one-meter pool fire. According to Equations (13) and (14), the relation between two turbulence kinetic energies, which one of them is calculated from an instantaneous velocity, and the other is calculated by the LES method, is given. In the LES method, SGS k is calculated using the assumption of the local equilibrium between production and energy dissipation in the SGS. The total turbulence kinetic energy is calculated according to (14). To calculate the sub-grid velocity fluctuations in equation (15), the isotropic sub-grid velocity fluctuations is assumed. An estimated value of sub-grid velocity fluctuations is obtained according to equation (16).

3-5. Comparison of mean squared fluctuation of velocity
In Figure 9, In Figure 11, uv  were compared in two height of Y=0.4 and 0.8 m (in Z=0 plane). uv  exhibit different behavior of the mean square of the vertical and horizontal fluctuations velocity; i.e. uv  is positive in one direction, which both of , uv  are marked and uv  is negative in the other direction that they are not marked. But uv  at two tails of the curve is close to zero due to the laminar flow.
As we approach the center (X=0), uv  first increases then decrease to zero at the center. Different models follow this trend, but the results of all these models are greater than predicted by the experimental results. Among these four cases, models that use an infinite fast chemistry combustion model has more over-prediction.
In general, flow is expected to be laminar in a region far from fire flames, and therefore the velocity fluctuations are minimal. By approaching the region of flames, the areas of turbulence and the fluctuations appear. When we get closer to the center of the flame, the amount of turbulence energy decreases.
In a summation, it can be seen that in predicting the values of turbulence parameters, the infinite fast chemistry model has an over-prediction of the turbulence parameters. The turbulence SGS models do not have much effect on the results, however, the one-equation SGS has more accurate results. In fact, the effect of the combustion model is more dominant than the turbulence SGS models; the other point is that nearly all models correctly predict the locations of maximal or minimal of the mean squared fluctuations velocity, but they have been modeled results by differences. The other point is that the far from the fuel bed (at heights of 0.4 and 0.8), the amount of fluctuations increases, due to the combustion intensification.

3-6. Comparison of turbulence kinetic energy results
In Figure  Due to the dominance of the perturbation in this direction, the mean turbulence kinetic energy follow the mean square of the vertical velocity fluctuation. The trend of variation in the mean turbulence kinetic energy at two tails of the curve goes to zero, due to the laminar flow; As we approach the center (X = 0), mean turbulence kinetic energy first increases, then decreases in the center due to the reduction of combustion energy. According to Figure 12, different models follow this trend, but the results of all these models are greater than the experimental results. Among these four cases, those using infinite fast chemistry combustion model, are farther from experimental results than others. But the results predicted by the EDM combustion model are closer to experimental values.
Prediction of the numerical method results depends on the nature of the two combustion models.
Infinite fast chemistry combustion model does the modeling with this assumption that combustion is implemented completely, as soon as the fuel and oxidizer are encounters. Then reaction accrues according to the irreversible stoichiometry reaction. In this model, combustion predictions become further than the actual values. Therefore the released energy increases and as a result, the turbulence becomes more apparent. However, the EDM combustion model accounts the characteristic time of the fuel and oxidizer mixing in the process of combustion, to improve the combustion results.

3-7. influence of combustion model constant
Considering the results of two combustion models and turbulence SGS in the previous sections, the EDM combustion model has much better results than the infinite fast chemistry combustion model. Therefore, in this section, by examining the coefficients of the EDM combustion model, the effect of these coefficients on the results accuracy has been examined. In general, two combustion models used, estimate the fuel consumption to be high, and the coefficients need to be corrected for using in the pool fire; one of the methods that can be used to correct the combustion model is the correction of the coefficients according to the flow regime. In [23], by examining small-scale methanol pool fire, the optimal EDM model constant is proposed to be one in a smallscale pool fire.
As indicated in equation (6), the default value of this coefficient is 4. As seen the results of the previous section, the values of turbulence parameter that carried out are somewhat higher than the experimental results, due to excessive combustion prediction. Therefore, in order to reduce the number of combustion parameters, the EDM combustion model must be reduced to decrees the combustion rate of the fuel; so one can reduce the amount of released energy and consequently, the turbulence parameters. In this modeling, by evaluating three coefficients 1, 2 and 4, the coefficient influence on the EDM combustion model is presented in the modeling of the turbulence parameters in large scale pool fire of methane.
In Figure 13 and Figure 14, the effect of the coefficients of the EDM combustion model on two cases O_E and S_E are investigated by changing the coefficients between 1-4. According to Figure   13 and Figure 14, the change in coefficients does not have much effect on the mean square of horizontal velocity fluctuations. In fact, the EDM model is not able to predict the exact mean square of the horizontal velocity fluctuation. If the coefficient one is used, the mean squared values of the vertical velocity, turbulence kinetic energy and uv  are in the acceptable range to the experimental results. So in general, the coefficient 1 is optimal EDM model constant.
A review of the results shows that although the coefficients of the combustion model reduced, the accuracy of the numerical results is closer to the experimental results; however, the use of other combustion models also can improve the numerical results. For example, the use of a steady laminar flamelet model (SLFM) or flamelet generated model (FGM) is one of these models [3].

Conclusion
The

Nomenclature and abbreviations
Greek symbol