7.3 Large-Eddy Simulation 大涡模拟[cfd-7-3]

LES is based on a spatial filtering operation, which decomposes any flow variable \(U\) into a filtered (large-scale, resolved) part \(\overline{U}\) and into a sub-filter (unresolved) part \(U'\), i.e.,

The filtered variable at the location \(\vec{r}_0\) in space is defined as

where \(\Omega\) denotes the entire flow domain, \(G\) represents the filter function, and \(\vec{r}\) is the position vector, respectively. The filter function determines the structure and size of the small scales. The filter function depends on the difference \(\vec{r}_0 - \vec{r}\) and on the filter width \(\Delta = \left(\Delta_1\,\Delta_2\,\Delta_3\right)^{1/3}\), with \(\Delta_i\) being the filter width in the \(i\)-th spatial coordinate. The following filter functions are the mostly used ones (see Fig. 7.3):

• the tophat filter:
• The sharp Fourier cut-off filter:
• The Gaussian filter:

The tophat and the Gaussian filter smooth the large-scale fluctuations as well as the small scales below the filter width. The cut-off filter affects only the scales below the cut-off wave-number. In practice, the Gaussian filter is always employed in conjunction with a sharp Fourier cut-off. Filters suitable for grids with varying cell sizes were proposed in Refs. [114], [115].

The spatial filtering, defined by Eq. (7.76) and Eq. (7.77), has to be applied to the Navier-Stokes equations in order to remove the small turbulent scales. The filter width \(\Delta\) as well as the filter function are considered as free parameters. In fact, the governing equations are usually not explicitly filtered. Instead, the grid as well as the discretisation errors are assumed to define the filter \(G\). For the discussion of explicit filtering see Refs. [116], [117].

Because of the differing treatment, we shall distinguish in the following between compressible (7.1) and incompressible (7.6) formulation of the Navier-Stokes equations.

For an incompressible flow of a Newtonian fluid, the filtered governing equations (7.6) take the form

where \(\nu\) denotes the kinematic viscosity coefficient. The equations (7.81) describe the temporal and spatial evolution of the large, energy-carrying scales of motion. The non-linearity of the convective term leads to the appearance of the so-called subgrid-scale stress (SGS) tensor

which describes the effects of the unresolved scales. The SGS tensor has to be modelled (see Subsection 7.3.3) in order to close the equations.

The SGS tensor can be decomposed into three parts [118], namely

The individual parts have the following physical meaning:

is the so-called Leonard stress term and represents the interactions between large-scale eddies which produce small-scale turbulence. This term only can be evaluated explicitly from the filtered velocity field \(v_i\). Further, the cross-stress term

describes interactions between large- and small-scale eddies. Finally,

is the so-called SGS Reynolds-stress tensor. It reflects interactions between the small-scale structures. The above decomposition (7.83) is no longer used mainly because \(L_{ij}\) and \(C_{ij}\) are not invariant with respect to Galilean transformation¹.

If LES is to be applied to compressible flows, we have to apply Favre averaging (Subsection 7.1.2) together with the spatial filtering to the Equations (7.1). Otherwise, the filtered Navier-Stokes equations would contain products between density and other variables like velocity or temperature. Thus, the velocity components, the energy and the temperature in Eq. (7.1) is decomposed as

The filtered variable at the location \(\vec{r}_0\) in space is given by

where the overbar denotes the filtering in Eq. (7.77). The Favre-filtered Navier-Stokes equations (7.1) read [111], [113]

with the terms

and

In the above equations (7.89)-(7.90), \(e\) denotes internal energy per unit mass, \(\tilde{S}_{ij}\) is the Favre-filtered strain-rate tensor, and \(\tau_{ij}^{SF} = \overline{\rho}\left(\widetilde{v_i v_j} - \tilde{v}_i\tilde{v}_j\right)\) represents the Favre-averaged subgrid-scale stress. Furthermore, \(\mu\), \(\mu_B\), and \(k\) stand for the molecular viscosity, the bulk viscosity, and for the thermal conductivity, respectively. Finally, \(\tilde{\mu}\), \(\tilde{\mu}_B\), and \(\tilde{k}\) are the corresponding values at the filtered temperature \(\tilde{T}\).

The right-hand side of Eq. (7.89) contains terms which have to be modelled. In the momentum equation, the SGS stresses \(\tau_{ij}^{SF}\) are approximated, but the second term, i.e., \(\left(\overline{\sigma}_{ij} - \hat{\sigma}_{ij}\right)\) is usually neglected. In the energy equation, term \(\mathcal{A}\) can be expressed through the SGS stresses [119], term \(\mathcal{B}\) can be neglected, and terms \(\mathcal{C}\), \(\mathcal{D}\) can be modelled as proposed in [120].

The main task of a subgrid-scale model is to simulate energy transfer between the large and the subgrid scales. On the average, the energy is transported from the large scales to the small ones (turbulent cascade process). Therefore, a subgrid-scale model has to provide means of adequate energy dissipation. However, in some instances the energy also flows from small to large scales - a process called backscatter. Thus, the model should account for this effect as well. Backscatter models are discussed, e.g., in [121].

Various subgrid-scale models were proposed in the past and the research still continues. The models can be divided into two basic classes. The first one consists of approaches which model the SGS tensor \(\tau_{ij}^S\) explicitly. A necessary condition is then that the numerical dissipation caused by the spatial discretisation scheme must be much lower than the subgrid-scale dissipation. The majority of explicit SGS models is based on the eddy-viscosity concept, which is explained next. Furthermore, we shall present the Smagorinsky model, which forms the basis of all subgrid-scale models. We shall also briefly discuss the basics of the so-called dynamic subgrid-scale models. A comparison of six different explicit subgrid-scale models was presented recently in [122].

The second class of subgrid-scale models consists of approaches where the SGS stresses are modeled implicitly by an appropriate discretisation of the convective fluxes (thus \(\tau_{ij}^S\) is omitted). These models are referred to as Monotonically Integrated LES (MILES). The methodology was first presented by Boris et al. [123] and recently advocated by Grinstein and Fureby [124], [125]. The condition in this case is that the numerical dissipation correctly models the subgrid-scale dissipation. This is not quite easy to achieve, as the investigations in Ref. [126] and [127] prove. Nevertheless, MILES was applied with some success to a variety of flow problems [128]-[134].

These explicit models are able to represent the global dissipative effects of the small scales, but they cannot reproduce the local details of the energy exchange. In the case of incompressible flows, the eddy-viscosity models relate the SGS stresses to the large-scale strain-rate \(\overline{S}_{ij}\) as follows

The strain-rate \(\overline{S}_{ij}\) is obtained from Eq. (7.3) by using filtered velocity components. The eddy viscosity \(\nu_T\) is in general evaluated from algebraic relations in order to save numerical costs. The isotropic part of the SGS stresses (\(\tau_{kk}^S\)) can either be added to the filtered pressure [135], modelled [136] or neglected.

Relation similar to Eq. (7.91) applies also in the case of compressible Navier-Stokes equations. The components of the Favre-averaged SGS stress tensor are approximated as

The components of the strain-rate tensor \(\tilde{S}_{ij}\) are given in Eq. (7.90).

The Smagorinsky model [75] is based on the equilibrium hypothesis which implies that the small scales dissipate entirely and instantaneously all the energy they receive from the large scales. The algebraic model assumes the form

where \(|\overline{S}| = \left(2\overline{S}_{ij}\overline{S}_{ij}\right)^{1/2}\) is the magnitude of the strain-rate tensor and \(C_s\) denotes the Smagorinsky constant. The theoretical value found by Lilly [137] is \(C_s \approx 0.18\). However, the Smagorinsky constant depends on the type of the flow. For example, in shear flows \(C_s\) has to be reduced to approximately 0.1. The filter width \(\Delta\) in Eq. (7.93) is usually chosen to be twice the average grid size, i.e., \(\Delta = 2\left(\Delta x_1\,\Delta x_2\,\Delta x_3\right)^{1/3}\).

In order to account for the reduced growth of the small scales near walls, the value of the eddy viscosity \(\nu_T\) has to be reduced. Thus, the Smagorinsky model Eq. (7.93) is modified according to Van Driest damping as

where \(y^{+}\) represents the dimensionless wall distance (for the computation of wall distances see, e.g., Ref. [39]).

The Smagorinsky model is numerically cheap and easy to implement. However, it has several serious disadvantages:

• it is too dissipative in laminar regions with mean shear;
• it requires special provisions near walls and at laminar-turbulent transition;
• the parameter \(C_s\) is not uniquely defined;
• the process of energy backscatter is not modelled.

Because of these shortcomings, various other approaches were proposed (see, e.g., [111]). Very popular are the following dynamic models.

The dynamic SGS models employ the same relation as the Smagorinsky model (Eq. (7.93)) for the evaluation of the eddy viscosity \(\nu_T\) in Eq. (7.91) or (7.92). The difference is that the Smagorinsky constant (adjusted a priori) is replaced by a parameter, which evolves dynamically in space and in time. Hence,

The parameter \(C_d\) is computed based on the energy content of the smallest scale of the turbulence. For this purpose, Germano et al. [138] proposed to employ a second filter - the so-called test filter \(\hat{\Delta}\). The width of the test filter has to be larger than that of the filter \(\Delta\) applied to the governing equations (usually \(\hat{\Delta} = 2\Delta\)). The application of the test filter to the filtered equations leads to the so-called subtest-scale stresses \(\tau_{ij}^{ST}\)

The subtest-scale stresses are related to the SGS stresses \(\tau_{ij}^S\) (Eq. (7.83)) via the Germano identity [139]

where \(\hat{L}_{ij}\) denotes the Leonard stresses associated with the test filter. It represents the contribution to the Reynolds stresses by the scales whose length is intermediate between the filter width \(\Delta\) and the test filter width \(\hat{\Delta}\).

If we express the subtest-scale and SGS stresses in Eq. (7.97) using the eddy-viscosity approach Eq. (7.91) together with Eq. (7.95), we obtain

with

The notation \([\,]^{\wedge}\) means that the whole term enclosed in the square brackets is test-filtered. The parameter \(C_d\) can be derived from Eq. (7.98) by using the least-squares minimisation of Lilly [140]. This leads to

The above formulation (7.100) is mathematically inconsistent since the parameter \(C_d\) was taken outside the test filter in Eq. (7.98). In practice, the numerator and denominator in Eq. (7.100) are therefore ensemble-averaged in the homogeneous directions, i.e.,

Improved dynamic SGS models were proposed, e.g., by Ghosal et al. [141], Carati et al. [142], Piomelli and Liu [143], and Held [90].

The computational costs of LES of wall-bounded flows at high Reynolds numbers (\(Re > 10^6\)) are still too high for engineering purposes. The reason is the excessively large number of grid points (cells) required to resolve the wall layer appropriately. In order to reduce the costs, it is possible to model the wall layer by specifying a correlation between the velocity in the outer flow and the stress at the wall. This approach is quite similar to using wall functions in RANS simulations. The basic assumption is that there is only a weak interaction between the near-wall and the outer region, which is supported by the investigations in [144] and [145].

Earlier implementations of the wall models were based on the assumption that the dynamics of the wall layer are universal and hence they can be approximated by a generalised law-of-the-wall. Basically, the models utilised the logarithmic law (see [77] and [146]-[148]). Balaras et al. [149] proposed recently a new zonal approach. Within the two-layer model, the filtered Navier-Stokes equations (7.81) are solved up to the first grid point above the wall. From this point to the wall 2-D boundary layer equation are solved on a refined embedded grid. The solution on the embedded grid is then used to prescribe the wall shear stress as a boundary condition for the LES. The zonal approach of Balaras et al. [149] allows it to place the first point in a region \(20 < y^{+} < 100\), which leads to significantly reduced grid size and hence computational time. The methodology was applied with success to turbulent flows in a plane channel, square duct and rotating channel. Later on, it was also employed for the LES of separated flows with encouraging results [150]-[153].

Even though the above wall models help to reduce the number of grid points (cells) considerably, LES still remains too costly for complex engineering configurations. For this reason, Spalart recently suggested another approach, the so-called Detached Eddy Simulation (DES), which is aimed at the simulation of high Reynolds-number massively separated flows [154], [155]. The methodology represents a hybrid between the RANS and LES. The idea is to employ highly stretched grids together with a RANS turbulence model (mostly the Spalart-Allmaras model from Subsection 7.2.1 or Menter's SST model from Subsection 7.2.3) to resolve the attached boundary layer(s), and to use LES outside the wall region together with an isotropic grid to capture the detached 3-D eddies. Thus, DES tries to combine the strengths of both methods in a single framework.

The implementation of DES for the Spalart-Allmaras turbulence model consists of replacing the wall distance \(d\) in the formulae (7.36) and (7.38)-(7.40) by the DES length scale

The length scale is dependent on the largest dimension of the control volume, i.e., \(\Delta = \max(\Delta x, \Delta y, \Delta z)\). The constant \(C_{DES}\) depends to some extent on the type of the flow. For a homogeneous turbulence, the value \(C_{DES} = 0.65\) was found optimal [156]. On the other hand, \(C_{DES} = 0.1\) was recommended for transonic and supersonic jets [157]. The definition of the length scale \(l\) in Eq. (7.102) makes sure that within the boundary layer, where \(d < C_{DES}\Delta\) and hence \(l = d\), the original RANS model is recovered. On the other hand, outside the boundary layer \(l = C_{DES}\Delta\) and the Spalart-Allmaras model serves as a one-equation SGS model for the LES (cf. Eq. (7.91) and (7.24)). When integrating the governing equations in time, the global time step has to be adjusted such as to resolve the time scales of the detached eddies. This usually means that the time step would by far exceed the stability margin of an explicit scheme for the boundary layer region. It is therefore more efficient to employ a time-accurate implicit scheme, such as the dual time-stepping approach described in Section 6.3. More details regarding the DES methodology and examples of simulations can be found in the above references or in [158]-[161].