Computers & Fluids 72 (2013) 30–45
Contents lists available at SciVerse ScienceDirect
Computers & Fluids j o u r n a l h o m e p a g e : w w w . e l s e v i e r . c o m / l o c a t e / c o m p fl u i d
Technical note
On the relevance of the type of contraction of the Germano identity in the new integral-based dynamic Smagorinsky model Filippo Maria Denaro ⇑ Dipartimento di Ingegneria Aerospaziale e Meccanica, Seconda Università di Napoli, via Roma 29, 81031 Aversa, Italy
a r t i c l e
i n f o
Article history: Received 20 February 2012 Received in revised form 7 December 2012 Accepted 10 December 2012 Available online 22 December 2012 Keywords: Large-eddy simulation Dynamic Smagorinsky model Finite volume methods Germano identity error
a b s t r a c t The 20-year old dynamic model is a well-established procedure in large eddy simulation that allowed us achieving several progresses in the field of numerical simulation of turbulence. The key of the procedure is in the derivation of the exact Germano identity, which is a tensor relation. After introducing the eddy viscosity model, in order to achieve a scalar relation and thereby determine a single value of the model constant, Germano et al. proposed contracting the tensor identity with the resolved strain rate tensor. Then, Lilly observed that this method can be efficient, but does raise the problem of indetermination when the resolved strain rate tensor cancels out. To remedy this problem, Lilly proposed calculating the model constant by a least-squares method, that is by minimizing the Germano identity error, observing that a negative value of the constant results to be consistent. Such a contraction is now a standard in the LES community, however, since the model function is arbitrarily extracted out from filtering, the necessity to adopt suitable averaging, clipping and some other cares to achieve numerically stable solutions are considered to be still open problems. Unlike the original Germano identity, which was developed in the framework of the differential formulation of the filtered equations, the recent integral-based filtered equations were shown to produce a new tensor identity wherein the model function is no longer acted on by filtering [8,9]. While the constant model preserves its full three-dimensional character, this new model was proved to produce stable and accurate solutions and is suitable also for complex flows. The aim of this study is to investigate the effects of three different contractions of the new integral-based tensor identity since each one of the projection represents a different minimization of the residual. The new expression of the Germano identity error is introduced to show the resulting accuracy in testing the turbulent flow in a plane channel. Unlike what is reported for the differential form of the dynamic model, the results highlight that the integral-based method is much less sensitive to the type of contraction than expected in the differential-based formulation. This fact confirms the better aspect of the integral-based Germano identity. 2012 Elsevier Ltd. All rights reserved.
1. Introduction Large Eddy Simulation (LES) was introduced more than 50-year ago [1–4] as a technique for simulating turbulent geophysical flows using the computational resources available at that time. LES results were obtained while using the eddy-viscosity concept in the static Smagorinsky model for the Sub-Grid Scale (SGS) flow fluctuations. Nevertheless, after many years of quite unsatisfactory results, in order to fit the static eddy viscosity model better to the local structure of the flow, the dynamic procedure for adapting the Smagorinsky model (DSM), by adjusting the function model at each point in space and at each time step, was introduced 20-year ago [5]. This model has represented a real cornerstone in LES, permitting us great advances in this field. ⇑ Fax: +39 081 5010264. E-mail address:
[email protected] 0045-7930/$ - see front matter 2012 Elsevier Ltd. All rights reserved. http://dx.doi.org/10.1016/j.compfluid.2012.12.002
The dynamic procedure requires using some explicit test-filtering to be applied on the primary (explicitly or implicitly) filtered variable. Because of the application of the test-filter on the filtered momentum equation, the exact Germano identity appears as a tensor relation incorporating the scale-invariance concept. The eddy viscosity approximation is introduced in the exact Germano identity, leading to a tensor relation in which the model function (i.e., the SGS viscosity) is the only unknown to be determined. In the dynamic procedure, the user must suitably prescribe as input only the ratio between the test and the primary filter lengths. As it is well known, in the DSM procedure the practical determination of the dynamic coefficient requires that the latter is arbitrarily extracted out from the test-filtering [5]. In order to overcome this mathematical inconsistence, it is advisable that the model function is made constant [1–4] in average sense, along suitable directions. However, doing this way the original goal of the local adaption of the model is somehow lost.
31
F.M. Denaro / Computers & Fluids 72 (2013) 30–45
Furthermore, as the modeled tensor Germano identity is not exactly satisfied by a single scalar function, a contraction of the tensor equation to a single scalar equation is required. This goal was initially obtained by means of a projection along the resolved strain tensor [5]. Successively, Lilly recognized [6] that such a projection had the risk of producing an inconsistence when the resolved strain tensor vanishes and therefore proposed a different projection corresponding to a minimization of the Germano identity error (GIE). At present, owing to the stable solutions, is such type of contraction to be accepted by the LES community as a standard. Nonetheless, it is worthwhile reminding that the original Germano identity was deduced in the framework of the differential form of the filtered momentum equation [5]. Owing to its multiscales foundation, the dynamic procedure appears perhaps much more general than it has been used so far and it is suitable to be applied in a more general formulation. Therefore, the new integral-based dynamic Smagorinsky model (IDSM) was successfully developed and illustrated, see Refs. [7–9]. In those papers, the main achievement was in developing the integral-based counterpart of the Germano identity in which the model function appears naturally out from filtering and no mathematical inconsistence is therefore introduced. In addition, the proper choice of the discrete filter parameters for a Finite Volume (FV)-based implicit filtering was analyzed [9]. On the other hand, even for the integral-based tensor Germano identity the contraction to a scalar equation remains necessary, requiring some suitable projection. A new expression of the GIE can be deduced in which, however, there is a priori no contribution due to the extraction of the model constant out of filtering. Therefore, as a measure of the modeling error was not tested in Refs. [7–9], the new element in this paper is the analysis of the effect of three different types of projections in terms of the GIE. Some studies on the reduction of the Germano identity error for the differential-based formulation were presented in literature, see Refs. [10,32]. Anderson and Meneveau studied the reduction of the GIE for the one and two-parameter dynamic mixed model and assessed that only the scale-similar part actively enters into the reduction of the GIE providing the 50% of reduction compared to the 3% of the DSM alone. Then, Park and Mahesh recently developed a complicated minimization of the GIE for the DSM. According to those authors ‘‘the pursuit of small GIE is at least a necessary condition for a good SGS model’’. However, the problem due to the eddy-viscosity hypotheses is present in both the DSM and the IDSM: the alignment between principal axes of the strain-rate tensor and SGS modeled one. The magnitude of the GIE is still remarkable as observed by Anderson and Meneveau even in simulations of isotropic turbulence. According to Park and Mahesh, this behavior may be attributed to the very low correlation coefficients, typically around 0.2, between true SGS stress and DSM in a priori tests. However, another possible cause for the large GIE is supposed to be the assumption of ‘‘frozen’’ velocity with respect to change in Cs. Therefore, the IDSM cannot change the intrinsic limitation of the eddy viscosity hypothesis but can avoid some mathematical inconsistence. Even if some conclusions on the reduction of the GIE were drawn by those papers, the expressions of the GIE in the DSM and IDSM are different. To the best of the authors’ knowledge, such differences, as well as the dependence of the GIE on the type of contraction, are quite unexplored and, in the authors’ opinion, that requires further testing and assessment. Therefore, this paper has the aim of illustrating some technical issues that complete the topics about the IDSM in Refs. [8,9]. The present study wants to assess the relevance on the LES solution of using three different types of projections for contracting the integral-based Germano identity to a single scalar equation. Specifically, one wonders what about
the GIE for each contraction of the new integral-based GIE expression. To this goal, for each contraction the eddy viscosity function is evaluated and used in the LES computation of turbulent channel flow at Res = 590, without any further spatial or time averaging but using its natural local form. Furthermore, the effect of possible negative values of the model function (that is the presence of some amount of energy backscatter) is tested. An investigation of the behavior of both the statistics and GIE for different discrete values of the grid-to-test filter ratio is also explored. Remembering that the Lilly contraction is strongly recommended to avoid problems for the DSM, the results obtained here for the IDSM are quite surprisingly: irrespective of the three different types of contractions, practically the same statistics are produced. Slight variations are obtained only when considering the presence of negative values of the eddy viscosity. Furthermore, also for the IDSM, the GIE profiles reproduce quite well the same behavior illustrated in Ref. [10] for the GIE in the DSM. Differences are highlighted because of diverse choice of the input test-to-grid filter ratio. Even though a reader could find excessive a full-length paper for concluding that results are insensitive to the type of contraction, it is authors’ opinion that such a result is a step ahead Ref. [8] as it assesses the successfully elimination of the commutation error in the IDSM. That fact gives the opportunity to propose the IDSM as a general dynamic procedure applicable for non-uniform filters distribution, such as that one on unstructured grids. The paper is organized as follows. Section 2 describes the integral-based filtered equations and the SGS dynamic modeling, Section 3 describes the three types of contractions and, finally, the results obtained in the solution of a turbulent channel flow are discussed in Section 4. In order to show the competitiveness of the IDSM, a sub-section illustrates the results of the LESinItaly project wherein several LES codes simulated the same channel flow problem on the same grid. 2. Dynamic SGS modeling via integral formulation In an LES approach, the scale separation is obtained via low-pass filtering [1–4]. Since FV-based LES are here considered [7–9], according to the historical [1] Schumann work dated 1975, the spatial filtering operation (no time-filtering is considered) is associated with a local average over an elementary box-volume X of linear homogeneous measure D. Thus, the volume-averaged velocity defined in the continuous form is
v ðx; t; DÞ ¼
1
D3
Z
zþD=2
zD=2
0
dz
Z
yþD=2
yD=2
dy
0
Z
xþD=2
v ðx0 ; tÞdx0 ;
ð1Þ
xD=2
the parameter D formally being a degree of freedom in the filtering procedure. If the elementary volume is homogeneous, it is possible to consider the volume average (1) in terms of the convolution product with the filter kernel 1/jXj, "x 2 X and vanishing elsewhere. The continuous filtered velocity (1) is known as the tophat filtered1 velocity. According to (1), the continuous governing equations for the are ‘‘formally’’ obtained by applying large-scale velocity field v the filtering operator to each term of the Navier–Stokes (NS) equations, commuting filtering and differential operators and supplying then some SGS modeling for the unresolved terms arising from the non-linear flux. Hence, for given boundary and initial conditions, a certain discretization in time and space allows us updating the LES solution. The term ‘‘formally’’ highlights that an explicit filtering 1 However, note that the velocity (1) still possesses an infinite number of Fourier components, its transfer function having an infinite number of zeros and vanishing asymptotically, Ref. [9]. Only the domain discretization introduces a scale separation due to the truncation of frequencies over the Nyquist grid cut-off frequency.
32
F.M. Denaro / Computers & Fluids 72 (2013) 30–45
operation is, usually, never applied in the solution algorithm but the numerical discretization acts to induce implicitly the effective filter parameters [9]. For the sake of clarity, let us start from the continuous form by applying the filtering to each term of the NS equations for incompressible flows, which results in the system
rv ¼0 @v þ r Fðv Þ ¼ 0; @t
ð2Þ
being v(x, t) the point-wise velocity field, F(v) = vv 2mrsv + Ip the total tensor flux (superscript s stands for the symmetric part having zero trace of the tensor and the gradient operator is underlined), sum of the nonlinear convective, the linear diffusive and the pressure (normalized by the homogeneous density) fluxes, respectively. In order to obtain the LES equations, the common procedure consists in commuting filtering and divergence operators in (2),
r v ¼ 0 @v Þ ¼ r FU ðv ; v Þ; þ r FR ðv @t Þ ¼ v v 2mrs v þ Ip being F R ðv
1 jXðxÞj
Z
v dx0 þ XðxÞ
1 jXðxÞj
0.2
0 0 2 0.5
1.5
1
2
2.5
3
2
2.5
3
x (2+cos(x)) (2+cos(y))/9 1
0.8
0.6
0.4
0.2
0 0
x
2 0
0.5
1.5
1
y
Z
ð4Þ Fig. 1. Exact (top) and numerical (bottom) transfer functions of the used 2D top-hat test filtering.
n F dSðx0 Þ ¼ 0;
@ XðxÞ
r v ¼ r ðv v Þ ð5Þ
This way, the commutation error is always eliminated from the filtered momentum equation, either the commutation property between filtering and divergence operator applies or not. If commu2
0.4
ð3Þ
where X(x) is an elementary volume, its frontier being denoted by oX, n being its unit vector, outward-oriented in normal direction. Owing to a proper decomposition of the flux F in the resolved Þ ; v Þ ¼ FR ðv Þ Fðv Þ ¼ ðv v vvÞ FR ðv and unresolved Tðv 2mrs v Þ, a fundamental consequence of using the inteð2mrs v gral-based formulation (4) instead of the differential (3) is that one does not demand anymore commuting filtering and differentiation. As a consequence, the unresolved flux2 T is very different from FU defined in the differential-based procedure (3), both in the non-linear part and for the presence of the residual diffusive flux. The integral-based filtered Eq. (4), governing the top-hat velocity v , are rewritten in compact form as
@v Þ ¼ r Tðv ; v Þ: þ r FR ðv @t
0.6
0
n v dSðx0 Þ ¼ 0
@ XðxÞ
@ 1 @t jXðxÞj
0.8
y
the resolved and ; vÞ ¼ v v vv the unresolved fluxes. Moreover, when the FU ðv commutation property does not apply, one can work using special commutative filters, Refs. [1–4]. Eq. (3) are appropriate to be discretized by means of either Finite Difference (FD) or Spectral Methods (SM) approaches. Nevertheless, the real shape of the filter kernel, which would be implicitly implied by the discretization, is well identified only for SM-based LES, changing significantly for FD methods of various accuracy order, Refs. [9,11]. On the other hand, numerical FV approximations, which born from the integral (weak) formulation of the conservation laws [12], are analyzed in this study. According to the integral formulation, the discretization of the divergence of the resolved flux is substituted by the discretization of the integral of the normal component of the flux [7–9]. Therefore, the symbolism in (2) actually implies that the overbar applied on the divergence operator stands for the integral of the normal component of the flux. Thus, the Eq. (2) are equivalent to write the integral equations
Z
sin(x) sin(y)/x/y
1
In the expression of T the isotropic residual pressure flux is not needed since computed within the total pressure in the elliptic equation.
¼ 0 and also the source term in the tation applies, r v ¼ r v continuity equation vanishes exactly; otherwise, depending on the accuracy order of the discretization, the second order magnitude source term can be disregarded and the filtered velocity field assumed divergence-free [7–9]. This is what used in the present paper because of the FV second order discretization. In order to close Eq. (5), the unresolved flux is modeled by introducing an ; v Þ ffi TSGS ðv Þ. The key-difference in the use of approximation Tðv either Eq. (3) or Eq. (5) is that the latter univocally imply only the volume average-based (top-hat) filtering since it is characterized by the presence of resolved and unresolved fluxes under explicit surface integrals. Conversely, the differential form (3) has not a univocally defined shape of the filtering, each discretization implying a different type of filtering [9,11]. The dynamic procedure, developed for the integral form of the equation (IDSM), has its foundation in the derivation of the integral-based counterpart of the exact Germano identity [7–9]. After e > D on the continuous applying the top-hat test-filter of width D filtered momentum Eq. (5) one has
e ge g e @v Þ ¼ r Mð v ; v Þ; þ r FM ð v @t
ð6Þ
F.M. Denaro / Computers & Fluids 72 (2013) 30–45
33
Fig. 2. LES computations of turbulent channel flow at Res = 590. Statistically averaged stream-wise velocity profiles for contractions P1, P2 and P3. Top figure: no negative SGS viscosity values allowed; Bottom figure: negative values allowed.
being
e e e e e vvÞ ð2mrs v 2mrs v Þ; ; v Þ ¼ FM ð v Þ Fðv Þ ¼ ð v v Mð v
ð7Þ
the integral-based exact sub-test scale tensor. Thus, by comparing (7) to the exact unresolved tensor T, the integral-based resolved tensor, denoted by L in analogy with the classic Leonard tensor, is obtained
e e e Þ ð2mrs v Þ: v v v 2mrs v M T L ¼ ðv
ð8Þ
Identity (8) represents the counterpart of the classical differentialbased Germano identity, extended to the integral form of the equation. It differs from the original for the presence of the diffusive flux as well as for the fact that the quadratic product is not acted upon f by the test-filter, the term v v being indeed absent in (8). This fact is relevant for determining the model coefficient in a consistent mathematical way. In fact, the eddy viscosity assumption for modeling the deviatoric part of the tensors (denoted by the subscript dev), e and Mdev ffi 2m0LES rs v , is introduced. Thus, that is Tdev ffi 2mLES rs v by assuming the same coefficient and the same rate of dissipation at primary filter and test-level (scale-invariance hypothesis), the
e 2 jrs v e j and m0LES ¼ C D j classic eddy viscosity model mLES ¼ C D2 jrs v can be introduced, resulting in a system of five independent equations [1–6] Mdev Tdev = Ldev with a single unknown
e e jrs v Þ ¼ Ldev : jrs v jrs v 2C D2 ða2 jrs v
ð9Þ
e =D, test to FV-based implicit filter In the expression (9), the ratio a ¼ D width, is introduced and it is the only input parameter the dynamic SGS model depends upon. Comparing expression (9) to the classical differential-based expression [1–4] derived from Eq. (3),
g e e e Þdev ; v fv 2 mLES v 2m0LES rs v rs v ¼ ð v
ð10Þ
the relevance of the new integral-based expression (9) appears more evident: the practical determination of the eddy viscosity mLES(x, t) does not require that the model function C(x, t) has to be arbitrarily extracted from the filtering operation. Clearly, while avoiding this mathematical inconsistency, no additional averaging on C is required, the (9) can be used to determine a real threedimensional eddy viscosity function. However, the model constant
34
F.M. Denaro / Computers & Fluids 72 (2013) 30–45
RMS of Streamwise velocity component (GNR) 4 DNS from Moser et al. No-model P1 P2 P3
3.5 3
RMS of u
2.5 2 1.5 1 0.5 0
0
0.2
0.4
0.6
0.8
1
Y RMS of Streamwise velocity component (GNR) 4 DNS from Moser et al. No-model P1 P2 P3
3.5 3
RMS of u
2.5 2 1.5 1 0.5 0 0
0.2
0.4
0.6
0.8
1
Y Fig. 3. LES computations of turbulent channel flow at Res = 590. Statistically averaged RMS of stream-wise velocity components for contractions P1, P2 and P3. Top figure: no negative SGS viscosity values allowed; Bottom figure: negative values allowed.
is still assumed to be the same at the grid and test-filter scale, such limitation causing well-known consequences, see Ref. [24] The last step, required for a practical computation of the eddy viscosity function, is the contraction of the tensor Eq. (9) to a scalar equation for determining the single unknown. 3. Three types of contractions for the modeled tensor equation In the original DSM procedure [5], after extracting the model constant out of filtering, the reduction of Eq. (10) to a scalar equation was obtained by contracting the tensor relation with the re . In the case of the IDSM relation solved strain rate tensor rs v produces now the contraction (9), the same projection along rs v
j 2mLES ¼ jrs v
Ldev : r
v
s
e e jrs v jrs v : rs v Þ jrs v : rs v ða2 jrs v
j 2mLES ¼ jrs v
ð11Þ
In principle, other projections can be chosen, for example the (9) e , producing the second contraction can be projected along rs v
se j
se
vrv
e Ldev : rs v : s e e jrs v : rs v Þ : r v jrs v
ð12Þ
In order to avoid indetermination of the eddy viscosity value when rs v vanishes, Lilly [6] proposed a modification of the original projection. The modification consists of projecting (10) along all its LHS, in such way Lilly recognized that one minimizes the square of the residual of (10). Owing to the good performances showed by this contraction, such procedure is a standard in the LES community that uses the DSM. According to Lilly modification extended to the IDSM, the projection of (9) produces the third type of contraction
j 2mLES ¼ jrs v :
ða r 2j
e e jrs v Þ jrs v jrs v Ldev : ða2 jrs v : e e e e jrs v Þ : ða2 jrs v jrs v Þ jrs v jrs v jrs v jrs v ða2 jrs v ð13Þ
Each one among (11)–(13) corresponds to an evaluation of the residual of (9), that is the integral-based expression of the GIE
35
F.M. Denaro / Computers & Fluids 72 (2013) 30–45
Averaged SGS viscosity Vertical Profile - GNR 2 Molecular kinematic viscosity P1 P2 P3
Y
1.5
1
0.5
0 0
0.002
0.004
0.006
0.008
0.01
SGS Viascosity Averaged SGS viscosity Vertical Profile - GNR 2 Molecular kinematic viscosity P1 P2 P3
Y
1.5
1
0.5
0 0
0.002
0.004
0.006
0.008
0.01
SGS Viscosity Fig. 4. LES computations of turbulent channel flow at Res = 590. Statistically averaged eddy viscosity function for contractions P1, P2 and P3. Top figure: no negative SGS viscosity values allowed; Bottom figure: negative values allowed. Black vertical line indicates the molecular viscosity value 1/590.
E ¼ 2mLES a2
e j s jrs v e rs v rv s j jr v
! Ldev ;
ð14Þ
which is clearly different from that one defined in the standard DSM [5,6,10]. In order to prevent numerical instabilities, several authors strongly advice for the DSM to use averaging of the model constant in homogeneous directions, as well as use some clipping of the eddy viscosity values. Some analyses of the GIE for the DSM were illustrated in Ref. [10], along with some proposals to reduce further the errors due to both the contraction and the mathematical inconsistence of the extraction of the model constant out of filtering. Actually, in the IDSM the latter error is naturally eliminated from the GIE (14). Hence, the goal of the present study is to assess the effective error introduced by each type of projection for the new IDSM. Therefore, a numerical experiment is designed for determining the eddy viscosity by means of the three different Eqs. (11)–(13). No supplementary averaging on the eddy viscosity functions are introduced in the computation, the local character mLES = mLES(x, t)
is therefore preserved. A series of practical numerical simulations of the turbulence in a plane channel is now illustrated. 4. Results: solutions of turbulent channel flow at Res = 590 While exploiting the IDSM, an assessment of the behavior of the three projections on a practical LES solution is illustrated for the classical problem of the turbulent flow in a plane channel at Res = 590 for which several LES solutions exist (see Refs. [10–15]) and a DNS database [16] can be used. Assuming x and z to be the directions of homogeneity of the flow and y the normal-to-wall direction, the channel height H = 2d, the stream-wise length Lx = 2pd, the span-wise length Lz = pd are chosen. The turbulent flow is therefore assumed to be periodic in x and z directions, sustained in the stream-wise direction by a driven force provided by a suitable constant pressure gradient. The corresponding non-dimensional lengths are obtained by means of the half-height d, the qffiffiffiffiffiffiffiffi non-dimensional velocity by means of us ¼ Dqp0Lxd, this way the 0
forcing non-dimensional pressure gradient is simply the unity. A
36
F.M. Denaro / Computers & Fluids 72 (2013) 30–45
Streamwise energy spectra of streamwise velocity component (GNR) 0.1
0.01
Euu
0.001
DNS from Moser et al. -5/3 slope No-model P1 P2 P3
0.0001
1e-005
1e-006 0.1
1
10
100
Kx Streamwise energy spectra of streamwise velocity component (GNR) 0.1
0.01
Euu
0.001
DNS from Moser et al. -5/3 slope No-model P1 P2 P3
0.0001
1e-005
1e-006 0.1
1
10
100
Kx Fig. 5. LES computations of turbulent channel flow at Res = 590. Statistically averaged stream-wise energy spectra along stream-wise wavenumber kx at y+ = 590 for contractions P1, P2 and P3. Top figure: no negative SGS viscosity values allowed; Bottom figure: negative values allowed.
second order accurate co-located FV-based discretization is used in the LES code, implicitly defining the type of primary filtering, see Refs. [7–9]. In order to separate velocity and pressure gradient vector fields without introducing spurious modes, a second order accurate double-projection method [17,18] is used. This procedure consists of a first step in which the intermediate non-solenoidal filtered velocity is integrated in time by means of a semi-implicit Adams-Bashforth/Crank-Nicolson scheme. Then, in a second step an approximate projection method allows computing a gradient of a scalar field that ensures that the cell-face velocities are divergence-free, last an exact projection step allows computing a gradient of a second scalar field that ensures the cell-centre velocities are divergence-free, too. This procedure has been proved to provide stable and spurious-free modes solutions [17]. A Cartesian structured computational grid, having Nx = Nz = 64 computational nodes uniformly distributed in the homogeneity plane (x, z) and corresponding to have Dx+ = 58 and Dz+ = 29 in wall
units, is used. Along the normal-to-wall direction y, the trigonometric stretching law
yj ¼ ð1 þ sin #j ÞLy =2;
#j ¼ p=2 þ ðj 1Þp=Ny
is introduced with Ny = 32 computational nodes, corresponding to have a nonfully resolved boundary layer, being the minimum cell size Dy+ = 1.42 close to walls. Such a quite coarse vertical grid is used for an assessment of the effect of the SGS model also along the non-homogeneous part of the flow. Unlike the DSM, in the IDSM the non-uniform filtering effect has no relevance in the GIE expression. This way the test appears particularly meaningful to highlight the possible different behavior of the LES solutions computed for the projections 12, 11 and 13, labeled in the following as P1, P2, P3, respectively. Furthermore, in order to isolate the effects of the SGS model, a solution obtained without any turbulence modeling (labeled no-model) is considered. Such solution can be considered
37
F.M. Denaro / Computers & Fluids 72 (2013) 30–45
Streamwise spectra of eddy viscosity function (GNR) -5/3 slope P1 P2 P3
0.1
E_nisgs
0.01
0.001
0.0001
1e-005
1e-006 1
10
Kx Streamwise spectra of eddy viscosity function (GNR) -5/3 slope P1 P2 P3
0.1
E_nisgs
0.01
0.001
0.0001
1e-005
1e-006 1
10
Kx Fig. 6. LES computations of turbulent channel flow at Res = 590. Statistically averaged spectra of eddy viscosity function along stream-wise wavenumber kx at y+ = 590 for contractions P1, P2 and P3. Top figure: no negative SGS viscosity values allowed; Bottom figure: negative values allowed.
as an under-resolved DNS. The computational time step is fixed for all cases to Dt = 5 105, the run being performed until an energy equilibrium is achieved after which the samples of the fields are memorized for computing the time-averaged statistics. As the test-filtering is concerned, a two-dimensional local averaging (top-hat filter) over the homogeneous plane with a computational stencil of uniform steps 2Dx, 2Dz is computed according to
Z zþDz Z xþDx 1 0 dz v ðx0 ; y; z0 ; tÞdx0 4DxDz zDz xDx 1 h di;j;k þ 4 v diþ1;j;k þ v di1;j;k þ v di;j;k1 þ v di;j;kþ1 16v ffi 36 i di1;j;kþ1 þ v diþ1;j;k1 þ v di1;j;k1 ; diþ1;j;kþ1 þ v þv ð15Þ
ve ðx; y; z; t; De Þ ¼
having used the fourth order accurate Simpson formula (superscript d indicate discrete grid-values). It is worthwhile reminding that the primary filter is implicitly applied by the second-order accurate FV discretization, which is supposed to be in effect as a 3D local average over a volume of
some characteristic length D (therefore, having a smooth transfer function [9]) plus a projective filter at the Nyquist frequencies, induced by the computational grid. As the test-filter (15) is concerned, the choice was to retain the smooth transfer function, which is truncated by the action of the computational grid, using a doubled computational area. According to Ref. [1], the use of the same value of the constant for the SGS model at the two filtering levels implicitly relies on two following self-similarity assumptions: (a) the two cutoff wavenumbers are located in the inertial range of the kinetic energy spectrum; and (b) the filter kernels associated to the two filtering levels are themselves self-similar. The use of a 2D test-filter in plane channel flows is in accord to what used by several authors but a 3D extension of (15) can be also supposed to be valid. More details can be found in Ref. [22] and future assessments on effects of a 3D shape of the test filtering are in progress. e is inObviously, the effective value of the test-filter length D duced [1,8,9,19–21] by the discretization (15). To clarify the discrepancy with the ideal top-hat filter, it is instructive to show
38
F.M. Denaro / Computers & Fluids 72 (2013) 30–45
Fig. 7. LES computations of turbulent channel flow at Res = 590. Statistically averaged Germano-identity error (GIE) profiles along the vertical direction for contractions P1, P2 and P3. Top figure: no negative SGS viscosity values allowed; Bottom figure: negative values allowed.
the behavior of the numerical transfer functions associated to (15). A simple Fourier analysis leads to
b d ðn; fÞ ¼ 1 ð2 þ cos nÞð2 þ cos fÞ; G 9
ð16Þ
where n = kxDx, f = kzDz, and the comparison to the exact one b fÞ ¼ sin n sin f is shown in Fig. 1. Gðn; nf Owing to the uncertain value of the primary filter length D which is only implicitly defined by the FV discretization [9,11], fixing a proper input value ad is not straightforward. In fact, consider that the only lengths exactly known are the computational mesh sizes. Defining a global mesh parameter h, the primary filter length is related to it by means of the sub-filter resolution parameter e =D ¼ ð D e =hÞ=ðD=hÞ ¼ [8,9,11] D = Qh. Therefore, since ad ¼ D e =hÞ=Q , both D e =h and the parameter Q are required for prescribðD ing a proper input value ad. An expression of the test-filter width based on the standard deviation of the corresponding discrete filter function [8,9,19,20] is
e qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi X D ¼ 12 I;K ðI2 þ K 2 ÞW IK ; h
ð17Þ
WIK being the coefficients of the discrete filtering expression, centered around the I, K = (0, 0) grid point in the homogenous plane x, z. An evaluation of (17) for the coefficients of the Simpson formula pffiffiffi e =h ¼ 2 2, and the mesh parameter h can (15) produces the value D pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi be evaluated by an Euclidean norm as h ¼ Dx2 þ Dz2 . More details about the analysis of the dependence of the LES solutions on ad can be found in Refs. [8,22]. 4.1. Testing the three types of projections In first instance, assuming a numerical value ad = 3 a series of IDSM-based simulations are performed for P1 to P3 cases. Furthermore, two clipping procedures on the eddy viscosity function are tested, the first one forcing the eddy viscosity to be exclusively positive, Clipping a), that is mLES(x, t) P 0, whilst the second admits some negative values such that the total viscosity is positive, Clipping (b), that is m + mLES(x, t) P 0. This way, the goal is to stress the three types of projections to show different behaviors and highlights potential numerical problems in the solution. The statistically averaged stream-wise velocity, stream-wise velocity fluctuations (RMS) and eddy viscosity profiles along the
39
F.M. Denaro / Computers & Fluids 72 (2013) 30–45
RMS of Streamwise velocity component (GNR) 4
DNS from Moser et al. No-model alpha = 2 alpha = 3 alpha = 4
3.5 3
RMS of u
2.5 2 1.5 1 0.5 0
0
0.2
0.4
0.6
0.8
1
Y Averaged SGS viscosity Vertical Profile - GNR 2
Molecular kinematic viscosity alpha = 2 alpha = 3 alpha = 4
Y
1.5
1
0.5
0 0
0.005
0.01
0.015
0.02
SGS Viscosity Fig. 8. LES computations of turbulent channel flow at Res = 590 for contractions P2 and ad = 2, 3, 4. Statistically averaged stream-wise velocity, RMS, Germano-identity error and eddy viscosity profiles. No negative SGS viscosity values were allowed.
vertical direction are shown in Figs. 2–4, the top figure for the Clipping (a) case, the bottom figure for the Clipping (b) case. In addition, the statistically averaged stream-wise velocity and eddy viscosity spectral distribution along the wavenumber kx at y+ = 590 are also shown in Figs. 5 and 6. The no-model and the unfiltered DNS solutions are superimposed. The results obtained with the IDSM appear quite surprising. We recall that for the DSM is commonly accepted that the Lilly contraction and the supplementary averaging of the model constant are very beneficial in the quality of the solution. Now, despite using the different contractions P1 to P3 the solutions are practically indistinguishable. Some differences appear only confronting the solutions for the two clipping procedures. Allowing the eddy viscosity to have some negative values seems producing some benefits in the solution but, in order to give a better indication, a quantitative analysis is required. According to Ref. [10], the Germano-identity error can provide a measure of the performances of the three projections. When considering the GIE (14), the IDSM-based scalar error expression
! ! se e j s s jrs v s 2 jr v j s e e r v r v a r v r v : j j jrs v jrs v ! e j jrs v a2 s rs ve rs v : Ldev þ Ldev : Ldev ð18Þ jr v j
eðx; tÞ ¼ E : E ¼ 4m2LES a2 þ 4mLES
is obtained. The profiles along the vertical direction of the statistically averaged GIE are shown, for each projection, in Fig. 7. Globally, all the curves show a peak of the errors at about y+ = 20 (the location where the velocity profiles in Fig. 2 show a transition to the logarithmic slope and the RMS in Fig. 3 show a peak) and a lesser magnitude of the GIE for the Clipping (a). The shape and the peaks location of the GIE appear in accord to those reported in Ref. [10] wherein, however, a Spectral Method (SM), which is much more accurate than the present second order FV-based, is used. As well, SM introduces the sharp cut-off, that is a different shape of the transfer function, see Refs. [9] and [11]. The location of the peaks in the GIE profiles seems to confirm that the isotropic form of the Smagorinsky model can behave poorly when approaching the wall, where different turbulence scales appear [10]. Despite the previous curves (Figs. 2–6) appeared almost coincident for all the three projections, now some differences are more visible in the GIE profiles: the P2 case (that is the projection along the resolved strain as originally proposed in Ref. [5]) appears producing the lower peak in the error both for Clipping (a) and (b). The Lilly projection P3 does not seem to produce a benefit over P2. A question appears about the use of Clipping (b) that seems to slightly improve the velocity profiles. The improvement appears in terms of a reduction of the overestimation of the unfiltered DNS velocity profile but the contemporary increasing in the RMS and GIE seems to drive us to a different conclusion. Again, as al-
40
F.M. Denaro / Computers & Fluids 72 (2013) 30–45 RMS of Streamwise velocity component (GNR)
4
DNS from Moser et al. No-model alpha = 2 alpha = 3 alpha variable along y
3.5 3
RMS of u
2.5 2 1.5 1 0.5 0 0
0.2
0.4
0.6
0.8
1
Y
2
Averaged SGS viscosity Vertical Profile - GNR Molecular kinematic viscosity alpha = 2 alpha = 3 alpha variable along y
Y
1.5
1
0.5
0 0.005
0.01
0.015
0.02
SGS Viscosity Fig. 9. LES computations of turbulent channel flow at Res = 590 for contractions P2 and ad = 2, 3 and variable depending on y. Statistically averaged stream-wise velocity, RMS, Germano-identity error and eddy viscosity profiles No negative SGS viscosity values were allowed.
ready observed in Ref. [9], it is questionable that if an LES velocity based on a smooth filter has a profile showing an overshoot of the DNS profile then the LES solution has to be considered poorly accurate. The friction velocity is filter-dependent therefore, the magnitude of the LES profile can also overestimate the DNS one. In actual fact, there are no theoretical reason to expect an exact concordance between the statistically averaged LES velocity and the logarithmic law produced by DNS, recalling Ref. [14]: ‘‘the analytical form of the law of the wall is not necessarily invariant to a filter operation since the friction velocity of filtered data can be different from the unfiltered data’’. In the authors’ opinion, this issue highlights a more general problem that is still debated in the LES community: depending on the filter shape, how much an LES solution should be expected close to a DNS one? Presumably, a good spectral-based (projective filter) LES solution should be quite similar to a DNS solution but is that necessarily true for an integral-based LES? Recently the shape of the transfer functions induced by some numerical discretization were analyzed, see Ref. [9]. The conclusions can be synthesized in this fact: increasing the accuracy of an integral-based FV scheme the shape remains smooth and close to the ideal top-hat unlike what happens for a FD scheme for which the shape tends to the spectral cut-off. The use of a smooth filter necessarily alters the spectral content at high resolved wavenumbers, this being not a merely numerical error but a response to the continuous integral-based formulation, independently from any kind of SGS
model. Ideally, owing to the smoothing of the top-hat filter, also an ‘‘exact’’ solution of the integral-based LES equations can differ from a DNS, at least if the filter width is not very small. The question that still requires a response is in what extent is acceptable a discrepancy of an integral-based LES solution from a DNS one. Such task is out of the aims of the present paper, however more details about such issue can be found in Refs. [8,9,25]. Furthermore, an assessment of the LES solutions to the filter shape can be found also in Refs. [27,28]. On the other hand, a further confirmation appears from a look to the no-model solution which is often considered in literature in order to isolate the effect of the subgrid scale model. However, this kind of solution has twofold meaning and to avoid misconceptions is here better addressed. On a side, in the framework of implicit LES (ILES – Ref. [29]), an explicit SGS model is not added as the use of asymmetric formula or the use of some limiter flux function produce a local truncation error of the scheme that can be shown to have a dissipative form and therefore mimics the same behavior of an eddy viscosity model at the length scale of the computational grid. On the other side, when the scheme is built in order to eliminate (or minimize) the numerical dissipation in the local truncation error, as happens for central formula as well as for spectral polynomial, the no-model solution should be considered as a way to assess the real contribution when the SGS model is explicitly added. For this latter case, the no-model simulation can be argued to produce inevitably unstable solutions due to
41
F.M. Denaro / Computers & Fluids 72 (2013) 30–45
RMS of Streamwise velocity component (with Models - GNR) 4 DNS from Moser et al. FD II order - Dynamic Anisotropic Smagorinsky Fluent II order FV - Dynamic Smagorinsky OpenFOAM II order FV - Dynamic Smagorinsky Present II order FV - Integral Dynamic Smagorinsky TransAT II order FV - Dynamic Smagorinsky Code_Saturne II order FV - Dynamic Smagorinsk IV order FD/Spectral - Smagorinsky
3.5 3
RMS of u
2.5 2 1.5 1 0.5 0 0
0.2
0.4
0.6
0.8
1
Y
Streamwise energy spectra of streamwise velocity component (with Models - GNR) 0.1
Euu
0.01 DNS from Moser et al. -5/3 slope II order FD - Dynamic Anisotropic Smagorisnky Fluent II order FV - Dynamic Smagorinsky OpenFOAM II order FV - Dynamic Smagorinsky II order FV - Integral-based Dynamic Smagorinsky TransAT II order FV - Dynamic Smagorinsky Code_Saturne II order FV - Dynamic Smagorinsky IV order FD/Spectral Smagorinsky
0.001
0.0001
1e-005
1e-006
0.1
1
10
100
Kx Fig. 10. Project LESinltaly: Comparison of LES codes in the solution of turbulent channel flow at Res = 590 for an unresolved boundary layer. From top to bottom: Statistically averaged stream-wise velocity, RMS and spectra profiles.
42
F.M. Denaro / Computers & Fluids 72 (2013) 30–45
Fig. 11. Project LESinltaly: Comparison of LES codes in the no-modeled solution of turbulent channel flow at Res = 590 for resolved (bottom-GR) and unresolved (top-GNR) boundary layer. Statistically averaged stream-wise velocity profiles.
the energy pile-up near the cut-off frequency. However, that not always happens but depends on the shape of the transfer function that is implicitly induced by the discretization. Indeed, spectral methods suffer of instability since the cut-off filter does not smooth the highest resolved. Often high order FD discretizations share the same problem. Actually, when using FV method, the numerical transfer function is always smooth Ref. [9] and, even if the scheme has not numerical dissipation (as the scheme adopted here, see Ref. [30]), the highest resolved frequencies do not pile-up their energy content and the solution remains stable. According to the above framework, owing to the lack in the SGS contribution, the no-model solution has to be assumed wrong. Nevertheless, the velocity profile is much closer to the DNS one than the LES ones. Therefore, adding the eddy viscosity-based contribution to the no-model case, the overshoot in the logarithmic region is somehow a fact to be expected. On the other hand, the spectral distribution of the eddy viscosity shows (Fig. 6) that the magnitude is prevalent at low wavenumbers and, as a result, the low wavenumbers of the velocity spectra (Fig. 5) are better in accord for the LES solutions in case P3.
As a further comment, increasing the grid resolution in vertical direction (for the sake of brevity the comparisons are not reported) does not change the fact that projection P1 to P3 produce similar results. However, some results on the fine grid are presented in Section 4.3. 4.2. Analysis of the GIE for locally variable test-to-grid filter ratio The use of a primary implicit filtering joined to the use of an explicit test-filtering performed only in the 2D homogeneous plane (having therefore a constant characteristic width), has some consequences on the input value ad that is worth to be addressed. Some analysis of the sensitivity to the parameter ad was described in Refs. [8,22]. Starting from those conclusions, the attention is here focused on a new aspect not analyzed there, which is the fact that such parameter depends locally on the non-uniform vertical grid. In principle, along the non-homogeneous direction two opposite situations can result when the computational grid is stretched. The first situation is that the vertically stretched grid is everywhere so fine to resolve all the scales (DNS-like): since the implicit filtering
F.M. Denaro / Computers & Fluids 72 (2013) 30–45
43
Fig. 12. Turbulent channel flow at Res = 590: statistically averaged integral-based Germano identity error profiles for the integral-based dynamic eddy viscosity model (IDSM), integral-based dynamic scale similar model (ISSM) and integral-based dynamic mixed model (TDMM) simulations.
length D can be assumed non-depending on y, the assumption e =D ¼ constant is justified. In a second situation (that is also ad ¼ D the present one), the vertical grid is stretched to improve the near-wall resolution but is unable to resolve all the scales, thus implicit non-uniform filtering effects are present such that D = D(y). In e =DðyÞ should be propthe latter case, a local dependence ad ðyÞ ¼ D erly considered in the dynamic procedure. Unlike DSM (10), this dependence does not introduce any mathematical inconsistence in the IDSM since in (9), ad appears naturally out of filtering. The law of variation ad(y) has to accord to the grid-stretching law. In fact, it is plausible to assume3 Dðyj Þ ¼ ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi r Dx2 þ Dz2 þ Dy2j where h(yj) represents a characterisQhðyj Þ ¼ Q tic grid length at position yj. On the other hand, considering the constant test-filtering length according to (17), one has pffiffiffi pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi e ¼ 2h 2 ¼ 2 2ðDx2 þ Dz2 Þ. Therefore, the following law appears D well suited:
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2ðDx2 þ Dz2 Þ r ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi Dx2 þ Dz2 þ Dy2j Q pffiffiffi sffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 2 2 1 : ¼ Q 1 þ Dy2j =ðDx2 þ Dz2 Þ
e D ad ðyj Þ ¼ ¼ Dðyj Þ
ð19Þ
e and Dyj approaching to zero close to a wall, For a quite large value D pffiffi one gets ad jwall ! 2Q 2 and the value decreases for y approaching the half-height channel. Hence, in order to explore the possibilities that both Q – 1 and ad depends on y, three more simulations were performed for the projection P2: two fixed values ad = 2 and 4 and ad variable along y according to (19), starting from the value 3 at walls. No negative SGS viscosity values were allowed, that is Clipping a) is used. The results are shown in Fig. 8 for the fixed values ad = 2, 3, 4 and Fig. 9 (fixed and variable values) for the statistically averaged stream-wise velocity, RMS, GIE and eddy viscosity profiles. From the curves, a someway conflicting result appears: the lower is the peak in the GIE (ad = 2), the higher is the overestimation in 3 As a matter of fact, a unique sub-filter resolution parameter Q is rigorously correct only for an homogeneous computational grid while a Qi set should be used for variable grid.
the velocity profile. The eddy viscosity profile at ad = 2 is also largely overestimated those obtained at ad = 3,4 and variable law. From the inspection of the GIE profile, it can be supposed that close to wall and in the incipient transitional region (y+ < 30) the value ad = 2 is better suited but in the logarithmic region a higher value such as ad = 4 is better since it makes the profile much closer to the DNS. Actually, this fact would be in contrast to the law supposed in (19) which provides the lowest value at half-channel. A response in this behavior can be in the fact that the GIE for the Smagorinsky model appeared large [11,21] due to a very low correlation between true and modeled SGS stresses. Hence, the reduction in the peak of the GIE profile at ad = 2 is associated to a ‘‘laminarization’’ of the flow due to the excess of the eddy viscosity, showing that the magnitude of the GIE is not, by itself, a measure of the quality of the LES solution. It is worthwhile reminding that in the IDSM there is no contribution to the GIE due to the arbitrary extraction of the model function out of filtering. On the other hand, deducing that an arbitrary high value of ad produces a solution according better to DNS is not correct, the eddy viscosity tends to vanish for high values (and the GIE magnitude increases) and the solution simply tends to the no-model one. Again, if one adds an eddy viscosity-type contribution for the unresolved scales, some amount of overestimation in the velocity profile seems a consequence of the LES approach based on implicit smooth filtering. 4.3. A comparison of the IDSM-based solution with other LES codes: the project LESinItaly The results previously illustrated can appear poor to a reader if compared to the unfiltered DNS data. Unfortunately, for this test, the full database of the DNS velocity fields for performing a proper post-filtering and a statistical averaging (as done in Ref. [8]) is not available. Even the fact that the no-model solution appears better in accord to DNS can be somehow surprising. However, it is worthwhile remarking that model-free solutions are notoriously surprisingly quite good in some conditions, e.g., see Refs. [13,25]. Furthermore, the computational vertical grid has a quite lower number of cells that cannot be compared to grids used by others authors in the literature for this Reynolds number. Nevertheless, the IDSM-based results can be now compared further; indeed, during the last year the writing author coordinated
44
F.M. Denaro / Computers & Fluids 72 (2013) 30–45
the project LESinItaly, involving several research groups, which used different computational LES codes with the aim of assessing a common running protocol and creating a comparative data-base, see Ref. [26]. More specifically, commercial (Ansys Fluent, TransAT), open-source (OpenFOAM, Code_Saturne) and own-made academic LES codes were tested by several groups in the simulation of the previous channel flow test-case, all using the same computational grids. Specifically, all the codes run on two grids, the first is exactly the same previously presented (does not resolve the boundary layer) whilst the second grid has a finer vertical grid (same horizontal resolution) built of Ny = 100 non-uniformly distributed computational nodes, allowing the boundary layer to be resolved. This way, the dependence of the results on the vertical filter size was assessed. The codes, based on FD, FV and pseudo-spectral methodologies, used the standard Smagorinsky/DSM model and the comparison of the results can further clarify the behavior of the present IDSMbased solution. The Fig. 10 shows the comparison for the unresolved grid in terms of stream-wise velocity, RMS and spectra, the IDSM being reported for the projection P2 and ad = 3. The comparison shows that the FV-based IDSM appears competitive with the classic Smagorinsky-based pseudo-spectral solution as well as performs better than other FV codes. Mainly it appears for all codes that the implicit smooth filter of second order discretization governs the spectral resolution in a similar way, only the pseudospectral solution shows a clear difference (mainly in the spectra) due to the sharp-cut off filter. All codes overestimate the DSN velocity profile, while being (quite surprisingly) stable and in good accord to the DNS for the no-model (here not reported) solutions [26]. Again, the typical behavior of adding an eddy viscosity model seems confirmed: a shift in the logarithmic-law magnitude owing to a reduced value of the stress velocity. A question arise about some compensation of errors that on a single grid can somehow drive to results better than they really are. Therefore, as reinforce to the comparison in Fig. 10 the solution on the grid resolving the boundary layer can clarify the picture. In Fig. 11 the no-model solutions obtained by all the codes on the coarse (upper) and refined (lower) grids are shown. For the sake of brevity, the solutions are limited to the statistically averaged velocity profiles. Fig. 11 demonstrates that Fluent, OpenFOAM, Code_Saturne and TransAT have a basic overestimation of the DNS profile assessing an intrinsic dissipative character of the discretization that is more relevant than the expected implicit filtering in effect for each scheme [9]. Thus, an eddy viscosity SGS model can only make worst such discrepancy. The other codes showed a different behavior of the no-model solutions. The present FV-based method obtains both resolved and unresolved solutions that are very similar, quite insensitive to the refinement. That confirms the potential beneficial effect of a good SGS model.
5. Conclusions The present study had the goal of continuing the analyses of the issues on the new integral-based dynamic SGS modeling recently presented [7–9,22,23]. As in the LES community is accepted that the Lilly contraction is necessary to ensure stable results for the DSM, the relevance of three types of contractions of the new expression of the modeled Germano tensor identity is studied. Furthermore, the effect of negative values of the eddy viscosity as well as the dependence on the discretization of the test-filtering is considered. The consequence of using a local value of the lengths filter ratio is explored. According to Ref. [10], the expression of the Germano identity error is used to assess the relevance of these parameters in the LES solutions. The IDSM has a key-advantage over the DSM in the absence of the mathematical inconsistence caused by
the arbitrary extraction of the model constant, indeed the latter appears naturally extracted out of filtering and, therefore, there is no need to further averaging the model function. However, the model constant is still assumed to be the same at the grid and test-filter scale and such limitation should be considered. The capability of using a real local eddy viscosity function is assessed in the simulation of the turbulence in a plane channel flow. The quality of the solution is shown to be governed mainly by the FV-based discretization inducing a smooth implicit filtering (no-model solution), the type of contraction of the Germano identity seems to have minor relevance. The dynamic eddy viscosity model shows that the DNS velocity profile is systematically overestimated in logarithmic region and, since the error due to the arbitrary extraction of the model function is eliminated in the IDSM, this behavior can be attributed to the isotropic eddy viscosity assumption. Therefore, in the future it appears suitable to extend the new integral-based dynamic procedure to more sophisticated SGS models such as the dynamic mixed model with an anisotropic assumption. It is worthwhile remarking that the IDSM is easily extendible to complex flow configurations in which a FV discretization over hybrid structured/non-structured grids is used. During the long revision process, the integral-based dynamic procedure was further tested and extended to the one- and twoparameters dynamic mixed model (IDMM) [31]. This new model is currently still under study and the preliminary results on the channel flow seem to indicate a significant improvement of the LES solutions and a drastic reduction of the GIE. The results are very promising as shown in Fig. 12 by the dramatic reduction of the GIE reduced. In the figure the label ISSM stands for the scalesimilar model alone (no eddy viscosity part), and C = 1 stands for the one-parameter mixed model. It is confirmed the improvement of the solution, therefore a full analysis will be presented in future. Acknowledgments The author acknowledges the CASPUR centre – Italy, to have supported this study providing the computational resources of Matrix. The author is also indebted to people of the LESinItaly group for the creation of the database. References [1] Sagaut P. Large eddy simulation for incompressible flows. An introduction. Springer; 2006. [2] Berselli L, Iliescu T, Layton W. Mathematics of large eddy simulation of turbulent flows. Springer; 2005. [3] Lesieur M, Métais O, Comte P. Large-eddy simulations of turbulence. Cambridge; 2005. [4] Geurts BJ. Elements of direct and large-eddy simulations. Edwards; 2004. [5] Germano M, Piomelli U, Moin P, Cabot W. A dynamic subgrid-scale eddy viscosity model. Phys Fluids A 1991;3:1760. [6] Lilly DK. A proposed modification of the Germano sub grid-scale closure method. Phys Fluids A 1992;4(4). [7] Denaro FM, De Stefano G, Iudicone D, Botte V. A finite volume dynamic largeeddy simulation method for buoyancy driven turbulent geophysical flows. Ocean Modell 2007;17(3):199–218. [8] Denaro FM, De Stefano G. A new development of the dynamic procedure in large-eddy simulation based on a finite volume integral approach. Appl Stratified Turbul Theor Comput Fluid Dyn 2011;25(5):315–55. [9] Denaro FM. What does finite volume-based implicit filtering really resolve in large-eddy simulations? J Comput Phys 2011;230(10). [10] Park N, Mahesh K. Reduction of the Germano-identity error in the dynamic Samgorinsky model. Phys Fluids 2009;21:065106. [11] Geurts BJ, van der Bos F. Numerically induced high-pass dynamics in largeeddy simulation. Phys Fluids 2005;17:125103. [12] LeVeque RJ. Finite volume methods for hyperbolic problems. Cambridge Press; 2002. [13] Sarghini F, Piomelli U, Balaras E. Scale-similar model for large-eddy simulations. Phys Fluids 1999;11(6):1596–607. [14] Stolz S, Adams NA, Kleiser L. An approximate deconvolution model for largeeddy simulation with application to incompressible wall-bounded flows. Phys Fluids 2001;13(4).
F.M. Denaro / Computers & Fluids 72 (2013) 30–45 [15] Meyers J, Sagaut P. Is plane-channel flow a friendly case for the testing of large-eddy simulation subgrid-scale models? Phys Fluids 2007;19:048105. [16] Moser RD, Kim J, Mansour NN. Direct numerical simulation of turbulent channel flow up to Res = 590. Phys Fluids 1999;11(4). [17] Denaro FM. A 3D second-order accurate projection-based finite volume code on non-staggered, non-uniform structured grids with continuity preserving properties: application to buoyancy-driven flows. Int J Numer Methods Fluids 2006;52(4):393–432. [18] Denaro FM. Time-accurate intermediate boundary conditions for large eddy simulations based on projection methods. Int J Numer Methods Fluids 2005;48. [19] Lund TS, Kaltembach H-J. Experiments with explicit filters for LES using a finite-difference method. Center for turbulence research. 1995: Ann. Research Briefs, Stanford University; 1995. pp. 91–105. [20] Lund TS. On the use of discrete filters for large eddy simulation. Center for turbulence research. Ann. Research Briefs, Stanford University; 1997. pp. 83–95. [21] Anderson R, Menevau C. Effects of the similarity model in finite-difference LES of isotropic turbulence using Lagrangian dynamic mixed model. Flow Turbul Combus 1999;62:201. [22] Denaro FM. On the relevance of discrete test-filtering in the integral-based dynamic modelling, on ERCOFTAC series. Direct Large-Eddy Simul VIII 2011;15(1):27–32. [23] Denaro FM, De Stefano G. A new development of the dynamic procedure for the integral-based implicit filtering in large-eddy simulation, on ERCOFTAC series. Qual Reliab Large-Eddy Simul II 2009;16(Part 1):79–90.
45
[24] Meneveau C, Lund TS. The dynamic Smagorinsky model and scale-dependent coefficients in the viscous range of turbulence. Phys Fluids 1997;9(12):3932–4. [25] Gullbrand J, Chow FK. The effect of numerical errors and turbulence models in large-eddy simulations of channel flow, with and without explicit filtering. J Fluid Mech 2003;495:323–41. [26] Denaro FM, Abbà A, Germano M, Icardi M, Marchisio D, Rolfo S, et al. LESinItaly group: a comparative test for assessing the performances of large-eddy simulation codes. In: Proc. of AIMETA conference, Bologna; September 2011. [27] De Stefano G, Vasyliev O. Sharp cut-off versus smooth filtering in large eddy simulation. Phys Fluids 2002;14(1):362–9. [28] P Lafon J Berland, Daude F, Crouzet F, Bogey C, Bailly C. Filter shape dependence and effective scale separation in large eddy simulations based on relaxation filtering. Comput Fluids 2011;47:65–74. [29] Grinstein FF, Margolin LG, Rider WJ. Implicit large eddy simulation: computing turbulent fluid dynamics. Cambridge University Press; 2007. [30] Aprovitola A, Denaro FM. A non-diffusive divergence-free, finite volume-based double projection method on non-staggered grids. Int J Numer Methods Fluids 2007;53(7):1127–72. [31] Denaro FM. A new integral-based two-coefficients mixed model for LES. In: European fluid mechanics conference, Roma; 10–14 September, 2012. [32] Anderson R, Meneveau C. Effects of the similarity model in finite-difference LES of isotropic turbulence using a Lagrangian dynamic mixed model. Flow Turbul Combust 1999;62(3):201–25.