**Arne Heinrich 1,\*, Guido Kuenne 1, Sebastian Ganter 1, Christian Hasse <sup>2</sup> and Johannes Janicka <sup>1</sup>**


Received: 25 July 2018; Accepted: 3 September 2018; Published: 6 September 2018

**Abstract:** Combustion will play a major part in fulfilling the world's energy demand in the next 20 years. Therefore, it is necessary to understand the fundamentals of the flame–wall interaction (FWI), which takes place in internal combustion engines or gas turbines. The FWI can increase heat losses, increase pollutant formations and lowers efficiencies. In this work, a Large Eddy Simulation combined with a tabulated chemistry approach is used to investigate the transient near wall behavior of a turbulent premixed stoichiometric methane flame. This sidewall quenching configuration is based on an experimental burner with non-homogeneous turbulence and an actively cooled wall. The burner was used in a previous study for validation purposes. The transient behavior of the movement of the flame tip is analyzed by categorizing it into three different scenarios: an upstream, a downstream and a jump-like upstream movement. The distributions of the wall heat flux, the quenching distance or the detachment of the maximum heat flux and the quenching point are strongly dependent on this movement. The highest heat fluxes appear mostly at the jump-like movement because the flame behaves locally like a head-on quenching flame.

**Keywords:** sidewall quenching; LES; premixed methane; flame–wall interaction; FGM

## **1. Introduction**

The world's primary energy demand will be increasing about 30% in the next 20 years, and 75% of this demand will be provided by combustion [1,2]. Since resources are limited and combustion can cause pollution, it is necessary to better understand the fundamentals. One important field, in which research is required, is the flame–wall interaction (FWI). Due to technical improvements, e.g., of internal combustion engines or gas turbines, its relevance rises [3,4]. Here, the reaction zone of the flame gets close to cold walls and quenches. This reduces efficiencies and increases pollution formation [5]. The FWI can be subdivided into two canonical configurations, namely head on quenching (HOQ) and sidewall quenching (SWQ) [5]. In the latter case, the flame moves perpendicular along the wall and only the flame tip is affected. In case of HOQ, the flame moves parallel towards the wall and the whole flame quenches. The dimensionless quantities for the quenching distance and heat flux differ clearly in both cases. Because of the importance of the FWI, an SWQ configuration is considered in this work.

In terms of premixed conditions, SWQ has been widely studied theoretically, experimentally and numerically. Primarily, fuels like hydrogen and gaseous hydrocarbons, such as methane, ethylene, propane or butane, were used. The most simple fuel, but experimentally difficult to handle, is hydrogen. It was investigated by Cheng et al. [6] for an experimental turbulent boundary layer over a heated

wall. They determined the limits of the reactions near the wall depending on the equivalence ratio. A more recent study regarding hydrogen combustion was conducted by Gruber et al. [7]. Within their three-dimensional direct numerical simulation (DNS), the authors used detailed chemistry (DC). Their data suggests that near-wall coherent turbulent structures are important for the wall heat flux and that the heat release rate at the wall is strongly controlled by exothermic radical recombination reactions. However, hydrogen behaves differently than hydrocarbon based fuels because it has a much higher consumption speed and it will not be investigated in this current work.

The effect on the quenching layer thickness of different hydrocarbons with variations in the equivalence ratio and different wall conditions was researched by Saffman [8]. The investigations of Lu et al. [9] and Ezekoye et al. [10] used different fuels and equivalence ratios as well, but in a constant volume chamber. They measured the unsteady heat transfer during SWQ. Boust et al. [11] developed a thermal formulation for single-wall quenching without any empirical coefficients. This formulation is based on experimental data and is applicable for lean and stoichiometric methane/air mixtures in a pressure range from 0.05–0.35 MPa. DNS studies combined with a single step mechanism, which mimics a lean methane flame with an inlet temperature of 600 K were carried out by Alshaalan and Rutland [12,13]. They investigated a V-shaped flame in a weakly turbulent Couette flow. Nevertheless, the number of simulations with a chemistry treatment beyond single step mechanisms for hydrocarbons is sparse and most turbulent premixed FWI investigations were focused on the HOQ configuration.

In the present numerical study, an SWQ burner configuration is investigated, which is based on the experimental design of Jainski et al. [14,15]. In the experimental setup, a premixed V-shaped flame interacts with a cooled wall under laminar and turbulent conditions. First numerical results for both cases have been already provided by Heinrich et al. [16,17] with a highly-resolved (large eddy simulation) LES with DNS-like spatial resolution. The turbulent and laminar cases were compared and analyzed with each other and the experiments. The global structure of the flame, the flame position, the temperature field and the velocities are given correctly by the methods used. Compared to other highly resolved simulations [7,12,13], which used simple configuration, this setup is more complex and the turbulent structures of the flow are not comparable with channel or pipe flows. In this work, a transient analysis of the flame–wall–turbulence interaction (FWTI) is conducted. This gives further insight into heat fluxes, quenching distances and quenching positions beyond the experimental findings. It will be shown that a strong dependency of the moving direction of the flame tip for characteristic properties exist. Moreover, it will be explained that the highest heat fluxes occur because the turbulent SWQ flame behaves locally like an HOQ flame.

The remainder of the paper is structured as follows: the numerical description, the chemistry treatment and the domain are given in Section 2. The transient analysis of the FWI is provided in Section 3. Finally, conclusions are given in Section 4.

#### **2. Methods**

## *2.1. Numerical Description*

The simulation in this work has been conducted with the Finite Volume code FASTEST. It is a block-structured code that uses hexahedral, boundary fitted grids to represent complex geometries. The Massage Passing Interface (MPI) is utilized for the inter-processor communication. The incompressible Favre-averaged Navier–Stokes equations with a low-Mach number formulation for the variable density *ρ* are discretized with a second order spacial scheme [18] and the time integration is done by an explicit low storage three-stage Runge–Kutta scheme with second order accuracy. In each Runge–Kutta step, the pressure correction procedure is applied. The corresponding mass and momentum equations are the following:

$$
\frac{\partial \overline{\rho}}{\partial t} + \frac{\partial}{\partial \mathbf{x}\_i} \left( \overline{\rho} \widetilde{u}\_i \right) = 0 \,, \tag{1}
$$

$$\frac{\partial}{\partial t} \left( \overline{\rho} \,\,\widetilde{u}\_{\dot{j}} \right) + \frac{\partial}{\partial \mathbf{x}\_{\dot{i}}} \left( \overline{\rho} \,\,\widetilde{u}\_{\dot{i}} \,\widetilde{u}\_{\dot{j}} \right) = \frac{\partial}{\partial \mathbf{x}\_{\dot{i}}} \left( \overline{\mathbf{r}}\_{\dot{i}\dot{j}} - \overline{\rho} \mathbf{r}\_{\mathbf{s}\mathbf{g}^{\sf s}} \right) - \frac{\partial \overline{p}}{\partial \mathbf{x}\_{\dot{j}}} + \overline{\rho} \mathbf{g}\_{\dot{j}} \,\,. \tag{2}$$

The quantity *u* is the velocity, *p* is the pressure, *τij* is the shear stress tensor and the operators and ¯ represent the Favre-averaging and filtering operation, respectively. The sub grid fluxes of the momentum *τsgs* are accounted by eddy viscosity approach proposed by Smagorinsky [19]. The model coefficient is determined by the dynamic procedure proposed by Germano et al. [20] with a modification by Lilly [21]. For the discretization of the convective terms in the scalar transport equations, a total variation diminishing (TVD) limiter as proposed by Zhou et al. [22] is used.

#### *2.2. Numerical Domain*

Figure 1 displays the numerical domain used in this work, which is based on the experimental setup of Jainski et al. [14]. The dimensions are given on the left side of this figure. On the right side, an instantaneous velocity field of the center plane with the flame front is shown in the three-dimensional domain. The turbulent inlet flow and acceleration of the fluid due to the flame front is clearly visible.

**Figure 1.** Numerical domain of the simulation. (**a**) dimensions of the domain; (**b**) snapshot of the 3D LES with velocity magnitude together with the flame front (orange iso-surface, *ω*˙ *CO*<sup>2</sup> = 60 kg/(m3 s)).

Like in the experiments [14], the wall temperature is kept constant at *Tw* = 350 K and the rod temperature is set to 300 K. The numerical rod is fully resolved with an O-grid. To stabilize the stoichiometric methane main flow, it is surrounded by a pure air co-flow (see Figure 1a). The inlet temperatures of the main and co-flow are 300 K. The turbulent inlet flow was calculated in a pre-processing step by simulating the experimental burner geometry. Each time step, a turbulent velocity field is set at the inlet. How this was calculated will be explained in the next paragraph. At the outlets, convective boundary conditions are used. A block structured grid with approximately 21.5 million control volumes is applied. Distant from the wall, the spatial resolution ranges from 200 μm to 400 μm. Towards the wall, the mesh is refined to 90 μm. As discussed in [17] this grid leads to an direct numerical simulation (DNS)-like spatial resolution and no additional wall modeling is needed.

The complete validation and verification of the simulation of the non-reacting turbulent inlet flow can be found in Heinrich et al. [17]. Here, the basic idea will be pointed out. Therefore, the body of the experimental burner was simulated (Figure 2a) and the experimental mass flows for the main and co-flow were set as inlet conditions [14].

**Figure 2.** Characterization of the simulated burner geometry used for the turbulent inlet profiles. (**a**) sketch of the burner geometry. The boundary conditions are color-coded; (**b**) velocity profiles 11 mm after the nozzle exit.

The simulated geometry includes the main and co-flow, the perforated plate, the flame stabilizing rod and the wall. The height of the nozzle exit is marked in Figure 2a. The velocity fields were extracted along this area for each time step and stored in a database. Within the burner, the turbulence was generated by the perforated plate. After the nozzle exit, a boundary layer developed at the wall. The Reynolds number based on the nozzle exit is *Re* = 5300. The resulting turbulent intensity at the nozzle exit section is approximately 8%. The analyses of the flow field revealed that the structures are inhomogeneous. In Figure 2b, the time averaged velocity and its rms-value are shown for the main component 11 mm downstream the nozzle exit. Both quantities are in very good agreement with the experimental values, since the maximum difference between the profiles is less than 0.18 m/s.

#### *2.3. Chemistry Treatment*

Depending on the quantities of interest, a very detailed description of the near-wall chemistry may be required. However, the application of full reaction mechanisms, typically covering tens of species and hundreds of elementary reactions, is not possible in complex configurations. Even for the laboratory scale system considered here, the computational cost given by the resolution requirement, the number of transported scalars and their stiff coupling prohibit its application. Accordingly, a chemistry reduction approach by means of tabulation is employed. In the framework of the Flamelet Generated Manifolds (FGM) approach introduced by van Oijen and de Goey [23], pre-computed detailed chemistry solutions of premixed flames enter the chemistry database. By that, vital information is included such that important quantities like the flames' propagation speed or its spatial structure including all species can be exactly reproduced by means of the controlling variables. The table generation, its assessment in a FWI situation, and its treatment in turbulent flow are briefly outlined in the following.

#### 2.3.1. Construction of the FGM Table

A three-dimensional chemistry table is used where the physical processes of this configuration being the reaction, the heat transfer to the wall and the mixing with the co-flowing air are accounted for by a reaction progress variable *YCO*<sup>2</sup> (being the mass fraction of CO2), the enthalpy *h*, and the mixture fraction *Z*, respectively. The table is illustrated in Figure 3b where the temperature and equivalence ratio have been chosen as the axis for illustrative purposes. As shown on the left, first a slice is generated by joining flamelets obtained with the GRI3.0 reaction mechanism [24,25] for different enthalpy levels. Hereby, region a) can be covered using different preheating temperatures while the burner stabilized method [26] is utilized in region b) where chemical reactions evolve under significant heat losses as present at the wall. Finally, in region c), the absence of activation energy suppresses most of the chemical activity and it mostly corresponds to a cooling of exhaust gases. Accordingly, the thermo-chemical states in this region can no longer be obtained by detailed chemistry simulations and Ketelheun et al. [27] approximated them by a thermo-dynamically consistent extrapolation. To cover mixing (i.e., different equivalence ratios), several of those slices are then combined to build the final three-dimensional table as given on the right of Figure 3b. Here, the source term *ω*˙ CO2 has been added representing one of the quantities extracted during the lookup in the simulation.

**Figure 3.** Representation of the table used. (**a**) cut at *φ* = 1 with three flamelet types, a) freely propagating, b) burner stabilized, c) extrapolated. (**b**) 3D table with a cut at *φ* = 1 and CO2 source term *ω*˙ *CO*<sup>2</sup> .

#### 2.3.2. Evaluation of Simplifying Assumptions for FWI Application

The above-mentioned tabulation approach has been successfully applied to turbulent laboratory scale combustors [28–30] where the FWI is not a crucial part of the physics. With respect to the latter, some works have been conducted to evaluate the accuracy of FGM near cold walls. Specifically, Meier et al. [31] and Ganter et al. [32] performed comparisons to detailed chemistry simulations. It turned out that FGM can generally satisfactorily predict the flame attachment position and the transported major species as well as the temperature profiles. Pollutants such as the accumulation of carbon-monoxide near the wall could not be predicted by FGM (Ganter et al. [33]). They showed that this is also possible by advanced or at least more specific tabulation techniques, which, however, require further development. Based on the analysis, we consider that the FGM approach is sufficiently accurate for the quantities analyzed in this study.

Specific for this work and the assumptions made, Figure 4 shows the evaluation of the employed FGM in a generic, two-dimensional SWQ configuration. Detailed chemistry simulations represent the verification basis for the tabulation in general and also to justify diffusivity assumptions. Regarding the latter, within the table generation, the unity Lewis number assumption is employed, since it simplifies the mapping process due to having a constant enthalpy and mixture fraction throughout a flamelet. The under prediction of the flame speed caused by this assumption (about 30% for the relevant conditions) is then corrected afterwards by consistently altering the chemical sources and diffusivities which accelerates the flame according to classical flame theories [34] while preserving its structure [35].

**Figure 4.** Comparison of different chemistry treatments: DC with mixture averaged diffusion (dashed lines), DC with *Le* = 1 assumption (solid lines) and FGM (dotted lines). (**a**) temperature field of the DC simulation. The black lines represent iso-lines of *T* = 1700 K and the white lines of *ω*˙ *CO*<sup>2</sup> = 60 kg/(m2 s); (**b**) the colors of the lines represent the extracted position from (**a**). The top figure shows the corresponding temperature lines and the bottom figure the mass faction of CO2.

First, Figure 4a shows the temperature field in the region around the quenching point for orientation. Added to this graphs are lines marking the flame by means of an isothermal as well as the chemical source term obtained by three different simulations. The first (dashed line) is a detailed chemistry simulation employing individual diffusivities for all species and can be considered as the most accurate physical description. The Second (solid line) is also a detailed chemistry simulation but with the *Le* = 1 assumption. Finally, the dotted line is the FGM solution. Due to the above-mentioned correction of the flame speed, the flame–wall attachment points and flame angles are in close agreement. Slightly upstream of the quenching point (coordinate zero) the individual diffusivities cause the flame to reach slightly closer to the wall while the FGM solution is close to the *Le* = 1 detailed chemistry indicating that these deviations are due to the remaining transport simplifications. Downstream of the quenching point, the isothermal of the detailed chemistry simulations then approach each other while the FGM slightly depart from them. Accordingly, and as also analyzed by Ganter et al. [33], this is due to tabulation assumption regarding diffusive fluxes orthogonal to the controlling variables. As mentioned, correcting this would require an immensely increased tabulation effort being currently not yet applicable to this configuration. Considering the rather small errors this seems acceptable. A quantification of the temperature and species evolution is furthermore provided on the right of Figure 4 where profiles are shown extracted along the horizontal lines added on to the left of Figure 4. Here, the blue, magenta and black lines represent the states at the quenching point and one millimeter up- and downstream of it. All of these positions being in the region of intense FWI confirm a very good approximation by the FGM approach for the purpose of this work.

#### 2.3.3. Treatment in the Turbulent FWI Flow

As mentioned above, the thermo-chemical state is described with three transported quantities, the mixture fraction *Z*, the reaction progress variable *Ypv* and the enthalpy *h*:

$$
\frac{\partial \overline{\rho} \mathcal{Z}}{\partial t} + \frac{\partial}{\partial \mathbf{x}\_{\dot{j}}} \left( \overline{\rho} \widetilde{\mu}\_{\dot{l}} \widetilde{\mathcal{Z}} \right) = \frac{\partial}{\partial \mathbf{x}\_{\dot{j}}} \left( \left( \mathcal{F} \mathcal{J} \frac{\overline{\mu}}{\text{Sc}} + (1 - \Omega) \frac{\mu\_{l}}{\text{Sc}\_{l}} \right) \frac{\partial \mathcal{Z}}{\partial \mathbf{x}\_{\dot{j}}} \right) \,, \tag{3}
$$

$$\frac{\partial \overline{\rho} \tilde{Y}\_{\text{PV}}}{\partial t} + \frac{\partial}{\partial \mathbf{x}\_{\dot{j}}} \left( \overline{\rho} \overline{u}\_{\dot{j}} \tilde{Y}\_{\text{PV}} \right) = \frac{\partial}{\partial \mathbf{x}\_{\dot{j}}} \left( \left( \mathcal{J} \mathcal{J} \frac{\overline{\mu}}{\text{Sc}} + (1 - \Omega) \frac{\mu\_{t}}{\text{Sc}\_{t}} \right) \frac{\partial \tilde{Y}\_{\text{PV}}}{\partial \mathbf{x}\_{\dot{j}}} \right) + \frac{\mathcal{J}}{\mathcal{J}} \dot{\omega}\_{\text{PV}} \tag{4}$$

$$\frac{\partial \overline{\rho} \overline{h}}{\partial t} + \frac{\partial}{\partial \mathbf{x}\_{\dot{j}}} \left( \overline{\rho} \overline{\mathbf{u}}\_{\dot{j}} \tilde{h} \right) = \frac{\partial}{\partial \mathbf{x}\_{\dot{j}}} \left( \left( \mathcal{J}^{\mathbf{x}} \mathcal{J} \frac{\overline{\mu}}{\text{Sc}} + (1 - \Omega) \frac{\mu\_{l}}{\text{Sc}\_{l}} \right) \frac{\partial \overline{h}}{\partial \mathbf{x}\_{\dot{j}}} \right) \tag{5}$$

where *μ* is dynamic viscosity, *Sc* the Schmidt number, F the thickening factor, E the efficiency function, Ω the flame sensor and *ω*˙ pv the reaction source term. To capture the flame–turbulence interaction correctly, the artificially thickened flame (ATF) approach is used. In this context, the flame has been thickened with F to resolve its structure on the LES grid. Therefore, the transformation of Butler and O'Rouke [36] is used, where the flame speed is preserved. A dynamic version of the procedure based on a sensor Ω is used in this work [37,38]:

$$
\mathbb{A}\mathcal{F} = 1 + (\mathcal{F}\_{\text{max}} - 1)\Omega \,, \tag{6}
$$

$$
\Omega = 16 \left( \chi\_{\rm n,pv} (1 - \chi\_{\rm n,pv}) \right)^2 \,\, \, \, \, \, \tag{7}
$$

where Fmax is the grid depending maximum thickening factor and *Y*n,pv is the normalized progress variable:

$$Y\_{\rm R,pv} = \frac{Y\_{\rm PV} - Y\_{\rm PV,min}}{Y\_{\rm PV,max} - Y\_{\rm PV,min}} \,. \tag{8}$$

Due to the small thickening factor Fmax of 1.4 in this work, the flame turbulence interaction is only slightly altered. Since ATF is used for the flame–turbulence interaction, the turbulent diffusion is switched of in the region of the flame. As explained in [16], the influence of the dynamic thickening on relevant quantities, like wall heat flux, the flame position and quenching point determination is negligible because the thickening can be retransformed. In general, the interaction with the turbulence is reduced if the flame is artificially thickened. To compensate these effects, Colin [35] used a so-called efficiency function E . The formulation of the implemented efficiency function E can be found in Charlette et al. [39]. Within this work, E is nearly negligible because the grid resolution is fine enough and Fmax is relatively small. Following the classical turbulent premixed combustion diagram from Peters [40], this setup operates in the wrinkled flamelet regime, meaning that the inner flame structure is close to a laminar flame and the thin flame front is wrinkled by the turbulent motion due to the fact that the flamelet concept is valid [40]. Additionally, the effects of curvature and stretch on the mass burning rate have been investigated for the thin reaction zone regime [41] and the corrugated flamelet regime [42]. The resulting errors have been found to be within 5% when using a single progress variable as done in this work. Therefore, FGM can be applied for the considered configuration.

#### *2.4. Determination of the Quenching Point*

Following Heinrich et al. [16], the quenching point *yqp* is determined by tracking the maximum value of the CO2 source term, which defines the flame front. Due to the enthalpy losses at the wall, this source term decreases and if it drops under a certain threshold the flame is quenched. The threshold is defined as the half of the maximum CO2 source term along the flame front. This procedure is analogous to the experimental determination where the OH gradient was used. For this FGM approach, using the CO2 source term is validated against detailed chemistry simulations and it appears to be the most reliable method [16].

#### *2.5. Comparison with Experimental Data*

As mentioned before, the models used lead to good results compared to the experimental laminar and turbulent SWQ burner. This was shown before by Heinrich et al. [16,17]. Within these references, the comparisons of the laminar and turbulent for the non-reacting and reacting configurations can be found. For the analyses conducted in Section 3, the velocity, flame position and quenching points are of importance. For the sake of completeness, these quantities are shown again in Figure 5.

**Figure 5.** Comparisons of the turbulent data with experiments at *z* = 60 mm in the conditioned coordinated system. (**a**) time averaged velocity profiles 5 mm upstream the quenching point; (**b**) averaged flame position; (**c**) quenching probability.

Firstly, a conditioned coordinated system has to be introduced (*x*q, *y*q) with respect to local averaged quenching point is introduced. It is centered at the wall with *x*<sup>q</sup> = 0 mm and *y*<sup>q</sup> = 0 at the averaged axial quenching point (exemplary marked with the black solid circles in Figure 6). The experimental data is taken from Jainski et al., 2016 [14]. The velocity profiles 5 mm before the time averaged quenching point are shown in Figure 5a. Both directions recover the experimental velocity very well. The difference in the main flow direction is less than 0.17 m/s, except for the point at the wall. Clearly, at the wall, an experimental inaccuracy occurred because the the velocity is not zero. In the wall normal direction, the maximum difference is slightly higher (0.3 m/s). The averaged flame positions (Figure 5b) are in excellent agreement. The maximum deviation is in order of the grid size, 0.2 mm. Figure 5c depicts a statistical quantity, the quenching probability. If it is zero, the flame is always burning and, if it is one, the flame is always quenched. Both profiles match perfectly. It shows that quenching events takes place in the region −8 mm < *yqp* < 8 mm.

#### *2.6. Dimensionless Parameters*

Two dimensionless parameters used for the analyses have to be explained. The Péclet number Pe describes the ratio of the wall normal distance of the flame and a typical flame thickness:

$$\mathrm{Pe} = \frac{\mathfrak{x}}{\delta}.$$

The diffusive flame thickness is used as reference to define the Péclet number at quenching Pe*q* [5]:

$$\mathrm{Pe\_q} = \frac{x\_{\mathrm{q}}}{\rho\_{\mathrm{u}} c\_{p\_{\mathrm{u}}} S\_L / \lambda\_{\mathrm{u}}} \frac{1}{f\_{\mathrm{atf}}},\tag{9}$$

where *ρ* is the density, *cp* is the specific heat capacity at constant pressure, S*<sup>L</sup>* is the laminar flame speed and *λ* refers to the thermal conductivity of the gas. The occurrence of *f*atf is due to the ATF transformation as detailed in Heinrich et al., 2017 [16]. The subscript *q* refers to values at the quenching point and the subscripts *u* and *b* refer to unburnt and burnt states. The wall heat flux Φ*w*,*<sup>q</sup>* is made dimensionless with the flame power *PQ* = *ρuSLcp*<sup>u</sup> (*T*<sup>b</sup> − *T*u) [5]:

$$F\_{\rm Q} = \frac{\Phi\_{w,q} \cdot f\_{\rm att}}{P\_{\rm Q}},$$

where *T* is its local temperature. The applied parameters are given in Table 1.

**Table 1.** Characteristic values from the CHEM1D calculation for a premixed freely propagating flamelet with 300 K inlet temperature, stoichiometric conditions and mixture averaged diffusion.


#### **3. Results**

To provide a proper understanding, this section is divided into three parts. First, a mostly descriptive insight into the transient FWI is given to provide a basic understanding on the typical physical behavior and the associated quantities like the heat flux. Therefore, first, Section 3.1 illustrates the basic evolution of the FWI depending on the type of the quenching point displacement by means of individual, representative events. The observations are then shown in Section 3.2 to also hold on a statistical basis where further observations are provided. Finally, Section 3.3 outlines the root causes by considering the full interaction of the flame, the wall and the turbulent field (FWTI). It is exemplified by following the mechanisms leading to very large heat fluxes associated with a departure from the SWQ towards the HOQ regime.

The averaged quenching points along the wall are shown as a red solid line in Figure 6 together with all quenching points from one snapshot (black solid line). These lines are called quenching lines. As it can be seen, the averaged quenching line is not constant and, for an instantaneous moment, the flame can quench far away from this line. In this study, three independent planes will be used for the transient analysis of the FWI, at the center (*z* = 60 mm), at *z* = 40 mm and at *z* = 80 mm (blue dashed lines). These planes are uncorrelated for the quenching point statistics since the correlation-coefficient drops under zero at a length of 10 mm. Since the averaged quenching line is not constant in the *z*-direction, the conditioned coordinate system (*x*q, *y*q) introduced before is used. With this transformation, the three extracted lines can be combined to one dataset.

**Figure 6.** Time averaged quenching points along the wall (red solid line) with one instantaneous quenching line (black dashed-dotted line) and the three extracted planes marked (blue dashed lines).

#### *3.1. Description of the Quenching Point Movement*

As pointed out by Heinrich et al. [17], the quenching positions, heat fluxes and wall angles are widely distributed. In this work, it will be shown that most of these statistics depend on the way the flame tip moves. The movement of the flame tip is tracked via the temporal quenching point position. A sequence from the center plane is shown in Figure 7, where the downstream and upstream movement are colored in green and red, respectively.

**Figure 7.** Temporal section of the axial quenching point location with colored downstream and upstream movement (green and red). The interval between the points is 0.1 ms.

Clearly, this movement seems to be chaotic and is not equally distributed. In 71% of the events, the flame tip is moving in the downstream and in 27% in the upstream direction. Consequently, the averaged upstream velocity of the flame tip must be faster than its averaged downstream velocity. The remaining 2% are a jump-like upstream movement, with multiple quenching points along the wall. This leads to the following categorization:


The following paragraphs together with the Figures 8–10 describe the scenarios in more detail. The earliest time step is shown at the top and the latest on the bottom. The individual time steps are organized as follows (from the left to the right):


The explanations and statistics for scenarios are provided in the next sections.

#### (Scenario A)

The moderate downstream movement of the flame tip appears in nearly 71% of the cases, which makes it the most frequent event. How the flame and the wall heat flux are connected is explained in conjunction with Figure 8.

**Figure 8.** (Scenario A): Consecutive snapshots (**a**–**c**) showing the time sequence of the downstream movement of the quenching point. From the left to the right are the wall heat flux with all quenching points shown, then Φ<sup>w</sup> along the wall at *z* = 0 mm, including the position of the maximum wall heat flux point *y*(Φw) (blue dashed line) and *yqp* (red dashed line). The next graph displays a 2D slice (*x*-*y*-plane) with T-isolines at z = 0 mm with the flame position (black line) overlaid with the velocity vectors and *y*(Φw) (blue dashed line) and *yqp* (red dashed line). The last graph shows the temperature field including the magnitude of the velocity and its vector for two positions at z = 0 mm (one at *yqp*, one at *yqp* − 4*δT*, where *δ<sup>T</sup>* is the thermal thickness of the flame).

The quenching point *yqp* and the location of the maximum wall heat flux *F*Q,max move together in an axial direction, while the heat flux constantly decreases. As visible, *F*Q,max is downstream of the quenching point. The difference between both locations can be up to 12 thermal flame thicknesses. Within this movement and upstream of *yqp*, mostly cold fresh gases flow parallel to the wall. This can be seen at the temperature field and at the extracted profiles. As visible on the right graphs, directly at the wall and 2 mm before the quenching point, the fluid is not heated by the flame. Compared to the laminar case, the heat fluxes are rather low and the quenching distances are large (Φw,l = 200 kW/m2, *Pe*q,l = 7.3).

#### (Scenario B)

The next scenario is depicted in Figure 9.

**Figure 9.** (Scenario B): Time sequence of the upstream movement of the quenching point. The explanation of the figure can be found in Figure 8.

This moderate upstream movement with increasing maximum wall heat flux occurs in 27% of the events. The maximum heat fluxes and the quenching point are closely together, like in the laminar case, while the heat flux is mostly higher than in (scenario A) or than the laminar reference value Φw,l = 200 kW/m2. Within this kind of movement, only high heat fluxes and low *Pe*<sup>q</sup> can be found. In contrast to (scenario A), no explicit temporal behavior of Φw,max can be distinguished; within this movement, it can increase or decrease. Additionally, it can be seen from the temperature profiles on the bottom at right side of Figure 9 that the fluid upstream of the quenching position is heated by the flame.

#### (Scenario C)

The last scenario is found only in 2% of the cases. However, it cannot be neglected because only the strongest heat fluxes arise there. It is also an upstream movement like (scenario B), but it appears

to look like a jump of the quenching point. Figure 10 will be used to explain the behavior. At the beginning (a), the flame front is nearly parallel to the wall and moves downstream like in (scenario A), but the wall heat flux upstream the quenching point increases at a large region. At the same moment, the flame gets weakened where the heat flux rises and the enthalpy near the wall decreases in the same area (not shown here). A widespread quenched region emerges (Figure 10b) and the heat flux reaches extreme values (up to *F*Q,max/*F*Q,max,*<sup>l</sup>* = 1.7). In this moment, three quenching points occur in the center plane and the moving direction of the flame tip flips. Since the quenching point most upstream is always used as the reference point, it seems like it jumps upstream. As it can be seen in the left graph, in the next moment (c), both quenching areas will merge together and only one quenching point will remain. During this process, the occurring values of *F*Q,max are only the highest values. Afterwards, the flame tip moves downstream again with a lower heat flux, being the transition to (scenario A).

**Figure 10.** (Scenario C): time sequence of the upstream jump of the quenching point. The explanation of the figure can be found in Figure 8.

#### *3.2. Statistics of the Axial Quenching Position and the Maximum Wall Heat Flux*

As shown in Figure 11a, the locations of the maximum wall heat flux *y*(*F*Q,max) (blue solid line) and the quenching point *yqp* (red dashed line) are at the same position in a laminar SWQ flame. The black crosses correspond to the dimensionless wall heat flux profile (top *x*-axis) and the magenta line displays the iso-line *ω*˙ *CO*<sup>2</sup> = 60 kg/(m2s), which envelops the flame front. These profiles are shown in Figure 11b for one instantaneous moment of the turbulent SWQ configuration. As already observed by Heinrich et al. [17], in a turbulent setup, *yqp* and *y*(*F*Q,max) do not always coincide. This phenomenon is visualized with the black arrow in Figure 11b. The quantity Δ(*y*(*F*Q,max), *yqp*) = *y*(*F*Q,max) − *yqp* expresses the distance between the axial locations. A positive distance implies that the maximum heat flux is farther downstream. *F*Q,max and *yqp* are spacial quantities that change during time.

**Figure 11.** Flame position with marked quenching point (*yqp*, red dashed line) together with the dimensionless wall heat flux (black crosses, top *x*-axis) and its maximum marked (*y*(*F*Q,max), blue solid line). The magenta line corresponds to *ω*˙ *CO*<sup>2</sup> = 60 kg/(m<sup>2</sup> s), which envelops the flame front. (**a**) laminar; (**b**) turbulent: the black arrow shows the difference between *yqp* and *y*(*F*Q,max).

First, the temporal distribution of *y*(*F*Q,max) and *yqp* are analyzed in Figure 12a,b.

For both moving directions, *yqp* is centered at zero and normal distributed (magenta dashed lines). However, the mean value of *y*(*F*Q,max) is zero for the upstream movement with equal number of counts the positive and negative side but with larger positive range. In contrast to that, the downstream direction has a mean value around 2 mm and more events occur on the positive side. The moving direction does not affect *yqp*, but *y*(*F*Q,max) because it is shifted for the downstream movement, which could be seen in Figure 8 in the graphs in the middle. The same observation can be made by looking at Figure 12c,d, where Δ(*y*(*F*Q,max), *yqp*) is shown for all movements. In the following, most of the figures include similar histograms. A solid black line represents the mean value of the shown data and the dashed black the laminar reference. A cumulative histogram in light blue is overlaid. As it can be seen, within the downstream movement, for more than 60% of all events, the maximum wall heat fluxes are farther downstream. As a result, the mean difference is 1.8 mm. Like expected from the opposite direction (Figure 12b) and also visible in Figure 12d, Δ(*y*(*F*Q,max), *yqp*) is very small (80% are in the range |Δ(*y*(*F*Q,max), *yqp*)| ≤ 1 mm). Only a few outliers occur. The jump-like movement was not considered here because it is a very fast process whereby the determination of the local detachment is error-prone. Figure 12 showed that the shift between *yqp* and *y*(*F*Q,max) is not constant and is mainly caused by the downstream movement of the flame tip.

**Figure 12.** Relation between *yqp* and *y*(*F*Q,max) in dependency of the flame tip movement. (**a**,**b**): *yqp* with corresponding normal distribution (magenta dashed line) and *y*(*F*Q,max) with corresponding mean value (black solid line). (**c**,**d**): Mean values (black solid lines) and laminar reference (black dashed line). (**a**) Downstream. (**b**) Upstream. (**c**) Downstream. (**d**) Upstream.

Figure 13 displays the probability distribution of the maximum dimensionless wall heat flux. All three scenarios show a different behavior. If the flame tip moves downstream, only low values appear: nearly 90% are lower than the laminar reference. On the contrary, in the upstream scenario, 80% are higher than the laminar value. Accordingly, the mean values for these two cases are 0.155 and 0.27, respectively. The jump-like movement is most extreme; here, in 99% of the events, *F*Q,max is higher than the laminar scenario; accordingly, mostly the highest observed values appear. The distribution of the mean value of *F*Q,max depends on the movement of the flame tip; this was also shown in the representative snapshots in Figures 8–10.

**Figure 13.** Dependency of *F*Q,max from the flame tip movement with its mean value (black solid line) and the laminar reference value (black dashed line). (**a**) downstream; (**b**) upstream; (**c**) jump-like.

Another quantity, which was not discussed before, is the dimensionless heat flux at the quenching point *F*Q,q. Since it has a similar qualitative behavior to *F*Q,max, its distribution will not be shown here. The difference is that the values for *F*Q,q are a bit lower and the maximum values only appear at the upstream and jump-like movement. At the upstream movement, nearly 70% of the events are higher than the laminar reference with a mean value of *F*Q,q = 0.255. The opposite applies for the downstream direction, almost 95% are smaller, with a mean value of 0.12. Both distributions (*F*Q,q and *F*Q,max) showed that high heat fluxes appear mostly while the flame tip moves upstream and the lowest values if it moves downstream.

As known from the literature (e.g., [5,17]), there is a negative correlation between *Pe*<sup>q</sup> and *F*Q,q. The cross-correlation between *Pe*<sup>q</sup> and *F*Q,q is very strong, the correlation coefficients are −0.87 and −0.97 for the downstream and the upstream movement, respectively. With this in mind, the distribution of *Pe*q Figure 14 is not unexpected.

**Figure 14.** Dependency of *Pe*q from the flame tip movement with its mean value (black solid line) and the laminar reference value (black dashed line). (**a**) downstream; (**b**) upstream; (**c**) jump-like.

For the downstream direction, almost 90% of the events are higher than the laminar reference and the mean value is about *Pe*<sup>q</sup> = 12. This means that the flame is much farther away from the wall. The smallest possible wall distances do not occur. In contrast to that, only small distances occur for the upstream movement and its mean value *Pe*<sup>q</sup> = 7.5 is slightly below the laminar case. *Pe*<sup>q</sup> has a narrow distribution centered around its mean value with over 50% of the events in its proximity. The distribution for the jump-like movement is qualitatively similar to the upstream direction and the mean value *Pe*<sup>q</sup> = 7.5 is slightly lower than the laminar value. Certainly, the most events occur around *Pe*<sup>q</sup> = 6.5.

#### *3.3. Physical Mechanism Governing the FWTI*

Figure 14 displayed an obvious dependency of *Pe*<sup>q</sup> from the moving direction and, with its help, the occurring wall heat fluxes can be explained, since the flame can get close to the wall and a strong temperature gradient occurs. The resulting heat fluxes are very high or the opposite can occur, the flame is far away and the heat fluxes are much smaller.

The behavior for the downstream and upstream distributions can be explained with the temperature field. For that reason, time averaged temperature fields with respect to the conditioned coordinate system (*x*q, *y*q) are shown in Figure 15b,c.

Figure 15a shows the laminar temperature field for comparing purpose. Compared to the laminar case, in the downstream movement, more cold unburnt fluid is enclosed between the flame and the wall, which pushes the flame farther away from the wall. The consequence of this was shown before, the resulting heat flux is smaller. While the flame tip moves in an upstream direction, the temperature field is very close to the laminar case as already shown in Figure 14b. In addition, the flame shapes are different, depending on the movement. In case of the downstream movement, the curvature is stronger.

For a laminar SWQ flame, Φ<sup>w</sup> is always smaller than in the corresponding HOQ flame [5], but, as shown before, in the turbulent case, the wall heat flux can clearly exceed the laminar reference values up to values of the HOQ case or even higher. Due to that, the turbulent case will be compared

qualitatively and quantitatively with an HOQ flame. Furthermore, a one-dimensional HOQ flame is calculated with FGM for this purpose. It has the same mixture properties as the main flow of the SWQ configuration and the same wall temperature. The one-dimensional grid used has the same resolution as the wall normal direction of the 3D SWQ case. As seen before in Figure 13, the jump-like movement has only the highest *F*Q,max values. This can be explained by analyzing this movement in a different conditioned coordinate system. Therefore, the location of the maximum wall heat flux *F*Q,max defines the zero on the *y*-axis. Figure 16 displays the time averaged temperature fields for all jump-like movements and the laminar reference.

**Figure 15.** Dependency of the temperature field from the flame tip movement within the conditioned coordinate system together with the flame (black line). The white dashed line in (**b**,**c**) corresponds to the laminar reference. (**a**) laminar case; (**b**) downstream; (**c**) upstream.

**Figure 16.** Temperature field in the conditioned coordinate system, where *y*(*F*Q,max) marks zero. (**a**) laminar case; (**b**) jump-like.

As it can be seen on the right side, for the jump-like movement, the flame is nearly parallel to the wall, like in an HOQ manner. In addition, the temperature field is much more compressed in the direction of the wall, which causes the high heat fluxes. For a one-dimensional HOQ flame, the velocities near the wall are nearly zero at the moment the flame quenches. The same can be seen for the jump-like scenario in Figure 17.

**Figure 17.** Wall normal velocity in the conditioned coordinate system, where *y*(*F*Q,max) marks zero. (**a**) laminar case; (**b**) jump-like.

Compared to the laminar case, where the flame is clearly an SWQ flame, the wall normal velocities are much smaller. Due to that, the flame can burn much closer to the wall because the convection away from it is much weaker.

Since the turbulent SWQ flame front is transient and three-dimensional, a point of interest has to be tracked in time and space. How it is done will be explained with Figures 18–20.

**Figure 18.** 3D flame front for the laminar case. (**a**) 3D flame front with wall (black surface); (**b**) zoom with monitored point (black dot).

A section of the 3D laminar flame is shown in Figure 18. On the left side, the wall is depicted for orientation purpose. The figure includes a temperature iso-surface of *T* = 1750 K and a slice of the wall normal velocity component at the height of the maximum wall heat flux at the center of the burner. On the right side, a close-up view of the same section without the wall is displayed and the monitored quenching point is marked (black dot). Additionally, the momentum vectors are projected onto the slice. As it can be seen, the momentum vectors are all aligned and are pointed perpendicular away from the wall. The fluid is accelerated away from the wall by the flame without any redirection. Due to this convection, the flame cannot get closer to the wall.

The temporal evolution of the tracked point in turbulent case is shown in Figures 19 and 20, (black dots). It follows the path of a point on the flame front, which leads to the maximum wall heat flux at *t* = 0 ms.

**Figure 19.** Temporal evolution of the 3D flame front near the wall, part 1.

**Figure 20.** Temporal evolution of the 3D flame front near the wall, part 2.

As visible in the first sub-figure (top, left), the flame front is shaped convexly towards the unburnt gases. In contrast to the laminar case, the momentum vectors, which are not in the proximity of the tracked point, are directed away from the point into positive and the negative *z*-direction. During this progress, the wall normal velocity decreases and the convex flame shape gets flattened. Due to the convection of the main flow in the axial direction, the tracked point moves upstream while it stays more or less at the same *z*-position (which does not look like it is due to the chosen angle of view). By following this point, the graphs in Figures 21 and 22 were created.

The temporal evolution of the tracked point from the turbulent SWQ setup, an HOQ flame and the laminar reference case are shown in Figure 21. It displays from the left to the right: the wall heat flux, the normal velocity and the position of the flame.

**Figure 21.** Temporal evolution of typical flame quantities for different flame types. (**a**) wall heat flux; (**b**) wall normal velocity at *T* = 1750 K; (**c**) position of the flame at *T* = 1750 K.

The time point *t* = 0 ms marks the moment with the highest wall heat flux. The position and the velocity of the flame front are extracted, where the temperature is equal to 1750 K. This location coincides with the position of the maximum source term of the flame and therefore it is associated with the center of the flame.

At first, the 1D HOQ case will be illustrated. At the beginning, the heat flux rises up to its maximum value (at *t* = 0 ms), while the wall normal velocity is positive and starts to decrease. At the same time, the flame moves towards the wall. It can been seen that the velocity is nearly zero when the highest heat flux appears, which is higher than the laminar SWQ case (black asterisk). After that, the velocity drops under zero, it stays negative, and the heat flux decreases again. The relation of the velocity and the temperature for the SWQ cases and the HOQ flame at characteristic points in time are plotted in Figure 22.

**Figure 22.** Temporal evolution of the flame. The solid line marks the HOQ and the dashed line the SWQ case.

This kind of representation decouples the dependency of the spacial coordinate. Before the flame interacts with the wall (*t* = −1.4 ms, solid red line), the flame has a linear dependency, as the temperature increases the velocities increases, too. Right before the flame quenches (*t* = −0.2 ms), the velocity over the whole flame is reduced. Like seen before in Figure 21, at the moment of quenching, the velocity is very small. The high heat fluxes appear because the velocities between the flame and

the wall decreases if the flame reaches the wall (Figure 21b,c), due to the fact that the flame can closely approach the wall, which leads to high temperature gradients. As already mentioned and known from the literature, the laminar wall heat flux is lower and the quenching distance is greater, which can be seen in Figure 21a,c. Due to the higher wall normal velocity in the laminar case, the convection away from the wall is greater, which explains why this flame is farther away from the wall (Figure 21b). The *T* − *ux* plot (Figure 22) reveals that the velocity in the laminar case is much greater than zero over most of the flame. Finally, the behavior of the turbulent SWQ flame will be discussed. The temporal evolution of the heat flux, the wall normal velocity and the flame position are also shown in Figure 21. Qualitatively, the shapes are similar to the one-dimensional HOQ flame. Certainly, the increase of the wall heat flux and the decrease of the velocity start earlier in the turbulent case because this flame is closer to the wall. The behavior right before the quenching occurs (*t* = −0.2 ms) is different too because the turbulent SWQ flame moves much slower. Furthermore, the linear dependency of the velocity from the temperature vanishes much faster for the SWQ flame; in the moment before quenching, it is gone. As it can be seen, at the moment of quenching (*t* = 0 ms), the velocity of the HOQ is nearly zero over the whole temperature range, while, for the SWQ flame, the velocity is very small in the range of *T* = 300 K–*T* = 1750 K. After quenching, both flames have the same behavior for the region between the wall and the center of the flame. The drop of the velocity is also visible in the 3D plots in Figures 19 and 20. As observed in these figures, the fluid moves away from the tracked point and, as a result, the wall normal velocity gets lower than in the laminar case. During this process, the convex shaped flame front gets flattened. As a result, the turbulent flame behaves locally like an HOQ flame because the effects of the velocity are getting negligible. The maximum wall heat flux is Φw,max = 311 kW/m<sup>2</sup> and the minimal wall distance is *x*min = 0.475 mm, which are very close to the HOQ flame (Φ*w*,*max*,*HOQ* = 313 kW/m2, *x*min,HOQ = 0.417 mm). After the turbulent flame reached the maximum Φw,max (*t* > 0 ms), the profiles from Figure 21 are nearly the same again. However, a difference can be seen in the *T* − *ux* plot (Figure 22), for the region far away from the wall, the velocity of the turbulent flame is much greater because the fluid of the SWQ flame can accelerate again.

With the knowledge gained above, the identified regimes are very natural and somehow obvious, but the characterization raises the question if also a jump-like downstream movement could appear. Therefore, pockets of fresh gases must appear downstream of the flame front inside the burnt region. In our configuration, this seems to be very unlikely and was not observed, but it could be possible e.g., if large eddies transport fresh gases into the burnt region. Due to the flame topology, these events would not be connected with high heat fluxes because the fresh gases do not preheat fluid near the walls.

#### **4. Conclusions**

This work was conducted to gain insight into the flame–wall–turbulence interaction. The analysis utilized results obtained from the simulation of a well-suited SWQ configuration where the statistical data was confirmed to be in good agreement with measurements. Building upon this predicting capability, the availability of the highly resolved, three-dimensional transient field gave access to the physical processes like the heat flux evolution that significantly depart from the laminar operation. First, several processes were identified that show a distinct phenomenological behavior associated with the quenching point displacement. Specifically, the three archetypes being the downstream, upstream and jump-like movements were shown to operate on individual time-scales and heat flux evolutions. This physical understanding was then supplemented with statistical evaluations confirming its generality. For this, the conditionally averaged data showed that the average downstream movement of the quenching point yields significantly lower heat fluxes than the laminar SWQ flame and vice versa for the upstream movement. Significantly larger heat fluxes were observed when the quenching point jumped. The root causes for this were subsequently identified by considering the full interaction of the turbulent field, the flame and the wall. It revealed that the near-wall flame wrinkling by turbulent structures locally creates a departure from the SWQ towards the HOQ regime. This latter is known for

conditions that enable the flame to approach the wall very rapidly causing a strong interaction with high heat fluxes. We showed that, by the curvature of the flame towards the wall, these conditions can even get surpassed by a flame-tangential divergence of the velocity field.

The identification of the mechanisms of FTWI was conducted based on simulations of a laboratory scale experimental SWQ burner. The configuration was chosen first, since its geometrical complexity is a reasonable compromise in between computational effort and practical relevance regarding the process of interest and, second, basic validation data are available. With respect to the generality of the results and their relevance for future combustor design, it should be mentioned that the main goal of this work was to show that simulations with the given approaches can handle such configurations and their results provide a sound basis for an understanding of the physical processes. Considering the phenomenological character of our analysis, it is certainly transferable for in terms of the fundamental FTWI mechanism. However, the parameter specific to a configuration will decide how they scale and accordingly quantitatively behave. In this regard, the Reynolds number is generally higher in real applications, intensifying the FTWI. Accordingly, the mechanism of regime alteration (SWQ-HOQ) will likewise be enhanced. Furthermore, for regimes of higher Karlovitz numbers, stretch might further modify the near wall flame structure, causing significantly different peaks of the heat release. Likewise, the specific fuel is of major importance e.g., with respect to the pollutants formed in the walls vicinity. Similarly, the wall treatment (e.g., effusion cooling) will have an influence. In this regard, we mentioned that the tabulation approach should be chosen with care depending on the quantities of interest.

**Author Contributions:** A.H. performed the simulations. A.H., G.K. and S.G. exploited the numerical data and analyzed the results. A.H. wrote the paper, while G.K. and C.H. further supported improving the manuscript. J.J. contributed by providing materials, computing resources and supervision.

**Funding:** The authors gratefully acknowledge the financial support by the Deutsche Forschungsgemeinschaft (DFG) through grants SFB/TRR 150 and in the framework of the Excellence Initiative, Darmstadt Graduate School of Energy Science and Engineering (GSC 1070).

**Acknowledgments:** All computations were performed on the Lichtenberg High Performance Computer (HHLR) at the TU Darmstadt.

**Conflicts of Interest:** The authors declare no conflict of interest.

#### **References**


c 2018 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (http://creativecommons.org/licenses/by/4.0/).

MDPI St. Alban-Anlage 66 4052 Basel Switzerland Tel. +41 61 683 77 34 Fax +41 61 302 89 18 www.mdpi.com

*Fluids* Editorial Office E-mail: fluids@mdpi.com www.mdpi.com/journal/fluids

MDPI St. Alban-Anlage 66 4052 Basel Switzerland

Tel: +41 61 683 77 34 Fax: +41 61 302 89 18