## 1. Introduction

An understanding of the behavior of turbulent flow and bottom sediment transport behind structures requires a series of experimental investigations. However, physical experiments demand an enormous time and include a similarity limitation between the prototype and the model (Ettema et al., 2004). Therefore, an alternative method is desirable. Many studies have been devoted to the development of numerical modeling in the past decades (Bosch and Rodi, 1998;Liou et al., 2002;Chen et al., 2003;Honzejk and Frana, 2008;Park et al., 1998;Park et al., 2009;Park et al., 2014).

Simple attempts have been made to model the nature of turbulence in the wake of the flow past the stationary objects engaged in 1-D or 2-D numerical simulations using the cross section averaged or depth averaged flows (Bosch and Rodi, 1998;Park et al., 1998;Honzejk and Frana, 2008). With computer hardware and numerical scheme developments, numerical modeling allows us to elucidate the comprehensive physical mechanism of the wake (Park et al., 2014). In previous studies, Olsen and Melaaen (1993) performed a 3-D numerical simulation with the k-ε model to estimate the turbulent flow field around the cylinder structure, and Ishino et al. (1993) examined wake characteristics behind the cylinder as the Reynolds number increased from small to large values. The turbulent motions near the structures were produced and compared with the experimental data from Richardson and Panchang (1998) who utilized the RNG (Renormalized Group) k-ε model. To analyze the flow around the square and circular structures, the Large Eddy Simulation (LES) turbulence model was adopted by Tseng et al. (2000), Liou et al. (2002) and Chen et al. (2003). Yang and Choi (2002) simulated the turbulent wake behind the bridge piers with the built-in RNG k-ε and LES models of the computational fluid dynamic commercial code (FLOW-3D). Park et al. (2014) recently simulated wakes behind a square cylinder using the zero- and one-equation turbulence models with FLOW-3D and investigated the applicable range of the major empirical constant in each model.

The wakes behind a structure depend on the drag coefficient and the flow characteristics. Accurate predictions of wakes are critical to the safety of structures in channel, rivers, harbors and coastal areas (Park et al., 2014). The objective of this study is to calculate the wakes behind a square cylinder using two types of two-equation turbulence models based on the eddy viscosity concept, the k-ε and RNG k-ε models. In addition, we evaluate each model performance by comparing its results with analytical solutions using three skill assessments: the correlation coefficient (*γ*) for the similarity of the wake shape, the error of maximum velocity difference (EMVD) for the accuracy of wake velocity, and the ratio of drag coefficient (RDC) for the pressure distribution around the structure (see Park et al., 2014). The study also investigates the applicable range of empirical constants and suitable values for each turbulence model in the same manner as Park et al. (2014).

## 2. Analytical Solutions and Turbulence Theories

### 2.1 Analytical Solutions

Schlichting (1930) investigated wakes behind a single body based on Prandtl’s mixing length hypothesis and suggested an analytical solution. Görtler (1942) and Reichardt (1951) introduced the eddy viscosity *ν _{t}* that is a function of the non-zero mean velocity gradient, and assumed that

*ν*is a constant

_{t}*ν*

_{0}. They provided an equation for a wake using

*ν*and its solution can be given by

_{t}

where *u _{D}* is the velocity difference of a wake,

*U*is the free stream velocity,

*C*is the drag coefficient of a resistance body,

_{D}*D*is the representative distance of the body,

*x*is the flowing distance,

*y*is the transverse distance and

*ν*=

_{t}*ν*

_{0}= 0.0222

*C*

_{D}*D U*. This theoretical solution Eq. (1) was in excellent agreement with Schlichting’s measurements (1979) and was used to evaluate, therefore, the wakes behind a square cylinder in this study (see Fig. 1 of Park et al., 2014).

### 2.2 Turbulence Theories

For the broadest range of applicability of the turbulent model, the characteristic length scale *L* and the velocity scale need to be specified through theoretical approaches. In the one-equation model, the characteristic velocity scale was expressed as the function of the turbulent kinetic energy k as ${\nu}_{t}=L\sqrt{2\text{k}/3}$. The turbulence energy dissipation ε scaled as *C*k^{3/2}/*L* where *C* was a constant (Chou, 1945). Thus k and ε could be combined to eliminate the empirical mixing length and to produce an eddy viscosity equation Eq. (2). The k-ε model, one of the two-equation models, solved two transport equations for k and ε and, after that, yielded the eddy viscosity.

where *C _{μ}* is an empirical constant and the typical

*C*is

_{μ}*C*

_{μ}_{0}(= 0.09) (Launder and Spalding, 1974).

Using a field-theoretical approach and statistical methods, the RNG k-ε model was used to calculate model constants. The basic idea was that the small turbulence scales could be removed through space-time Fourier decomposition of the velocity field and repeated substitution in the Navier-Stokes equations. Also, the model constants were evaluated from the large scale field and the modified viscosity induced by the elimination of the small scales (Yakhot and Orszag, 1986;Yakhot et al., 1992). The model introduced the modified eddy viscosity *ν _{eddy}* in lieu of an eddy viscosity in the k-ε model (Yakhot et al., 1992).

where *ν* is the kinematic viscosity and *C _{μR}* is an empirical constant and its typical value is known as

*C*

_{μR}_{0}(= 0.085) (Yakhot and Smith, 1992;Yakhot et al., 1992).

## 3. Numerical Models and Simulation Conditions

### 3.1 Numerical Models

FLOW-3D models (Flow Science, 1999) based on the Marker and Cell (MAC) and the Solution Algorithm Volume Of Fluid (SOLA-VOF) methodologies were used in this study. The velocity and pressures were computed using the continuity and momentum equations, and the water elevations were solved using the VOF (Volume of Fluid) method (Flow Science, 1999). To calculate the turbulence flow around a square cylinder, we utilized built-in turbulence modules such as the k-ε and RNG k-ε models in FLOW-3D. A detailed description of these models is given in the FLOW-3D User’s Manual (Flow Science, 1999).

The principal goal of any turbulence model is to provide a mechanism for estimating the influence of turbulent fluctuations on mean flow quantities. This influence is usually expressed by the additional diffusion terms in the equations for mean mass and momentum. In the one-equation or k-ε models, k is governed by Eq. (4).

where *P* is a shear production and *D _{k}* is a diffusion term. The transport equation for the ε is

where *C _{∈}*

_{1}and

*C*

_{∈}_{2}are the user defined constants and

*D*is a diffusion term. In the RNG k-ε model (Yakhot and Smith, 1992;Yakhot et al., 1992), the additional term

_{∈}*R*on the righthand side of Eq. (5) is an extra strain term given by

Where $\mu =\sqrt{P}\text{k}/\in $, *η*_{0} = 4.38 and *β _{R}* = 0.012.

### 3.2 Drag Coefficients

As described in the previous study by Park et al. (2014), the drag coefficient in the *x*-direction is defined as

where *F _{x}* is the drag force in the

*x*-direction,

*A*′ is the projected area in the

_{x}*x*-direction and

*U*is the upstream velocity in the

_{x}*x*-direction. Park et al. (2009) and Park et al. (2014) calculated the drag coefficient of an object in the numerical basin using FLOW-3D, and the same approaches are applied to this study.

### 3.3 Simulation Conditions

For the adequate calculation of wakes behind a square cylinder, we used the numerical basin that was 280 m long, 60 m wide and 5 m deep. Also a non-uniform grid containing 190 × 44 horizontal cells and 11 vertical layers was used as shown in Fig. 1. A constant flow propagates toward a square cylinder that is 30 m from the boundary inlet in which the drag coefficient of the square cylinder is 2.2 at a Reynolds number Re ~ 1.0 × 10^{5} (Blevins, 1984). The constant velocity, ** V** = (

*u*, 0, 0) = 0.1 m/s, was specified at the boundary inlet. And the continuative boundary condition was applied at the flow outlet, while the smooth and the no-slip condition were used at the bottom boundary conditions. Here, all the calculation conditions were the same as in the previous study by Park et al. (2014).

FLOW-3D with the inclusion of the two different turbulence models, the k-ε and RNG k-ε models, was used to simulate wakes behind the square cylinder. Launder and Spalding (1974) and Yakhot et al. (1992) suggested the model constant *C _{μ}* = 0.09 (=

*C*

_{μ}_{0}) and

*C*= 0.085 (=

_{μR}*C*

_{μR}_{0}) for the k-ε and RNG k-ε models, respectively. According to the experiment of Rodi (1972), the eddy viscosity for the simulation of wakes behind the cylinder must be 3 to 13 times higher than that for the round jet, plane jet and pipe flow. Thus we calculated wakes using constant values like 4

*C*

_{μ}_{0}, 5

*C*

_{μ}_{0}, 6

*C*

_{μ}_{0}, 7

*C*

_{μ}_{0}, 8

*C*

_{μ}_{0}and 9

*C*

_{μ}_{0}or 4

*C*

_{μR}_{0}, 5

*C*

_{μR}_{0}, 6

*C*

_{μR}_{0}, 7

*C*

_{μR}_{0}, 8

*C*

_{μR}_{0}and 9

*C*

_{μR}_{0}in both the turbulence two-equation models. In addition, the simulation was performed for over 5000 seconds with varying time steps controlled by the program. After the flow reached a steady state condition, the simulation results were analyzed (Park et al., 2014).

## 4. Model Results and Discussion

### 4.1 General Horizontal and Vertical Distribution of Wakes

Using the RNG k-ε model, we calculated the velocity fields, and the longitudinal (*x*-direction) velocity distribution is especially shown in Fig. 2. The wakes behind the square cylinder were similar to the long tail of a shooting star in the sky and the width of the wake increased as the flow moved downstream. The magnitude of the horizontal velocity at the surface layer was larger than that of the middle and bottom layers (Park et al., 2014).

From the longitudinal velocity along the centerline, we could find that a reverse flow occurred near the square cylinder due to the low pressure. After passing the cylinder, the flow rapidly recovered around *x*/ (*C _{D}*

*D*) = 10. The velocity gradient slowly varied till

*x*/(

*C*

_{D}*D*) = 50 and became stabilized after

*x*/(

*C*

_{D}*D*) = 50. These results were quite close to those obtained by Schlichting (1930). Resulting from the strong vertical mixing and interaction of the vortex, the differences in the horizontal velocity at the three layers were insignificant just behind the structure. These velocity characteristics looked very similar to the previous study by Park et al. (2014). On the other hand, the results of the k-ε model were similar to those of the RNG k-ε model, and so are not presented here.

### 4.2 The k-ε Model’s Results

The k-ε model, which solves two transport equations for the k and the ε, does not need the information of the mixing length. Although the eddy viscosity in the model depends on the k and ε, an empirical constant *C _{μ}* in eddy viscosity equation (Eq. 2) has to be specified. The standard

*C*(

_{μ}*C*

_{μ}_{0}= 0.09) from the empirical observation was used everywhere (Launder and Spalding, 1974). Rodi (1972) noted that the

*C*for the wake experiments behind a sphere must be 3 to 13 times higher than that of for the jet or pipe flow experiments. Due to the fact that there were insufficient studies to suggest reasonable values of

_{μ}*C*for the simulation of wakes behind a square cylinder, we conducted several scenarios using the k-ε model with different

_{μ}*C*that ranged from 4

_{μ}*C*

_{μ}_{0}to 9

*C*

_{μ}_{0}, found the best model constant and evaluated the applicability of the model constant for the simulation of wakes.

For all scenarios, the k-ε model results underestimated the wake width (Fig. 3) and the values of *γ* ranged from 0.843 to 0.877 (Table 1), which were worse than that of the one-equation model and mixing length model in the previous study by Park et al. (2014). Even if the EMVD was similar to that of the one-equation model, each velocity difference at three vertical layers was significant which meant that the model created weak vertical mixing. The RDC ranged from 56.2 % to 56.4 % in Table 1 which was relatively high when compared to the other models, so that the model had a limitation to require decent RDC at any *C _{μ}* . Despite the fact that the k-ε model employed two transport equations, the model could not successfully calculate wakes behind the square cylinder. The reason stemmed from the eddy viscosity hypothesis and ε-equation. The ε-equation included two empirical constants with suggested

*C*

_{∈}_{1}= 1.44 and

*C*

_{∈}_{2}= 1.92 in Eq. (5). It performed well for simple flows but not for complex flows as mentioned by Stephen (2000). To improve the accuracy of model calculation for wakes behind the square cylinder, the model constants must be adjusted and close investigations for their ranges are necessary.

### 4.3 The RNG k-ε Model’s Results

Many studies have contributed towards improving the performance of the k-ε model (Pope, 1978;Hanjalic and Launder, 1980;Bardina et al., 1983). Yakhot and Smith (1992) developed the RNG k-ε model based on theoretical grounds. The RNG k-ε model included an additional term in the ε-equation which could not be neglected and was important in the highly strained flows and low Reynolds number flows (Yakhot et al., 1992;Maurizi, 2000;Kim and Baik, 2004). Since the model is based on a theoretical ground, all the constants in the ε-equation derived from the theory itself and the values are *C _{∈}*

_{1}= 1.42 and

*C*

_{∈}_{2}= 1.68. Therefore, the empirical constant

*C*in the modified eddy viscosity equation plays a major role in calculating the eddy viscosity. The typical

_{μR}*C*(=

_{μR}*C*

_{μR}_{0}) in the RNG k-ε model is 0.085, which is similar to the value

*C*

_{μ}_{0}(= 0.09) used in the k-ε model. We calculated the wakes behind the square cylinder using the RNG k-ε model with various ranges of

*C*. The

_{μR}*C*ranges 4

_{μR}*C*

_{μR}_{0}to 9

*C*

_{μR}_{0}were the same as those used in the k-ε model because of the closeness of two constant values as well as for the comparisons of two models.

As the *C _{μR}* increased, the wake width increased and the increasing trend was similar to that of the results of the one-equation model in the previous study by Park et al. (2014). As shown Fig. 4, the range of

*γ*was 0.963 to 0.995 (Table 2) which was the best when compared to other models. The EMVD showed large variations from -28.1 % to 22.4 % in Table 2 corresponding to the

*C*. But the variations of EMVD become relatively insignificant like 5 % to 7 % when the

_{μR}*C*changed from 6

_{μR}*C*

_{μR}_{0}to 8

*C*

_{μR}_{0}. At these ranges, the model predictions agreed well with the analytic solutions. The RDC generally increased along with the increase of

*C*and its range was 0.662 to 1.315 (Table 2). The best fit of the RDC was 1.032 which showed 3.2 % error when it occurred at

_{μR}*C*= 6

_{μR}*C*

_{μR}_{0}as shown in Fig. 4 and Table 2. The AME (Absolute Mean Error) with the

*γ*, EMVD and RDC ranged from 3.8 % to 19.3 % (Table 2) and these values were much better than those of the k-ε model (Fig. 3 and Fig. 4). This confirmed that the additional term and the model constants derived from the theory in the ε-equation were key components in simulating the wakes behind the square cylinder. From the RNG k-ε model results, the agreement of model predictions and analytic solutions was exceptionally good when the

*C*was 6 times higher than the suggested value. In addition, on the basis of the experiment of Rodi (1972) at P/ε < 0.5, the eddy viscosity for wakes behind the cylinder must be 3 to 13 times higher than that for the round jet, plane jet and pipe flow simulations. As a result, the high value of

_{μR}*C*might be appropriate for the simulation of wakes behind the square cylinder.

_{μR}## 5. Conclusions

The study simulated wakes behind a square cylinder using two-equation turbulence models based on the eddy viscosity concept, the k-ε and RNG k-ε models, and investigated the applicable range of the major empirical constant in each model. For the comparisons between the model predictions and analytical solutions, we employed three skill assessments:, the correlation coefficient for the similarity of the wake shape, the error of maximum velocity difference for the wake velocity distribution, and the ratio of drag coefficient for the flow patterns. Based on the above results, the following conclusions were arrived at.

The k-ε model results underestimated the wake width, and its performance was worse than the one-equation and mixing length models in the author’s previous study. The reason stemmed from the eddy viscosity hypothesis and ε-equation. The use of the ε -equation, including empirical constants with suggested values, might be inappropriate for the simulation of wakes behind the square cylinder. Therefore, adjusting the model constants was necessary for the accuracy of the model calculations.

Excellent results were obtained with the RNG k-ε model, resulting from an additional term and the model constants derived from the theory in the ε-equation. The AME averaging the *γ*, EMVD and RDC were less than 3.8 % at *C _{μR}* = 6

*C*

_{μR}_{0}, which was 6 times higher than the suggested value. It indicated that the high value of

*C*might be appropriate for the simulation of wakes behind the square cylinder.

_{μR}