j /r Center for TurbulenceResearch Annual Research Briefs 199f2 A normal stress subgrid-scale eddy viscosity model in large eddy simulation By K. Horiuti, 1 N. N. Mansour 2 AND J. Kim 2 1. Motivation and objectives The---'Smagorinsky subgrid-scale eddy viscosity model (SGS-EVM) is commonly used in large eddy simulations (LES) to represent the effects of the unresolved scales on the resolved scales. This model is known to be limited because its constant must be optimized in different flows, and it must be modified with a damping function to account for near-wall effects. The recent dynamic model (Germano et al. 1991) is designed to overcome these limitations but is compositionally intensive as compared to the traditional SGS-EVM. In a recent study using direct numerical simulation data, Horiuti (1993) has shown that these drawbacks are due mainly to the use of an improper velocity scale in the SGS-EVM. He also proposed the use of the subgrid-scale normal stress as a new velocity scale that was inspired by a high-order anisotropic representation model (Horiuti 1990). The testing of Horiuti (1993), however, was conducted using DNS data from a low Reynolds number channel flow simulation. It was felt that further testing at higher Reynolds numbers and also using different flows (other than wall-bounded shear flows) were necessary steps needed to establish the validity of the new model. This is the primary motivation of the present study. The objective is to test the new model using DNS databases of high Reynolds number channel and fully developed turbulent mixing layer flows. The use of both channel (wall-bounded) and mixing layer flows is important for the development of accurate LES models because these two flows encompass many characteristic features of complex turbulent flows. 2. Accomplishments The subgrid-scale stress tensor, rij, that results from filtering the Navier-Stokes equations consists of three terms (Bardina 1983): rij = Lij + Cij + Rij, (1) Lij = uiuj -- uiuj, Cij : -ffiu_ h- ui_u--j, Rij : ui_u j where overlineui denotes the filtered velocity component and uIi = ui - _i denotes the SGS component of ui. Lij is the Leonard term, Cij is the cross term, and Ri3 is the SGS Reynolds stress. The indices i = 1,2,3 correspond to the directions 1 Institute of Industrial Science, University of Tokyo 2 NASA Ames Research Center 62 K. Horiuti, N. N. Mansour gJJ. Kim x, Y, and z, respectively, with x the streamwise (ul = u), Y the "major-gradient" (wall-normal or cross-stream) (u2 = v), and z the spanwise (us = w) directions. The Leonard term in eq. (1) is not modeled but is treated explicitly by applying the filter, while the other two terms (Cij and Rij) need to be modeled. A successful model for the cross term is a model suggested by Bardina (1983) where _w w-- Cij = u_i_j + _iutj This model has been tested by Bardina (1983) for homogeneous flows and by Horiuti (1989) for the channel flow and was found to be a good model for the cross terms. This model will not be tested further in this work. For the Rij terms, the eddy viscosity model by Smagorinsky (Smagorinsky 1963): 2- ve(0 , 0 i. R,i ~ _Ec_i - 0xj + 3-_-_)' (2) v. = (csa)2[ si_s,j] 1/_, s_, = ox--_+ ox---7' and the Bardina model nit _ C(:_i - _i)(-_j - _j). (3) are two of several models which are used in LES computations. In these models, Cs and C are model constants, and Ec = u_u_/2 and ve are, respectively, the SGS turbulent kinetic energy and SGS eddy viscosity coefficient. A is the characteristic SGS length scale whose value is defined as (AxAyAz)l/s; Ax, Ay, and Az are the grid intervals in the x, y, and z directions, respectively. The Smagorinski model is a "Prandtl-type" mixing length model that can be derived by starting with the eddy viscosity approximation to the subgrid-scale Reynolds stresses and assuming production and dissipation are in balance. In an eddy viscosity approximation, v, is written as the product of a characteristic time scale r and a velocity scale E I/2, v,=C_rE (4) where C_ is a model constant, r is then expressed as (Horiuti 1993) -EG Ou_Ou_ -#3/2_a r = _ e (5) e ' =UOzlOx_=C_ _ , where e is the dissipation rate of Ea and C, is a model constant. The Smagorinsky model assumes that E = Ea in (4). In the present study, we make use of the direct numerical simulation flow fields available at CTR to directly test the various approximations. The fields we con- sider are homogeneous in two-directions. To compute the large-eddy flow fields, we filter the DNS fields by applying a two-dimensional Gaussian filter in the i = 1,3 directions. In the inhomogeneous direction (i = 2), a top-hat filter is applied to the SGSnormal stress model 63 0.8 0,6' 0.4, ¢,tt- O.2 A 0, _"_ -02, V -111.4 Model _.0 DNS -11.8 o:1 o:2 0:3 0:4 o:e o:o 1.o FIGURE 1. y-distribution of the SGS-Reynolds shear stress in turbulent channel at Re,. = 790. (Model with E = EG/u_.) 0."14 0.12 0.10 0.08 tN_l,- 0.06 0.04 A 0.02 I:::0.° -0.04 V -°.°6 -0,08 Model -0.10 .... DNS 41.12 -0.14 o o:1 0:2 o:3 0:4 0:5 o:s 0:7 o:e o:o 1.o y/26 FIGURE 2. y-distribution of the SGS-Reynolds shear stress in turbulent channel at Re,. = 790. (Model with E = u2, u,2/u2,.) channel flow fields. No filter was applied in this direction (i = 2) to the mixing layer flow field. This is due to the fact that occasionally the doubly filtered (_) grid-scale variables were larger than the singly filtered ones (i]), owing to the inaccuracy of a top-hat filter in regions where grid spacing is coarse. The DNS databases we used were the fully developed incompressible channel flows at Re_ (Reynolds number based on the wall-friction velocity, u,., and the channel height, 26)= 360 (Kim et aL 1987) and 790 (Kim 1990), and the incompressible mixing layer at Reo (the Reynolds number based on the momentum thickness, 6m, 64 K. Horiuti, N. N. Mansour _J J. Kim and the velocity difference, AU)= 2400 (Moser and Rogers 1992). We started with the low Reynolds number channel flow data as a confidence test. We found that the results obtained in this case are consistent with the previous work of Horiuti (1993), who used a different set of DNS data but at the same Reynolds number. The details of this testing are not shown in the present report. The high Reynolds number channel flow field (with 256 x 193 ×192 grid points) was filtered to 64 ×97 ×48 grid points. The mixing layer flow field (with 512 x 210 × 192 grid points) was filtered to 64 x 210 x 48 grid points. These LES grid point numbers were chosen so that the turbulent kinetic energy retained in the SGS components is large. This is needed to make a fair assessment of the SGS models. SGS model evaluations were conducted by comparing the y-distribution of the mean values averaged in the x - z plane (denoted by < 0 >), and also by comparing the y-distribution of the root-mean-square (rms) values of the exact terms with the model predictions. Only the y-distribution of the mean values are shown in the present report because the rms values were found to give similar results. 2.1 A proper eddy viscosity velocity scale 2.1.1 Channel flow The y-distribution of the SGS Reynolds shear stress < u_u_ > obtained with E = Ec and Cv = 0.1 in (4) is compared with the DNS data in Fig. 1. While the agreement between the model and the term is good in the central portion of the channel, the agreement deteriorates near the wall where the model predicts a very large peak compared to the actual data. This overprediction of the shear stress near the wall when Ec is used for E in (4) implies that a damping function is needed to account for the presence of the wall. This near-wall overprediction of the stress is similar to the near-wall behavior of one-point closure models (see Rodi & Mansour 1991). This behavior of one-point closure models is attributed to the rapid reduction of the Reynolds shear stress (as the wall is approached) due to the preferential damping of the normal stress (Launder 1987, and Durbin 1992). Horiuti (1993) reasoned that the same wall damping effects should hold true for the SGS field. Indeed, when the SGS normal stress u2tu t2 is used for E (with C_ = 0.23, see Fig. 2), the model agrees well with SGS Reynolds shear stress near the wall without an additional damping function. The model is, however, less effective as compared to using the total energy in the core region of the channel. The main deficiency in the core region is attributed to excessive grid stretching in the y-direction because of the mapping used in conjunction with Chebyshev expansions. In an actual LES computation, finite differences with a more uniform grid are used in the y direction and, therefore, a more isotropic energy distribution can be expected in this case. The effects of the anisotropic grid can be evidenced by the y-distribution profile of the 'flatness parameter' A (Lumley 1978) averaged in the x-z plane. In this case A is defined as A=[1 -9 (6) _{A2 - A3}], .42 = aijaij, A3 = ai.iajkaki, SGS normal stress model 65 0.5 ¸ 0.4, A 0.3 ¸ V 02, 0.1 Van Driest o:, o:2 o:_ o:4 o.s FIGURE 3. y-distribution of the flatness parameter A and the Van Driest function (channel at Rer = 790). 1 1 aij = {u_u_ - = ! I _UkUk We find (see Fig. 3) that in the core region of the channel, A -_ 0.35, which is much smaller than the expected A = 1 when the small scale turbulence is isotropic. In the region around y ,,_0.1, A peaks around A -,, 0.5, and then gradually decreases to 0.35 at the channel center. The y-distribution of A for the unfiltered DNS data does not show this overshoot and is close to A = 1 around the centerline. The grid spacing in the central region of channel seems to be too coarse; therefore, a considerable anisotropy exists in the SGS turbulence fluctuations. In fact, when the SGS-EVM model with E = u22u' ' in (4) was used in an actual LES channel flow calculations using a more uniform grid at high Reynolds number (Rer = 1280) (Horiuti 1993), a good agreement with experimental data was found. The present comparisons for the high Re channel flow confirm the conclusions of Horiuti (1993) based on the low Re channel flow fields. For the record, the y-distribution of the conventional Van Driest damping func- tion ((1 - exp(-y+/26.0)) (normalized with value of A at the channel center) is included in Fig. 3. It should be noted that the 'flatness parameter' A has a similar distribution across the channel as the Van Driest function, suggesting that A may be used as an alternative method to damp the eddy viscosity near the wall (Horiuti 1992). 2.1.2 Mizing layer The y-distribution of ux2u' ' obtained using E = EG (C_ = 0.20) and E = u2u2r, (C_ = 0.26) in (4) are compared with the DNS data in Fig. 4 and 5, respectively. Both cases show a good agreement of the model with the DNS data, indicating that the two models are equivalent in this case. It should be noted that the optimized C_ values obtained for the u2'u '2model in the clmnnel flow at lower Re (0.22), at high Re 66 K. Horiuti, N. N. Mansour _ J. Kim 0 ¸ `0._8 ¸ `0.OOLO A a ,_lt :: -0.0015 v Model -0.0020 t .... DNS i1 le -0.00_ -20 -15 .;o ._ o _ _o _'s 2o FIGURE 4. y-distribution of the SGS-Reynolds shear stress in turbulent mixing layer. (Model with E = Ea//(AU)2.) 0.0005 0 J <1 A -0.0010 i J_r° i :_ `0.o¢1s v ,, :,, Model `0.0020 DNS `0.0025 -- 30 -'6 _ 's _'o ;s -20 -1'5 20 Y/ m FIGURE 5. y-distribution of the SGS-Reynolds shear stress in turbulent mixing layer. (model with E = u_u_/(AV) 2 ) (0.23), and in the mixing layer (0.26) are very close to each other. This implies that the model constant of the SGS normal stress model is rather universal independent of the type of flow, whereas the optimized Cv values for the E = Ea model were 0.11, 0.10, 0.20, respectively. This is further indication of the potential strength of the normal-stress model. It is interesting to note that, particularly in the outer edge region of the mixing layer (y -_ 10 or ,,_-10), the magnitude of normal component was the largest among the three components of the SGS turbulence fluctuations. A possible relationship of this phenomenon with the significant intermittency in these regions will be investigated in future work. SGS normal stress model 67 0.14 ' 0.12 ' 0.10 , ]I 0°.06o, :: ..." 0.02, ..o''" /_c'4 0 ..- ..... - ........ "° _: -0.02 ' . -0.04 k ,,'J ::1°_1/[ -Model ...... DNs 0 0.1 0.2 0.3 0.4 0.5 0.6 0.7 0.8 0.0 1.0 y12_ FIGURE 6. y-distribution of the SGS-Reynolds shear stress using Generalized nor- mal stress model (Channel at Rer = 790) 0 s <1 i *_s A -0.0010 :_ -0.001S V Model DNS iJ -0.0025 -2o .is .;o :s _) "s ;o 1'5 y/6_ FIGURE 7. y-distribution of the SGS-Reynolds shear stress using Generalized nor- mal stress model (Mixing layer) P.._,Generalized SGS normal-stress model Although it is shown that the SGS-EVM based on the 'major-gradient stress' model (Eq. (4) with E = u'2u'2) shows a high correlation with the DNS data, this model does not preserve the tensoral invariance of the SGS Reynolds stresses. This drawback can be circumvented by generalizing this model as follows (Durbin 1991 and Horiuti 1993): (7) " 'ui1u'= 6,j( 3EG + _P) - v,,, --cgx,- v,j, 68 K. Horiuti, N. N. Mansour g4J. Kim where A ,, 0-_,,_ v_,i = Cv-_-i_uiuj, P = v_.._ Oxt " EG This model is not identical to the conventional eddy viscosity model. In this model, the velocity scale is chosen as the component normal to the principal shear plane, for example, for ' ' the most dominant term in the RHS becomes nu_O-ffa/x2. Ul_t 2, The generalized normal stress model is tested in the same manner as in the last section. The results for the high-Re channel flow (C_ = 0.23) and the mixing layer (C_ = 0.26) are shown in Figs.6 and 7, respectively. It is evident that the generalized normal stress model shows a high correlation with the previous normal stress model. It should be pointed out that the generalization of the normal stress model with- out losing the tensorial invariance is not unique, e.g. terms such as O-a_ 0-at ax--S- v j' (8) can be added to Eq. (7). When these terms are included, however, the term u_ x O-o.1/Ox2 causes a large peak in the model for the U_lU2t profile in the channel flow. For this reason, the terms in Eq. (8) were excluded when generalizing the normal stress model. _.3 Approzima_ion method of the SGS _urbulent energy To effectively use the model advocated in the previous sections, a model for the normal stresses is needed. We can either carry equations for the normal stresses or estimate the energy in the subgrid-scales from the energy in the large scales. In testing the scale-similarity model of Bardina (1983), Horiuti (1993) found good correlation between the model and the data. The model reads, Ec = CK( t - - (9) 2Iu2! = CN( 2 - - where a constant different from unity was needed. It was pointed out that the optimized model constants Ch" and CN were not equal to unity because the scale- similarity model provides a partial estimate of the whole SGS fluctuations which resides in the vicinity of the cutoff-wave number (= r/A). The poor performance of the model when these coefficients are set equal to unity can be evidenced by the fact that in this case the SGS flatness parameter A becomes identically zero (purely two-dimensional state). We have optimized CI,'/CN for the low-Re, the high-Re channel, and the mixing layer flows and found 7.0/12.0, 7.0/9.0, and 9.0/ 12.0 to be representative values for these flows. We note that they are slightly (but tolerably) sensitive to the type of flow fiel, and that they are generally close to each other. A representative comparison of the model prediction with the DNS data for u_u_ is shown for the high-Re channel flow and the mixing layer in Figs. 8 and 9, respectively. We find a good agreement with the DNS data for both flows. SGS normal stress model 69 1.2 1.0 gAt Model ¢'_J" 0.8 A 0.6 v O.2 0 0 011 0'.2 013 01, o'.s ore o:7 0:e o'.o 1.0 y/26 FIGURE 8. y-distribution of ul'u]/'u 2r using the Bardina model (Channel at Rer = 790) 0.012 Model 0.010 DNS A :: 0.004 v 0.002 0 -20 -I5 -10 -5 0 5 10 1'5 20 yl m FIGURE 9. y-distribution of u'_u_/(AU) 2 using the Bardina model (Mixing layer) 3. Future plans An ultimate goal of the present study is to develop an SGS model which yields good predictions of turbulent flows in a complex geometry. Of particular interest is the flow over a backward-facing step. In this case, while the flow is bounded by the walls, the internal mixing layer present in this flow plays a major role in setting the turbulence levels. In this work, a proper velocity scale for the SGS-EVM viscosity was determined for the fully developed channel and the mixing layer flows. In the channel, a clear advantage over more conventional treatments was shown by using the normal stress. It was also shown that the SGS normal stress is equally 70 K. ttoriuti, iV. iV. Mansour _ J. Kim useful as the total SGS turbulent energy for modeling in mixing layers. The model constant in the normal-stress model was found to be fairly independent of the type of flow. It was demonstrated that a generalized normal stress model can be used as an alternative method for the normal-stress model, and the tensorial invariance which is violated in the normal-stress model can be recovered. This generalized normal-stress model will be tested in a backward-facing step flow in both 'a priori' and 'a posteriori' manner in the future. Although the Bardina model constants CK and CN in Eq. (9) are rather consistent in three different flow fields, some variance was noticed. An attempt to determine these coefficients more accurately using the Dynamic scale model approach is currently underway. Acknowledgements K. H. is grateful to the Center for Turbulence Research for providing the oppor- tunity to access the DNS database located at CTR. We are grateful to Drs. P. Durbin, P. Moan, and J. H. Ferziger for stinmlating discussions. Gratitude is also extended to Drs. R. D. Moser and M. M. Rogers for their kindness in providing the recent DNS data. This work was partially supported by the International Exchange Program, University of Tokyo. REFERENCES BARDINA, J. 1983 Improved turbulence models based on large eddy simulation of homogeneous, incompressible turbulent flows. Ph.D. dissertation. Stanford University, Stanford, California. DURBIN, P. A. 1991 Near wall turbulence closure modelling without damping functions. Theor. Comp. Fluid Dynamics. 3, 1. GERMANO, M., PIOMELLI, U., MOIN, P. _: CABOT, W.H. 1991 A dynamic subgrid-scale eddy viscosity model. Phys. Fluids A. 3, 1760. HORIUTI, K. 1989 The role of the Bardina model in large eddy simulation of turbulent channel flow. Phys. Fluids A. 1,426. HORIUTI, K. 1990 Higher-order terms in the anisotropic representation of Reynolds stresses. Phys. Fluids A. 2, 10. HORIUTI, K. 1992 Damping effects in large eddy simulation subgrid-scale model- ing. In preparation. HORIUTI, K. 1993 A proper velocity scale for modeling subgrid-scale eddy viscosity in large eddy simulation. Phys. Fluids A. (1). KIM, J., MOlN, P. & MOSER, R. D. 1987 Turbulent statistics in fully developed channel flow at low Reynolds nmnbcr. J. Fluid Mech. 177, 133. KIM, J. 1990 unpublished data. LAUNDER, B. E., REECE, G. J. _ RODI, W. 1975 Progress in the development of a Reynolds-stress turbulence closure. J. Fluid Mech. 68, 537.