Capability Assessment of Five Different RANS-Based Turbulence Models to Simulate the Various Regions of Slot Turbulent Impingement Jet Flow

In this paper, impingement of a turbulent rectangular flow to a fixed wall is investigated. The jet flows from bottom-to-top and the output jet Reynolds is 16000. The nozzle-to-plate distance is equal to 10 ( H/e =10). Five turbulence models, including k-ε, RNG k-ε, k-ω SST, RSM and v 2 f model have been used for two-dimensional numerical simulation of the turbulent flow. Because of the complexities of the impingement flow, such as curved streamlines, flow separation, normal strains and sudden deceleration in different areas, different turbulence models are proposed to simulate different regions of the flow. To investigate the capability of these turbulence models in simulating different regions of the impinging jet, the mean flow velocity field and turbulent kinetic energy are extracted and compared with the experimental data of a two-dimensional particle image velocimetry (PIV). The calculated error of these five turbulence models was presented for the various flow regions, while it have not been clearly investigated earlier. Results indicate the highest conformity of the v 2 f model with the experimental data at the jet centerline. However, this model does not predict well the flow at the shear layer and wall-jet areas. RSM Gibson and Lander model has the highest conformity with the experimental data in these regions.


INTRODUCTION
Impinging of a jet flow to a wall is an effective and efficient way for the heat and mass transfer in the industrial applications. A high speed jet flow impinges a plate directly and increases the heat and mass transfer to a large extent. Many applications have been mentioned for this phenomenon in the literature. For example it is possible to mention metal industry, cooling of electronic components, defogging by heating optical surfaces, cooling turbine components and cooling sensitive parts in machining [1][2][3]. In this study, the impinging jet flow is divided into four regions. They are potential core, developing, impact and wall-jet regions (see Fig. 1). The characteristics of the mentioned regions are comprehensively described by Zuckerman and Lior [4]. If the distance between the nozzle and plate is high, the jet behaves like a free submerged jet when it leaves the nozzle. In this case, the velocity variation in the jet edges creates a shear layer and consequently, momentum transfer happens. Therefore, more flow is drawn from the surrounding towards the jet and the mass flow increases. Hence, the energy dissipates and the velocity profiles become wider and its maximum reduces gradually. Due to the existence of different regimes as well as complications in the impingement of the turbulent jet to the walls, this problem has always used to test the validity and accuracy of the turbulence models. Since this phenomenon is widely used in simulating industrial processes, and the simulation of industrial processes by LES model is not cost effective, so investigating the strengths and weaknesses of the widely used RANS-based models is necessary. In the following paragraphs some more investigation found in the literature are presented.
Heyerichs and Pollard [5] simulated the turbulence jet impinging to the fixed wall using three different wall functions with five different damping functions. They concluded that the previously developed damping functions, which were dependent to distance from the wall (y+), are sufficient to simulate the parallel and transient flows from the wall. However, they are not accurate for the impact area. Craft et al. [6] compared the numerical solution data of their model with experimental data of an impinging two dimensional jet to a wall. They observed that velocity at the centerline is estimated four times bigger than that of experimental data. Kubacki and Dick [7] simulated the impingement of a rectangular jet to the wall using k-ω model and three different hybrid models of k-ω/LES. They investigated the shear stress of the wall and the heat transfer rate and reported a large difference between k-ω model and experimental data. Recently, Ortega-Casanova and Granados [6] studied the heat transfer between a circular air jet and a curved plate. They used k-ε, k-ω SST and RSM model to simulate flow at different Reynolds numbers (between 7000 and 19000). They also examined the influence of the distance ratio on Reynolds number by varying the distance between the nozzle and wall to find the desired Reynolds number and distance ratio for optimum heat transfer. However the flow structure is not considered at all. Zuckerman and Lior [4], compared different investigation results conducted in the field of simulating impinging jet in their review article.
Although considerable research have been conducted in the field of simulating jet impingement using RANS models, the majority of them were conducted to heat transfer of the circular jets and short impingement distances (Nozzle to plate distance). While, in the most of studies, the working fluid is air, the novelty of this work is the use of a fine space resolution of 2d-PIV measurements in water instead of air, in order to well check the validity of the chosen RANS-based turbulent models. In addition, in the most of previous works, characteristics of flow such as mean velocity field and turbulent kinetic energy at different regions of the jet impinging flow have not been clearly investigated.
In this study, the impingement of a rectangular turbulent jet to a wall is investigated both experimentally and numerically using five turbulence models k-ε, RNG k-ε, k-ω SST, RSM and v 2 f for Reynolds number of 16,000. The distance of nozzle from the impingement wall is 10 times the nozzle width. Particle image velocimetry (PIV) data of Koched et al. [8] are extracted to evaluate the numerical simulation results.

GOVERNING EQUATIONS
To simulate the impingement of a turbulent isothermal jet to a wall, the time averaged continuity and momentum equations are used. Equations (1) and (2) represent the dimensionless averaged time averaged continuity and momentum equations for an incompressible fluid.
where i u and P are the time average values of velocity and pressure, respectively, and ij  is the Reynolds stress tensor. These equations are in dimensionless form using the output jet velocity and nozzle width. Reynolds stress tensor is given by equation (3).

 
Many models have been proposed for modeling Reynolds stress term. In this study, five RANSbased turbulence models are used, which are illustrated with a brief explanation in the following sections.

2-1 k-ε standard model
The first two-equation model used in this study is k-ε standard model. This is the most widely used and economical model among the different turbulence models. The k-ε standard model is a semiempirical model, solving the transfer equations for the turbulent kinetic energy k and its dissipation rate ε, respectively.
where σk and σε are kinetic turbulent Prandtl number and dissipation turbulent Prandtl number and are equal to 1 and 1.3, respectively. The model constants are Cμ=0.09, C1ε=1.44 and C2ε=1.92. Also, ρ is the fluid density and P k is the production rate of turbulent kinetic energy. It is defined as follows.
where Sij is strain rate tensor and μt is the dynamic turbulent viscosity.
The RNG k-ε model includes a modification to the transport equation for ε stemming from renormalization group theory. The formulation of RNG k-ε model is consisted of the transfer equations for the turbulent kinetic energy k and dissipation rate ε, similar to k-ε model formulation (equations (4) and (5)), substituting C2ε with C * 2ε illustrated in equation (9).
where η=Sk/ε and the turbulent viscosity being calculated in the same manner as with the standard k-ε model (equation (8)). The coefficients in equations (4), (5) and (9)

2-3 k-ω SST model
The standard k-ω SST was first introduced by Menter [9] in 1994 and then, some of its coefficients were developed by Menter et al. [10] in 2003. In this study, the modified SST model is used as follows.
where F1 is the blending function and is defined as follows.
and y is the distance from the nearest wall. The value of F1 in the far from the surfaces and inside the boundary layer is equal to zero (k-ε model) and equal to one (k-ω model), respectively. In this model, the turbulent viscosity is illustrated by equation (13).
In which, S is the invariant of strain rate and F2 is the second blending function illustrated by equation (14).
It should be noted that in SST model, a restrictive term is used for the production term in order to not overestimate the turbulence in stagnation regions (equation 15).
A review of second-moment computations for engineering flows has been provided by Leschziner [11] and Launder [12]. Results of these computations demonstrate the superiority of RSMs over eddy-viscosity models for curved flows, swirling flows, buoyant flows and recirculating flows.
In the present study, the model proposed by Gibson and Launder [13] has been chosen mainly because it takes into account the influence of a wall on the pressure field. The mentioned characteristic causes some advantage to simulate the decelerating and stagnation regions with high fluctuations in the pressure field.

2-5 v 2 f model
This model was first introduced in 1991 by Durbin [14]. In this model, two equations are added to the k   model and therefore, a 4-equation model is achieved. One equation for turbulent stresses perpendicular to the flow stream lines (v 2 ) and another as a function f to estimate the effects of the wall on the changes of v 2 . In this study, a modified model by Davidson et al. [15] is used. Equations 17 and 18 are used to calculate v 2 and f.
where Lscale and Tscale shows length and time scales respectively, and

Numerical method
The numerical platform is an open source CFD code based on the field operation and manipulation C++ class library for continuum mechanics (OpenFOAM 2.3.0) and is used to simulate the flow and represent the mean and instantaneous flow field 1 semi-implicit method for pressure-linked equations characteristics. The program uses collocated grid arrangement and the pressure-velocity coupling is handled with SIMPLE 1 [16]. The convective terms in the momentum equations are discretized using the second-order, bounded van Leer scheme [17]. The convective terms in the equations for turbulent quantities are discretized with hybrid upwind/central differencing. Fig. 2 represents the schematic of two dimensional computational domain and boundary conditions of the problem. As shown in Fig. 2, the nozzle width (e) is equal to 20 mm and the distance of nozzle from the plate (H) is 10e based on the Maurel and Solliec [18] study. They stated that a maximum turbulence intensity of the vertical component of the velocity field is occurred at the distance H/e=10. Therefore, in this distance one can be sure that the jet flow is fully turbulent when it hits the plane.

Computational domain
The channel length (Lch) is 15e and the total domain length (Lx) is equal to 85e. There are two outlet boundary conditions in both sides of the impingement plate. The length of these two outlet domains (Lout) are 5e. The aforementioned sizes of the computational domain are chosen from experimental setup data [8].

Boundary conditions
All boundary conditions are shown in Fig. 2. The inlet bulk velocity is uin=0.75 m/s. A parametric study of the inlet velocity is performed to approximate the inlet profile with conformity the experimental data used. As shown in Fig. 3, comparison to the experimental inlet values, it is observed that Lch/e =15 is the value that better fits our requisites, Re=16000. Turbulent kinetic energy in the nozzle exit is equal to that obtained from the PIV data [8]. The output boundary condition is considered as constant pressure and its value is equal to the atmospheric pressure.
No slip boundary condition is considered as at the channel wall, chamber wall and impingement plate. A wall function has been used for all simulation and the first grid point near the wall has been ensured to fall in the logarithmic region, i.e., 30<Y + <100 where Y + =yuτ/ν, uτ being the friction velocity.
The rate of dissipation (ε) and the rate of specific dissipation (ω) are taken equal to 3/4 3/2 / 0.07 C k e    and 1/2 1/4 / 0.07 kC    at the inlet for k-ε model (also RNG k-ε) and k-ω SST. In Reynolds stress model turbulence is assumed to be isotropic at the inlet and all normal components of Reynolds stresses are given the value 2k/3, whereas all shear Reynolds stresses are set to zero. Similarly, v 2 and f are equal to 2k/3 with zero gradient at the inlet for v 2 f model. At the walls v 2 and f are assumed to be 0 and (20v 2 v 2 )/(εy 4 ), respectively.

3-4 Grid study
To ensure the independency of results from computational grid size, four different grids are studied. These grids are named as case 0 to case 3 and contain 19,280, 34,280, 80,480 and 101,980 cells, respectively. In all cases, the mesh sizes are fine near the impingement plate and near the jet center. Fig. 4 indicates the longitudinal components of the mean flow velocity field near the impingement plate (y/H=0.9) and also on the parallel line with the impingement plate for k-ε and k-ω SST models. As seen in these figures, all four grids predict very close results for the mean velocity components, while according to the turbulent flow characteristics (such as turbulent kinetic energy) this condition is not acceptable. From results obtained for TKE in Figs. 4 a-d, case 2 is used in this study.

RESULTS
In this section, results obtained for five turbulence models, k-ε, RNG k-ε, k-ω SST, RSM and v 2 f, are presented for impingement of the water jet to a wall with Re=16000 in comparison with the experimental data at different impinging jet regions. Fig. 5 shows the longitudinal velocity component profile obtained from aforementioned RANS-based turbulence models at jet centerline in comparison with k-ω model (Kubacki and Dick [7]) and PIV data. It is seen that the k-ω model, has more error compared to the other models and estimates the potential core length about 0.79H (where the mean velocity is almost constant). As it can be seen in Fig.  5, the potential core lengths are 0.72H, 0.78H, 0.72H, 0.49H and 0.31H for k-ε, RNG k-ε, k-ω SST, RSM and v 2 f models, respectively. The potential core length obtained from PIV data is about 0.3H. In the next sections, flow simulation is done using different turbulence models at three main impinging jet regions.

4-1 Potential core region
The main characteristic of this region is uniformity of the velocity at jet centerline. When the nozzle-to-plate distance is 10 times the nozzle width, the potential core region continues up to y/e=3 (Fig.  5). As it is seen from Fig. 6, by moving away from the nozzle, the jet grows and the impact region on the surrounding fluid extends to x/e=1 (Fig. 6-c). However, the maximum velocity on the jet centerline is almost constant. In this region, the best agreement between the numerical results and experimental data is observed. By moving away from the nozzle, only the RNG k-ε model indicates the value of the longitudinal mean velocity component slightly less than the experimental one. As shown in Fig. 7, k-ε, RNG k-ε and k-ω SST models predict the dimensionless turbulent kinetic energy in the potential core close to zero, which is much lower than that of experimental data. While RSM and v 2 f predict higher values (about 0.005 for y/e=3). However, the value of the dimensionless turbulent kinetic energy from the experimental data is higher and is equal to 0.004, 0.008 and 0.0014 for y/e=1, y/e=2 and y/e=3, respectively. Moving away from the jet centerline to the either sides, the turbulent kinetic energy is increased due to the impingement of the jet to the surrounding fluid and reaches to its highest value. In the potential core, all models but RNG k-ε, over predict the maximum value of turbulent kinetic energy.  Fig. 8 represents the longitudinal mean velocity profiles in the developing region, calculated by different turbulence models in comparison with experimental data for different sections. By approaching the stagnation region, the difference between the simulated velocity profiles at the centerline with the experimental data gradually increases. The longitudinal velocity component calculated by k-ε and k-ω SST models in the developing region (up to y/e=7) is almost constant and is equal to the nozzle exit velocity. It then decreases slightly up to y/e=8. These two models overestimate the value of the longitudinal velocity component by about 70% at the centerline. Accuracy of the mean longitudinal velocity predicted by RNG k-ε model is slightly better than the other models.

4-2 Developing region
However, among the mentioned models, the maximum under prediction of the velocity in the shear layer region is obtained by k-ε model. The RSM model predicts the average longitudinal velocity much better than the other three models. By approaching to the decelerating region, the pressurestrain rate increases compared to the other terms in the momentum equation. So, this model is able to predict the velocity in this region better than previous models. However, the constant coefficients in the presented models are not accurate for strainpressure term in impinging flows. Therefore, the model is faced with restrictions to simulate the flow in this area. The vertical velocity component results obtained by v 2 f model have the best agreement with PIV data. However, the accuracy decreases by approaching the impingement wall (Figs. 8 a-e). Figures. 9a to 9e represent the non-dimensional TKE achieved from various turbulence models and also from experimental data. Unlike the potential core region in which most of the models over-estimate the turbulent kinetic energy compared to the experimental data, in this region the maximum energy (in the shear layer) is under-estimated by most of the models. Only the v 2 f model overestimates the value of the maximum energy. As it is seen from Fig. 9, similar to the potential core region, three models of k-ε, RNG k-ε and k-ω SST, estimate the TKE value in the jet centerline near zero. However, by moving away from the jet center line and shear layer, TKE value is over predicted. Among the five models, RSM and v 2 f models estimate the kinetic energy in the shear layer region more accurately. However, The RSM model estimates the value of TKE less than the experimental data and v 2 f model estimates the value higher than experimental data.

4-3 Impact region
Figs. 10-a to 10-c represent the velocity profiles in the impact region. As it can be seen, the RSM and v 2 f models provide a good estimation of the mean longitudinal velocity near the jet centerline.
RSM not only allows transport and development of the individual Reynolds stresses. They also have the advantage that terms accounting for anisotropic effects are introduced automatically into the stress transport equations. These non-isotropic characteristics of the turbulence play very important role in the flows with streamline curvature, swirl or strong recirculation. However, despite the reasonable estimation of RSM model, it does not have good accuracy in calculating the velocity field in the shear layer. In fact, all of the studied models under-estimate the value of velocity in the shear layer region (x/e>1). The largest difference between the numerical and experimental results is observed for turbulent kinetic energy in the impact region. The dimensionless turbulent kinetic energy of the jet is increased in the centerline of the jet in this region and reaches to 0.035 in y/e=9.5 ( Fig. 11-b). The flow tends to the either sides and a boundary layer near the wall-jet region forms. Therefore, the value of TKE is not decreased to zero by moving away from the jet centerline. As it is seen from Figs. 11-a to 11-c, RNG k-ε and k-ω SST models estimate the value of TKE in the jet centerline very poor. However agreement between the models results and experimental ones in the wall-jet region is good. Estimation of TKE by k-ε model is higher than the two previously mentioned models and still a large difference between these results and the experimental data is observed. TKE values obtained by RSM and v 2 f models are slightly over-estimated on the jet centerline at y/e=9 and 9.5. However TKE at y/e=9.8 obtained from PIV data is much higher than simulated using different turbulence models.
The maximum errors percentage of velocity stated by equation (19)

5-CONCLUSIONS
In this paper, the turbulent slot impinging flow is simulated using five RANS-based turbulence models at Re=16000 and the results are compared with the experimental data. In order to have fully turbulent flow, the nozzle-to-wall distance was considered sufficiently large (H/e=10). In order to investigate the impingement of a jet to a wall numerically and experimentally, the mean velocity field and the turbulent kinetic energy are considered. The numerical simulations were conducted by five most used models in the industrial processes, including k-ε, RNG k-ε, k-ω SST, RSM and v 2 f and the obtained results were compared with each other. The main novelty of this work is the use of a fine space resolution of 2D-PIV measurements (resolved in time) in a turbulent impinging water jet configuration (low velocity compared to that of air using similarity, so good accuracy) in order to well check the validity of the tested RANS turbulent models in different flow regions.
 All checked models simulate correctly the mean velocity field with rather great accuracy in the potential core region so that the maximum error is about 7% belonging to the RNG k-ε model. However the error estimated for turbulent kinetic energy in this region is large and reaches its highest value for the shear layer region.  The k-ε, RNG k-ε and k-ω SST models estimate the potential core length about 2.5 times more than the experimental values while the RSM and v 2 f models estimate this length 1.6 and 1.03 times bigger, respectively.  Comparison of the average velocity field obtained from the numerical and experimental results indicates that the estimated value by twoequation turbulence models (RNG k-ε, k-ω SST and k-ε) involved with errors up to 100%, which is maximum in the decelerating regions.  Precision of calculated velocity profiles at the jet centerline is very high in the RSM model. However, this model shows a high error (71%) in predicting the kinetic energy in the developing region at y/e=6. In this region, the results of v 2 f present the highest consistency with the experimental data, which indicates the validity of considering the pressure perturbation and wall effect in calculating v 2 and f functions, respectively.  The maximum error in calculating TKE was observed in the impact region, in which the sudden change in the direction of the velocity vector in addition to non-isotropic nature of the flow was resulted in an inaccurate modeling of the turbulence viscosity (νt) and the turbulent kinetic energy. The RSM and v 2 f results show the highest agreement with the experimental data. It seems that by improving the constant in these models, better agreement is achieved.
 Although the v 2 f model is successfully estimates the time averaged velocity field in three potential core, developing and impact regions, especially on jet centerline, however, it showed the highest error in estimating the turbulent kinetic energy in the wall-jet region. This reflects the failure of this models to predict the vertex deformation in the impact area and formation of the boundary layer.