**2. Large eddy simulation (LES)**

Spatial scale separation is done by applying a low‐pass filter to the governing balance equations. This set of filtered equations governs the dynamics of the large scales. Spatial fluctuations, lower than a defined filter cut‐off length ‐Δ‐, are smoothed or removed. The spatial filtering operation of a flow variable *Q,* being a function of time and space, is defined by the convolution integral:

Transport and Mixing in Liquid Phase Using Large Eddy Simulation: A Review http://dx.doi.org/10.5772/63993 403

$$\{\mathcal{Q}\}\_{\iota} \equiv \frac{1}{\{\rho\}\_{\iota}} \int\_{\mathcal{D}} \rho G(\mathbf{x} - \mathbf{x}') \mathcal{Q}(\mathbf{x}) d\mathbf{x}' = \frac{\langle \rho \mathcal{Q} \rangle\_{\iota}}{\langle \rho \rangle\_{\iota}} \tag{1}$$

where *Q* = *Q <sup>L</sup>* + *Q*". Favre filtering has been applied to Eqs. (1) and (2). *G* is the filter kernel defined over the entire domain. The filter kernel depends on the filter width, Δ, and must satisfy a set of conditions (see Ref. [16]). Favre‐filtered mass, momentum, and scalars' balance equations are as follows:

constitutive LES equations are grid‐dependent, that is, LES is an incomplete model. Referenc‐

Compared to DNS, LES saves the computational cost of solving the lowerturbulent scales (99% of DNS calculations), and can be used for the simulation of systems of practical importance.

In shearflows, LES prediction capabilities are diminished by the fact that in the viscous region, the energy‐containing structures are of the same order of magnitude as the viscous length‐ scales (*δ<sup>v</sup>* =*v ρ* / *τw*), and LES cannot solve them. To face this issue, it is possible to refine the computational mesh in near‐wall zones in order to consider the small coherent structures developed on those zones, but this refinement implies more computational effort. There are three additional possibilities: incorporating the fluid behavior in the wall region by wall models, simulating those zones with RANS (hybrid RANS‐LES), or incorporating wall‐

Most research efforts in LES for reactive flows are concentrated on the SGS fluxes and Favre‐ filtered source terms [30]. In diffusion flames, chemical reactions take place after mixing of reactive species is achieved at the smallest scales of turbulence (unresolved in LES); so, in LES

According to Ref. [31], the uncertainty initially contained in the nonresolved scales is propa‐ gated to the resolved ones; hence, it is said that SGS modeling is not a well‐posed problem. LES results can be interpreted as a different realization of the flow. They could have the same statistical properties of the flow and may predict the same spatially organized structures but

LES has been successfully applied to turbulent flows, in particular to complex geometries in liquid phase such as stirred tank reactors (see Refs. [10, 24, 32–34]). Because LES equations are unclosed, SGS models must be supplied in order to specify the SGS stresses, SGS scalar fluxes,

There are still unanswered questions about the behavior of high Schmidt turbulent flows. LES could be an alternative method for studying these flows, as long as the flow's physical behavior can be captured by the subgrid‐scale models [15]. The LES equations are presented in the

Spatial scale separation is done by applying a low‐pass filter to the governing balance equations. This set of filtered equations governs the dynamics of the large scales. Spatial fluctuations, lower than a defined filter cut‐off length ‐Δ‐, are smoothed or removed. The spatial filtering operation of a flow variable *Q,* being a function of time and space, is defined

However, it is more expensive than a *k‐ε* model for a given grid resolution.

es [27–29] review the LES literature.

402 Numerical Simulation - From Brain Imaging to Turbulent Flows

adapted SGS closures.

at a different location.

and filtered reaction terms.

by the convolution integral:

**2. Large eddy simulation (LES)**

following section.

the combustion process must be modeled.

$$\frac{\partial \left\langle \rho \right\rangle\_{\iota}}{\partial t} + \frac{\partial}{\partial x\_{\wedge}} \left( \left\langle \rho \right\rangle\_{\iota} \left\langle u\_{\wedge} \right\rangle\_{\iota} \right) = 0 \tag{2}$$

$$\frac{\left\|\hat{c}\left\{\rho\right\}\_{i}\left\{u\_{i}\right\}\_{l\_{L}}+\frac{\left\|\left<\rho\right\rangle\_{i}\left\{u\_{i}\right\}\_{l\_{L}}\left\{u\_{j}\right\}\_{l\_{L}}\right\}}{\left\|\mathbf{x}\_{j}\right\|} = -\frac{\left\|\left<\rho\right\rangle\_{i}\left\{\left\{\mathbf{r}\_{\overline{\boldsymbol{\eta}}}\right\}\_{l} + \mathbf{r}\_{\overline{\boldsymbol{\eta}}}^{\mathrm{asy}}\right\}\right\| - \frac{\left\|\left<\rho\right\rangle\_{l}\right\|}{\left\|\mathbf{x}\_{i}\right\|} + \left\{\rho\right\}\_{l}F\_{i} \tag{3}$$

$$\frac{\left\|\left<\rho\right>\_{\cdot}\left<\phi\_{a}\right>\_{L}}{\left<\gamma\right>\_{\cdot}}+\frac{\left<\left<\rho\right>\_{\cdot}\left<\left<\phi\_{a}\right>\_{L}\left<\mu\_{\cdot}\right>\_{\cdot}+\left\_{L}+J\_{a,\cdot}^{\ast\ast}\right>\right>}{\left<\gamma\right>\_{\cdot}}\right>\_{L}\tag{4}$$

where the symbol ∂ denotes the partial differential operator, and the summation convention is used for repeated indices. Time (*t*) and spatial coordinates (*xj* ) are the independent varia‐ bles; *ρ* is the fluid mass density; *ϕα* is the mass‐weighted value of the scalar field (e.g., internal energy, mass fraction); *uj* is the mass‐averaged velocity in the *j*‐direction; *p* denotes pressure; *τij* is the shear stress tensor; and *J* is the diffusive flux of the scalar.

New terms appear in LES equations, as it occurs in RANS. These terms account for resolved and subgrid scales' interactions. The SGS stresses (*τij SGS* ) and SGS scalar fluxes ( *j α*, *j SGS* ) are defined by

$$
\pi\_{\circ}^{\text{sys}} = \left\langle \boldsymbol{\mu}\_{i} \boldsymbol{\mu}\_{j} \right\rangle\_{L} - \left\langle \boldsymbol{\mu}\_{i} \right\rangle\_{L} \left\langle \boldsymbol{\mu}\_{j} \right\rangle\_{L}; \quad j\_{\alpha,j}^{\text{sys}} = \left\langle \boldsymbol{\mu}\_{j} \boldsymbol{\phi}\_{\alpha} \right\rangle\_{L} - \left\langle \boldsymbol{\mu}\_{j} \right\rangle\_{L} \left\langle \boldsymbol{\phi}\_{\alpha} \right\rangle\_{L} \tag{5}
$$

In addition, a model is needed for filtered source terms. Subgrid‐scale (SGS) models must be supplied in order to specify the SGS stresses, SGS scalar fluxes, and filtered reaction terms.

The main task of SGS stress tensor model is to dissipate the energy transferred from turbu‐ lent large scales. The dynamics of the subgrid scales affect that of the resolved flow field through the subgrid‐scale stress tensor.

There are several SGS stress tensor models (for review, see Refs. [17, 35]). The simplest models are based on the gradient assumption. More complex models are of first and second order. As the model complexity increases, so does the required number of equations. Those models based on eddy viscosity show good balance between accuracy and numerical computation and present good prediction abilities in turbulent combustion and other highly interacting processes. It is also possible to close the LES approximation by PDF methods, such as the velocity‐filtered density function (VFDF) [36] and the velocity‐scalar filtered mass density function (VSFMDF) [37], but this approach increases considerably the number of equations to be solved. This procedure increases the computational cost 15–30 times longer than the Smagorinsky [38] and dynamic Smagorinsky [39, 40] models, respectively.

#### **2.1. The Smagorinsky model**

The Smagorinsky model [38] is one of the pioneer SGS models in the development of LES. Its simple formulation and good performance have made it very popular, and it is used in this work. The Smagorinsky model is an eddy‐viscosity based model, which assumes equilibri‐ um between the turbulent kinetic energy dissipation and production rates to obtain a relation between the characteristic velocity and the resolved strain rate. A Boussinesq approximation is applied to the deviatoric part of the SGS stress.

$$
\tau\_{\vec{\nu}}{}^{\rm sys} - \frac{1}{3} \tau\_{kk}{}^{\rm sys} \delta\_{\vec{\nu}} \approx -\nu\_T \left[ 2 \left\{ S\_{\vec{\nu}} \right\}\_L - \frac{2}{3} \left\{ S\_{kk} \right\} \delta\_{\vec{\nu}} \right] \tag{6}
$$

where *νT* is the subgrid eddy viscosity of residual motions and *δ* is the Dirac delta function. The filtered strain rate is given by

$$\left\{ S\_y \right\}\_L = \frac{1}{2} \left( \frac{\partial \left\langle u\_i \right\rangle\_L}{\partial x\_j} + \frac{\partial \left\langle u\_j \right\rangle\_L}{\partial x\_i} \right) \tag{7}$$

The subgrid eddy viscosity is modeled in a similar way as the mixing length:

$$\boldsymbol{\nu}\_{T} = \left(\boldsymbol{C}\_{s}\boldsymbol{\Delta}\right)^{2} \left|\boldsymbol{S}\right|\_{L} \quad , \quad \left|\boldsymbol{S}\right|\_{L} = \left(\boldsymbol{2}\left\{\boldsymbol{S}\_{\boldsymbol{\vartheta}}\right\}\_{L} \left\{\boldsymbol{S}\_{\boldsymbol{\vartheta}}\right\}\_{L}\right)^{\underline{\mathsf{X}}} \tag{8}$$

where *Cs* is the Smagorinsky coefficient. This Model parameter is not universal and depends on the flow configurations. While for isotropic turbulence, the Smagorinsky coefficient is about 0.17 [41], for other configurations this value may not be correct; for a channel flow, the Smagorinsky coefficient is about 0.1 [42]. In addition, the spatial variation of this coefficient makes it difficult to find a proper value. This model parameter can be dynamically calculat‐ ed by using the dynamic procedure [39, 40], in which, by assuming scale invariance, a test filter greater than the filter width, Δ, can be used for calculating *Cs* from the resolved flow field.

#### **2.2. The dynamic Smagorinsky model**

The major drawback of the Smagorinsky model is that a single universal constant cannot correctly represent different turbulent flows. Germano et al. [39] proposed a dynamic procedure for computing the Smagorinsky model coefficient, based on the instantaneous

information given by the filtered velocity field. The coefficient is locally recalculated during the simulations; so, it is no longer necessary to specify its value as an input parameter.

The dynamic procedure is based on the idea that the information given by the smallest resolved scales can be used to model the largest unresolved scales, as they have a similar behavior. The latter can be done by employing a test filter, with the test filter width larger than the filter width (Δ), usually taken as *Δ* ^ =2*Δ*.

A model for the SGS stress is

(6)

(7)

(8)

processes. It is also possible to close the LES approximation by PDF methods, such as the velocity‐filtered density function (VFDF) [36] and the velocity‐scalar filtered mass density function (VSFMDF) [37], but this approach increases considerably the number of equations to be solved. This procedure increases the computational cost 15–30 times longer than the

The Smagorinsky model [38] is one of the pioneer SGS models in the development of LES. Its simple formulation and good performance have made it very popular, and it is used in this work. The Smagorinsky model is an eddy‐viscosity based model, which assumes equilibri‐ um between the turbulent kinetic energy dissipation and production rates to obtain a relation between the characteristic velocity and the resolved strain rate. A Boussinesq approximation

where *νT* is the subgrid eddy viscosity of residual motions and *δ* is the Dirac delta function.

where *Cs* is the Smagorinsky coefficient. This Model parameter is not universal and depends on the flow configurations. While for isotropic turbulence, the Smagorinsky coefficient is about 0.17 [41], for other configurations this value may not be correct; for a channel flow, the Smagorinsky coefficient is about 0.1 [42]. In addition, the spatial variation of this coefficient makes it difficult to find a proper value. This model parameter can be dynamically calculat‐ ed by using the dynamic procedure [39, 40], in which, by assuming scale invariance, a test filter greater than the filter width, Δ, can be used for calculating *Cs* from the resolved flow field.

The major drawback of the Smagorinsky model is that a single universal constant cannot correctly represent different turbulent flows. Germano et al. [39] proposed a dynamic procedure for computing the Smagorinsky model coefficient, based on the instantaneous

The subgrid eddy viscosity is modeled in a similar way as the mixing length:

Smagorinsky [38] and dynamic Smagorinsky [39, 40] models, respectively.

**2.1. The Smagorinsky model**

The filtered strain rate is given by

**2.2. The dynamic Smagorinsky model**

is applied to the deviatoric part of the SGS stress.

404 Numerical Simulation - From Brain Imaging to Turbulent Flows

$$
\tau\_y^{\ast \text{gs}} - \frac{1}{3} \tau\_{kk}^{\ast \text{gs}} \mathcal{S}\_y = -2C\_{\text{dyn}} \Delta^2 \left| \left\{ S\_y \right\}\_L \right| \left( \left\{ S\_y \right\}\_L - \frac{1}{3} \left\{ S\_{kk} \right\}\_L \delta\_y \right) = -2C\_{\text{dyn}} \alpha\_y \tag{9}
$$

where *Cdyn* is a Model parameter, and the term should be modeled. If the test filter is applied to the filtered momentum equation, the filter stress has the form:

$$T\_{\psi} = \left\langle u\mu\_{\rangle} \right\rangle\_{L} - \left\langle u\_{i} \right\rangle\_{L} \left\langle u\_{j} \right\rangle\_{L} \tag{10}$$

As it was the case for the SGS stress, a model for the filter stress is

$$T\_y - \frac{1}{3} T\_{ik} \delta\_y = -2 C\_{dyn} \Delta^2 \left| \left\langle S\_y \right\rangle\_L \right| \left( \left\langle S\_y \right\rangle\_L - \frac{1}{3} \left\langle S\_{kk} \right\rangle\_L \delta\_y \right) = -2 C\_{dyn} \beta\_y \tag{11}$$

The coefficient has been assumed to be independent of the filtering process. By applying the test filter on the unresolved SGS scalar flux (*τij SGS* ) and subtracting it from the test filter stress (*Tij*), the resolved part of the SGS stress tensor is obtained. This is the Germano identity [40]:

$$L\_y = T\_y \ -\hat{\tau}\_{\vec{v}} = \left\langle u\_i \right\rangle\_L \left\langle u\_j \right\rangle\_L - \left\langle u\_i \right\rangle\_L \left\langle u\_j \right\rangle\_L \tag{12}$$

The Germano identity can be partially satisfied by replacing model equations for the SGS and filter stresses [Eqs. (9) and (11), respectively] into Eq. (12), giving an approximation for the model coefficient (*Cdyn*):

$$L\_{\circ} - \frac{1}{3}L\_{kk}\mathcal{S}\_{\circ} = -2C\_{\phi\circ} \left(\mathcal{B}\_{\circ} - \alpha\_{\circ}\right) = -2C\_{\phi\circ}M\_{\circ} \tag{13}$$

Equation (13) is overspecified, because there are five independent equations and one un‐ known value. There are different ways to calculate the coefficient. By applying a least‐squares methodology, Lilly [40] minimized the error incurred in the calculation of the coefficient. The error is

$$E\_y = L\_y - T\_y + \tau\_y = L\_y - \frac{1}{3}L\_{kk}\delta\_y + 2C\_{d\phi n}\mathcal{M}\_y \tag{14}$$

The minimum corresponds to

$$\frac{\partial \left(E\_{\psi} E\_{\psi}\right)}{\partial C\_{\psi\eta}} = 0\tag{15}$$

leading to

$$C\_{\phi\nu} = -\frac{1}{2} \frac{L\_{\tilde{\nu}} E\_{\tilde{\nu}}}{M\_{\tilde{\nu}} M\_{\tilde{\nu}}} = 0 \tag{16}$$

The calculation of the parameter generates spatial and temporal fluctuations, as well as negative values [43]. Regions with a negative coefficient may be interpreted as regions where backscatter takes place. According to Carati et al. [44], the dynamic Smagorinsky model does not have information about the amount of energy that is available in the subgrid scales.
