Large eddy simulation of turbulent channel flow with transverse roughness elements on one wall

International Journal of Heat and Fluid Flow xxx (2014) xxx–xxx

Contents lists available at ScienceDirect

International Journal of Heat and Fluid Flow journal homepage: www.elsevier.com/locate/ijhff

Large eddy simulation of turbulent channel flow with transverse roughness elements on one wall C.D. Dritselis ⇑ Mechanical Engineering Department, University of Thessaly, Pedion Areos, 38334 Volos, Greece

a r t i c l e

i n f o

Article history: Received 19 February 2014 Received in revised form 4 August 2014 Accepted 18 August 2014 Available online xxxx Keywords: LES Subgrid-scale model Immersed boundary method Roughness elements Ribs

a b s t r a c t This study examines the feasibility of large eddy simulation for predicting turbulent channel flows with two-dimensional roughness elements of square, circular and triangular shapes transversely placed on the bottom wall. Results are obtained for several values of the cavity width to the roughness height ratio using various subgrid-scale turbulence models. The present large eddy simulation predictions of mean streamwise velocity, root-mean-square velocity fluctuations, and skin frictional and form drags agree reasonably well with previously published results of direct numerical simulations at a low Reynolds number. All the subgrid-scale models examined here are capable of reproducing the relevant physics associated with the effect of the rough surface on the turbulent flow, exhibiting similar performances. Moreover, the use of the turbulence models leads to an improvement in the predictions of several turbulence statistics as compared with the case when no model is considered. Large eddy simulation can be combined easily with an immersed boundary method yielding satisfactory results based on a coarser grid resolution than in direct numerical simulation and, thus, it is suitable for the investigation of turbulent channel flows with riblets of various shapes. Ó 2014 Elsevier Inc. All rights reserved.

1. Introduction Roughness elements are often incorporated on the walls of a duct in order to promote turbulence and to enhance heat transfer (see, for example, Nagano et al., 2004; Saha and Acharya, 2004; Luo et al., 2005; Saidi and Sunden, 2000; Kim et al., 2004; Yucek and Dinler, 2006). These elements may have random shape and structure or well-defined dimensions and ordered pattern. Recently, several studies of direct numerical simulation (DNS) considered the turbulent flow over riblets with square, circular, or triangular cross-sections (Leonardi et al., 2003a,b, 2004, 2007; Ashrafian et al., 2004; Orlandi et al., 2006; Burattini et al., 2008; Leonardi and Castro, 2010). While restricted to flows with small Reynolds number, the use of DNS has shown itself to be a very suitable tool in obtaining detailed information of the relevant physics in turbulent flows with rough walls. Larger Reynolds numbers can be achieved potentially by using large eddy simulation (LES). Based on the standard and the dynamic Smagorinsky subgrid-scale turbulence models, flows with or without heat transfer in ribbed channels or pipes have also been investigated numerically

⇑ Tel.: +30 24210 74337; fax: +30 24210 74085. E-mail address: [email protected]

(Ciofalo and Collins, 1992; Cui et al., 2003; Iacono et al., 2005; Leonardi et al., 2006; Vijiapuparu and Cui, 2007; Khan et al., 2010). In the standard Smagorinsky model (Smagorinsky, 1963), the value of the model parameter is determined a priori and it changes depending on the flow. A wall function is always necessary to damp the eddy viscosity at the walls of the turbulent flow. Germano et al. (1991) proposed a dynamic procedure to calculate the model parameter during the numerical simulation. The subgrid-scale eddy viscosity of the initial version may acquire negative values causing subsequently numerical instabilities. Thus, averaging of the model coefficient in homogeneous directions is usually performed in order to overcome this difficulty. Some other alternatives include averaging locally over a stencil of three or more grid points in each direction (Zang et al., 1993; Gullbrand, 2004), averaging based on the dynamic localization model of Ghosal et al. (1995), time averaging of Meneveau et al. (1996) and ensemble averaging over many simultaneous computations proposed by Carati et al. (1996). Several others models suitable for engineering applications have also been proposed and used frequently (see, for example, Shimomura, 1994; Nicoud and Ducros, 1999; Yoshizawa et al., 2000; Vreman, 2004). In these models, the turbulent eddy viscosity is locally determined with some fixed values of the model parameters. Kobayashi (2005) proposed a local subgrid-scale model based on the coherent structures of the turbulent flow field. The model

http://dx.doi.org/10.1016/j.ijheatfluidflow.2014.08.008 0142-727X/Ó 2014 Elsevier Inc. All rights reserved.

Please cite this article in press as: Dritselis, C.D. Large eddy simulation of turbulent channel flow with transverse roughness elements on one wall. Int. J. Heat Fluid Flow (2014), http://dx.doi.org/10.1016/j.ijheatfluidflow.2014.08.008

2

C.D. Dritselis / International Journal of Heat and Fluid Flow xxx (2014) xxx–xxx

parameter is composed of a predefined, fixed value and a function describing the coherent structures that plays the role of walldamping. The coherent structure model has been tested in a series of canonical turbulent flows, such as rotating and non-rotating channel flows (Kobayashi, 2005) and magnetohydrodynamic (MHD) flows in channels (Kobayashi, 2006) and in ducts (Kobayashi, 2008). The results based on the coherent structure model were in good agreement with DNS and/or experimental data, exhibiting a level of accuracy at least similar to that obtained by using the standard dynamic model. The applicability of the coherent structure model was also assessed in large eddy simulations of a flow over a backward-facing step, a flow in an asymmetric plane diffuser, and staggered jets in cross-flow (Kobayashi et al., 2008). For all configurations, it gave almost the same performance as the dynamic localization model and the dynamic Smagorinsky model with averaging in homogeneous directions. However, the coherent structure model was inexpensive relative to the dynamic model and it was numerically stable without requiring averaging. The objective of the present study is to assess the applicability and the predictive performance of LES in turbulent channel flows with roughness elements placed on one wall. For this reason, the LES results produced for several cases with square, circular and triangular bars on the bottom wall of the channel based on the standard dynamic Smagorinsky model, the Lagrangian dynamic Smagorinsky model, and the coherent structure model are compared with those obtained by using DNS at a relatively low Reynolds number (see, for example, Leonardi et al., 2003a,b, 2004, 2007; Orlandi et al., 2006), which has not been done previously. The main interest here is to stably apply the aforementioned subgrid-scale turbulence models to the present flows with complex geometries and to obtain reasonable results with respect to the aforementioned DNS studies, and not to conclude which model is better. An immersed boundary technique is adopted to treat the two-dimensional (2d) wall disturbances of various cross sections, which allows the numerical solution of flows over complex geometries without the need of computationally intensive body-fitted grids (see, for example, Fadlun et al., 2000; Iaccarino and Verzicco, 2003; Mittal and Iaccarino, 2005; Orlandi and Leonardi, 2006). 2. Simulations overview

A schematic of the flow configuration is shown in Fig. 1. The velocity satisfies periodic boundary conditions in the x- and zdirections, and no-slip conditions at the walls. The 2d roughness elements are numerically treated by an immersed boundary method, which consists of imposing zero values to all fluid velocity components on the stationary boundary surface that does not necessarily coincide with the computational grid (see, for details, Fadlun et al., 2000; Orlandi and Leonardi, 2006). In accordance with the aforementioned studies, zero velocities are imposed in the grid points within the roughness elements. At the first grid point outside each roughness element, all the viscous derivates and the terms of the extra subgrid-scale stress tensor in the filtered flow equations are discretized by using the distance between the grid point and the boundary of the wall disturbance and not the actual mesh size. This is done in order to avoid describing the geometry in a stepwise way and to perform numerical simulations by maintaining a constant fluid mass flow rate in the channel. 2.2. Subgrid-scale turbulence models In the present study, the extra subgrid-scale stress tensor srij appearing in the LES equations due to the filtering procedure is modeled by the standard and Lagrangian dynamic Smagorinsky models and the coherent structure model. 2.2.1. Standard dynamic Smagorinsky model In the standard dynamic Smagorinsky model, the subgrid-scale stress srij is modeled as

srij ¼ 2cD 2 jsjsij ;

ð4Þ

where sij is the resolved strain-rate tensor and |s|=(2sijsij )1/2 is its  = (DxDyDz)1/3, where magnitude. The characteristic length is D Dx, Dy, and Dz are the grid spacings in the x-, y-, and z-directions, respectively. The model parameter c is determined by the dynamic procedure proposed by Germano et al. (1991) based on Lilly’s (1992) modification

c¼

hlij mij ixi ;hom hmkl mkl ixi ;hom

;

ð5Þ

with lij and mij given as

2.1. Flow governing equations

j  u  ^i u ^i u ^j ; lij ¼ u

ð6Þ

The governing equations of the present flow are the threedimensional (3d), unsteady, incompressible filtered continuity and Navier–Stokes equations

^ 2 js^j s^ ; mij ¼ 2D2 jsj^sij  2D ij

ð7Þ

i @u ¼ 0; @xi

ð1Þ

@ srij i @ u i u j  i @u 1 @p @u þ ¼ þm  þ Pdi1 ; @t @xj qf @xi @xj @xj @xj

ð2Þ

 i is the component of the large-scale fluid velocity in the where u  is the fluid pressure, qf and m are the fluid density i-direction, p and kinematic viscosity, respectively, P is the extra pressure gradient in order to keep the fluid mass flow rate constant, and dij is the x , u  y , and u  z are the fluid velocity comKronecker delta. Note that u ponents at the grid-scale in the streamwise (x), wall-normal (y), and spanwise (z) directions, respectively, and the indices i,j denote x, y, or z. The term srij represents the effect of the subgrid-scale motions on the resolved grid-scale velocities of the fluid. It is the traceless part of the subgrid-scale stress tensor, defined as

1 3

srij ¼ sij  smm dij with sij ¼ ui uj  ui uj :

ð3Þ

where (^) denotes variables calculated on the test filter and h ixi ;hom averaging over homogeneous directions. The present large eddy simulations of the turbulent channel flow with roughness elements on the bottom wall are performed by using a box filter in the physical space based on the trapezoidal rule. No filtering is applied in the wall-normal direction, while the width of the test filter is twice the size of the uniform grid spacing in the x- and z-directions. A cut-off is set to ensure positive values of the total viscosity in order to prevent any numerical instability. 2.2.2. Lagrangian dynamic Smagorinsky model The model parameter c in Eq. (4) is also calculated based on the approach developed by Meneveau et al. (1996) as

c¼

iLM : iMM

ð8Þ

In principle, the quantities iLM and iMM are obtained from the solution of separate transport equations. However, the numerator and denominator of Eq. (8) are calculated by using a simple time discretization resulting in

Please cite this article in press as: Dritselis, C.D. Large eddy simulation of turbulent channel flow with transverse roughness elements on one wall. Int. J. Heat Fluid Flow (2014), http://dx.doi.org/10.1016/j.ijheatfluidflow.2014.08.008

C.D. Dritselis / International Journal of Heat and Fluid Flow xxx (2014) xxx–xxx

3

Fig. 1. Flow configuration.

nþ1

nþ1

iLM ðxÞ ¼ Hfe½lij mij 

n

 n DtÞg; ðxÞ þ ð1  eÞiLM ðx  u

nþ1

n

 n DtÞ; iMM ðxÞ ¼ e½mij mij nþ1 ðxÞ þ ð1  eÞiMM ðx  u

ð9Þ ð10Þ

where

e¼

Dt=t ns 1=8 n n and t ns ¼ 1:5D ðiLM iMM Þ ; 1 þ Dt=tns

ð11Þ

and where H{x} = max(x, 0) is the ramp function and bold symbols  nDt) and inMM(x  u  nDt) are indicate vectors. The values of inLM(x  u obtained through linear interpolation. In accordance with Meneveau et al. (1996), the initial values of iLM and iMM are set equal to iMM(x, 0) = mijmij(x, 0) and iLM(x,0) = 0.0256 mijmij(x, 0), respectively. The wall-normal derivates of both iLM and iMM are set to zero at the walls. 2.2.3. Coherent structure model Details about the coherent structure model may be found in Kobayashi (2005) and only the main features of the model are summarized here. In the coherent structure model, the model parameter c is ‘‘locally’’ defined as

c ¼ cCSM jf CS j3=2 f X ;

ð12Þ

with

C CSM ¼ 1=22; f CS ¼ q=e; f X ¼ 1  f Cs ; j @ u i 1 1 @u  ij x  ij  sijsij Þ ¼  ; q ¼ ðx 2 2 @xi @xj  2 j 1 1 @u  ij x  ij þ sij u  ij Þ ¼ e ¼ ðx ; 2 2 @xi   j @ u i 1 @u  ij ¼ ; x  2 @xi @xj

ð13Þ ð14Þ ð15Þ ð16Þ

where cCSM is a fixed model constant, fCS is the coherent structure function defined as the second invariant q normalized by the magnitude of the velocity gradient tensor e, fX is the energy-decay suppression function, which becomes approximately 1.1 in homogeneous isotropic turbulences and at the center of turbulent  ij is the vorticity tensor in a resolved flow field. channel flows, and x Moreover, fCS and fX have definite upper and lower limits:

1 6 f CS 6 1; 0 6 f X 6 2:

ð17Þ

As a result, the parameter of the coherent structure model has smaller variance than that of the standard dynamic model. In addition, it is locally determined since no averaging is performed. 2.3. Solution method and numerical details The numerical algorithm used to solve Eqs. (1) and (2) is based on a semi-implicit, fractional step method (see Orlandi, 2000). The

time integration incorporates an implicit second-order Crank– Nicolson method for the diffusion terms and a low-storage explicit third-order Runge–Kutta method for the convection terms and the extra subgrid-scale stress tensor. All spatial derivates are discretized by using a second-order central differencing scheme. A Poisson equation is solved for the pseudo-pressure to enforce mass continuity. The present implementation of the immersed boundary method is similar to that used in the DNS studies of Orlandi, Verzicco and coworkers (see, for example, Fadlun et al., 2000). It satisfies the continuity and preserves the correct behavior near the wall, while the second-order accuracy of the LES solver is maintained. Therefore, the numerical error introduced by the immersed boundary method is similar to that associated with the numerical schemes of the solver. The projection step of the numerical algorithm ensures that continuity is globally and locally satisfied with a maximum error of the order of 1015 in the present double precision calculations of LES combined with the immersed boundary method. In addition, for all the cases shown here it was verified that the largest change in the imposed values at the immersed boundaries during the projection step was less than 1011. The values of the cavity width w to the roughness height ratio k (=0.2h) examined here are w/k = 0, 1, 3, and 7 for square bars, and w/k = 0 and 1 for circular rods and triangular elements. These cases consist of a sufficient number of characteristic flows to examine the feasibility of LES coupled with the immersed boundary method for predicting the turbulent flows in channels with a rough lower wall. Note that the distance between two circular or triangular elements varies, while the streamwise wavelength k (=w + k) is constant, and that here w denotes the minimum distance between two successive roughness elements. The case with w/k = 0 corresponds to the channel flow with smooth walls, in which the lower smooth wall is also described by the immersed boundary method. All simulations are performed at the same Reynolds number of Reb (=ub2h/m) = 5600 based on the bulk velocity ub and the distance between the smooth wall and the crest of the rough wall 2h. This relatively low Reynolds number is adopted to facilitate the comparison with the available low Re DNS results, and it corresponds to a DNS value of Res  178 based on the friction velocity and the half distance h in the case with both smooth walls (Leonardi et al., 2006). The predictions of the LES models and without any subgrid-scale model are in the range 171 6 Res 6 175 for the case with both smooth walls, while Res varies in the range of 168 6 Res 6 317 in the cases with square roughness elements, and 225 6 Res 6 246 in the cases with the circular rods and triangles. For all the simulations, a fixed fluid mass flow rate was considered by imposing a mean pressure gradient P in the streamwise momentum equation, while the CFL number was fixed at a constant value of 0.5. The dimensions of the channel in the streamwise and spanwise directions are lx = 8h and lz = ph, respectively, and the

Please cite this article in press as: Dritselis, C.D. Large eddy simulation of turbulent channel flow with transverse roughness elements on one wall. Int. J. Heat Fluid Flow (2014), http://dx.doi.org/10.1016/j.ijheatfluidflow.2014.08.008

4

C.D. Dritselis / International Journal of Heat and Fluid Flow xxx (2014) xxx–xxx

computational mesh is 201  94  65 in the x-, y-, and z-directions, respectively. Each roughness element is discretized by using 30 equidistant grid points over the range of 1.2 < y/h < 1, while the rest 64 of the grid points are unevenly distributed in 1 < y/h < 1. A maximum of five grid points are used in the x-direction. It should be noticed that away the roughness elements, a mesh with 64  65 grid points is used in the flow domain of 2h  p. For the case with both smooth walls, the grid resolution in wall units (denoted by the symbol +) is 0.48 6 Dymin+ 6 0.49, 10.6 6 Dymax+ 6 10.8, and 7.7 6 Dz+ 6 7.9, depending on the subgrid-scale model, which is similar to that adopted by the majority of the LES studies at the same Reynolds number and consistent with the requirements for obtaining accurate results (see Piomelli and Balaras, 2002). The grid resolution in the flows with wall roughness is 0.48 6 Dymin+ 6 0.89, 10.5 6 Dymax+ 6 19.6, and 7.7 6 Dz+ 6 14.3. The present LES mesh size is at least 8 times smaller than in the corresponding DNS studies, which is a significant gain and consistent with the general spirit of LES. Consequently, it results in a significant effect of the subgrid-scale models, as it will be discussed and shown in the Results. Although several other turbulence models might also have been used, the standard and the Lagrangian dynamic Smagorinsky models and the coherent structure model were chosen in order to conduct the present study within a manageable parameter range, considering also the influence of the averaging procedure of the model parameter. In the standard dynamic model, the averaging of the model parameter is performed only in the spanwise direction, the Lagrangian dynamic model incorporates an alternative averaging over time, while the coherent structure model does not implement an averaging procedure of the model parameter. The present numerical code has been verified successfully for both DNS and LES of turbulent channel flow with smooth walls (see, for example, Orlandi, 2000; Dritselis et al., 2011; Dritselis and Vlachos, 2011). The accuracy of the present numerical treatment of the roughness elements was verified in two test cases: (a) the flow around a sphere and (b) the experiment of Furuya et al. (1976) for a boundary layer over 2d transverse circular rods. More specifically, several numerical simulations were carried out for those cases of axis-symmetric flow around a sphere examined by Fornberg (1988) who used body-fitted meshes and by Fadlun et al. (2000) who used a similar immersed boundary method. Fig. 2a shows that the present values of the drag coefficient Cd agree very well with those of Fornberg (1988) and Fadlun et al. (2000) for various values of the Reynolds number. Regarding the second test case, additional simulations were performed in a channel with circular rods placed on one wall, while the other was smooth at a lower Reynolds number of Reb = 5600 relative to that of the experiment of Furuya et al. (1976). For these simulations, the computational domain was 8h  2h  ph and the grid was 401  140  97 in accordance with Orlandi et al. (2006). In the x-direction, each roughness element was discretized with ten grid points, while 40 grid points with uniform spacing were used over the range 1.2 < y/h < 1. The ensemble averaged results shown in Fig. 2b were obtained by post-processing a smaller number of instantaneous flow fields as compared to Orlandi et al. (2006). However, Fig. 2b shows that the agreement for the pressure around the circular rods with w/k = 3 and 7 is rather good relative to the experimental results of Furuya et al. (1976) and the DNS results of Orlandi et al. (2006). The successful verification implies that the present numerical method is accurate and it can be used to reproduce the flow past the surface of the wall disturbances studied here. In addition, the good agreement between the numerical and experimental results in Fig. 2b indicates also that the influence of the Reynolds number is not so important in flows with rough walls as it is in flows with smooth walls (see also, Leonardi et al., 2006; Hanjalic and Launder, 1972).

(a)

(b)

Fig. 2. (a) Comparison of Cd based on the present model against the results of Fornberg (1988) and Fadlun et al. (2000) for the axis-symmetric flow around a sphere and, (b) distribution of averaged pressure around a circular roughness element normalized by qub2. Lines: present results, open symbols: DNS results of Orlandi et al. (2006), filled symbols: experimental results of Furuya et al. (1976). —, h, j: w/k = 3; —, r, .: w/k = 7. h is oriented clockwise and its origin is at the bottom wall.

3. Results and discussion 3.1. Square roughness elements Fig. 3 shows the mean streamlines averaged in time and in the spanwise direction for various values of w/k predicted by the coherent structure model, while qualitatively similar distributions were obtained by the standard and the Lagrangian dynamic Smagorinsky models. For w/k = 1, the whole cavity is occupied by a single recirculation region. The size of the vortex is increased for w/k = 3 and a smaller secondary vortex with opposite circulation is observed close to the left vertical wall of the cavity. For these two cases of w/k = 1 and 3, the separation of the flow occurs at the edge of each element, while the flow reattaches at the vertical wall of the following roughness element. For w/k = 7, a large recirculation zone is located downstream of each element, which is accompanied by a smaller vortex. The flow reattaches on the bottom wall of the cavity and separates again as the vertical wall of the next element is approached. The present LES results are in good qualitative agreement with those of Leonardi et al. (2003a) who used a similar immersed boundary technique together with DNS, and with the results of Cui et al. (2003) who used a body-fitted numerical code coupled with LES in order to study the turbulent

Please cite this article in press as: Dritselis, C.D. Large eddy simulation of turbulent channel flow with transverse roughness elements on one wall. Int. J. Heat Fluid Flow (2014), http://dx.doi.org/10.1016/j.ijheatfluidflow.2014.08.008

C.D. Dritselis / International Journal of Heat and Fluid Flow xxx (2014) xxx–xxx

5

(a)

Fig. 3. Mean streamlines for square bars with w/k = 1, 3, and 7. The flow is from left to right.

flow in similar geometries. All subgrid-scale models predict a reattachment length of about 5k for w/k = 7, which is close to the value of 4.8k found in the DNS study of Leonardi et al. (2003b). Fig. 4 shows the distribution of the non-dimensional mean  x ix;z;t =ub Þ as a function of y/h, in streamwise velocity hU x ix;z;t ð¼ hu the range of –1 6 y/h 6 1, for turbulent flows in a channel with smooth walls (w/k = 0) and with square bars placed on the lower wall and w/k = 1, 7. The brackets hix,z,t denote that averaging has been performed with respect to the x and z-directions and time t, in order to assess the predictability of the LES models regarding the overall effect of the roughness elements on the overlying flow (y/h P –1). The results shown in Fig. 4 correspond to the cases with the standard dynamic Smagorinsky model (S. Dyn-SD), the Lagrangian dynamic Smagorinsky model (L. Dyn-LD), the coherent structure model (CS), and without using subgrid-scale turbulence model (No model-NM). The main effect of the transverse roughness elements on the mean streamwise velocity hU x ix;z;t at the outer flow region is the shift of its maximum value toward the upper smooth wall with increasing w/k, which is well captured in the present LES study. Small differences between the distributions of hU x ix;z;t can be seen and all profiles almost coincide with the DNS results of Orlandi et al. (2006). This indicates a good performance of all the subgrid-scale models used in the present well-resolved large eddy simulations with either smooth or rough walls. The DNS studies (Leonardi et al., 2003; Orlandi et al., 2006) at this low Reynolds number have shown that, despite the large shift of the maximum value of hU x ix;z;t , the effect of w/k on the friction velocity of the upper smooth wall Us,u is smaller than that of the lower rough wall Us,l, where the subscripts u and l refer to the upper and lower walls, respectively. This is also confirmed here; more specifically, the values of the friction velocity of the upper smooth wall normalized by ub for w/k = 0 are (Us,u)SD = 0.061, (Us,u)LD = 0.061, (Us,u)CS = 0.0625, and (Us,u)NM = 0.0625, relative to the value of (Us,u)DNS = 0.063 yielded in Orlandi et al. (2006). For w/k = 7, these values are (Us,u)SD = 0.0728, (Us,u)LD = 0.0737, (Us,u)CS = 0.0738, (Us,u)NM = 0.0738, and (Us,u)DNS = 0.073, exhibiting small differences from those corresponding to the case with w/k = 0. On the other hand, the friction velocity at the rough wall Us,l increases significantly, which is also verified in the present LES study. For w/k = 7, it is found that (Us,l)SD = 0.1424, (Us,l) LD = 0.1423, (Us,l)CS = 0.1401, (Us,l)NM = 0.1344, and (Us,l)DNS  0.16, as compared with the values of (Us,l)SD = 0.0609, (Us,l)LD = 0.0612, (Us,l)CS = 0.0625, (Us,l)NM = 0.0625, and (Us,l)DNS = 0.0634 in the w/k = 0 case. As it can be seen, Us,l and Us,u are underestimated by all the subgrid-scale models with respect to the DNS predictions. However, their variations with w/k can be reproduced by all the LES models. It should be noted that the friction velocity

(b)

(c)

Fig. 4. Distribution of mean streamwise velocity normalized by ub for square roughness elements with w/k = 0 (a), w/k = 1 (b), and w/k = 7 (c). Lines: present LES results, symbols: DNS results of Orlandi et al. (2006).

for rough walls is given by the square root of the sum of the skin frictional drag and the form drag, whereas for a smooth wall only the skin friction contributes to Us. In low Reynolds number turbulent flows in channels with smooth walls, a relatively large value of the model parameter

Please cite this article in press as: Dritselis, C.D. Large eddy simulation of turbulent channel flow with transverse roughness elements on one wall. Int. J. Heat Fluid Flow (2014), http://dx.doi.org/10.1016/j.ijheatfluidflow.2014.08.008

6

C.D. Dritselis / International Journal of Heat and Fluid Flow xxx (2014) xxx–xxx

increases the subgrid-scale dissipation, which in turn attenuates the fluid vorticity of the turbulent coherent structures, while their size is enlarged. The logarithmic region is then shifted at larger distances from the wall and the friction velocity is decreased (see, for example, Dritselis and Vlachos, 2011). The above scenario is observed in the flow case with w/k = 0 examined here for both the standard and the Lagrangian versions of the dynamic model, and at a lesser extent for the coherent structure model. For a turbulent channel flow with transverse roughness elements, the coherence is decreased with increasing w/k, enhancing the isotropy of the stresses and decreasing the skin friction drag (Leonardi et al., 2003b, 2004; Orlandi et al., 2006). The present LES study indicates that the viscous drag for w/k = 7 is smaller than that for w/k = 0 and 1, as shown in Fig. 4 by the smaller velocity gradient at y/h = –1 observed in the case with w/k = 7 for all models. This is further corroborated by the non-dimensional viscous shear stress at the wall, hC f iz;t ¼ Re1 ð@hU x iz;t =@YÞY¼1 , where Y = y/h and the averaging is performed over the z-direction and time. In Fig. 5, hC f iz;t is shown for w/k = 1 and 7 along the horizontal wall, with x/k = 0 at the leading edge of a square roughness element. The LES predictions of hC f iz;t are compared with the DNS results of Leonardi et al. (2003b) and (2007). It can be seen that the trends of the DNS results are well captured by the LES models from a qualitative point of view. More specifically, hC f iz;t attains negative values over a range of x/k values within the cavity, exhibiting one negative peak for w/k = 1 and two for w/k = 7. On the crests of the square element (0 < x/k < 1) a positive peak is observed at the leading edge, while hC f iz;t is reduced with increasing w/k. However, there are certain quantitative discrepancies. For example, hC f iz;t is underestimated by all the LES models for w/k = 1 and 0 < x/k < 1, while its positive peak at the leading edge of the element cannot be reproduced correctly, especially in the case with w/k = 7. The main reasons for these differences are the smaller grid resolution and the slightly different distribution of the grid points used in the present LES analysis with respect to the fully-resolved DNS studies. Fig. 6 shows the pressure contours around a square roughness element predicted by the LES models and without subgrid-scale model in the representative case with w/k = 7. Maximum (positive) and minimum (negative) values of pressure are observed on the upstream face and crest of the square bar in qualitative agreement with the DNS findings (see, for example, Leonardi et al., 2003b). The pressure field in Fig. 6 indicates the existence of a strong pressure gradient that drives the flow outward around the vicinity of the leading edge of the element. For w/k = 7, the pressure gradient

(a)

is not strong enough to observe flow separation above the crests, but it results in an enhancement of the recirculation zone and a significant decrease of the shear friction as compared to the cases with smaller values of w/k. The total shear stress on a rough wall is the sum of the skin frictional drag Fd and the form drag Pd. These are obtained by integrating hC f iz;t and hPiz;t over a streamwise wavelength k, Rk Rk F d ¼ k1 0 hC f iz;t s  x ds and P d ¼ k1 0 hPiz;t n  x ds; where s is a coordinate that follows the wall contour, s and n are the unit vectors parallel and normal to the wall, respectively, and x is the unit vector in the x-direction. Both DNS and LES analyses indicate that Fd is reduced with increasing w/k, while the opposite is observed for Pd as shown in Fig. 7. It can be seen that the values of Fd produced by all the subgrid-scale models are close to each other, while larger differences are found for Pd. More specifically, Fig. 7 shows that the differences between the DNS and LES predictions of Fd in the cases with w/k = 3 and 7 are smaller than those with w/k = 0 and 1. In contrast, the discrepancies of the LES values of Pd from the DNS ones are increased with increasing w/k. This different behavior can be explained as follows: the skin frictional drag depends on the mean velocity, which at the present grid resolution is well reproduced as shown in Figs. 3 and 4. Moreover, the differences in the mean velocity produced by the LES models with respect to the DNS results near the wall become smaller with increasing w/k. For example, the mean velocity predicted by the LES models are 34–37% smaller than in DNS for w/k = 0 at y/h = 0.092, 13–14% for w/k = 1 at y/h = 0.962, and 10–11% for w/k = 7 at y/h = 0.985. Consequently, the changes in the mean velocity are also reflected on the LES predictions of the skin frictional drag. On the other hand, by integrating the Navier–Stokes equation for the streamwise velocity component over the cavity region (see Fig. 1) as shown in the analysis of Leonardi et al. (2003b), it turns out that the form drag Pd depends also on the Reynolds shear stress at y/h = 1. The latter quantity is generally under-estimated by the LES models, similarly to the diagonal Reynolds stresses that are discussed below, leading to smaller LES values of Pd as compared with the DNS results, depending on the value of w/k. For all the LES models and without subgrid-scale model, the largest discrepancies from the DNS results are observed in the skin frictional and form drags and to a lesser extent in the friction velocity due to its square root dependence on the sum of Fd and Pd. Such quantities depend on the wall contour resolution, which is not as good as it is in the DNS studies. For example, the LES models, including no model, cannot reproduce the large peak of hC f iz;t at

(b)

Fig. 5. hC f iz;t for w/k = 1 (a) and w/k = 7 (b). Lines: present LES results, symbols: DNS results of Leonardi et al. (2007) for w/k = 1 and Leonardi et al. (2003b) for w/k = 7.

Please cite this article in press as: Dritselis, C.D. Large eddy simulation of turbulent channel flow with transverse roughness elements on one wall. Int. J. Heat Fluid Flow (2014), http://dx.doi.org/10.1016/j.ijheatfluidflow.2014.08.008

C.D. Dritselis / International Journal of Heat and Fluid Flow xxx (2014) xxx–xxx

(a)

(b)

(c)

(d)

7

Fig. 6. Contours of pressure around a square roughness element with w/k = 7 predicted by the standard dynamic Smagorinsky model (a), the Lagrangian dynamic Smagorinsky model (b), the coherent structure model (c), and without any model (d). Lines: positive values, dashed lines: negative values. The interval between successive contour levels is 0.01. The flow direction is from left to right.

(a)

(b)

Fig. 7. Variation of the skin frictional drag Fd (a) and of the form drag Pd (b) with w/k for square roughness elements. Lines: present LES results, symbols: DNS results of Leonardi et al. (2003b).

the leading edge of the roughness element in Fig. 5, influencing the integrated quantity Fd. As discussed above, Pd depends indirectly on the effect of the subgrid-scale flow motions. Considering also that the flow physics imposed by the rough wall leads to an increase of the local Reynolds number, the overall effect of the subgrid-scales will be more significant on the prediction of Pd. For w/k = 3, the LES predictions of Fd differ from that obtained without subgrid-scale model about 4.5–5.5%, while, on the other hand, Pd differs about 2–15%. The small variance of 1% in the former case of Fd indicates that it does not depend explicitly on the subgridscale motions, while the opposite is true for Pd which is more sensitive to the modeling effects of the small unresolved scales. It

should be noted also that the LES predictions of Fd and Pd differ from those without subgrid-scale model. The improvement associated with the subgrid-scale models with respect to the case without turbulence model is encouraging in the flows with the rough wall. For example, at w/k = 7, the friction velocity at the lower wall based on the standard dynamic model is 11% smaller than the DNS value, while this is 16% when no model is used. The LES predictions of Pd are also better than those without a subgrid-scale model as w/k increases, as it is evident in Fig. 7. Moreover, all the LES models predict reasonably the dependence of Fd, Pd, and Us with w/k consistent with the DNS trends, ensuring that the present results will converge to the DNS results by adopting a better grid resolution. Indeed, results from grid refinement tests for the case of square roughness elements with w/k = 7 toward a DNS based on a computational mesh with 401  160  129 in the x-, y-, and z-directions, respectively, reveal that approximately the same friction velocity at the lower rough wall of Us,l  0.16 as in DNS is obtained, while both Fd and Pd differ less than 2% as compared with the DNS values of Leonardi et al. (2003b) and Orlandi et al. (2006). Fig. 8 shows the LES predictions of the root-mean-square (rms) velocity fluctuations Ux,rms, Uy,rms, Uz,rms normalized by the bulk velocity ub as a function of y0 /h (=y/h + 1) in the range of 0 6 y0 /h 6 1 for w/k = 0, 1, and 7 in comparison with those yielded without using a subgrid-scale model and the DNS results of Orlandi et al. (2006). In the case with smooth walls w/k = 0, a good agreement is observed for Ux,rms and Uy,rms between the LES and the DNS results, while Uz,rms is under-predicted by all the subgrid-scale models. At the grid resolution of the present study, the molecular viscosity provides sufficient dissipation to avoid numerical instabilities, even without any turbulence model. If finite differences simulations are not fully resolved, the small scales generated by the insufficient resolution act as a subgrid-scale model resulting in satisfactory statistics (Orlandi and Leonardi, 2004; Dritselis and Vlachos, 2011). Consequently, more dissipation is given when a LES model is used, reducing the velocity fluctuations and increasing the mean velocity. For roughness with w/k = 1, it is seen that Ux,rms is underestimated near the wall and in the central region of the channel, while it is little overestimated in the rest of the flow region. The LES models yield smaller values of Uy,rms and Uz,rms throughout the channel relative to the DNS results. The larger discrepancies are observed for the spanwise component of the velocity fluctuations similarly to the case with smooth walls. For w/k = 7, it is seen that Ux,rms obtained without turbulence model is close to the DNS distribution, while the agreement is reasonable for the LES models, with the coherent structure model producing a distribution of Ux,rms closer to the DNS one. In contrast, Uy,rms is under-predicted as compared with the DNS distribution near the wall by all the LES models, but the agreement is rather good in the central flow region. The Uz,rms profile from the Lagrangian dynamic model agrees closely with the DNS one, whereas the standard dynamic

Please cite this article in press as: Dritselis, C.D. Large eddy simulation of turbulent channel flow with transverse roughness elements on one wall. Int. J. Heat Fluid Flow (2014), http://dx.doi.org/10.1016/j.ijheatfluidflow.2014.08.008

8

C.D. Dritselis / International Journal of Heat and Fluid Flow xxx (2014) xxx–xxx

(a)

(b)

(c)

(d)

(e)

(f)

Fig. 8. Distributions of the rms velocity fluctuations for square bars with w/k = 0, 1 (a, c, e), and w/k = 7 (b, d, f): (a, b) streamwise, (c, d) wall-normal, and (e, f) spanwise component. Lines: present LES results, symbols: DNS results of Orlandi et al. (2006).

model and the coherent structure model predict similar distributions with slightly smaller values of Uz,rms near the wall and in the central region of the channel relative to the DNS results of Orlandi et al. (2006). The differences between the LES and the DNS predictions of Uz,rms for w/k = 7 are smaller than those found in the case without using subgrid-scale model, as well as the other components of the velocity fluctuations, in contrast to what is found for w/k = 1. The above discussion indicates that all the LES

models are capable of reproducing the large increase of the velocity fluctuations in the region near the rough wall with increasing the value of w/k, as observed previously in the DNS studies. For the comparison between the LES and DNS results, the reduced (deviatoric) turbulence intensities are more proper quantities than the turbulence intensities shown in Fig. 8 (see Winckelmans et al., 2002). Since the DNS reduced turbulence intensities are not available, the rms velocities based on the

Please cite this article in press as: Dritselis, C.D. Large eddy simulation of turbulent channel flow with transverse roughness elements on one wall. Int. J. Heat Fluid Flow (2014), http://dx.doi.org/10.1016/j.ijheatfluidflow.2014.08.008

C.D. Dritselis / International Journal of Heat and Fluid Flow xxx (2014) xxx–xxx

large-scale components were shown in the present study and, therefore, the discussion was focused mostly on highlighting qualitatively the agreement or the differences between the LES and DNS results with respect to the relevant physics due to the effect of the rough surface. In this sense, Fig. 8 indicates rather the differences in the results between the LES models and without subgrid-scale model, as a consequence of their different impact. For example, for w/k = 0, Uy,rms differs from the DNS value about 28%, 26%, 13%, and 8% at y0 /h = 0.103 based on the standard dynamic model, the Lagrangian dynamic model, the coherent structure model, and without subgrid-scale model, respectively. For w/k = 1, the differences in Uy,rms are 46% (SD), 42% (LD), 34% (CS) and 31% (NM) at y0 /h = 0.087, while for w/k = 7 these are 7% (SD), 5% (LD), 6% (CS) and 8% (NM) at y0 /h = 0.114. The differences between the LES and DNS results generally depend on the wallnormal location, while similar observations can be made for the other two components of the rms velocity fluctuations. Fig. 9(a) shows the ratio of the subgrid-scale eddy viscosity mt and the molecular viscosity m for the standard dynamic model, the Lagrangian dynamic model, and the coherent structure model in the cases with w/k = 0, 1, and 7. The quantity mt serves as an indicator of the effect of the subgrid-scale flow motions that are

(a)

(b)

9

actually modeled during the simulations. The values of vt yielded by all the LES models are smaller than those of v for w/k 6 3, while mt is comparable to m in the region adjacent the rough wall for w/k = 7, as shown in Fig. 9(a). Consequently, although the grid size used in the present LES study is adequate to resolve a significant part of the energy of the turbulent flow scales in the cases with either smooth or rough walls for Reb = 5600, the effect of the subgrid-scale models is significant. It should be noted that similar distributions of mt/m have been reported in previous LES studies (see, for example, Najjar and Tafti, 1996; Meneveau et al., 1996; Gullbrand, 2004; Fröhlich et al., 2005; Kobayashi et al., 2008). Fig. 9(a) shows that, for w/k = 1, the ratio mt/m is slightly increased relative to the case with w/k = 0 mainly near the rough wall, while small differences are observed for y/h > 0 between these two cases. On the other hand, mt/m is significantly increased throughout the channel for w/k = 7. This denotes that the subgrid-scale turbulence models have more influence on the resolved flow field with increasing w/k, not only near the region of the rough wall, but for the whole channel flow domain. This behavior can be attributed to the fact that wall roughness generally improves the flow isotropy near walls, the local Reynolds number is increased and smaller flow scales appear with increasing w/k (Leonardi et al., 2003b, 2004). Thus, a larger impact of the subgrid-scale models can be anticipated in the cases with roughness elements, since more unresolved fluid motions have to be modeled. The different distributions of mt are due primarily to differences in the calculations of the model parameter as opposed to changes in the strain-rate magnitudes. Regarding the standard and the Lagrangian dynamic models, the discrete profiles of mt can be considered an outcome of the different averaging procedure of the model parameter. The values of mt predicted by the coherent structure model are smaller than those of the other two LES models, which is more obvious in the case with w/k = 7. This can be attributed to the smaller variance of the model parameter and its dependence on the computational grid which is stretched near the channel walls, in contrast to the standard and Lagrangian dynamic models (see Kobayashi et al., 2008). Despite the differences in the values of vt, the statistics of the first and second order moments produced by all the LES models are close each other. Fig. 9(b) shows that the differences in the eddy viscosity are reflected on the subgrid-scale dissipation rate. The subgrid-scale dissipation rate is calculated as Tij Sij and the molecular dissipation rate as – (2/Reb) Sij Sij : Fig. 9(b) reveals that the dynamic models contribute more than 10% of the total dissipation rate, while for the coherent structure model it is about 5%, over much of the channel for w/k = 0. The contribution of the total dissipation rate is substantially increased to 30–45% for the dynamic models, and 20–27% for the coherent structure model covering a significant part of the channel in the case with a lower rough wall and w/k = 7. This is significant since it indicates that all the subgrid-scale models play a key role in the resolved-scale turbulence energy budget for all the cases with smooth and rough walls examined here. Similar ratios of the subgrid-scale dissipation rate to the total dissipation rate as in the case with w/k = 0, have also been reported in previous LES studies in turbulent flows with smooth walls (see, for example, Meneveau et al., 1996). 3.2. Roughness elements with circular or triangular cross-section

Fig. 9. Distributions of the ratio of the subgrid-scale eddy viscosity mt and the molecular viscosity m (a) and of the ratio of the subgrid-scale dissipation rate to the total (subgrid-scale plus molecular) dissipation rate (b) for the standard dynamic model, the Lagrangian dynamic Smagorinsky model, and the coherent structure model in the cases with square bars and w/k = 0, 1, and 7. The distributions for w/k = 1 and 7 are shifted upwards by 0.3 and 0.6, respectively.

In this section, the performance of the subgrid-scale turbulence models is also examined in turbulent channel flows with 2d roughness elements of circular or triangular cross-section placed on the lower wall. Fig. 10 shows the mean streamlines averaged in time and in the spanwise direction for circular rods and triangular elements with w/k = 1 predicted by the coherent structure model,

Please cite this article in press as: Dritselis, C.D. Large eddy simulation of turbulent channel flow with transverse roughness elements on one wall. Int. J. Heat Fluid Flow (2014), http://dx.doi.org/10.1016/j.ijheatfluidflow.2014.08.008

10

C.D. Dritselis / International Journal of Heat and Fluid Flow xxx (2014) xxx–xxx

while similar results are obtained by the other models. In accordance with the case of square bars, the cavity between two rods or triangles with w/k = 1 is occupied by a single recirculation, as shown in Fig. 10. The streamlines above the roughness elements exhibit a weak undulation. For circular rods, the separation point is a function of the Reynolds number, occurring slightly downstream of the stagnation point. This is in contrast to the cases with square bars and triangles for which separation is due to the surface discontinuity. Fig. 11 shows the distributions of the non-dimensional mean streamwise velocity hU x ix;y;t as a function of y/h for both circular rods and triangular elements with w/k = 1. Once again, it can be seen that all the subgrid-scale models are in good agreement with the DNS results of Orlandi et al. (2006). The LES models reproduce well all the relevant DSN trends, namely, the stronger effect of the rods and triangles on the mean streamwise velocity relative to the square bars, as that is indicated by the generally higher values of hU x ix;z;t and the larger shift of its maximum value toward the upper smooth. The values of the friction velocity of the upper smooth wall Us,u normalized by ub for w/k = 1 vary between 0.0631 (SD) and 0.0658 (NM) for circular rods, while for triangles these vary between 0.0668 (SD) and 0.0687 (NM, CS). The predictions of the friction velocity at the upper smooth wall are close to those in the cases with square bars, indicating a generally small impact of the rough surface with w/k = 1 on the overlying turbulent flow, independent of the shape of the 2d wall disturbances. On the other hand, the values of the friction velocity of the lower rough wall are (Us,l)SD = 0.095, (Us,l)LD = 0.095, (Us,l)CS = 0.096, and (Us,l)NM = 0.093 for circular rods and (Us,l)SD = 0.103, (Us,l)LD = 0.104, (Us,l)CS = 0.101, and (Us,l)NM = 0.099 for triangular elements. It can be seen that Us,l is significantly larger than that in the corresponding cases with square bars. This can be attributed to the significant increase of the form drag Pd caused by the rods or the triangles, which counteracts the reduction in the skin frictional drag Fd, that may also obtain negative values. More specifically, Fd varies between 8.18  104 (LD) and 103 (NM), and Pd between 4.79  103 (NM) and 5.22  103 (LD) for circular rods, while 4.25  104 (SD) 6 Fd 6 2.87  104 (NM) and 1.01  102 (NM) 6 Pd 6 1.12  102 (LD) for triangular elements. The LES predictions of Fd and Pd are in reasonable agreement with the DNS values Fd = 1.1  103, Pd = 3.9  103 for circular rods and Fd = 6  104, Pd = 4.5  103 for triangles (Leonardi et al., 2003a,b; Orlandi et al., 2006). Based on the above observations, it is evident that the influence exerted by the roughness elements with circular or triangular cross-section on the mean streamwise velocity and friction velocities is reproduced satisfactorily by the subgrid-scale models. In agreement with the present findings for the square roughness element, the largest discrepancies from the DNS results in the cases with circular rods or triangles are also observed in the skin frictional drag Fd and form drag Pd. These discrepancies can be similarly attributed to the lower resolution of the wall contour and the stronger impact of the subgrid-scale fluid motions. For example, the coherent structure model predicts 10% smaller Fd than the DNS value, while it is 25% for the Lagrangian model and

Fig. 10. Mean streamlines for circular rods (left) and triangles (right) with w/k = 1. The flow is from left to right.

(a)

(b)

Fig. 11. Distribution of mean streamwise velocity normalized by ub for circular (a) and triangular roughness elements (b) with w/k = 1. Lines: present LES results, symbols: DNS results of Orlandi et al. (2006).

9% for no model. The coherent structure model produces 25% smaller Pd than the DNS value, while for the Lagrangian model it is 34% and for the no model it is 22%. Similar differences can be observed for the triangular roughness elements, with the exception that the LES predictions of Fd are now closer to the DNS results and, therefore, better as compared with that obtained without subgrid-scale model, which is very encouraging for the use of the LES models. Results from grid refinement test simulations toward a DNS of the present flow configurations with triangles, reveal that Fd = 4.2  104, Pd = 7.9  103 for a mesh with 201  129 grid points in the x-, z-directions, respectively, and Fd = 5.1  104, Pd = 5.98  103 for a mesh with 401  129, indicating that the values of Fd, Pd converge to those of the DNS study as the computational mesh is refined. Fig. 12 shows the LES distributions of Ux,rms, Uy,rms, and Uz,rms based on the resolved grid-scale velocity components as a function of y0 /h for circular rods and triangles with w/k = 1; it also shows the corresponding DNS results of Orlandi et al. (2006) for comparison purposes. Fig. 12(a) and (b) shows that Ux,rms is adequately predicted by all the LES models with respect to the corresponding DNS results for both types of wall disturbances. It is also seen that the LES predictions of Ux,rms are better than that without using any subgrid-scale turbulence model (NM). It is clear that the roughness elements with circular and triangular cross-section have a stronger impact on Ux,rms, which is significantly larger in the region near the

Please cite this article in press as: Dritselis, C.D. Large eddy simulation of turbulent channel flow with transverse roughness elements on one wall. Int. J. Heat Fluid Flow (2014), http://dx.doi.org/10.1016/j.ijheatfluidflow.2014.08.008

C.D. Dritselis / International Journal of Heat and Fluid Flow xxx (2014) xxx–xxx

(a)

(b)

(c)

(d)

(e)

(f)

11

Fig. 12. Distributions of the rms velocity fluctuations for circular (a, c, e) and triangular roughness elements (b, d, f) with w/k = 1: (a and b) streamwise, (c and d) wall-normal, and (e and f) spanwise component. Lines: present LES results, symbols: DNS results of Orlandi et al. (2006).

rough wall, as compared to the case with square bars and w/k = 1. This roughness shape effect, which is pronounced even for small values of w/k, is accurately reproduced by both versions of the dynamic Smagorinsky model and the coherent structure model. Fig. 12(c–f) show the wall-normal distributions of Uy,rms, Uz,rms, respectively, for both circular rods and triangular elements with w/k = 1. It can be seen that Uy,rms is underestimated relative to the DNS distribution mainly in the region of 0.05 6 y0 /h 6 0.5 for

rods, and in a larger region of 0.04 6 y0 /h 6 0.7 for triangles, while a good agreement is observed near the rough wall and the central channel region for both types of roughness elements. Fig. 12(c)–(f) indicate that the level of under-prediction for all the LES models is larger for Uz,rms than Uy,rms, while their cross-examination with Fig. 8(c)–(f) reveals that the values of Uy,rms and Uz,rms yielded by the subgrid-scale models are more under-estimated in the case of square bars. Once again, the DNS trends are well captured by

Please cite this article in press as: Dritselis, C.D. Large eddy simulation of turbulent channel flow with transverse roughness elements on one wall. Int. J. Heat Fluid Flow (2014), http://dx.doi.org/10.1016/j.ijheatfluidflow.2014.08.008

12

C.D. Dritselis / International Journal of Heat and Fluid Flow xxx (2014) xxx–xxx

the various subgrid-scale turbulence models which exhibit similar performances, as indicated by the small differences in the distributions of Ux,rms, Uy,rms, and Uz,rms. Fig. 13(a) shows the ratio of the subgrid-scale eddy viscosity to the molecular viscosity mt/m yielded by the standard and the Lagrangian dynamic models, and the coherent structure model for rods and triangles with w/k = 1. Similar observations to those of the square bars discussed above can also be made. However, the stronger influence of the rods and triangular elements is apparent, as indicated by the larger values of mt/m in the vicinity of the rough wall. It should be noted that the triangles have a stronger impact than the circular rods. Fig. 13(a) also shows that the distributions of mt/m are similar for the standard and the Lagrangian dynamic models, while the coherent structure model yields smaller values than the former two LES models. However, this has no effect on the quality of numerical predictions of the turbulence statistics, in accordance with the findings for the square bars. Fig. 13(b) shows a relatively large contribution of all the models on the total dissipation rate. It should be noted that larger ratios of the subgrid-scale dissipation rate to the total dissipation rate are found in the cases with triangle or circular rods, as compared with the square bars and w/k = 1, denoting the stronger impact of the LES models in the former types of roughness elements.

(a)

4. Conclusions Results from large eddy simulations of turbulent flows in a channel with 2d roughness elements of square, circular and triangular cross-sections transversely placed on its bottom wall were obtained in order to assess the feasibility of several subgrid-scale models. For all flow cases with smooth and rough walls examined here, the effect of the LES models was important as indicated by the subgrid-scale eddy viscosity. In addition, all models contributed from 5% up to 50% of the total dissipation rate over much of the channel, playing thus a key role in the grid-scale turbulence energy budget. It was found that the standard dynamic Smagorinsky model, the Lagrangian dynamic model and the coherent structure model were capable of reproducing the relevant physics and the main effects of the rough surface on the overlying turbulent flow, exhibiting similar performances. The mean streamwise velocity, the rms velocity fluctuations, the skin frictional drag and the form drag were predicted reasonably by all the LES models with respect to the DNS results, based on a coarser grid resolution. In addition, an improvement associated with the use of the subgrid-scale turbulence models as compared with the case when no model is adopted was observed in several statistic results. This improvement was rather evident in the case with square roughness elements and w/k = 7, in which the influence of the subgrid-scale models was the strongest, while the differences between the LES and DNS results were generally smaller than those when no model was used. Consequently, LES coupled with an immersed boundary method is a suitable, efficient and accurate option for the numerical investigation of turbulent channel flows with riblets of various shapes. Acknowledgments The author would like to thank Professor Emeritus Nicholas Vlachos for his insightful remarks. Financial support from the Association EURATOM–Hellenic Republic is gratefully acknowledged. The content of this paper is the sole responsibility of its author and it does not necessarily represent the views of the European Commission or its services. Appendix A. Effect of grid resolution

(b)

Fig. 13. Distributions of the ratio of the subgrid-scale eddy viscosity mt and the molecular viscosity m (a) and of the ratio of the subgrid-scale dissipation rate to the total (subgrid-scale plus molecular) dissipation rate (b) for the standard dynamic Smagorinsky model, the Lagrangian dynamic Smagorinsky model, and the coherent structure model in the cases with circular and triangular roughness elements and w/k = 1.

The choice of the grid resolution to conduct accurate LES requires a careful consideration of several numerical and modeling issues, the physics of the flow under investigation and the complexity of the flow geometry. In this Appendix, the effect of grid resolution on the LES results is shown. The grid resolution adopted in the present LES study is much coarser than those used in the corresponding DNS studies of turbulent channel flow with roughness elements on one wall. Based on this grid resolution, a different and significant LES effect is modeled in the calculations by the subgrid scale models, and the LES predictions exhibit differences from those obtained without using any subgrid scale model. Figs. 14 and 15 show results in a representative case with square roughness elements and w/k = 3, produced using the standard dynamic Smagorinsky model on grids appropriately coarsened in the streamwise, wall-normal or spanwise direction. Similar results were obtained by the other LES models. In order to denote the LES impact on the numerical predictions, the results obtained without any model are also shown in Figs. 14 and 15. For all cases, the channel size was 8h  (0.2 + 2)h  ph. The results in Figs. 14(a) and 15(a) show that the refinement in the streamwise direction is less important; in fact, the profiles using 81  (30 + 64)  65 points do not differ largely from those using 201  (30 + 64)  65 points. This can potentially explain

Please cite this article in press as: Dritselis, C.D. Large eddy simulation of turbulent channel flow with transverse roughness elements on one wall. Int. J. Heat Fluid Flow (2014), http://dx.doi.org/10.1016/j.ijheatfluidflow.2014.08.008

C.D. Dritselis / International Journal of Heat and Fluid Flow xxx (2014) xxx–xxx

(a)

(b)

(c)

(d)

13

Fig. 14. Distributions of Ux,rms for square roughness elements with w/k = 3 predicted using a grid coarsened in the streamwise (a), wall-normal (b and c), or spanwise direction (d). Lines: present LES results based on standard dynamic Smagorinsky model, symbols: DNS results of Leonardi et al. (2003b).

the significant subgrid-scale effect that is modeled when using the initial grid resolution, despite the apparent large number of points in the streamwise direction. It is recalled here, that the use of 201 points in x is an outcome of the decision to discretize each roughness element with 5 points in the x direction, as opposed to the 10 points used in the DNS studies. The grid with 81 points uses only two points in the streamwise direction for each element. Obviously, more complex wall contours cannot be reproduced correctly by such grid resolution. This has a strong effect on the predictions of the skin friction drag Fd and the pressure form drag Pd. For example, the 81  (30 + 64)  65 simulations predict 2.3  104 6 Fd 6 3.7  104 and 4.9  103 6 Pd 6 5.7  103, which are lower than 3.8  104 6 Fd 6 4.9  104 and 6.7  104 6 Pd 6 7.7  104 obtained in the 201  (30 + 64)  65 simulations. The latter values are closer to the DNS predictions. The results from the 201  (30 + 32)  65 simulations, which use fewer points in the outer flow region in the wall-normal direction than in the initial grid, are in generally good agreement with the results obtained by the initial grid, as shown in Figs. 14(b) and 15(b). This is due to the non-uniform grid, which allows 32 points to be sufficient to locate the first point at the same location close to the wall as in the initial grid arrangement. However, a more uniform grid spacing, or locating the first point further away from the wall, leads to significant discrepancies from the DNS results. Small changes are also observed by using fewer points to discretize the

roughness element in the y direction as shown in Figs. 14(c) and 15(c). The results produced based on 201  (15 + 64)  65 points are close to those of the simulations with 201  (30 + 64)  65 points. The values of Fd and Pd produced by the 201  (30 + 32)  65 and 201  (15 + 64)  65 simulations are relatively smaller than those found for the 201  (30 + 64)  65 simulations. For the aforementioned grid arrangements, all the LES models produce better predictions as compared with those without any model. However, the flow reattachment on the bottom wall of the cavity for w/k = 7 was barely reproduced when using 15 points to discretize each roughness element in the y direction. A poor performance is indicated for all LES models in the 201  (30 + 64)  33 simulations. Figs. 14(d) and 15(d) show that the refinement in z has a greater influence than that in x. A better resolution in the spanwise direction is generally needed to account properly for the changes of the flow structures. For example, the values of Pd predicted are 4.8  103 (SD), 5.7  103 (LD), and 5.3  103 (CS) that are smaller than those obtained without using a subgrid-scale model Pd = 6  103, which in turn is closer to the DNS value. At the initial grid of 201  (30 + 64)  65, it is found that Pd = 7.1  103 (SD), 7.7  103 (LD), 6.8  103 (CS), and 6.6  103 (NM), indicating better predictive performances of all the subgrid-scale turbulence models. The above analysis supports the results and conclusions based on the initial grid size; however, larger discrepancies are observed

Please cite this article in press as: Dritselis, C.D. Large eddy simulation of turbulent channel flow with transverse roughness elements on one wall. Int. J. Heat Fluid Flow (2014), http://dx.doi.org/10.1016/j.ijheatfluidflow.2014.08.008

14

C.D. Dritselis / International Journal of Heat and Fluid Flow xxx (2014) xxx–xxx

(a)

(b)

(c)

(d)

Fig. 15. Distributions of Uy,rms for square roughness elements with w/k = 3 predicted using a grid coarsened in the streamwise (a), wall-normal (b and c), or spanwise direction (d). Lines: present LES results based on standard dynamic Smagorinsky model, symbols: DNS results of Leonardi et al. (2003b).

with respect to the DNS results, which confirm the need of using the initial grid in order to obtain accurate LES results in the cases of turbulent channel flows with rough walls studied here.

References Ashrafian, A., Andersson, H.I., Manhart, M., 2004. DNS of turbulent flow in a rodroughened channel. Int. J. Heat Fluid Flow 25, 373–383. Burattini, P., Leonardi, S., Orlandi, P., Antonia, R.A., 2008. Comparison between experiments and direct numerical simulations in a channel flow with roughness on one wall. J. Fluid Mech. 600, 403–426. Carati, D., Wray, A., Cabot, W., 1996. Ensemble averaged dynamic modeling. In: Proc. Summer Program, 1994, Cent. Turbl. Res., Stanford Univ., CA V, pp. 237– 38. Ciofalo, M., Collins, M.W., 1992. Large-eddy simulation of turbulent flow and heat transfer in plane and rib-roughened channels. Int. J. Numer. Methods Fluids 15, 453–489. Cui, J., Patel, V.C., Lin, C.-L., 2003. Large-eddy simulation of turbulent flow in a channel with rib roughness. Int. J. Heat Fluid Flow 24, 372–388. Dritselis, C.D., Vlachos, N.S., 2011. Large eddy simulation of gas-particle turbulent channel flow with momentum exchange between the phases. Int. J. Multiphase Flow 37, 706–721. Dritselis, C.D., Sarris, I.E., Fidaros, D.K., Vlachos, N.S., 2011. Transport and deposition of neutral particles in magnetohydrodynamic turbulent channel flows at low magnetic Reynolds numbers. Int. J. Heat Fluid Flow 32, 365–377. Fadlun, E.A., Verzicco, R., Orlandi, P., Mohd-Yusof, J., 2000. Combined immersed boundary finite difference methods for three-dimensional complex flow simulations. J. Comput. Phys. 161, 35–60. Fornberg, B., 1988. Steady viscous flow past a sphere at high Reynolds numbers. J. Fluid Mech. 190, 471–489.

Fröhlich, J., Mellen, C.P., Rodi, W., Temmerman, L., Leschziner, A., 2005. Highly resolved large-eddy simulation of separated flow in a channel with streamwise periodic constrictions. J. Fluid Mech. 526, 19–66. Furuya, Y., Miyata, M., Fujita, H., 1976. Turbulent boundary layer and flow resistance on plates roughened by wires. J. Fluids Eng. 98, 635–644. Germano, M., Piomelli, U., Moin, P., Cabot, W.H., 1991. A dynamic subgrid-scale eddy viscosity model. Phys. Fluids A 3, 1760–1765. Ghosal, S., Lund, T.S., Moin, P., Akselvoll, K., 1995. A dynamic localization model for large-eddy simulation of turbulent flows. J. Fluid Mech. 286, 229–255. Gullbrand, J., 2004. Dynamic modeling in large-eddy simulation of turbulent channel flow. Investigation of two-dimensional versus three-dimensional test filtering. Int. J. Numer. Methods Heat Fluid Flow 14 (4), 467–492. Hanjalic, K., Launder, B.E., 1972. Fully developed asymmetric flow in plane channel. J. Fluid Mech. 51, 301–335. Iaccarino, G., Verzicco, R., 2003. Immersed boundary technique for turbulent flow simulations. Appl. Mech. Rev. 56, 331–347. Iacono, G.L., Tucker, P.G., Reynolds, A.M., 2005. Predictions for particle deposition from LES of ribbed channel flow. Int. J. Heat Fluid Flow 26, 558–568. Khan, M.A.I., Luo, X.Y., Nicolleau, F.C.G.A., Tucker, P.G., Iacono, G.L., 2010. Effects of LES sub-grid flow structure on particle deposition in a plane channel with ribbed wall. Int. J. Numer. Method Biomed. Eng. 26, 999–1015. Kim, J.H., Jansen, K.E., Jensen, M.K., 2004. Analysis of heat transfer characteristics in internally finned tubes. Numer. Heat Transfer A 46, 1–24. Kobayashi, H., 2005. The subgrid-scale models based on coherent structures for rotating homogeneous turbulence and turbulent channel flow. Phys. Fluids 17, 045104. Kobayashi, H., 2006. Large eddy simulation of magnetohydrodynamic turbulent channel flows with local subgrid-scale model based on coherent structures. Phys. Fluids 18, 045107. Kobayashi, H., 2008. Large eddy simulation of magnetohydrodynamic turbulent duct flows. Phys. Fluids 20, 015102. Kobayashi, H., Ham, F., Wu, X., 2008. Application of a local SGS model based on coherent structures to complex geometries. Int. J. Heat Fluid Flow 29, 640–653.

Please cite this article in press as: Dritselis, C.D. Large eddy simulation of turbulent channel flow with transverse roughness elements on one wall. Int. J. Heat Fluid Flow (2014), http://dx.doi.org/10.1016/j.ijheatfluidflow.2014.08.008

C.D. Dritselis / International Journal of Heat and Fluid Flow xxx (2014) xxx–xxx Leonardi, S., Castro, I.P., 2010. Channel flow over large cube roughness: a direct numerical simulation study. J. Fluid Mech. 651, 519–539. Leonardi, S., Mostarda, E., Orlandi, P., Antonia, R.A., 2003a. Turbulent channel flow with rod and square bar roughness on one wall. 16th AIMETA Congress of Theoretical and Applied Mechanics, 1–10. Leonardi, S., Orlandi, P., Smalley, R.J., Djenidi, L., Antonia, R.A., 2003b. Direct numerical simulations of turbulent channel flow with transverse square bars on the wall. J. Fluid Mech. 491, 229–238. Leonardi, S., Orlandi, P., Djenidi, L., Antonia, R.A., 2004. Structure of turbulent channel flow with square bars on one wall. Int. J. Heat Fluid Flow 25, 384–392. Leonardi, S., Tessicini, F., Orlandi, P., Antonia, R., 2006. Direct numerical and largeeddy simulations of turbulent flows over rough surfaces. AIAA J. 44 (11), 2482– 2487. Leonardi, S., Orlandi, P., Antonia, R.A., 2007. Properties of d- and k-type in a turbulent channel flow. Phys. Fluids 19, 125101. Lilly, D.K., 1992. A proposed modification of the Germano subgrid-scale closure method. Phys. Fluids A 4, 633–635. Luo, D.D., Leung, C.W., Chan, T.L., Wong, W.O., 2005. Flow and forced convection characteristics of turbulent flow through parallel plates with periodic transverse ribs. Numer. Heat Transfer A 48, 43–58. Meneveau, C., Lund, T.S., Cabot, W.H., 1996. A Lagrangian dynamic subgrid-scale model of turbulence. J. Fluid Mech. 319, 353–385. Mittal, R., Iaccarino, G., 2005. Immersed boundary methods. Annu. Rev. Fluid Mech. 37, 239–261. Nagano, Y., Hattori, H., Houra, T., 2004. DNS of velocity and thermal fields in turbulent channel flow with transverse-rib roughness. Int. J. Heat Fluid Flow 25, 393–403. Najjar, F.M., Tafti, D.K., 1996. Study of discrete test filters and finite difference approximations for the dynamic subgrid-scale stress model. Phys. Fluids 8, 1076–1088. Nicoud, F., Ducros, F., 1999. Subgrid-scale stress modelling based on the square of the velocity gradient tensor. Flow, Turbulence Combust. 62, 183–200.

15

Orlandi, P., 2000. Fluid Flow Phenomena: A Numerical Toolkit. Kluwer, Dordrecht. Orlandi, P., Leonardi, S., 2004. Passive scalar in a turbulent channel flow with wall velocity disturbances. Flow, Turbul. Combust. 72, 181–197. Orlandi, P., Leonardi, S., 2006. DNS of turbulent channel flow with two- and threedimensional roughness. J. Turbul. 7, 1–22. Orlandi, P., Leonardi, S., Antonia, R.A., 2006. Turbulent channel flow with either transverse or longitudinal elements on one wall. J. Fluid Mech. 561, 279–305. Piomelli, U., Balaras, E., 2002. Wall-layer models for large-eddy simulation. Annu. Rev. Fluid Mech. 34, 349–374. Saha, A.K., Acharya, S., 2004. Unsteady simulation of turbulent flow and heat transfer in a channel with periodic array of cubic pin-fins. Numer. Heat Transfer A 46, 731–763. Saidi, A., Sunden, B., 2000. Numerical simulation of turbulent heat transfer in square ribbed ducts. Numer. Heat Transfer A 38, 67–88. Shimomura, Y., 1994. Subgrid-scale algebraic stress model of turbulence. J. Phys. Soc. Jpn. 63, 5–9. Smagorinsky, J., 1963. General circulation experiments with the primitive equations. I. The basic experiment. Month. Weather Rev. 91, 99–152. Vijiapuparu, S., Cui, J., 2007. Simulation of turbulent flow in ribbed pipe using large eddy simulation. Numer. Heat Transfer, Part A 51, 1137–1165. Vreman, A.W., 2004. An eddy-viscosity subgrid-scale model for turbulent shear flow: algebraic theory and applications. Phys. Fluids 16, 3670–3681. Winckelmans, G.S., Jeanmart, H., Carati, D., 2002. On the comparison of turbulence intensities from large-eddy simulation with those from experiment or direct numerical simulation. Phys. Fluids 14, 1809–1811. Yoshizawa, A., Kobayashi, K., Kobayashi, T., Taniguchi, N., 2000. A nonequilibrium fixed-parameter subgrid-scale model obeying the near-wall asymptotic constraint. Phys. Fluids 12, 2338–2344. Yucek, N., Dinler, N., 2006. Numerical study of laminar and turbulent flow through a pipe with fins attached. Numer. Heat Transfer A 49, 195–214. Zang, Y., Street, R.L., Koseff, J.R., 1993. A dynamic mixed subgrid-scale model and its application to turbulent recirculating flows. Phys. Fluids A 5, 3186–3196.

Please cite this article in press as: Dritselis, C.D. Large eddy simulation of turbulent channel flow with transverse roughness elements on one wall. Int. J. Heat Fluid Flow (2014), http://dx.doi.org/10.1016/j.ijheatfluidflow.2014.08.008