Large eddy simulation of flow over a wall-mounted cube : Comparison of different semi dynamic subgrid scale models

ABSTRACT In this paper the ability of different semi dynamic subgrid scale models for large eddy simulation was studied in a challenging test case. The semi dynamic subgrid scale models were examined in this investigation is Selective Structure model, Coherent structure model, Wall Adaptive Large Eddy model. The test case is a simulation of flow over a wall-mounted cube in a channel. The results of these models were compared to structure function model, dynamic models and experimental data at Reynolds number 40000. Results show that these semi dynamic models could improve the ability of numerical simulation in comparison with other models which use a constant coefficient for simulation of subgrid scale viscosity. In addition, these models don’t have the instability problems of dynamic models.


INTRODUCTION
In Large-Eddy Simulation, the flow is divided to resolved scales and subgrid scales.The resolved scales are simulated directly and subgrid scales are modeled.There are several methods for simulation subgrid scales such as spectral models, physical space models, deconvolution models and etc. Physical space models, particularly, viscose models are more practical than other models because of their usage in engineering problems.These models are based on this hypothesis "The action of the subgrid scales on the resolved scales is essentially an energetic action, so that the balance of the energy transfers alone between the two scale ranges is sufficient to describe the action of the subgrid scales.The energy transfer mechanism from the resolved to the subgrid scales is analogous to the molecular mechanisms represented by the diffusion term, in which the viscosity appears" [1].The Smagorinsky model is one of famous models in this category.In this model the subgrid scale viscosity is defined as: Results show that C s is not constant in different problems.In addition the constant coefficient makes the subgrid scale viscosity acts on whole of domain which is not according to reality, for instance v sgs should be zero in vicinity of wall or in laminar regions.Different numerical simulations showed that the ability of simulation can improve by adapting subgrid scale models to local state of flow.Many new models have been introduced for this purpose in recent years.Dynamic procedures for computing subgrid scale models constant are the more famous than others.Dynamic models usually are based on Germano [2] identity; in this procedure constant coefficient is computed locally in space and time by reducing error and a second filter.Germano dynamic procedure does not change the prior form of model; therefore it can be used for different models.Lagrangian dynamic procedure [3], the constrained localized dynamic procedure [4], approximate localized dynamic procedure [5] are some dynamic models which are based on Germano procedure.It should be noted that there are some dynamic procedure that are not based on Germano procedure such as multilevel dynamic procedure by Terracol and Sagaut [6].Numerical results show that the dynamic models usually are unstable and need numerical contrivances such as clipping or averaging in the directions of statistical homogeneity.Consequently, the ability of dynamic models can be restricted in complex geometries.
Many researchers have been involving to introduce easier model to adapt v sgs to local state of the flow, for example Selective Structure model (SSF) [7], Coherent structure model (CSM) [8], Wall Adaptive Large Eddy model [9].
The main purpose of this paper is to investigate the ability of different subgrid scale models which are semi dynamic (the constant coefficient is computed locally in space and time) in simulation of complex phenomena such as horseshoe vortices, flow separation, arc vorticity and etc. for turbulent flow over a wall-mounted cube confined in a channel.It should be mentioned that these phenomena make this test case to a challenging test case for numerical simulation.There is a vast amount of literature about experiments undertaken for this geometry [10][11][12][13][14]; among them being the comprehensive work of Martinuzzi and Tropea [13] at Reynolds number equal 40000.

MATHEMATICAL FORMULATION
As mentioned in previous section, Large Eddy Simulation is based on this idea that is not necessary to model the total of cascade of energy (transform of energy from large scale to small scale) such as RANS models, RSM models and etc.In large eddy simulation, the resolved scales are simulated directly and subgrid scales (SGS) are modeled.For this purpose a filter has been introduced. ( Where f and f ′ are the resolved and subgrid scale components, respectively.G ∆ x (xi) is the filtering function.The famous one of filtering function which is usually used in large eddy simulation is top hat filter.By applied this filter on Navier-Stokes equations, these equations can be represented as follow: (2) Where ut , P -and τ ij are the velocity components, pressure and stress tensor of subgrid scales, respectively.In above equations the stress tensor of subgrid scales should be modeled.In Functional modeling for large eddy simulation, one is supposed that the effect of subgrid scales on resolved scales is same as viscous effect.Therefore the stress tensor of subgrid scales is computed as follow: ( Where sgs is the viscosity of subgrid scales and modeled by different subgrid scale models.

SELECTIVE STRUCTURE FUNCTION MODEL
M′etais and Lesieur [15] introduced Structure Function Model (SF) by transferring their constant effective viscosity model into the physical space.They supposed that the energy at cutoff can be estimated by the second-order velocity structure function.In this model, the eddy viscosity is evaluated according to; (5) Where ∆c = (∆x 1 × ∆x 2 × ∆x 3 ) 1 ⁄ 3 is the geometric mean of the meshes in the three spatial directions.C k is Kolmogrov constant and F 2 is the local structure function constructed with the filtered velocity field u -(x, t): (6) Results show that this model is very dissipative.In order to improve the prediction of SF and remove its shortcomings, Selective Structure Function (SSF) was introduced and developed by David [16].SSF switches off the subgrid scales eddy viscosity when the flow is not sufficiently three-dimensional.The criterion is the angle (α) between the vorticity at a given grid point and the average vorticity at the six closest neighboring points.When α less than 20, the flow is supposed laminar and calculating eddy viscosity is not necessary.
Where Φ 20°( x, t) is the indicating function based on the value of (ϰ): Suksangpanomrung et al. [17] used a smoothly varying function rather than an abrupt cutoff, Φ′ 20°( x, t) instead of Φ 20°( x, t), which is evaluated using a smoothly varying function defined as: Ackermann and Metais [18] proposed Modified Selective Structure Function.The used a function for critical angle based on the ratio of the cut-off wave number and the spectrum peaks.The results of MSSF were found to be very close to those obtained with Smagorinsky's dynamic model of Germano.It should be mentioned that the utilization of MSSF is difficult for complex geometries.Farhadi and Rahnama [19][20] studied different SSF models for turbulent flow over obstacles.They represented that the SSF model with smooth Function [17] can achieve better results than others.Therefore, this model was used in this study.

WALL-ADAPTIVE LOCAL EDDY VISCOSITY (WALE)
Nicoud and Ducros [9] stated that a better model can be obtained by using both strain and rotation rates in computing the subgrid scale viscosity.However, this requirement is included in some models such as SF but another major problem in LES subgrid scale modeling is behavior in vicinity of the solid walls.Most subgrid viscosity model (SM, SF and etc) do not exhibit the correct behavior in the vicinity of solid walls.Results show that the behavior of sgs should be O(y +3 ).SM gives non-zero value in vicinity of wall and SF behaves O(1) in this region.A common way for solving this problem is using damping function.But these damping function usually don't exhibit correct behavior in vicinity of wall, for example Van Driest [21] damping function which is used widely for this purpose varies O(y 2 ) instead of O(y +3 ).Moreover these functions usually require to the distance to the wall and the skin friction as input parameters, consequently their application are limited in complex geometries.Nicoud and Ducros [9] tried to find a combination of resolved velocity spatial derivatives that use both of strain and rotation rates and exhibits the correct behavior in vicinity of wall.They introduced Wall Adaptive Large Eddy model (WALE).In this model, the viscosity of subgrid scale is defined as [9]: Where S ij is the strain-rate tensor: (11) Is the traceless symmetric part of the tensor, and C w = 0.55 − 0.60.

COHERENT STRUCTURE MODEL (CSM)
Kobayashi [8] proposed an adaptive model based on the coherent structures.The model parameter is composed of a fixed model-parameter and a coherent structure function.Consequently the SGS model coefficient is calculated locally.Moreover Kobayashi [22] applied this model to complex geometry such as diffuser, jet and hill and showed that the result of CSM is same as dynamic model.The CSM is an eddy viscosity model that the coefficient of model is not constant and calculated locally as follow [8]: (12) Where: ) Where C CSM is a fixed model constant, F CS is the coherent structure function defined as the second invariant normalized by the magnitude of a velocity gradient tensor E, F Ω is the energy-decay suppression function and W ij is the vorticity tensor in a resolved flow field.Moreover, F CS and F Ω have definite upper and lower limits: (17) Kobayashi [22] mentioned that the C CSM has smaller variance than the model parameter of the Dynamic Smagorinsky Model without averaging and the numerical simulation with the CSM is more stable.

COMPUTATIONAL METHOD AND BOUNDARY CONDITION
In this investigation, the computational domain, grid spacing and the average time were selected based on Farhadi et al. [19] (Figure . 1).The minimum grid spacing used for this test case was 0.03 H in x, y and z directions with grid expansion ratio of 1.06.The number of grid points is 113 × 51 × 100.The average time in the simulation was 150 H/U inlet where H is the cube height and U inlet is the velocity at the inlet.
The governing equations presented in the preceding section were discretized using a finite volume method with the staggered grid.The convective terms were discretized using QUICK scheme.As described by Farhadi and Rahnama.[19 and 20], However the QUICK scheme has some deficiencies, such as large numerical dissipation, as compared with the Central Difference (CD) scheme, but QUICK scheme could be better in simulation of turbulent flow over obstacles.In addition instability problem decreases by using QUICK scheme instead of central difference scheme.A third order Runge-Kutta algorithm is used for the time integration in conjunction with the classical correction method at each sub-step.The semi-implicit fractional step method used which provides an approach that does not use pressure in the predictor step as in the pressure corrector method (such as the well-known SIMPLE family of algorithms).The no-slip boundary condition was used for all the walls and surfaces of obstacles.Outlet boundary condition was of convective type given by: (18) Where U c is mean velocity or bulk velocity of inflow.At inlet a fully developed turbulent flow was considered.

RESULTS
Large eddy simulation was employed for the complex problem of flow over a wallmounted cube.The study was conducted for three adaptive and a non-adaptive subgrid scale models.The results were compared with Martinuzzi and Tropea [13] experimental data at Re = 40000.The time and space averaged streamwise velocity at plane z = 0 is presented for different models in Figure 2. By computing the subgrid scale viscosity locally, discrepancy between numerical results and experimental data decrease.In fact the semi dynamic models could predict the horse shoe vortex and separation regions around and behind the wall mounted cube better than SF model.Figures 3 and 4 show the streamlines in vicinity of channel floor and at plane z = 0 for different models used in the present simulations and other experiments, respectively.As shown in these Figures, the SSF and CSM models predicted largest and smallest recirculation zone behind of obstacle, respectively.Moreover, all semi dynamic models show a horseshoe vortex around the cube and the separation regions on the roof, lateral sides and behind of cube.
The main point about horseshoe vortex is its converging-diverging behavior which is shown in the simulation of these semi dynamic models correctly.
The flow pattern for upstream of obstacle was very similar for different semi-dynamic models (Figure 4).As discussed by Martinuzzi and Tropea [13], two recirculation regions exist upstream of the cube.All simulations predicted the primary recirculation.The size of this recirculation zone is very small for these semi dynamic models.The verity of velocity at the vicinity of channel floor among the channel reveals this recirculation zone (Figure 5).It is noted that the saddle point in front of obstacle was calculated same as experimental data by these semi dynamic models (X f 2 ≅ 0.6h).No author has detected the secondary recirculation zone, which is very small and close to the front side of cube, through numerical computation except DNS [23].
The effect of discretizing scheme for convective term in momentum equations was studied for a low Reynolds number too.The streamlines for Re = 3200 represented in Figure (6).Both QUICK and Central schemes could simulate the general pattern of flow over obstacle same as experimental data.However, QUICK scheme displayed some discrepancy in the prediction of the wake center behind the obstacle and recirculation zone above the obstacle.
Generally, by computing model's coefficient locally in space and times better results were achieved.For instance, the SF model estimated is far distance for the center of recirculation zone behind of obstacle, but the results of other dynamic and these semi dynamic models are Figure 4 Comparison of time-average streamlines plots at plane z = 0 obtained from different models and experimental results [13].

Modeling Method
Grid point %Er X F X R %Er Ref. [13] 1 compares various lengths of separation regions defined in Figure 4 and 5.
It is observed that none of the used models could predict both upstream and downstream recirculation lengths correctly.Time-averaged resolved velocities and turbulent statistics are computed and compared with the experiments in Figure 7.The locations are selected from x = H = 0:5 to 4.0 at z = 0.All presented variables in this Figure are in non-dimensional form.The obtained Results for locations upstream of the cube were not shown, because no significant flow feature exists in that region except for a small recirculation region.
The computed mean velocity for this region shows reasonable correspondence with the experiment, in spite of the differences that exist in the upstream recirculation length, as mentioned in the preceding discussion.All semi-dynamic models results follows the trend of experimental variable correctly; however there are some discrepancies with the experimental data and these discrepancies become less along the channel.Variables distribution near the top surface of the cube shows discrepancies with the experimental data for most of the region between the top side of the cube and the channel wall.The main reason for such a prediction could be a result of low grid resolution and the QUICK scheme used for the discretization of convective terms.SSF model shows better agreement with experimental data especially for cross and streamwise velocity and maximum point of turbulent statistics.The results of WALE and CSM are very similar to each other.
Verity of turbulent kinetic energy at plane z = 0 is illustrated in Figure 8.The maximum points in these graphs represent the center of recirculation zone behind the obstacle.As shown in this figure the SSF model have minimum kinetic energy in this region.In fact this model dissipated the momentum of fluid more than other models in this area.different semi dynamic subgrid scale models x/H = 1.5 x/H = 0.5 x/H = 4.0 1 0.5  x/H = 0.5 x/H = 1.5 x/H = 0.For comparison the behavior of different models to denotation the turbulent regions in the domain, the instant contour of C (Eq. 12), and selective sensor (Φ) for CSM, WALE and SSF models was presented in Figure (9), respectively.As shown in this figure, the regions should be modeled by CSM and WALE similar to each other approximately.This is the reason that the results of these models are similar to each other.
The sensor of SSF model affected on more regions in comparison to other models and most part of domain were modeled for SSF simulation.Second conclusion which can be achieved from this figure is that CSM and WALE models are more sensitive to grid resolution than SSF model because most regions in these models were simulated directly and no model effect on this parts.Our computation for simulation of a fully develop turbulent flow in a channel showed this one.The WALE and CSM simulations diverged for a coarse grid (∆z + > 40) but SSF simulation converged for this coarse grid.
It should be mentioned that the time step computed dynamically in our simulations.Results show that CSM and WALE models need to small time step rater that SSF model and the simulation time for these models is more than SSF model.

CONCLUSION
In this investigation the ability of Selective Structure model, Coherent structure model, Wall Adaptive Large Eddy model for large eddy simulation of flow over a wall-mounted cube confined in a channel was studied and compared to other numerical results and experimental data [13,14].Results show that by computing model's coefficient locally in space and times better results were achieved.In addition, it was observed that the results of CSM and WALE models were very similar to each other and these models were more sensitive to grid resolution than SSF model because most regions in these models were simulated directly.As a result, the CSM and WALE need to small time step rater that SSF model and the simulation time for these models is more than SSF model.

Figure 1
Figure 1 a) Geometry of problem for flow over wall-mounted cube.b)general pattern of flow over wall-mounted cube.

Figure 2 Figure 3
Figure2The time and space averaged streamwise velocity at plane z = 0.

Figure 6 Figure 5
Figure 6 Time-averaged computations of streamlines plots at plane z = 0 (a) SSF with CD (b) SSF with QUICK and vicinity of channel floor (c) SSF with CD and (d) SSF with QUICK.

Figure 7
Figure 7 Time-averaged resolved velocities and turbulent statistics from x = H = 0:5 to 4.0 at z = 0.

Figure 8
Figure 8 Verity of average turbulent kinetic energy at plane z = 0 for WALE, CSM and SSF.

Figure 9
Figure 9 Contour of turbulent region which defined by different models.The blue region is laminar.