This preview shows page 1. Sign up to view the full content.
Unformatted text preview: Introduction to Large Eddy Simulation of Turbulent Flows 1
J. Frohlich, W. Rodi Institute for Hydromechanics, University of Karlsruhe, Kaiserstra e 12, D76128 Karlsruhe, Germany Contents
1 Introduction
1.1 Resolution requirements of DNS . . . . . . . . . . . . . . . . . . . . . . . . . . . . 1.2 The basic idea of LES . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2.1 2.2 2.3 2.4 3.1 3.2 3.3 3.4 3.5 Schumann's approach Filtering . . . . . . . . Variable lter size . . Implicit versus explicit ..... ..... ..... ltering . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 2 2 2 3 4 5 6 2 Governing equations and Filtering 3 3 Subgrid{scale modelling Introduction . . . . . . . . . . . . . . . . . . Smagorinsky model . . . . . . . . . . . . . . Dynamic procedure . . . . . . . . . . . . . . Scale similarity models . . . . . . . . . . . . Further models and comparative discussion .6 .7 .8 . 10 . 10 6 4 Numerical methods 4.1 Discretization schemes in space and time . . . . . . . . . . . . . . . . . . . . . . . . 12 4.2 Analysis of numerical schemes for LES . . . . . . . . . . . . . . . . . . . . . . . . . 12 4.3 Further developments . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 14 5.1 5.2 5.3 5.4 5.5 Resolution of the near{wall region Wall functions . . . . . . . . . . . Other approaches . . . . . . . . . . In ow and out ow conditions . . . Sample computations . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . 12 5 Boundary conditions 14 15 15 16 17 17 6 Concluding remarks
1 19 to appear in B.E. Launder, N.D. Sandham (eds.): Closure Stratigies for Turbulent and Transitional Flows, Cambridge University Press 1 Introduction
This chapter is meant as an introduction to Large{Eddy Simulation (LES) for readers not familiar with it. It therefore presents some classical material in a concise way and supplements it with pointers to recent trends and literature. For the same reason we shall focus on issues of methodology rather than applications. The latter are covered elsewhere in this volume. Furthermore, LES is closely related to direct numerical simulation (DNS) which is also discussed in several chapters of this volume. Hence, we concentrate as much as possible on those features which are particular to LES and which distinguish it from other computational methods. For the present text we have assembled material from research papers, earlier introductions and reviews (Ferziger(1996), Hartel (1996), Piomelli (1998) ), and our own results. The selection and presentation is of course biased by the authors' own point of view. Supplementary material is available in the cited references. 1.1 Resolution requirements of DNS
The principal di culty of computing and modelling turbulent ows resides in the dominance of non{linear e ects and the continuous and wide spectrum of observed scales. Without going into details (the reader might consult classical text books such as Tennekes & Lumley (1972)) we just recall here that the ratio of the size of the largest turbulent eddies in a=4 ow, L, to that of the smallest ones determined by viscosity, , behaves like L= Re3 . Here, Reu = u L= with u being a characteristic velocity uctuation and u the kinematic viscosity. Let us consider as an example a plane channel, a prototype of an internal ow. Reynolds (1989) estimated Reu Re0:9 from u u c1=2, cf Re 0:2, where f Re is based on the center line velocity and the channel height. In a DNS no turbulence model is applied so that motions of all size have to be resolved numerically by a grid which is su ciently ne. Hence, the computational requirements increase rapidly with Re. According to this estimate a DNS of channel ow at Re = 106 for example would take around hundred years on a computer running at several GFLOPS. This is obviously not feasible. Moreover, in an expensive DNS a huge amount of information would be generated which is mostly not required by the practical user. He or she would mostly be content with knowing the average ow and some lower moments to a precision of a few percent. Hence, for many applications a DNS which is of great value for theoretical investigations and model testing is not only una ordable but would also result in computational overkill.
0 0 0 0 0 0 ; 1.2 The basic idea of LES
Suppose somebody wants to perform a DNS but the grid that would be required exceeds the capacity of the available computer so a coarser grid is used. This coarser grid is able to resolve the larger eddies in the ow but not the ones which are smaller than one or two cells. From a physical point of view, however, there is an interaction between the motions on all scales so that the result for the large scales would generally be wrong without taking into account the in uence of the ne scales on the large ones. This requires a so{called subgrid{scale model as discussed below. Hence, LES can be viewed as a `poor man's DNS'. The poor man, however, has to compensate by cleverness in that a model for the unresolved motion has to be devised and an intricate coupling between physical and numerical modelling is generated. On the other hand, the resolution of the large scales of the ow while modelling only the small ones { not the entire spectrum { is an advantage of the LES approach compared to methods based on the Reynolds averaged Navier{Stokes equations (RANS). The latter methods often have di culties when applied to complex ows with pronounced vortex shedding or special in uences of buoyancy, curvature, rotation or compression. Finally, LES gives access to the dominant unsteady motion so that it can, for example, be used to study aero{acoustics, uid{structure coupling or the control of turbulence by an appropriate unsteady forcing. 2 Governing equations and Filtering
The Navier{Stokes equations (NSE) constitute the starting point for any turbulence simulation. Here, we consider incompressible, constant{density uids for which these equations read @ui = 0 @xi @ui + @ (uiuj ) + @ = @ ( 2Sij ) @t @xj @ xi @xj where Sij = (@uj [email protected] + @[email protected] )=2 is the strain{rate tensor and
!1 (1) (2) = p= . For later reference we introduce Reynolds averaging which isR used in statistical turbulence mod1 elling (RANS) as time averaging: hui = limT T 0T u dt. Reynolds averaging has the properties hhuii = hui huhvii = huihvi: (3) According to the idea of LES a means is required to distinguish between small, unresolved, and larger, resolved structures. This is accomplished by the operation u 7! u de ned below. Unlike the above Reynolds time averaging it is an operation in space. The fact that RANS and LES methods employ averaging in di erent dimensions inhibits an easy link between them. Several attempts have been made to put both in a common framework (Speziale 1998, Germano 1999) but will not be discussed here. We now turn to the ways of de ning u and illustrate them in the one{dimensional case. 2.1 Schumann's approach
The `volume{balance approach' of Schumann (1975) starts from a given nite volume mesh. The integral of a continuous unknown u(x) in (1),(2) over one cell is denoted V u = 1 R u(x)dx as illustrated in Fig. 1 (indices referring to cells are dropped). Integrating xV the NSE over a cell and using Gauss' theorem relates these values to surface{averaged quantities denoted j , such as j uv. These need to be expressed in terms of the cell{ averages, which is done in two steps. If the discretization is su ciently ne replacing j uv by j uj v as is usual in nite volume methods is possible with only a minor approximation error. This is done in DNS. If the grid is not ne enough, however,the di erence can be signi cant and the unresolved momentum ux j uv ; j uj v has to be accounted for by a model, the so{called subgrid{scale (SGS) model. Subsequently, j u are related to the V u either by setting them equal to cell averaged quantities if a staggered arrangement is used or by interpolating from neighbouring values. The nal SGS contribution to be modelled therefore also depends on the expressions used for j u , i.e. on the discretization scheme. To sum up, the equations are discretized and with discretizing them the splitting into large and small scales is performed since the latter cannot be resolved by the discrete system. Observe that the operations u 7! V u and u 7! j u map an integrable function onto discrete values, a continuous function u(x) is not constructed. Thus, with Schumann's approach, scale separation, discretization, and SGS model are not separated conceptually but intimately tied together. This has advantages in that anisotropies and inhomogeneities of the grid can easily be incorporated. However, it renders the analysis of the various contributions to the solution relatively di cult and hence is considered too restrictive by many workers in the eld. 2.2 Filtering
Leonard (1974) proposed to de ne u by u(x) = Z+ 1 ;1 G (x ; x ) u (x ) dx :
0 0 0 (4) An integral of this kind is called a convolution. Here, G is a compactly supported or at least R rapidly decaying lter function with qG(x) dx = 1 and width . The latter can be de ned by theq second moment of G as = 12 R x2G(x) dx. Fig. 2 displays the Gaussian Filter GG = 6= 1= exp(;6x2= 2) and the box lter de ned by GB = 1= if jxj =2 and GB = 0 elsewhere. In fact, already Deardor (1966) used (4) in the special case G = GB . Figs. 3a and 3b illustrate the ltering with smaller or larger lter width: the larger , the smoother is u. According to (4) u is a continuous smooth function as displayed in Fig. 3 which can subsequently be discretized by any numerical method. This has the advantage that one can separate conceptually the ltering from the discretization issue. It is helpful to transfer eq. (4) to Fourier space by means of the de nition u(!) = ^ R i!x dx, since in Fourier space, where the spatial frequency ! is the independent u(x) e variable, a convolution integral turns into a simple product. Eq. (4) then reads
; b b b u(!) = G(!) u(!) : (5) Figure 4 illustrates the ltering in Fourier space. Equation (5) allows the de nition of b another lter, the Fourier cuto lter with GF (!) = 1 if j!j = and 0 elsewhere. From b b (5) it is obvious that only this lter yields u = u since (GF )2 = GF . In all other cases the identity is not ful lled. This can be appreciated by comparing u and u for the box lter in Figs. 3 and 4. The second relation in (3) is never ful lled except in trivial cases, so that for general ltering we have u 6= u uv 6= u v (6) which distinguishes clearly the ltering in LES from Reynolds averaging (see Germano (1992) for a detailed discussion). The vertical line in Fig. 4 represents the nominal cuto d at = related to the grid. The Fourier cuto lter GF would yield a spectrum of u which is equal to the one of u left of this line and zero right of it. Eq. (5) and Fig. 4 therefore demonstrate that when a general lter is applied, such as the box lter, this does not yield a neat cut through the energy spectrum but rather some smoother decay to zero. This is important since SGS modelling often assumes that the spectrum of the resolved scales near the cuto follows an inertial spectrum with a particular slope and a particular amount of energy transported from the coarse to the ne scales on the average. We see that even if u ful lls this property this can be altered by the ltering (for further remarks see Section 4.3). Nevertheless, it is convenient and common to use the notion of a simple cuto as a model in qualitative discussions. Eq. (5) is also helpful to illustrate that derivative and lter commute, i.e. @ [email protected] = (@[email protected]). Any convolution lter (4) can be written as in (5) regardless of the choice of G. Di erentiation appears as multiplication by i! in Fourier space, eq. (28) below, which is commutative. Applying the three{dimensional equivalent of the lter (4) to the NSE (1) and (2), the following equations for the ltered velocity components ui result where Sij and @ ui = 0 @xi @ ui + @ (ui uj ) + @ = @ ( 2S ij ) ; @ ij @t @xj @ xi @xj @xj (7) (8) (9) are de ned analogous to the un ltered case. The term ij = ui uj ; ui uj represents the impact of the unresolved velocity components on the resolved ones and has to be modelled. In mathematical terms it arises from the nonlinearity of the convection term which does not commute with the linear ltering operation. An important property of ui is that it depends on time. Hence, an LES necessarily is an unsteady computation. Furthermore, ui always depends on all three space{dimensions (except for very special cases). Symmetries of the boundary conditions generally produce the same symmetries for the RANS variable hui i, e.g. vanishing dependence on a homogeneous direction. However, due to the very nature of turbulence this does not hold for ui since the instantaneous turbulent motion is always three{dimensional. The fact that a three{ dimensional unsteady ow is to be computed makes LES a computationally demanding approach. We nally note that for any lter the term in (9) vanishes in the limit ! 0, since then u ! u according to (4), and all scales are resolved so that the LES turns into a DNS. 2.3 Variable lter size
It should be mentioned here that ltering as de ned by (4) is not easily compatible with boundary conditions. Applying a box lter of constant size for instance yields u 6= 0 within a distance =2 outside of the computational domain and raises the question of how to impose boundary conditions for u. This issue is removed by supposing G to be x{dependent and locally asymmetric. However, if G(x ; x ) is generalized to some G(x x ) or if the prolongation of u from a nite domain to the real axis induces discontinuities, the commutation property is lost and additional commutator terms arise in (7), (8) (Ghosal & Moin 1995). In contrast to the usual SGS term ij , which is generated by the nonlinearity of the convection term, the commutator also appears for linear expressions (see the discussion by Geurts (1999) and Section 4.3). This issue is relevant for pronounced grid stretching in the interiour of the domain and close to walls but has been disregarded until recently. Studies for a channel ow are reported in (Frohlich et al. 1998, 2000).
0 0 2.4 Implicit versus explicit ltering
The ltering approach relaxes the link between the size of the computed scales and the size of the grid since the lter can be coarser than the employed grid. Consequently, the modelled motion should be called sub lter{ rather than subgrid{scale motion. The latter labelling results from the Schumann{type approach and is frequently used for historical reasons to designate the former. In practice, however, the lter G does not appear explicitly at all in many LES codes 2 so that in fact the Schumann approach is followed. Due to the conceptual advantages of the ltering approach reconciliation of both is generally attempted in two ways. The rst observation is that a nite di erence method for (7), (8) with a box lter employs the same discrete unknowns as Schumann's approach, such as u(xk ) = V u with k referring to a grid point. Choosing appropriate nite di erence formulae the same or very similar discretization matrices are obtained in both cases. Another argument is that the de nition of discrete unknowns amounts to an `implicit ltering' { i.e. ltering with some unknown lter ( but one that in principle exists ) { since any scale smaller than the grid is automatically discarded. In this way the lter is more or less used symbolically only to make the e ect of a later discretization appear in the continuous equations. This is easier in terms of notation and stimulates physical reasoning for the subsequent SGS modelling. In contrast to implicit ltering one can use a computational grid ner than the width of G and only retain the largest scales by some (explicit) ltering operation. This explicit ltering is recently being advocated by several authors such as Moin (1997) since it considerably reduces numerical discretization errors as the retained motion is always well resolved. On the other hand it increases the modelling demands since for the same number of grid points more scales of turbulent motion have to be modelled and it is up to now not fully clear which approach is more advantageous (Lund & Kaltenbach 1995). The ltering approach of Leonard is almost exclusively introduced today in papers on LES and has triggered substantial development, e.g. in subgrid{scale modelling. In practice, however, it is most often used rather as a concept than as a precise algorithmic construction.
k 3 Subgrid{scale modelling
3.1 Introduction
Subgridscale modelling is a particular feature of LES and distinguishes it from all other approaches. It is well{known that in three{dimensional turbulent ows energy cascades in the mean from large to small scales. The primary task of the SGS model therefore is to ensure that the energy drain in the LES is the same as obtained with the cascade fully resolved as in a DNS. The cascading, however, is an average process. Locally and instantaneously the transfer of energy can be much larger or much smaller than the average and can also occur in the opposite direction ("backscatter") (Piomelli et al. 1996). Hence, ideally, the SGS model should also account for this local, instantaneous transfer. If the grid scale is much ner than the dominant scales of the ow, even a crude model will su ce to yield the right behaviour of the dominant scales. This is due to two reasons. First, the larger the distance in wavenumber space between di erent contributions the looser is
2 Apart from some ltering operations for the dynamic model, discussed below, which is of a somewhat di erent nature their coupling. Second, as a consequence of this as well as of the energy cascading, the ner scales exhibit a more universal character which is more amenable to modelling. On the other hand, if the grid scale is coarse and close to the most energetic, anisotropic, and inhomogeneous scales the SGS model should be of better quality. Obviously, there exist two possible approaches, one is to improve the SGS model and the other is to re ne the grid. In the limit, the SGS contribution vanishes and the LES turns into a DNS. Re ning the grid, however is restricted due to the rapid increase of the computational cost. The alternative strategy, for example solving an additional transport equation in a more elaborate SGS model, can be comparatively inexpensive. Another aspect results from the numerical discretization scheme which introduces a difference between the continuous di erential operators and their discrete equivalents. This di erence is particularly large close to the cuto scale. For DNS this is not so disturbing, but with LES precisely these scales have a substantial in uence on the modelled SGS contribution as will be illustrated below. Hence, in LES the discretization scheme and the SGS model have to be viewed together. Indeed, some schemes such as low order upwind discretizations generate a considerable amount of numerical dissipation as discussed in Section 5.2 below. Therefore certain authors perform LES without any explicit SGS model (Tamura, Ohta & Kuwahara 1990, Meinke et al. 1998). The grid is re ned as much as possible to decrease the importance of the SGS terms, and the energy drain is in one way or another accomplished by the numerical scheme. Although yielding valuable results in some cases, this kind of modelling can hardly be evaluated or controlled. Hence, in most LES central or spectral schemes are used and the SGS term is represented by an explicit model. We shall now turn to the description of some basic SGS models before giving a summarizing view at the end of this section. 3.2 Smagorinsky model
The Smagorinsky model (SM) (Smagorinsky 1963) was the rst SGS model and is still widely used. As with most of the current SGS models, it employs the concept of an eddy a viscosity, relating the traceless part of the SGS stresses, ij , to the strain rate Sij of the resolved velocity eld: 1 (10) 3 ij rk = ;2 t Sij The advantage of (10) is that the resulting equation for ui to be solved looks like (2) with ui instead of ui, + 1 ij kk instead of , and + t instead of . Hence, it is very easy to 3 incorporate this into an existing solver for the unsteady NSE. The second part of this model is the determination of the eddy viscosity t. Dimensional analysis yields (11) t / l qSGS where l is the length scale of the unresolved motion and qSGS its velocity scale. From the above discussion it is natural to use the lter size as the length scale, hence setting l = Cs . Similar to Prandtl's mixing length model, the velocity scale is related to the gradients of ui expressed by q (12) qSGS = l jS j jS j = 2S ij S ij
a ij = ij ; = (Cs )2jSj (13) This amounts to assuming local equilibrium between the production of the SGS kinetic 3 a energy, P = ; ij S ij and dissipation " expressed by qSGS =l. Introducing (10) and (12) in P = gives (13). The constant Cs can be determined assuming an inertial{range Kolmogorov spectrum for isotropic turbulence which yields Cs = 0:18. This value has turned out to be too large for most ows so that often Cs = 0:1 or even lower values are used. Close to walls t has to be reduced to account for the anisotropy of the turbulence. This is generally accomplished by replacing Cs in (13) with CsD(y+ ). Most often the van Driest damping
t which yields (14) is used which is known from statistical models. However, this yields t (y+)2 for small y+ while t should behave like y+3. The correct behaviour is achieved by the alternative damping function (Piomelli, Moin & Ferziger 1993)
; D (y + ) = 1 ; e y+ =A+ A+ = 25 (15) The main reason for the frequent use of the SM is its simplicity. Its drawbacks are that the parameter Cs has to be calibrated and its optimal value may vary with the type of ow, the Reynolds number, or the discretization scheme. The kind of damping to be applied near a wall is a further point of uncertainty. Also, the SM, like any other model based on (10) with t 0, is strictly dissipative and does not allow for backscatter. It is furthermore not appropriate for simulating transition since it yields t 0 even in laminar ows.
; D (y + ) = 1 ; e (y+ =A+ )3 1=2 3.3 Dynamic procedure
From the previous section it is apparent that for physical reasons one would prefer to replace the constant value Cs by a value changing in space and time. The dynamic procedure has been developed by Germano et al. (1991) in order to determine such a value from the information provided by the resolved scales, in particular the ones close to the cuto scale. m a In fact this procedure can be applied with any model ij od(C u) for ij or ij containing a parameter C . 3 The basic idea is to employ the model chosen not only on the grid scale, or lter scale, but also on a coarser scale b as illustrated in Fig. 5. This is the so{called test scale with, e.g., b = 2 : sub{grid scale stresses ( {level) : sub{test scale stresses ( b {level) : = uiuj ; ui uj Tij = udj ; ci uj uc iu
ij mod u) ij (C mod b u) b ij (C (16) (17) u From the known resolved velocities ui the velocities ci can be computed by applying the b lter :c: to ui using an appropriate function G. Similarly, the term Lij = udj ; ci uj can : uc iu be evaluated. It is this part of the sub{test stresses Tij which is resolved on the grid as sketched in Fig. 5 : The total stresses uiuj in the expression for Tij can be decomposed
3 In this subsection we distinguish between exact and modelled SGS stresses for clarity. into the contribution ui uj resolved on the grid and the remainder ij . Inserting this in (17) gives Tij = Lij + cj (18) i known as Germano's identity. Hence, on one hand Lij can be computed, on the other hand the SGS model yields a model expression when inserting (16),(17) in (18) : Lmod = ij mod b u) ; mod (C b ij (C ij u) (19) (20) Ideally, C would be chosen to yield Lij ; Lmod = 0 ij but this is a tensor equation and can only be ful lled in some average sense, minimizing e.g. the root mean square of the left hand side as proposed by Lilly (1991). Principally, the b b consecutive application of G and G to obtain u yields an e ective lter, of width b 6= b , b and G (e.g. when the box lter is used). This which generally is even of di erent type as G issue is generally neglected in the literature. For that reason and since b = b with the Fourier cuto lter presently used for illustration we write b instead of b in this section. E ectively, it is the ratio b = which is required by the dynamic models. We now apply the dynamic procedure to the SM (10),(13) and get (21) with C = Cs2 for convenience. Classically, the model is developed by extracting C from the ltered expression in the second term although in fact C will vary in space. The right hand side of (21) can then be written as ;2CMij so that inserting into (20) with the least{squares minimization mentioned above yields 1L C = ; 2 Mij Mij (22) ij Mij The advantage of (22) or a similar equation is that now the parameter of the SM is no longer required from the user but is determined by the model itself. In fact, it is automatically reduced close to walls and vanishes for well{resolved laminar ows. Negative values of C are possible and can be viewed as a way to model backscatter. The resulting \backward di usion" can however generate numerical instability so that often + t 0 is imposed. Furthermore, C determined by (22) as it is, exhibits very large oscillations which generally need to be regularized in some way. Most often Lij and Mij are averaged in spatially homogeneous directions space before being used in (22). However, this requires the ow to have at least one homogeneous direction. Another way is to relax the value in time according to C n+1 = C + (1 ; )C n using C n from the previous time step (Breuer & Rodi 1994). Yet another way is to use the known value C n in the rightmost term of (21) so that it need not be extracted from the test lter (Piomelli & Liu 1995). This yields smoothing in space without any homogeneous direction required.
bb Lmod = ;2 C b 2jS jSij + 2 C 2jSjS ij ij 3.4 Scale similarity models
Scale similarity models (SSM) have been created to overcome the drawbacks of eddy viscosity{type models. Filtering the decomposition ui = ui + ui yields the (exact) relation
0 ui = ui ; ui:
0 0 0 0 0 0 0 0 (23) This can be interpreted as equality between the largest contributions of ui and the smallest contributions of ui (see Fig. 4). Furthermore, it is computable from ui. Introducing the decomposition of ui into (9) and modelling uiuj ui uj and uiuj ui uj , respectively, a yields the model ij = Lm a with ij Lm = ui uj ; ui uj ij (24) and :::a indicating the traceless part of a tensor. A model constant is not introduced as this would destroy the Galilean invariance of the expression. For a spectral cuto lter u p b is replaced by u with b = 2 since for this lter u = u as discussed above. The SSM allows backscatter, i.e. transfer of energy from ne to coarse scales, and does not impose alignment between the SGS stress tensor ij and the strain rate S ij . On the other hand, (24) turns out to be not dissipative enough so that it is generally combined with a Smagorinsky model. Horiuti (1997) subsumes some current SSMs in the model
a ij
0 0 0 0 = CLLm a + CB LR a ; 2C ij ij 2 j S j S ij (25) with LR = ui uj ; ui uj , evaluated using (23). A further step is to combine (25) with the ij dynamic procedure for the determination of the constants: (a) C with CL = 1 CB = 0 (Zang, Street & Kose 1993) (b) CL and C with CB = 0 (Salvetti & Banerjee 1995) (c) CB and C with CL = 1 (Horiuti 1997). Di erent tests in the cited references as well as by Piomelli, Yu & Adrian (1996) show that SSMs in conjunction with the dynamic procedure perform quite well for low order nite di erence or nite volume methods. Apart from the ability to represent backscatter this may also be due to the fact that no spatial derivatives are involved in the SSM which reduces the impact of numerical discretization errors. 3.5 Further models and comparative discussion
Let us sum up a few strategies or concepts which are currently followed in SGS modelling. One, already mentioned in the beginning of this section, is to employ a crude model and to compensate by grid re nement which decreases the impact of the model. Another strategy is to employ the same approaches as in RANS modelling. The Smagorinsky model based on an eddy{viscosity and an algebraic mixing{length expression is the most prominent example. But as with RANS, more elaborate methods can be used to compute the turbulent viscosity such as a model employing a transport equation for the SGS kinetic energy kSGS = 1=2 kk which furnishes a velocity scale qSGS = kSGS 1=2 (Schumann 1975, Davidson 1997). Obviously, the lter width constitutes an adequate reference length 1=2 so that according to (11) t = C kSGS is a reasonable model and no second length{scale determining transport equation is required. Spalart et al. (1997) have developed an approach called Detached Eddy Simulation (DES). They start from a one{equation RANS turbulence model (Spalart & Allmaras 1994) based on a transport equation for t. In this equation the distance from the wall is introduced as a length scale in the destruction term. Replacing this physical length scale by a resolution{based scale CD (where CD is a parameter) turns the model into a SGS model. This method furthermore o ers a particular way of wall modelling which is discussed below. Still more complex approaches have been carried over from RANS. Fureby et al. (1997) employ the SGS equivalent of a Reynolds{stress model and obtain satisfactory results in some tests. The cost increase is claimed to be moderate as solving the pressure equation requires most of the work. A third strategy, applied with SSM and the dynamic procedure, is based on the multiscale nature of turbulence. It could only be developed with scale separation de ned independently from the discretization according to (4) since ltering is used as an individual operation. By analyzing experimental ow elds along these lines Liu, Meneveau & Katz (1994) propose (26) ij = CL Lij with Lij de ned in (18). This di ers from (24) since two lters of di erent size are used. A fourth strategy is to relate SGS models to classical theories of turbulence. An elementary example is the determination of the Smagorinsky constant assuming a Kolmogorov spectrum. This strategy is also pursued when de ning a wave{number dependent eddy viscosity to be employed with a spectral Fourier discretization and using EDQNM theory to determine t(k) (Chollet & Lesieur 1981). The spectral eddy viscosity model has also been reformulated in physical space for application in complex ows yielding the structure function model (Metais & Lesieur 1992)
t = 0:063 q F2 ( ) F2(r) = (ui(x + r) ; ui(x))2 (27) with F2 spatially averaged in an appropriate way. Di erent variants have been developed (Lesieur & Metais 1996). It can be shown that when implemented in a Finite Di erence context, this yields a Smagorinsky{type model with jS j in (13) replaced by j@ [email protected] j. The last strategy we mention concerns the testing of SGS models. Of course, as with other turbulence models, prototype ows can be computed and the results then compared with experimental data or DNS results. This is still the ultimate test to pass. However, another kind of testing particular to LES has been developed, namely the so{called a priori test: a fully resolved velocity eld from a DNS is used to explicitly compute the terms which have to be modelled in an LES on a coarser grid. The large{scale velocity on that grid is extracted to determine the SGS stresses by means of a SGS model. The di erence between exact and modelled stresses re ects the quality of the model. This information should however be taken with some caution as the test involves discretization e ects in a substantially di erent way than in the actual LES. Finally, one has to bear in mind that a perfect SGS model is impossible. Assume the exact grid{scale velocity u being known at all points. The perfect SGS model would then amount to inferring from u on the exact instantaneous SGS velocity u to deduce the exact instantaneous SGS stresses at all points. Since, however, in nitely many velocity elds u are compatible with the same u, even the best SGS model cannot decide which of them is realized in the actual DNS. In fact, the error introduced by missing SGS information propagates in an inverse cascade to larger scales (Lesieur 1997).
0 0 4 Numerical methods
4.1 Discretization schemes in space and time
With the ltering approach discussed in Section 2, physical modelling and numerical discretization are conceptually independent. Hence any available numerical method can in principle be used to discretize the ltered equations. A minimal requirement for precision and cost{e ectiveness is that the discretization scheme is at least of second order in space and time. Classically, spectral methods were frequently used for LES and are still employed for problems with simple geometry. Derivatives are discretized most accurately and ltering and de ltering as discussed below is naturally applied in this framework. For more complex boundary conditions, nite{di erence or nite{volume methods are prefered. Here, one current trend goes to unstructured meshes, another to grids with special local treatment at the boundary if this has an irregular shape. Some numerical methods favour particular modelling ideas such as spectral methods which allow to use a spectral eddy viscosity (Chollet & Lesieur 1991) and explicit ltering by means of (5). Others need particular care in certain points. For example if implicit ltering is used together with a type of nite element that has a di erent number of degrees of freedom for velocity compared to pressure (which is classically the case for stability reasons), this results in a di erent amount of ltering for these quantities and can deteriorate the result (Rollet{Miet, Laurence & Ferziger 1999). Discretizations in space can be selected according to relevant properties such as ability to treat complex geometries, cost per grid point, etc. If possible, however, equispaced grids are used since the in uence of grid{inhomogeneity and grid{anisotropy on SGS modelling is not yet fully mastered. A comparative study of a structured and an unstructured method for the same problem was undertaken by Frohlich et al. (1998) where, for the particular case considered, adaptivity of the former method was roughly compensated by higher cost per node. For more complex geometries unstructured methods are certainly favourable. Concerning the time scheme we already noted that temporal resolution has to be compatible with resolution in space so that C = u t= x O(1). Since this type of limit is equivalent to the stability limit of explicit methods, in many cases the latter are typically used for LES. Adams{Bashforth, Runge{Kutta or leap{frog schemes are the most popular ones. If the di usion limit is stricter in a computation, semi{implicit time stepping can be more e cient. 4.2 Analysis of numerical schemes for LES
The essential feature that distinguishes LES from DNS is that the smallest resolved grid{ scale components, which are just a little larger than the cuto scale, typically carry more energy. Hence, without explicit ltering which employs a lter coarser than the mesh size, the smallest resolved scales are by de nition substantially a ected by the employed numerical scheme. These scales however in uence most strongly the contribution determined by the SGS model. In fact a complex discrete model for the SGS e ect on the resolved ow is created which results from physical as well as numerical modelling. Consequently, the order of a method is not necessarily an appropriate notion in the context of LES. It rather has to be supplemented with a re ned analysis like the modi ed wavenumber concept as e.g. discussed by Ferziger (1996). Let us illustrate these statements by means of Fig. 6. Refering to eq. (5), the exact spatial derivative of u formulated in Fourier space is
b b (28) @x (!) = i! G(!) u(!): The numerical evaluation of @ [email protected] by a nite{di erence formula corresponds to replacing the factor ! in (28) by a modi ed wavenumber !eff (!) which depends on the particular d! @u scheme employed. Derivation and formulas are given, e.g., by Ferziger & Peric (1996). For symmetric schemes this is a real quantity, otherwise it is complex. Starting from !eff (0) = 0, j!eff j increases and then drops down to zero again. The point ! where this takes place is determined by the numerical grid employed it is the highest frequency resolved by the grid. In a DNS this point would be pushed as far as possible to the right (Fig. 6a). The order of a discretization scheme can be reformulated in terms of the exponent p in lim! 0 j! ; !eff (!)j / !p. Obviously, the information about the order of a scheme is su cient only if ! ! 1 so that the solution to be computed is located entirely at !=! 0. If however ! approaches the relevant scales of the solution to be discretized, the behaviour of the whole curve !eff is decisive, not only the limit ! ! 0. It is rather b visible that computing a derivative with !eff 6= ! amounts to replacing G(!) in (28) by b b Geff = !eff =! G(!). Hence the nite di erence formula results in additional ltering applied to the derivative of u (Salvetti & Beux 1998). Fig. 6b furthermore shows that the decay of the error obtained with grid re nement in an LES depends on the behaviour of the solution itself, e.g. on the decay rate of its spectrum. This information is indispensable when aiming to assess the numerical error in an LES and to compare it with the size of the SGS term (Ghosal 1996). So far we have discussed the discrete derivative operator which is a building block when discretizing the whole system of equations. Qualitatively, the real part of !eff =! in the convection term introduces spurious or numerical dispersion while its imaginary part results in additional numerical dissipation. Analyzing the fully discrete system is much more complicated but can be achieved when disregarding boundary conditions etc. by the modi ed equation approach (Hirt 1968). This has been applied by Werner (1991) to a staggered nite volume discretization with Adams{Bashforth time scheme, central di erences for the viscous term, and the QUICK convection scheme. Recall that the QUICK scheme (Leonard 1979) is a third order upwind interpolation scheme for the ux over the surface of a control volume. Werner observed that this combination results in a spurious 4th order dissipative term proportional to the cell Reynolds number Recell = u x= . The same analysis for a leap{frog time scheme with second order central di erencing yields a 4th order error term which is independent of Recell. Such an analysis nicely shows that the upwind scheme produces excessive damping for large Recell . This, however, is precisely the working range of LES with typically Recell = O(10000), even if is replaced by + t. To demonstrate the e ect in a real LES, Werner also computed a plane channel with Re = 1954 employing the Smagorinsky model. With a modied leap{frog scheme t= = 1:2 C = 0:14 Recell = 1350 on the centerline. Analysis yields j num = j 0:2. With the QUICK upwind scheme the corresponding numbers are 180. Hence the numerical dissipation introduced by the upt = = 0:9 and j num = j winding exceeds the one by the SGS model by two orders of magnitude. Similar though mostly less detailed experiences have been reported in several papers by comparing the solution obtained with di erent schemes. The QUICK scheme and lower order upwind schemes gave worse results than a second order central scheme in LES of a circular cylinder (Breuer 1998). This was, though with decreasing impact, also observed for 5th and 7th
! order upwinding (Beaudan & Moin 1994). Further studies of the numerical error in LES were performed by Vreman, Geurts & Kuerten (1994), Kravchenko & Moin (1997) and others. Hence, on the one hand upwind schemes can spoil the result by excessive damping. On the other hand, some researchers omit explicit SGS modelling and let the numerical dissipation of the employed scheme remove the energy. The MILES approach (Boris et al. 1992), e.g., falls into this class. Comte and Lesieur (1998) however found such schemes to be inferiour to explicit SGS modelling. In contrast to the numerical dissipation, the dispersion of a scheme is of lower importance as it has no e ect on the energy drain which has been claimed to be the principal task of the SGS modelling. Dispersion, however, is related to the generation of spurious wiggles which in some LES of blu bodies pose problems (Rodi et al. 1997). 4.3 Further developments
In order to improve the current status, attempts are made to separate more clearly the di erent ingredients in an LES. The aim is to study and improve each of them in a separate and controlled way. One of the directions pursued is explicit ltering as mentioned above. A similar approach is used by Vreman, Geurts & Kuerten (1994,1997) who use a value of larger than the mesh size of the grid, e.g. by a factor of two in the SGS model, which leads to increased SGS dissipation e ectively damping the solution in a similar way as explicit ltering. A second direction is the use and improvement of higher order energy{conserving discretization schemes (Morinishi et al. 1998). They ensure that the total dissipation is entirely controlled by the SGS model and not by the discretization. Bearing in mind the uncertainty in SGS modelling, when for example determining the parameter C in the DSM, the practical importance of an energy{conserving scheme is presently not clear. A third direction is to use higher order methods as they narrow the range of scales which are in uenced by the discretization of the ltered equations. Furthermore, lters recently have been devised so as to commute with discrete derivatives (Vasilyev, Lund & Moin 1998). This ideally ensures that apart from the SGS term ij in eq. (8) no commutator term arises which would require modelling. Finally, \de ltering" has been applied to invert the attenuation of resolved scales by the implicit ltering related to the discretization (Stolz, b b Adams & Kleiser 1999). It uses an operation like multiplication of u(!) with G(!) 1 to devise an estimate u for the true velocity u based on the resolved velocity u, as may be illustrated with eq. (5) and Fig.4. This procedure requires higher order methods and principal adjustments such as restriction to certain scales since inverse ltering is ill{ conditioned (imagine obtaining u from u in Fig.3a or 3b by backward di usion). First results look promising.
; 5 Boundary conditions
We have mentioned already the mathematical problems when boundary conditions for ltered quantities have to be de ned. From a physical point of view, the ow near a solid wall exhibits substantially di erent structures than away from it. In this region the \large scales" { in the sense that they signi cantly determine the overall properties { are of the order of the boundary layer thickness and hence typically much smaller than in the core of the ow, in particular if the Reynolds number is large. In addition, the small scales in this area exhibit substantial anisotropy and energy transfer mechanisms are di erent compared to the core ow (Hartel 1996, Piomelli et al. 1996). This makes subgrid{scale modelling in the vicinity of walls a di cult task. 5.1 Resolution of the near{wall region
The most natural boundary condition at a wall is the no{slip condition. It requires however that the energy{carrying motion is resolved down to the wall. In an attached boundary layer this motion is mainly constituted by the well{known streaks of spanwise distance z+ = 15 ; 40. (Piomelli z 100 resolution of which requires y + < 2 x+ = 50 ; 150 & Chasnov 1996). The resulting simulation is in fact a hybrid one between an almost{DNS near the wall and an LES in the main part of the ow. If locally a ne grid is required an e cient discretization calls for a block{structured or an unstructured method. Care has to be taken however, since for example a low order FV method that locally splits each cell into a number of smaller ones introduces a sudden decrease in the size of the implicit lter by a factor of two at least in one direction. This may lead to problems with the SGS modelling. Kravchenko, Moin & Moser (1997) have developed a discretization that employs overlapping B{splines and hence results in a smoother transition of the e ective resolution. It was successfully applied to channel ow up to Re = 109410. Another possibility is an unstructured nite element method as used by RolletMiet et al. (1999). Due to particular discretization issues discussed in this reference only very few codes of this kind are capable of LES up to now. With an unstructured code, one is still left with the task of generating a grid that ful lls the needs, in particular with respect to its in uence on SGS modelling. Regardless which method is used to discretize and compute the near{wall region, the wall{resolving approach can result in substantial complexity and computational e ort. Spalart et al. (1997) stressed that re nement needs to be performed not only in the wall{normal but also in the streamwise and spanwise directions and estimated that O(1011) grid points would be necessary for a wing at Re = 6:5 106 while "108 is impressive today". Also, resolving the ow in space is worthless if it is not also resolved properly in time, hence the CFL number has to be of order unity, a fact that even further increases the computational burden. To conclude: although a wall{resolving LES is appropriate for lower Re and transitional ows, a di erent approach is needed for higher Re, particularly when the interest of a simulation focuses on features away from the wall. 5.2 Wall functions
When for higher Re a wall{resolving LES is not possible, the way out is to use a near{wall model approximating the overall dynamic e ects of the streaks on the larger outer scales which are resolvable by the LES. The most commonly used models are wall functions for bridging a region very close to the wall, often the viscous sublayer. Such wall functions are classically used in RANS methods, where they take the form: h ! i = W (hu1i y1) (29) Here, y1 is the distance of the rst grid point from the wall, h ! i the average wall shear stress, hu1i the average tangential velocity at y1, and W a functional dependence. Note that any relation hui+ = f (y+ ) can be converted to the form (29). This can be the log{ law, the 1/7{power law, a linear viscous law or even a numerical t to DNS data. An + appropriate blending is generally used so that W is de ned from y1 0 to y+ of several hundreds. Even if many physical properties such as the low{order moments and hence W are well{ know for a certain ow { this is the case for the developed ow in a plane channel, for example { it is a delicate task to introduce this knowledge in the context of LES. The available information is of a statistical nature whereas the ltered velocity is an instantaneous, uctuating quantity. On the other hand it has been demonstrated that the inner and outer regions of a turbulent boundary layer are only loosely coupled (Brooke & Hanratty (1993)) so that an arti cal boundary condition bridging the inner layer has a chance to be successful (Piomelli 1998). One of these wall{function methods (Schumann 1975) employs the mean velocity hu(y1)i, which is successively computed during the LES, to determine the average wall shear stress h ! i from (29) with W being the logarithmic law of the wall. The same proportionality as between hu(y1)i and h ! i is then assumed to hold also between the instantaneous quantities u(y1) and ! , in particular they are supposed to be in phase. This yields the instantaneous wall stress as
! h! = hu(y i)i u(y1)
1 (30) which is used as a boundary condition in the LES. Werner & Wengle (1993) employed the 1=7{power law instead of the log law to avoid an iterative evaluation of W . Furthermore, they replaced (29) by ! = W (u1 y1) so that the average velocity is no longer needed. Other combinations and variants are possible as well. A generalization of the approach for the subcritical ow around a cylinder is described by Frohlich (2000). In the technically relevant case of a rough surface, using the wall{function approach is unavoidable since resolving the ow around each roughness element is impossible. The roughness e ect is brought in by the roughness parameter in the log{law (Grotzbach 1977). Wall function boundary conditions work reasonably well in simple ows and save substantial CPU time due to the reduced resolution requirements. They have also been applied to some complex ows around obstacles (Rodi et al. 1997) which, however, were found not to be very sensitive to variations in the boundary conditions. 5.3 Other approaches
Wall functions establish a relation between the local wall shear stress and the velocity at the wall{adjacent grid point. This can be generalized to the case where information on a line or a whole plane at some distance parallel to the wall is used to generate the wall{shear stress at a certain point. Such information can be introduced as a boundary condition in unsteady turbulent boundary layer equations which are solved along the wall within the wall cell using an embedded grid (Balaras, Benocci & Piomelli 1996) (cf. Fig. 7d). In these equations turbulence is modelled with an eddy viscosity depending on the wall distance. Similar work has been done by Cabot (1995,1996) where di erent models of this type were devised and applied to ow in a plane channel and over a backward facing step. Although yielding better results in some computations compared to wall functions, these methods have not found wide application yet due to their implementational complexity. Another approach that can be introduced more easily in an LES is based on using the no{slip condition which in turn requires re nement in the wall{normal direction. Parallel to the wall, however, the step size of the outer region is maintained leading to substantial savings. The idea then is to replace the unresolved near{wall structures by elements from RANS simulations. Schumann (1975) decomposed ij into an isotropic part for which the Smagorinsky model is used and an anisotropic part resulting from the mean ow gradient. The latter part is modelled with an eddy viscosity t an = min(c x z d)dhui=dy This is a RANS{like model in which close to the wall the size of the grid x z is replaced by the distance from the wall d as a length scale. A similar switch is used in the DES approach of Spalart et al. (1997) mentioned above so that close to a wall the original RANS model is employed (see Fig. 7b). DES, although conceived for di erent applications, has been tested for channel ow by Nikitin et al. (2000). The authors observe a spurious bu er layer re ecting di culties in connecting the quasi{steady RANS layer close to the wall to the outer unsteady computation. Further adjustments need to be introduced to apply DES in such cases. Bagett (1998) discussed the issue of blending RANS with LES turbulence models and points out the requirements for adequate spanwise resolution in the near{wall region in order to capture the energetically dominant features. If these are not captured, unphysical structures are generated which degrade the result. 5.4 In ow and out ow conditions
After discussing the boundary conditions at solid walls we brie y mention the conditions at arti cial boundaries, an issue shared with DNS. Turbulent out ow boundaries are relatively uncritical. Here, damping zones or convective conditions are generally applied which allow vortices to leave the computational domain with only small perturbations of the ow in its interior. A convective condition for a quantity reads @ + U @ =0 (31) conv @t @n applied on the outlet boundary with n the outward normal coordinate and Uconv an appro priate convection velocity such as the bulk velocity. The di culty posed by turbulent in ow conditions stems from the fact that LES computes a substantial part of the spectrum and hence requires speci cation of the in ow conditions in all this spectral range, not just the mean ow. The need for this information can be avoided by imposing streamwise periodicity with a su cient periodic length, but this is inapplicable in many practical ows. Imposing the mean ow plus random perturbations is generally not successful since these perturbations are unphysical so that a large upstream distance must be computed to produce the correct turbulence statistics. With more sophisticated perturbations the distance can be shortened. This is a subject of current research. If feasible the best solution is to impose some fully developed ow at the inlet. A separate companion LES, e.g. with streamwise periodicity, can then be performed to generate velocity signals at the grid points in the in ow plane of the main LES. An example is the ow around a single cube investigated in Rodi et al. (1997). 5.5 Sample computations
In order to illustrate the above discussion, we present results from a standard test case for LES calculations, namely fully developed plane channel ow. For this ow between two in nitely extended plates, periodic conditions can be imposed in the streamwise direction x and the spanwise direction z, with typical domain sizes of Lx = 2 and Lz = , respectively. Reference quantities are the channel half{width and the bulk velocity Ub . DNS of this ow has been performed for low and medium Reynolds number of which the currently highest is Reb = Ub = = 10935 (Moser, Kim & Mansour 1999) employing a high{precision spectral method with 384 257 384 points. This Reynolds number has been used in the computations below. The results have been obtained with the structured collocated nite volume code LESOCC developed by Breuer & Rodi (1994). The dynamic Smagorinsky model was used with test{ ltering and averaging in planes parallel to the walls. The bulk Reynolds number was xed and an external pressure gradient adjusted so as to yield the desired ow rate. Fig.8 shows a computation where resolving the near{wall ow has been attempted, i.e. no wall{function was used. The number of points in the y direction is 65 and a stretching of 11% has been applied to cluster them close to the walls. The gure shows the average streamwise velocity and the rms{ uctuations. The computed shear stress w yields Re = q 504 which is much below the Re = 590 in the DNS. The value of u = w = determines the scaling of the axes, and in particular the hui+ = f (y+){curve is quite sensitive to it. If the v{ and w{ uctuations are plotted in outer scaling, i.e. not by u , they are even further below the DNS curves. The observed failure occurs because resolving the ow near the wall requires adequate discretization in all three directions, not just normal to the wall. Here in particular the spanwise resolution is too coarse. In Fig.9 the wall{normal resolution has been improved using 159 points in y with clustering in the bu er layer accompanied by a substantially better resolution in spanwise direction. The computed shear stress yields Re = 598:5 and the Reynolds stresses compare quite satisfactorily with the DNS data. It is obvious that with a structured discretization the grid in the interiour of the channel is ner than it really needs to be, due to the requirements near the wall. To avoid this, a method with local re nement is bene cial as discussed above. Fig. 10 presents a computation using the wall function of Schumann (1975). In the wall{normal direction 39 equidistant volumes are used. In this case the rst cell center is still located in the bu er layer, but the viscous sublayer is well bridged. The computed wall shear stress gives Re = 589:5 which compares very well with Re = 590 in the DNS. The Reynolds stresses are well predicted, the v{ uctuations being somewhat too small. With this approach it is of course not possible to reproduce the peak in the u{ uctuations close to the wall. For an LES using a wall function the present Reynolds number is relatively low. With higher Re the rst point usually lies beyond the bu er layer at y+ 100. In Fig. 11 we nally show cuts of the instantaneous u{ and w{velocity of the third case. Straight lines have been inserted connecting the data points. The angles they form show that, as discussed in Section 4, on the grid level the discrete solution is not smooth, i.e. the velocity scales close to the cuto are hardly resolved. Hence, any gradient computed from these values by, e.g., a second order scheme can only be a crude approximation to the "true" gradient. Recall that gradients enter in the contribution of the SGS model. This gure illustrates the close interplay between the numerical discretization and the subgrid{scale modelling. The amount of SGS dissipation in eddy{viscosity models can be monitored by the ratio t= . It varies locally and in the above computations attains values up to 7 in the last case which shows the dominance of the SGS dissipation with respect to the resolved dissipation. Further applications of LES, in particular to blu body ows, are discussed in other chapters of this volume. 6 Concluding remarks
We have described the kconcept of the Large{Eddy Simulation technique which in fact is extremely simple and makes it appealing. It turns out, however, that several issues are not simple for numerical or physical reasons. We have aimed at making the reader aware of these points and at clarifying related concepts. In practice, LES is characterized by a large number of decisions concerning the numerical and physical modelling which have to be taken and which all in uence the nal result. Thorough testing is still a major occupation of the community, and this will presumably not change in the near future. On the other hand, LES has a potential on several levels. The rst is the determination of statistical quantities such as the average ow eld with a higher accuracy than obtained by statistical models. This is based on interchanging the order \ rst averaging { then computing" (RANS) to \ rst computing { then averaging" (LES). To pay o , the drastic increase in cost has to be justi ed by an improved quality of the results. The next level is the determination of statistical quantities which are inaccessible to RANS such as two{point correlations. The third level is to use the instantaneous information on the structure of the ow in order to improve the understanding of vortex dynamics, transition phenomena, etc. or to determine dynamic loading. Finally, this information can be coupled to other physical processes either within the ow eld such as the generation of sound, the transport of scalars (temperature, sediment, ...), chemical reactions, etc., or to external processes such as the dynamical response of a solid structure. This type of feature is required for important elds of research such as uid{structure aerodynamic coupling and turbulence control. In the future, the applications of LES will turn from the present academic cases to more applied con gurations. References
J.S. Baggett. On the feasibility of merging LES with RANS for the near{wall region of attached turbulent ows. In Annual Research Briefs 1998, pages 267{277. Center for Turbulence Research, 1998. E. Balaras, C. Benocci, and U. Piomelli. Twolayer approximate boundary conditions for large{eddy simulations. AIAA Journal, 34:1111{1119, 1996. P. Beaudan and P. Moin. Numerical experiments on the ow past a circular cylinder at sub{critical Reynolds number. Technical Report TF62, Stanford University, 1994. J.P. Boris, F.F. Grinstein, E.S. Oran, and R.L. Kolbe. New insights into large eddy simulation. Fluid Dyn. Res., 10:199, 1992. M. Breuer. Large eddy simulations for the ow past a circular cylinder, numerical and modelling aspects. Int. J. Numerical Methods in Fluids, 28, 1998. M. Breuer and W. Rodi. Large eddy simulation of turbulent ow through a straight square duct and a 180o bend. In P.R. Voke, R. Kleiser, and J.P. Chollet, editors, Fluid Mech. and its Appl., volume 26. Kluwer Academic, 1994. J.W. Brooke and T. J. Hanratty. Origin of turbulenceproducing eddies in a channel ow. Phys. Fluids, 5:1011{1022, 1993. W. Cabot. Large{eddy simulations with wall models. In Annual Research Briefs 1995, pages 41{50. Center for Turbulence Research, 1995. W. Cabot. Near{wall models in large eddy simulations of ow behind a backward facing step. In Annual Research Briefs 1996, pages 199{210. Center for Turbulence Research, 1996. J.P. Chollet and M. Lesieur. Parameterization of small scales of three dimensional isotropic turbulence. J. Atmos. Sci., 38:2747{2757, 1981. P. Comte and Marcel Lesieur. Large eddy simulation of compressible turbulence. In Advances in Turbulence Modelling, number 1998{05 in Lecture Series. Von Karman Institute for Fluid Dynamics, Rhode Saint Genese, Belgium, 1998. L. Davidson. Large Eddy Simulation: A dynamic one{equation subgrid model for three{dimensional recirculating ow. In 11th Symp. on Turbulent Shear Flows, volume 3, pages 26.1{26.6, Grenoble, 1997. J. Ferziger and M. Peric. Computational Methods for Fluid Dynamics. Springer Verlag, 1996. J.H. Ferziger. Large eddy simulation. In T.B. Gatski, M.Y. Hussaini, and J.L. Lumley, editors, Simulation and modelling of turbulent ows, ICASE/LaRC Series in Comp. Sci. and Eng., pages 109{154. Oxford University Press, New York, 1996. J. Frohlich. LES of vortex shedding past circular cylinders. to appear in Proceedings of ECCOMAS 2000, Barcelona, 11{14 September, 2000. J. Frohlich, W. Rodi, Ph. Kessler, S. Parpais, J.P. Bertoglio, and D. Laurence. Large eddy simulation of ow around circular cylinders on structured and unstructured grids. In E.H. Hirschel, editor, Notes on Numerical Fluid Mechanics, volume 66, pages 319{338. Vieweg, 1998. J. Frohlich, W. Rodi, J.P. Bertoglio, U. Bieder and H. Touil. Large eddy simulation of ow around circular cylinders on structured and unstructured grids, II. In E.H. Hirschel, editor, Notes on Numerical Fluid Mechanics, to appear. C. Fureby, G. Tabor, H.G. Weller, and A.D. Gosman. Di erential subgrid stress models in large eddy simulations. Phys. Fluids, 9:3578{3580, 1997. M. Germano. Turbulence: the ltering approach. J. Fluid Mech., 238:325{336, 1992. M. Germano. From RANS to DNS: Towards a bridging model. In P.R. Voke, N.D. Sandham, and L. Kleiser, editors, Direct and Large{Eddy Simulation III, volume 7 of ERCOFTAC Series, pages 225{236. Kluwer Academic, 1999. M. Germano, U. Piomelli, P. Moin, and W.H. Cabot. A dynamic subgrid{scale eddy viscosity model. Phys. Fluids A, 3:1760{1765, 1991. B. J. Geurts. Balancing errors in large{eddy simulation. In P.R. Voke, N.D. Sandham, and L. Kleiser, editors, Direct and Large{Eddy Simulation III, volume 7 of ERCOFTAC Series, pages 1{12. Kluwer Academic, 1999. S. Ghosal. An analysis of numerical errors in large{eddy simulations of turbulence. J. Comput. Phys., 125:187{206, 1996. G. Grotzbach. Direct numerical simulation of secondary currents in turbulent channel ow. In H. Fiedler, editor, Lecture Notes in Physics, volume 76, pages 308{319. Springer{Verlag, 1977. C. Hartel. Turbulent ows: direct numerical simulation and large{eddy simulation. In R. Peyret, editor, Handbook of Computational Fluid Mechanics, pages 283{338. Academic Press, 1996. C. W. Hirt. Heuristic stability theory for nite{di erence equations. J. Comp. Phys., 2:339{355, 1968. K. Horiuti. A new dynamic two{parameter mixed model for largeeddy simulation. Phys. Fluids, 9:3443{3464, 1997. A. G. Kravchenko and P. Moin. On the e ect of numerical errors in large eddy simulation of turbulent ows. J. Comp. Phys., 131:310{322, 1997. A.G. Kravchenko, P. Moin, and R. Moser. Zonal embedded grids for numerical simulations of wall{bounded turbulent ows. J. Comp. Phys., 127:412{423, 1997. A. Leonard. Energy cascade in large eddy simulations of turbulent uid ows. Adv. Geophys., 18A:237, 1974. M. Lesieur. Turbulence in Fluids, volume 1 of Fluid mechanics and its applications. Kluwer Academic Publisher, P.O. Box 322,3300 AH Dordrecht, The Netherlands, 3 edition, 1997. M. Lesieur and O. Metais. New trends in large{eddy simulations of turbulence. Ann. Rev. Fluid. Mech., 28:45{82, 1996. D.K. Lilly. A proposed modi cation of the Germano subgridscale closure method. Phys. Fluids A, 4:633{635, 1991. S. Liu, C. Meneveau, and J. Katz. On the properties of similarity subgrid{scale models as deduced from measurements in a turbulent jet. J. Fluid Mech., 275:83{119, 1994. T.S. Lund and H.J. Kaltenbach. Experiments with explicit ltering for LES using a nite{di erence method. In Annual Research Briefs 1995, pages 91{105. Center for Turbulent Research, 1995. M. Meinke, Th. Rister, F. Rutten, and A. Schvorak. Simulation of internal and free turbulent ows. In H.J. Bungartz, F. Durst, and C. Zenger, editors, High Performance Scienti c and Engineering Computing. Springer, 1998. O. Metais and M. Lesieur. Spectral large{eddy simulation of isotropic and stable{strati ed turbulence. J. Fluid Mech., 239:157{194, 1992. P. Moin. Numerical and physical issues in large eddy simulation of turbulent ows. In Proceedings of the International Conference on Fluid Engineering, Tokyo, July 1316, 1997, volume 1, pages 91{100. Japan Society of Mechanical Engineers, 1997. R.D. Moser, J. Kim, and N.N. Mansour. Direct numerical simulation of turbulent channel ow up to Re = 590. Phys. Fluids, 11:943{946, 1999. N.V. Nikitin, F. Nicoud, B. Wasisto, K.D. Squires, and P.R. Spalart. An approachto wall modellingin large{eddysimulations. submitted, 2000. U. Piomelli. Large eddy simulation turbulent ows. In Advances in Turbulence Modelling, number 1998{05 in Lecture Series. Von Karman Institute for Fluid Dynamics, Rhode Saint Genese, Belgium, 1998. U. Piomelli, W.H. Cabot, P. Moin, and S. Lee. Subgrid{scale backscatter in turbulent and transitional ows. Phys. Fluids A, 3:1766{1771, 1991. U. Piomelli and J. Liu. Large eddy simulation of rotating channel ows using a localized dynamic model. Phys. Fluids, 7:839{848, 1995. U. Piomelli, P. Moin, and J. H. Ferziger. Model consistency in large eddy simulation of turbulent channel ows. Phys. Fluids, 31:1884{1891, 1988. U. Piomelli, Y. Yu, and R. J. Adrian. Subgrid{scale energy transfer and near{wall turbulence structure. Phys. Fluids, 8:215{224, 1996. W. C. Reynolds. The potential and limitations of direct and large eddy simulations. In J.L. Lumley, editor, Lecture Notes in Physics, volume 357, pages 313{343. Cornell University Ithaca, Springer{Verlag, 1989. W. Rodi, J.H. Ferziger, M. Breuer, and M. Pourquie. Status of large eddy simulation: Results of a workshop. J. Fluid Eng., 119:248{262, 1997. P. RolletMiet, D. Laurence, and J.H. Ferziger. LES and RANS of turbulent ow in tube bundles. Int. J. Heat Fluid Flow, 20:241{254, 1999. M. V. Salvetti and S. Banerjee. A priori tests of a new dynamic subgrid{scale model for nite{di erence large{eddy simulations. Phys. Fluids, 7:2831, 1995. M.V. Salvetti and F. Beux. The e ect of the numerical scheme on the subgrid scale term in largeeddy simulation. Physics of Fluids, 10:3020{3022, 1998. U. Schumann. Subgrid scale model for nite di erence simulations of turbulent ows in plane channels and annuli. J. Comput. Phys., 18:376{404, 1975. J.S. Smagorinsky. General circulation experiments with the primitive equations, I, the basic experiment. Mon. Weather Rev., 91:99{164, 1963. P.R. Spalart and S.R. Allmaras. A one{equation turbulence model for aerodynamic ows. La Recherche Aerospatiale, 1994. P.R. Spalart, W.H. Jou, M. Strelets, and S.R. Allmaras. Comments on the feasibility of LES for wings, and on a hybrid RANS/LES approach. In C. Liu and Z. Liu, editors, Advances in DNS/LES. Greyden Press, 1997. C. G. Speziale. A combined large{eddy simulation and time{dependent rans capability for high{speed compressible ows. J. Sci. Comput., 13:253{274, 1998. S. Stolz, N.A. Adams, and L. Kleiser. The approximate deconvolution procedure applied to turbulent channel ow. In P. Voke, N.D. Sandham, and L. Kleiser, editors, Direct and Large{Eddy Simulation III. Kluwer Academic, 1999. T. Tamura, I. Ohta, and K Kuwahara. On the reliability of two{dimensional simulation for unsteady ows around a cylider{type. J. of Wind Eng. and Indust. Aerodyn., 35:275{298, 1990. H. Tennekes and J.L. Lumley. A rst course in turbulence. MIT Press, 1972. U.Piomelli and J.R. Chasnov. Large{Eddy Simulations: theory and applications. In M. Hallback et al., editor, Turbulence and Transition Modelling, pages 269{331. Kluwer Academic, 1996. O.V. Vasilyev, T.S. Lund, and P. Moin. A general class of commutative lters for LES in complex geometries. J. Comput. Phys., 146:82{104, 1998. B. Vreman, B. Geurts, and H. Kuerten. Discretization error dominance over subgrid terms in Large Eddy Simulation of compressible shear layers in 2d. Communications in Numerical Methods in Engineering, 10:785{790, 1994. B. Vreman, B. Geurts, and H. Kuerten. On the formulation of the dynamic mixed subgrid{scale model. Phys. Fluids, 6:4057{4059, 1994. B. Vreman, B. Geurts, and H. Kuerten. Large{Eddy simulation of the turbulent mixing layer. J. Fluid Mech., 339:357{390, 1997. H. Werner. Grobstruktursimulation der turbulenten Stromung uber eine querliegende Rippe in einem Plattenkanal bei hoher Reynoldszahl. PhD thesis, Technische Universitat Munchen, 1991. H. Werner and H. Wengle. Large{Eddy Simulation of turbulent ow over and around a cube in a plane channel. In U. Schumann et al., editor, 8th Symp. on Turb. Shear Flows, 1993. Y. Zang, R.L. Street, and J.R. Kose . A dynamic mixed subgrid{scale model and its application to turbulent recirculating ows. Phys. Fluids, 5(12):3186{3196, 1993. uH j H
Vu XXX z x
Figure 1: Illustration of Schumann's approach to LES as discussed in the text. 6
G(x) GG GB  x Figure 2: Gaussian lter GG and box lter GB as de ned in the text, both plotted for the same lter width . a u PPP q u u 6 b x u u HH j H u 6  x Figure 3: Filtered functions u and u obtained from u(x) by applying a box lter, a) narrow lter, b) wide lter, log E (!) u u u log !
Figure 4: E ect of ltering on the spectrum. Here the box lter is used corresponding to Figure 3, but the curves are similar for other lters such as the Gauss Filter. u and u are illustrated by the area between the curves for u and u, and u and u, respectively. The vertical line is related to the Fourier cuto lter on the same grid.
0 0 resolved scales log E (!)
resolved turbulent stresses HHH  unresolved scales Tij ij  H H Lij =c = log ! Figure 5: Illustration of the dynamic modelling idea as discussed in the text. a b ! ! log ! Figure 6: Discetization of derivatives (Sketch) in case of a DNS (a) and an LES !eff =!, the additional lter when (b) { spectrum of u or u,    !eff numerically computing a derivative as discussed in the text, here for a second{ order central formula. The vertical axis has an arbitrary scale. In the LES case we consider u to be obtained by an ideal low pass lter for illustration. Observe that due to the logarithmic frequency scale 88 % of the discretization points correspond to the range between the maximum of !eff and ! for a second order scheme in a three dimensional computation. a b
d< ∆ c
u1 d
u1 y 1 τw τw Figure 7: Schematical pictures for the di erent approaches close to solid walls: a) resolving the near{wall structure, b) blending with a RANS model, c) application of a wall function, d) determination of wall stress by boundary layer equation solved along the wall on an imbedded grid. 4
D NS : uurm+ s D NS : vvrm+ s + D NS : wwrms + D NS : < uv> + LE S : urms + LE S : vrms + LE S : wrms + LE S : < uv> 25 3
20 15 2 10
loglaw 1
+ 5 D NS : < u> + LE S : < u> 0 10 0 10 1 y + 10 2 10 3 0 0 0.25 0.5 0.75 1 y + Figure 8: Computation without wall function using x+ = 62 z+ = 30 y1 = 1:8. Left: + + u+, right: u+ vrms wrms ;huvi+. Continuous lines are DNS data, symbols LES. rms 3
20 D NS : uurm+ s D NS : vvrm+ s + D NS : wwrms + D NS : < uv> + LE S : urms + LE S : vrms + LE S : wrms + LE S : < uv> 15 2 10 1
loglaw 5 D NS : < u> + LE S : < u> + 0 10 0 10 1 y + 10 2 10 3 0 0 0.25 0.5 0.75 1 y + Figure 9: Computation without wall function using x+ = 50 z+ = 16 y1 = 1. Labels as in Fig. 8. 3
20 D NS : uurm+ s D NS : vvrm+ s + D NS : wwrms + D NS : < uv> + LE S : urms + LE S : vrms + LE S : wrms + LE S : < uv> 15 2 10 1
loglaw 5 D NS : < u> + LE S : < u> + 0 10 0 10 1 y + 10 2 10 3 0 0 0.25 0.5 0.75 1 y + Figure 10: Wall function computation with x+ = 62 z+ = 30 y1 = 31. Labels as in Fig. 8. 25 umean u w wmean 1 20 0.75 15 0.5 10 0.25 5 0 0
0 1 2 3 loglaw D NS : < u> + 0 10 0.5 10 1 y y + 10 1.5 12 0 Figure 11: Velocities u (upper curves) and w (lower curves) at three arbitrary cuts x = const:, z = const: in the computation of Fig.10. Thin lines connect the instantaneous values, thick lines show the corresponding averages. ...
View
Full
Document
This note was uploaded on 09/16/2011 for the course ME 563 taught by Professor Staff during the Spring '11 term at Auburn University.
 Spring '11
 Staff

Click to edit the document details