Turbulent Combustion

Turbulent Combustion models...

Under construction

ATF

The Artificial Thickening Method (ATF) approach essentially involves "thickening" the flame, allowing the mesh to resolve the scalar gradients across it. This serves a dual purpose: it enables capturing the flame dynamics and reduces numerical diffusion errors associated with resolving the sharp scalar jumps in the flame. The aim of ATF methods is to achieve adequate resolution within the flame front (i.e., the region with the sharpest gradients) on the transformed mesh, even when the coarse LES mesh lacks sufficient resolution. The geometry transformation ratio, denoted as $\mathcal{F}$, is referred to as the thickening factor.

$$\frac{\partial \rho Y_k}{\partial \tau} +\frac{\partial \rho u_i Y_k}{\partial \xi_i}= \frac{\partial}{\partial \xi_i }\left(\rho D \mathcal{F}\frac{\partial Y_k}{\partial \xi_i}\right) + \frac{ \dot{\omega}}{\mathcal{F}}$$

The transformation "thickens" the flame but maintains the correct flame speed by increasing diffusion accordingly.

The sub-grid wrinkling requires modelling and this is usually done in the ATF context by the inclusion of a so called efficiency function. The efficiency function, $E$, is defined by a dimensionless wrinkling factor $\mathcal{E}$, and its ratio between a laminar flame compared with its thickened counterpart. Inclusion into the ATF model results in a modification to the reactive scalar transport equation:

$$\frac{\partial \rho Y_k}{\partial \tau} +\frac{\partial \rho u_i Y_k}{\partial \xi_i}= \frac{\partial}{\partial \xi_i }\left(\rho D E\mathcal{F}\frac{\partial Y_k}{\partial \xi_i}\right) + \frac{ E\dot{\omega}}{\mathcal{F}}$$

$$E=\frac{\mathcal{E}\vert_{\delta_f=\delta_f^0}}{\mathcal{E}\vert_{\delta_f=\hat{\delta}_f^0}} \geq 1$$

The factor $\mathcal{E}$ essentially relates the the total flame front wrinkling with its resolved component. It can be approximated as:

$$\mathcal{E} \approx 1+\beta \Delta \vert <\nabla \vec{n} >_{sgs}\vert$$

where $\beta$ is a constant and $\vert <\nabla \vec{n} >_{sgs}\vert$ is the sub-grid surface curvature.

PasR

Broadly speaking, the Partially Stirred Reactor Model (PaSR) model assumes that within a computational cell (or sub-grid region), the reacting mixture is represented by a statistical ensemble of partially stirred reactors. The core idea is that chemical reactions and mixing occur over separate timescales: $\tau_c$ A chemical reaction timescale, $\tau_c$ , and a turbulent mixing timescale, $\tau_{sgs}$. The filtered reaction rate of i-species scales with

$$\dot{\omega}_i^{\text{PaSR}} = \frac{\tau_{c}}{\tau_{sgs} + \tau_{\text{c}}} \cdot \dot{\omega}_i (\overline{\phi}) = \gamma^\ast \dot{\omega}_i$$

Or introducing the sub-grid Damkholer number $\text{Da}_{sgs} \equiv \tau_{sgs}/\tau_c$ then

$$\gamma^\ast = \frac{1}{1 + \text{Da}_{sgs}}$$

If sub-grid turbulent mixing is fast relative to chemistry, reactions take place in a well-mixed environment. $\gamma^\ast \rightarrow 1$ If chemistry is fast, chemical equilibrium is sough $\gamma^\ast \rightarrow 0$. This is similar to the Eddy Dissipation Concept (EDC) where turbulence causes reactions to only take place in localized pockets (fine structure) and $\gamma^\ast$ represents the volume fraction of these pockets

To write

Eulerian Stochastic Fields

The Stochastic fields equations for the joint-velocity-scalar energy PDF equations to be solved are

Continuity

$$\frac{d \varrho^n}{d t} + \frac{\partial \varrho^n \mathscr{U}^n_j }{\partial x_j} = 0$$

Momentum

$$\frac{d \varrho^n \mathscr{U}^n_i}{d t} + \frac{\partial \varrho^n \mathscr{U}^n_j \mathscr{U}^n_i}{\partial x_j} = -\frac{\partial \mathscr{P}^n}{\partial x_i} + \frac{\varrho^n}{\overline{\rho}} \frac{\partial \widetilde{\tau}_{ij}}{\partial x_i} + \varrho^n G_{ij} \left( \mathscr{U}^n_j - \widetilde{u}_j \right) + \varrho^n \sqrt{C_0 \frac{\epsilon_{sgs}}{\overline{\rho}}} \frac{d W^n_i}{d t}$$

Mass Fraction of specie k

$$\frac{\partial \varrho^n \mathscr{Y}^n_k}{\partial t} +\frac{\partial \varrho^n \mathscr{U}^n_i \mathscr{Y}^n_k}{\partial x_i} = \frac{\varrho^n}{\overline{\rho}}\frac{\partial \widetilde{J}_{k,i} }{\partial x_i} + \varrho^n \dot{\omega}_k - \frac{1}{2} C_{Y} \frac{\epsilon_{sgs}}{k_{sgs}} \varrho^n \left( \mathscr{Y}^n_k - \widetilde{Y}_k \right)$$

Energy

$$\frac{\partial \varrho^n \mathscr{E}^n_t}{\partial t} + \frac{\partial \varrho^n \mathscr{U}^n_i \mathscr{E}^n_t}{\partial x_i} = \frac{\varrho^n}{\overline{\rho}}\frac{\partial \widetilde{q}_i}{\partial x_i} - \frac{\varrho^n}{\overline{\rho}}\frac{\partial \overline{p} \widetilde{u}_i}{\partial x_i} + \frac{\varrho^n}{\overline{\rho}}\frac{\partial \widetilde{\tau}_{ij}\widetilde{u}_j}{\partial x_i} - \frac{1}{2} C_{E} \frac{\epsilon_{sgs}}{k_{sgs}} \varrho^n \left( \mathscr{E}^n_t - \widetilde{e_t} \right)$$

The employed closure relation for the dissipation of the sub-grid kinetic energy is

$$\epsilon_{sgs} = C_\epsilon k_{sgs}^{3/2}/\Delta$$

where the constant $C_\epsilon = 1.05$. The micro-mixing constants $C_Y = 2$ and the Langevin constant is set to 2.1,

$$k_{sgs}= \frac{1}{2} \left( \frac{1}{N_f} \sum_{n=1}^{N_f} (\mathscr{U}^n_i-\widetilde{u}_i)^2 \right)$$

The filtered variables can be obtained from the average of the Eulerian stochastic fields.

$$\overline{\phi} = \frac{1}{N_f}\sum_{n=1}^{N_f} \phi^n; \ \ \ \ \widetilde{\phi} = \frac{\sum_{n=1}^{N_f} \varrho^n \phi^n}{\sum_{n=1}^{N_f} \varrho^n}$$

References

[1] Omer Rathore, "Numerical simulation of combustion instability: flame thickening and boundary conditions" PhD Thesis, Imperial College London (2022)

[2] Yuri Almeida "Large Eddy Simulation of Supersonic Combustion using a Probability Density Function method" PhD Thesis, Imperial College London (2019)

[3] Tin-Hang Un and Salvador Navarro-Martinez, “Stochastic fields with adaptive mesh refinement for high-speed turbulent combustion”, Comb. Flame, 272, 113897 (2025)