Large Eddy Simulation of Bluff Body Stabilized Turbulent Premixed Flame

: A turbulent reacting flow in a channel with an obstacle was simulated computationally with large eddy simulation turbulence modeling and the Xi turbulent combustion model for premixed flame. The numerical model was implemented in the open source software OpenFoam. Both inert flow and reactive flow simulations were performed. In the inert flow, comparisons with velocity profile and recirculation vortex zone were performed as well as an analysis of the energy spectrum obtained numerically. The simulation with reacting flow considered a pre-mixture of propane ( C 3 H 8 ) and air such that the equivalence ratio was equal to 0.65, with a theoretical adiabatic flame temperature of 1,800 K. The computational results were compared to experimental ones available in the literature. The equivalence ratio, inlet flow velocity, pressure, flame-holder shape and size, fuel type and turbulence intensity were taken from an experimental set up. The results shown in the present simulations are in good agreement with the experimental data.


INTRODUCTION
In turbulent reactive flow simulations, computational models have achieved a great development in recent years with the improvement of computer power.This development allowed more accurate solution of problems such as the instability caused by the turbulence in combustion chambers of rocket engines, gas turbines, turbojet afterburners, ramjets and scramjets.These models have been used to study bluff body stabilized flames, which allow combustion devices to operate at very high free stream velocities.Advanced afterburner design methods have been discussed by Lovett et al. (2004), who outlined the fundamental combustion sciences and engineering challenges that need to be addressed.Among other problems, Lovett et al. (2004) highlighted the requirements for flame stabilization and combustion dynamics.These authors discuss the need for advanced design methodologies and tools, and they stress the limitations of existing computational models to capture the physics of those phenomena.
In turbulent combustion, the behavior of the turbulent flame front is predominantly dictated by the turbulence (Peters, 2000).Combustion instability is also directly related to turbulence (Weller, Marooney and Gosman, 1990).Therefore, it becomes mandatory in numerical simulations to use turbulence models that reproduce these dynamic processes, which are mainly produced by large scale turbulence.Many researchers have used Reynolds Averaged Navier Stokes (RANS) and eddy viscosity turbulence models to simulate reactive flows behind bluff bodies.However, important discrepancies were observed due to shortcomings in the RANS methodology, especially in complex flows with circular obstacle (Saghafian et al., 2003;Frendi, Skarath and Tosh, 2004) and triangular one (Bai and Fuchs, 1994;Dourado, 2003;Eriksson, 2007).Significant differences in results were found for the turbulent velocity, integral length scale and turbulent viscosity distribution between RANS models.These differences affect the computation of flame front diffusion, which is underpredicted.The Kelvin-Helmholtz effects behind the obstacle are not captured well and the length of the recirculation zone as well as the turbulent flame speed are not recovered.
Experimental studies were conducted by Sjunesson, Henriksson and Lofstrom (1992) and Sjunnesson, Olovsson and Sjoblom (1991) at Volvo (Sweden).The Volvo test rig had an equilateral triangular bluff body with a blockage of 33% and a Reynolds number (Re) based on the inlet velocity and two times the channel height of 204,000.For this Re, the Strouhal number (St) observed was 0.417.The premixed gases were air and propane at equivalence ratios of 0.65 and 0.85.The Damköhler number (Da) was 10 and the Karlovitz number (Ka) was 4. The resulting flame for these conditions was on the thickened-wrinkled flame regime of the Borghi diagram rather than in the wrinkled flame regime.Sanquer (1998) also presented experimental results for premixed flame stabilized by a triangular prismatic bluff body at the University of Poitier (France).The blockage in his experiments was also 33% for a Re in the range of 6,690 to 23,150, much lower than the Re on the Volvo experiments.The St for the lower Re experiment was 0.276.Therefore, the vortex shedding characteristics and the turbulent scales (integral and Kolmogorov) were significantly different from the Volvo experiment.Sanquer's results will be discussed and compared to the present numerical simulations on the following sections.Specifically, the present simulations correspond to Sanquer's experiment that falls on the Borghi diagram where Da<1 and Ka>1 corresponding to the thin-wrinkled flame region.
These and other experimental results (Cheng, 1984, Cheng, Shepherd and Gokalp, 1989, Cheng and Shepherd, 1991, Kiel et al., 2007, Chaudhuri et al., 2011) describe the interaction between the turbulent structures and the flame.Numerical simulations of turbulent combustion must rely on models that are able to capture the complex turbulent vortical structures and this explains why results obtained with RANS models are less accurate.In order to improve numerical simulations, turbulent combustion models can be ported from RANS models to Large Eddy Simulation (LES) models, which are known to capture the large scale turbulent structures.This approach does not always result in more accurate simulations.It is necessary to identify those models that have the best results and show greater potential for future improvements when used in LES.
A number of investigators have used LES to simulate reactive flows with bluff body flame holders.Porumbel and Menon (2006), Akula, Sadiki, and Janicka (2006), Ge et al. (2007), Park and Ko (2011) and Manickam et al. (2012) simulated the Volvo experiment (Sjunnesson, Olovsson and Sjoblom, 1991).Porumbel and Menon (2006) have simulated bluff body flows with the Linear Eddy Mixing (LEM) LES based on the model proposed by Kerstein (1989) and developed into a sub-grid model by Menon et al. (1993) for premixed combustion.They discussed the differences between their model and the Eddy Breakup model (LES-EBU) and which of the two is best to represent the physics of the reacting flow behind an obstacle besides the ability to resolve the turbulent eddies that wrinkle the flame front.This model is specially adequate for high turbulent intensities where the chemical time is considered infinitely small, corresponding to a Da much greater than 1.Their conclusions were confirmed by experimental results presented by Chaudhuri et al. (2011).Akula, Sadiki, and Janicka (2006) also performed LES simulations with a flame surface density formulation adapted from Boger and Veynante (2000).This model includes the resolved progress variable in the transport equations and the flame wrinkling model avoids unrealistic detachment of the flame from the flame holder.The compressible model presented by Tabor and Weller (2004), along with a similar progress variable approach, is also able to capture the flame wrinkling and thus capture the physics as proposed by Porumbel and Menon (2006).The model of Tabor and Weller (2004) is used in the present simulations.
Park and Ko (2011) used a dynamic sub-grid scale combustion model based on the G-equation that describes the flame front propagation (Peters, 2000).The G-equation is a model for the chemical species conservation equations.In this model, a new LES dynamic sub-grid combustion model is introduced along with a new turbulent flame speed model based on the sub-grid turbulent diffusivity.These models are intended to better represent the flame characteristics.Park and Ko (2011) present comparisons between non-reacting and reacting flows and discuss the effect of chemical reaction on the development of the wake behind the triangular bluff body.They show that their model is able to capture the stabilization of the Karman vortex street behind the obstacle due to combustion.The results for the dynamic model show better agreement with experimental results when compared with the Smagorinsky model.The unsteady flow field was captured with good accuracy and the temperature and reaction rate profiles are well capture by both models, showing that the two combustion models tested are reasonable for the simulation of pre-mixed combustion behind a bluff body.Manickam et al. (2012) used an algebraic flame surface wrinkling (AFSW) reaction model based on the progress variable approach.Inert and reacting cases were analysed and compared with experimental results and with another well validated turbulent flame speed closure model denominated Turbulent Flame Speed Closure (TFC) (Zimont and Lipatnikov, 1995).For the non-reacting test case, they found a shedding frequency equal to 122Hz, while the experimental measured frequency was 110Hz.The comparisons for the mean flow variables and root mean square axial and normal velocity components were in good agreement with the experiments.In spite the fact that the recirculation zone length was well captured, the rms velocity distributions were underpredicted in this zone.Manickam et al. (2012) presented a detailed discussion on the influence of grid refinement and three subgrid scale models on the results and concluded that a coarse grid is too dissipative.Contrary to what might be expected, the grid refinement has a weak influence on the computation of the St, showing that at least the large scale structures were captured by a course grid LES simulation.The St dependence on the grid was also not observed for the reactive case, independently of the reaction model.The three sub-grid scale models tested were the Smagorinsky (SM), the dynamic Smagorinsky (DS) and the sub-grid scale kinetic energy (KSGS).They showed that the performances of the three different models are more or less equivalent.
Flame Surface Density (FSD) models rely on geometric parameters of the flame front to evaluate its progress.In this case, the laminar flamelet model used in RANS has been adapted to LES by considering a locally laminar flame wrinkled by turbulence.The amount of wrinkling (Σ) is measured by the flame surface area per unit volume (Boger and Veynante, 2000).Transport equations for the progress variable and for the wrinkling variable are solved to describe the evolution of the flame since the wrinkling increases the burning rate.
Similar to the Σ variable, Weller (1993) proposed a model based on the density of wrinkling (Ξ), which is the flame area per unit area resolved in the mean direction of propagation.Weller (1993) originally developed this model for RANS and later Tabor and Weller (2004) adapted the model for LES.The advantage of using Ξ is that it should be easier to model the transport terms as discussed in Weller (1993), Weller et al. (1998) and Tabor and Weller (2004).This model is used in the present study and will be discussed in detail on the following sections.
In the present study, the SM and dynamic sub-grid models for turbulence were used.For the combustion model, the flame surface wrinkling (Ξ) formulation, developed by Tabor and Weller (2004), was used.The objective of the study was to investigate the performance of a turbulent combustion model when applied to the simulation of pre-mixed turbulent flames behind a triangular obstacle.The review of the literature shows that such type of performance investigation has been conducted previously, with results compared to the Volvo experiment, which has a Da in the range of 10 and a Ka of 4, corresponding to a thickened, wrinkle flame.The present investigation considered an experiment which has a Da of 4.5 and a Ka of 1.3; on the range of thin, wrinkled flame, previous investigations found in the open literature did not consider LES simulations in this combustion regime, and the thinner flame thickness poses a more severe test for the turbulent combustion model.The Re for the simulation was 6,690, lower than the Re on the Volvo experiment, which was 204,000 and thus have a significant different flow dynamics.The chosen reactive LES model was evaluated by comparing results with experimental ones obtained by Sanquer (1998), which have not been analyzed before with other LES numerical models.

PROBLEM FORMULATION
In order to account for turbulence and combustion, the reactive flow governing equations are filtered using the LES concept, and the combustion process is accounted for by following the flame front.Therefore, it is necessary to define the filtering and a variable to account for the regions of burned and unburned gases.In this section, the formulation derived by Weller (1993) and Tabor and Weller (2004) is presented.

PreliminAry definitions
In premixed flames, a reaction wave propagates from the burned to fresh gases.The progress variable c, that identifies this propagation, varies between 0 for fresh gas and 1 for burned gas.The transitions between these values describe the flame front.A progress variable c can be defined based on the normalized temperature (T) or on the reactant mass fraction (Y).Using the temperature, it results in: Where b subscript stands for burned gas and u subscript, for fresh unburned gas.The flame front propagation is modelled by solving a transport equation for the density-weighted mean reaction regress variable denoted by b, where b = 1 -c.
In LES, it is assumed that the dependent variables can be divided into grid scale (GS) and sub-grid scale (SGS) components, such that, for any given dependent variable, it results in: (2) Here, D is the computational domain with boundaries ∂D, t is the time and x, the coordinate directions.The kernel G=G(x, ∆) is any function of x and of the filter width ∆.G has the properties Introducing a conditional filter (Tabor and Weller, 2004), with an indicator function l, results in: is in the unburned gas region otherwise.For a tensor ψ of any rank, one may define the phaseweighted value of ψ at any point.
Introducing the combustion progress variable b as a GS indicator function, we can obtain: where b (x, t) is the probability of the point (x, t) being in the unburned gas. .
In compressible flow, there is density (ρ) variation, and the product b ρψ ρψ = .
Where the subscript u indicates the unburned gases.
Defining a density-weight average u u ψ ρ ρψ = .in the unburned phase, we can obtain: From Eq. 8 in Eq. 9, results in: . ~u

filtered continuity equAtion
The governing equations will be written in a coordinate system placed at the flame surface, such that ⊥ n and n || are the unity vectors pointing the normal x ⊥ and parallel || x directions to the flame surface.The metrics of this coordinate system are h ⊥ and || h .This coordinate system is used in order to include the conditional filter based on the progress variable b.
The filtered continuity equation reads (Tabor and Weller, 2004) Where ʹ is the area element on the surface interface.Σ is interpreted as the amount of interface for the filtered component, the flame surface density.
The surface filtering operation ⎩ ⎨ ⎧ is defined as: From Eq. 10 in Eq. 11, results in: This surface filtering operation applied to ⊥ n results in: can be related to the GS with the n f direction of the interface: Where Ξ represents the total sub-grid surface area by the smoothed surface area in the f n direction: Where represents the area of the GS surface.
Combining the burned and unburned gases into the weighted total density results in: The final continuity equation reads: where Su is the laminar flame speed.

filtered momentum And energy equAtions
Conditionally filtering the momentum equation gives: Where p is the pressure, S u is the laminar flame speed, Σ is the flame surface density, μ λ is the stress tensor and is the symmetric part of the strain tensor.The terms in brackets represent the effect of the interface on the momentum balance.
B u represents the SGS stress tensor.
This term requires modeling.
The filtered energy equation is: and , ) ( represent the SGS pressure dilatation π and dissipation , ) ( (Tabor and Weller, 2004).The total energy at the interface is presented in brackets.

turbulence model
In the present study, three turbulence models available in OpenFoam were used.These turbulence models are described by Fureby et al. (1997).First, the Smagorinsky classical model is presented: Where S C 2 s sgs is the sub-grid viscocity and The constant C s is equal to 0.18.
The second turbulence model is the one-equation model with sub-grid kinetic energy ( ) proposed by Yoshizawa and Horiuti (1985) and given by: ( ) where the operation a stands for a ~in the product The sub-filter stress is: and turbulent viscosity is: To close the system of equations, the transport equation for the kinetic energy is used. where and The third turbulence model is the dynamic one-equation model that used the Germano identity, following Ghosal et al. 1995, considering a scaling law and homogeneous flow.The kinetic energy equation is: where L is the Germano identity (Ghosal et al., 1995): subgrid combustion model Models for the SGS stress tensor, flux vectors, dissipation and filtered reaction rates are used to close the governing equations.The models for the SGS stress tensor and flux vectors are the standard ones used in LES, since they do not depend on reacting flow data.
A flamelet model with conditional filtering for LES is used to derive the transport equations.This model considers conditions of the Klimov-Williams criterion (Ka=1) in the Borghi diagram (Borghi and Destriau, 1998).Instead of using the flame propagation speed in terms of the laminar flame area per unit volume Σ and the degree of wrinkling of the flame at a point in the domain, the present model (Weller, 1993) uses the flame surface density and a wrinkling surface function Ξ.The function Ξ is the average flame area per unit volume divided by the projected area in the mean direction of propagation.
In the flamelet regime, turbulent motions are slow and do not affect the flame structure.The disturbance velocity u´ considers how a rotation speed of the larger turbulent motion wrinkles the flame surface front (Veynante, 2006).In the reaction zone, the characteristic scales for the reaction processes are below the filter, so that, in reacting LES, a proper treatment (modeling) of the reaction zone is needed.In other words, equations for the geometric variables b and Ξ are needed.

Sub-grid scale model for b
In Eq. 20, the right hand side needs a SGS model, which is based on the conditional filter of unburned gas velocity . This term is modeled using: where is the slip velocity of the unburned minus burned gas .~c This is similar to the properties of laminar flame for LES: where is the diffusion coefficient of the sub-grid, and is the flame normal direction.
From Eq. 20, one arrives at the final equation for b (Tabor and Weller, 2004): Sub-grid scale modeled equation for Ξ From the transport equation for the sub-grid flame area density Σ, proposed by Weller (1993), an equation for Ξ is obtained from the relation and the resolved unburned gas volume fraction b (Tabor and Weller, 2004): where is the surface-filtered effective velocity of the fla-me and is the local instantaneous velocity of flame surface.
The proposed model for Eq.38 considers wrinkling generation , such that (Tabor and Weller, 2004): and The models for G and R, which accounts for the interaction between turbulence and flame, are: and where Finally, an equation for the laminar flame speed (Tabor and Weller, 2004) is needed.A proposed transport equation is: The convective velocity of the laminar flame is the filtered surface speed

NUMERICAL MODEL
Th e present investigation uses OpenFoam, an open source C++ collection of libraries for transport equations.Th e model equations are solved numerically based on a cell centered unstructured fi nite volume scheme.Th e solution is based on a segregated approach.For time integration, a fi rst order explicit Euler method is used.For spatial discretization, second and third order Total Variation Diminishing (TVD) schemes are used.Th e PISO algorithm (pressure-implicit split-operator) is used to solve the pressure-velocity coupling.Th e discretization method is the standard Gauss fi nite volume integration.Th e combustion solver called Xifoam, available on OpenFoam, was used.Th is premixed turbulent combustion model is described in the previous section -Preliminary defi nitions (Weller et al., 1998;Tabor and Weller, 2004).
Both reactive and non-reactive cases were run on a SGI Altix XE 1300 cluster with Intel x86 64 processors, running SUSE Linux Enterprise Server version 10.The cluster has 144 cores with 432 GB, DDR3, 1066 MHz RAM.For the non-reactive cases running with 64 cores, the average seconds per iteration computational time for a simulation using SM was 1.47 seconds.The clock time was about 12h52m.For the reactive case running with 136 cores, the average seconds per iteration for a simulation using the dynamic model was 8.04 seconds.) is 6,690.Th e obstacle used as fl ame holder is an equilateral triangular cross section obstacle whose backside is located 0.160 m from the entrance.Th e obstacle blockage is 33% of the total area and corresponds to the r-65 test case in Sanquer's experiment (Sanquer, 1998).
First, only inert simulations are considered with an initial temperature of 300K.At the outfl ow boundary, a pressure wave transmissivity boundary condition is used for pressure (Candel, 1992).A uniform velocity profi le is imposed at the entrance with a speed of U axe =3.1 m/s to match the experimental value.On top of the uniform average velocity profi le at the inlet, a fl uctuation is added using a routine available on OpenFoam that mimics the turbulent statistical properties.Th is procedure is necessary since the inlet turbulence determines the fl ow behavior on the domain, as discussed by Tabor and Baba-Ahmadi (2010).
Th e fl ow is considered periodic in the spanwise direction.A wall function is used in the channel walls, as shown schematically in Fig. 2 (Jayatilleke, 1969).
For reacting fl ows, propane (C 3 H 8 ) premixed with air is considered.Th e mixture is ignited behind the obstacle in the recirculation zone to achieve a proper performance and avoid fl ame blow off .Th e ignition point is located at 0.05m behind the obstacle in the center of the recirculation zone.A combustion time of 3ms was used before collecting data to avoid numerical transients and allow time to achieve stable combustion behavior.Th e imposed initial fl ame speed was 0.256 m/s and the initial condition for the regress variable was b=1.

Impose velocity pro le Wall
Wall channel    Th e grid used in the present study is composed of 388,355 volumes for two-dimensional simulations and 2,300,000 volumes for three-dimensional simulations.Figure 3 shows the grid topology around the obstacle.Lower grid densities were tested but it was not possible to recover the recirculation zone size obtained by Sanquer (1998).Th e adequacy of the grid was also tested through the turbulent decay rate, which should follow the Kolmogorov -5/3 decay, as discussed in the next section.

Wall function Obstacle
A detailed study of grid requirement using grid quality assessment techniques was presented by Manickam et al. (2012).Th ey simulated the Volvo experiment using LES and a fl ame surface wrinkling reaction model based on the progress variable, similar to the model used in this investigation.Th eir fi nest grid size has 2.4 million cells, close to the grid size used in the present investigation.Th eir intermediate and fi ne grids were able to recover all relevant experimental results for mean and turbulent quantities as well as the St.Nevertheless, even the fi ne grid presented lower quality near the fl amelet region.Considering that the Volvo experiment has a Re much larger than the Sanquer's experiment, it is safe to say that the fi ne grid for the Volvo higher Re should capture all the relevant structures for the Sanquer's lower Re experiment and will allow the necessary level of accuracy for the sub-grid scale model.
Th e instantaneous velocity was monitored in order to establish the time for which the fl ow can be considered periodic stationary.Figure 4 shows the history of longitudinal velocity.As can be seen, the transition lasts 0.18 seconds, so the fl ow can be considered periodic stationary beyond that.

RESULTS
In this section, the performance of the turbulent combustion models presented in problem formulation section and available on OpenFoam is assessed for the fl ow conditions considered by Sanquer (1998).
First, results are presented in terms of frequencies f q and Strouhal numbers (St=f q d/U axe ) of the vortex emission frequency behind the obstacle, where d is the obstacle height and U axe is the maxim velocity at the channel entrance.Th e size of the average recirculation zone X r is also compared.Next, the energy spectrum is analysed in order to verify the adequacy of the computational grid and to compare the inertial range against the -5/3 turbulence decay rate for both inert and reactive cases (Pope, 2000).Th en, streamwise and normal velocity profi les are compared to the experimental results in the recirculation zone and downstream of the recirculation zone.Finally, the progress variable and temperature profi les are presented.Th e recirculation zone and the energy spectra are determined through the velocity data aft er the simulation reaches the stationary periodic state.Once the length of the recirculation zone X r is determined, one can calculate the St to compare it with experimental results.
To determine the values of the energy spectrum, the Fast Fourier Transform (FFT) of the oscillatory instantaneous velocity is taken.Th is information is used to determine if the simulation captures the large turbulent scales and models the sub-grid scales correctly.

inert flow cAse
First, inert fl ow results are presented.Table 1 presents   LES simulation and for RANS simulations performed by Dourado (2003).Th e results show good agreement with those obtained experimentally.
Th e length of the recirculation zone based on the evolution of the average velocity Ux was determined taking the average values along the X axis downstream of the obstacle.Figure 5 shows the profi le of the average velocity used to defi ne the recirculation zone as shown schematically in Fig. 6.Th e length of the recirculation zone found in the present simulation is 0.023m, while the length found by Sanquer (1998) is 0.0204m, and a simulation based on a RANS model (Dourado, 2003) gives 0.0222m.Th e recirculation zone length diff ers approximately 11.5% from the experimental value.
To determine the vortex emission frequency, a temporal Fourier analysis was used (Dourado, 2003).For the numerical simulation using LES, a value of 93.5Hz was found, as shown in Fig. 7, while the experimental result was f q = 89Hz.Th ese values correspond to a St=0.284 for the simulation and St=0.276 from the experiment with a diff erence of only 3%.Th e RANS model gives fq=87 Hz and St=0.2694, with 2.39% diff erence in the St to the experimental result in this inert case.
For the three-dimensional simulations, the energy spectrum shown in Figs. 8 through 11 has the -5/3 energy decay rate expected for large Re turbulent fl ows.Figures 8 and 9 are the spectrum for the streamwise velocity component at X/X r =1.4, in two diff erent distances from the channel centerline, y/h=0.1 and y/h=0.41,respectively.Figures 10 and 11 are spectrum at those two same positions, but for the normal velocity component.On the energy spectrum in Figs. 9 and 11, peaks at about 80Hz can be identifi ed.Th ese peaks correspond to the large scale vortex shedding frequency measured above the channel centerline.At the centerline on Figs. 8 and 10, the peaks seem to be at 100Hz but the spectrum is not conclusive.
Filtered longitudinal Ux and normal Uy component velocity profi les at X/X r =0.8 and X/X r =1.4 are shown in Figs. 12 to 15. Figure 12 shows the longitudinal velocity profi les in the normal direction obtained in the simulation at X/X r =0.8.Th e results are in good agreement with the experimental results, but the dynamic model overpredicts the velocities in the region behind the obstacle and does not show the reverse fl ow in the recirculation zone.Downstream of the recirculation zone at X/X r =1.4 (Fig. 13), the longitudinal velocity distributions are also close to the experimental results, even for the dynamic model, which now shows a moderate underprediction.Th e results obtained with the SM and one equation models are more accurate.
Th e normal velocity component is underpredicted by all models in the regions behind the obstacle (y/h<0.6;Fig. 14), with the dynamic model showing the worse results.Approaching the upper wall (y/h greater than 0.6) the turbulence models give better predictions, but the Smagorinsky model gives results a little higher than the experiments.Never the less, the results are within the limits of the experimental error.the -5/3 decay rate.A peak on the energy spectrum can be observed at 94Hz, which is close to the experimental value of 89Hz observed by Sanquer (1998).
Streamwise and normal velocity profiles at two different positions in the streamwise direction (X/X r =0.8 and X/X r =1.4) are presented in Figs. 17 through 20.Again, dynamic, one equation and Smagorinsky SGS model results are presented.The streamwise velocity distribution at X/X r =0.8 is very close to the experimental observed velocity distribution for all the three models, as shown in Fig. 17.At X/X r =1.4, the results for the streamwise velocity also compare very well with the experiments, as shown in Fig. 18, but the one equation model overpredicts the velocity distribution.Th e values of the normal velocity component at X/X r =0.8 are presented in Fig. 19.Th e one equation and the Smagorinsky models are in good agreement with the experimental results in the region behind the obstacle (y/h<0.4),while the dynamic model underpredicts the experimental results.Away from the central region, approaching the upper wall, the trend is reversed and the dynamic model shows better results than the Smagorinsky and one equation models.Th ese results are not conclusive to which turbulence model is more suitable for the reactive case simulation.Th is same behavior is also observed for the normal velocity component at X/X r =1.4 (Fig. 20), with diff erent turbulent models performing better or worse at diff erent regions of the fl ow domain, but all models capturing the general trends of the experimentally observed velocity distributions.
Th e progress variable computed by the dynamic model is presented in Figs.21 and 22 and compared to the experimental distribution.Th e simulation captures the general behavior of the progress variable quite well.Despite the fact that the simulation overpredicts the progress variable inside the fl ame zone, the fl ame front at c=0.05 around y/h=0.6 matches the experimental value in X/X r =0.35 and underpredicts for X/ X r =1.4.Also, Sanquer (1998) states that the experimental results seem to be displaced to the left .
For the 3-D reactive case, a value of 1,750K was obtained for the flame average temperature with the dynamic model.The theoretical adiabatic flame temperature for a equivalence ratio of 0.65 corresponds to 1,750K for this fuel.Figures 23 and 24 show the temperature distribution in the normal direction with a profile corresponding to a premixed flame.
The vorticity distribution in the inert and reactive cases can be compared with the help of Figs. 25 and 26.These figures show the flow structure in the (x, y) cut at the center of the channel.As in Park's discussion (Park and Ko 2011), the chemical reaction has a stabilizing effect on the vortex shedding.Th e inert case clearly shows the Karman vortex street behind the bluff body, while in the reactive case the characteristic alternating vortex structure starts further downstream with a much weaker strength.Th e vorticity near the bluff body has a stretched topology on the reactive case and the wake spreading is lower than the wake spreading of the inert case, which is in agreement with the behavior described by Park and Ko (2011).

CONCLUSIONS
Results from numerical simulations of inert and reactive fl ows were compared to experimental results obtained by Sanquer    (1998).Smagorinsky, one equation and dynamic LES models were used and the reaction was simulated using the Ξ density of wrinkling model.The experiments are in the thin, wrinkled flame regime with a Re of 6,690, a Da of 4.5 and a Ka of 1.3.For the inert flow case, the simulations using LES good results were obtained for the St, with an error of the order of 3%.The size of the recirculation bubble is close to the size observed experimentally but the margin of error is larger, around 12%.But these values are not absolute values since the experimental results also have a margin of error.The energy decay on the inertial range was correctly captured at -5/3.As far as the velocity distribution is concerned, the streamwise velocity component shows good agreement with the experimental values, but the normal velocity components are somewhat off, in part as a consequence of the over-predicted recirculation zone length.
The results for the reactive flow case show similar performance with the inert case.The energy decay rate in the inertial range was correctly captured and the vortex sheading frequency was close to the experimental value.The progress variable and the correct behavior considering the stabilization of the vortex sheading behind the obstacle for the reactive case were captured by the simulation.Again, on the reactive case, the velocity distribution was somewhat off from the experimental measurements.Regarding the velocity profile, it is still interesting to observe that the SM gives a better result than the dynamic model used in these simulations.The results may be improved by tailoring the dynamic turbulent model coefficients.In general, the results for all three models are relatively similar and the differences between these models are small.Therefore, since the dynamic model has a higher computational cost, it may be worthwhile to use the Smagorinsky or the one equation model.
transformed coordinates for , || x x ⊥ ) ( : time, u´ is the intensity of turbulence in the sub-grid and And the superscripts 0 and ∞ stand for the value of unstrained flame speed and the value at equilibrium, respectively.

initiAl
And boundAry conditions Th e geometry of the channel with a bluff body fl ame holder is presented in Fig. 1.It consists of a channel 0.600 m long, 0.160 m wide and 0.0288 m tall.Th e Re based on the inlet velocity and twice the channel height (

Figure 2 .
Figure 2. Schematic boundary conditions.Turbulent velocity profi le at the channel entrance, walls and obstacle.

Figure 3 .
Figure 3. Grid structure around the obstacle.
the frequencies, St and recirculation zone length for the present t

Figure 5 .
Figure 5. Profi le of the average velocity u in the recirculation zone.

Figure 24 .
Figure 24.Temperature contours behind the bluff body.

Table 1 .
Comparison of numerical Large Eddy Simulation (LES) and Reynolds Averaged Navier Stokes (RANS) results with experimental results.