Turbulence Model Sensitivity and Scour Gap Effect of Unsteady Flow around Pipe: A CFD Study

A numerical investigation of incompressible and transient flow around circular pipe has been carried out at different five gap phases. Flow equations such as Navier-Stokes and continuity equations have been solved using finite volume method. Unsteady horizontal velocity and kinetic energy square root profiles are plotted using different turbulence models and their sensitivity is checked against published experimental results. Flow parameters such as horizontal velocity under pipe, pressure coefficient, wall shear stress, drag coefficient, and lift coefficient are studied and presented graphically to investigate the flow behavior around an immovable pipe and scoured bed.


Introduction
Scouring is a phenomenon caused by erosion of sediment of sand bed around an obstruction, that is, bridge piers and abutments in a flow field [1][2][3]. In other words, scouring basically happens due to the movement of the foundation of the bed under the flow filed conditions in which the flow surrounding the obstruction gets accelerated and induces high shear stress over the seabed surface [4]. The sand bed reduced around the obstruction under the flow level is named the scour depth. A scour hole is a pit or void that forms as a sequence of the sand bed sediment removal from the river bed [5].
Prediction of scour around bridge piers and submarine pipelines attract the hydraulic and ocean engineers. Cylindrical bridge piers are the most commonly used structures in coastal, offshore, and river engineering. Local scouring around the bridge piers is considered to be one of the most common causes of bridge failure [6,7]. The local scour around river hydraulic structures is a disaster mitigation of the engineering structure [8]. It leaves them in unsafe conditions requiring maintenance and occasionally results in loss of life. Damage of hydraulic structure because of local scouring is a global concern and it has been studied by many researchers experimentally and numerically for several decades. For bridge engineering practice, accurate prediction of local scour, such as the maximum depth of scour around the bridge piers, is critical for bridge design, maintenance, and evaluation. Scour surrounding obstructions, that is, bridge piers and abutments, are considered as a common reason causes bridge failures compared to other causes in the history [9,10]. The flow in the horizontal channels e.g. river obstructed by the vertical column gets separated, and when viewed from the top, looks like the shoe of a horse. Because of the occurrence of such shape, it is known as horseshoe vortex. This separated layer rolls up along the bridge piers to form a vortex which is known as horseshoe vortex because of its shape.
Mao [11] studied the interaction between a pipeline and erodible bed. Author observed the scour around horizontal cylinders in steady current and wave conditions, as well as with different Reynolds numbers (Re), Shields parameters, and pipeline gaps. These experiments investigated scour features such as shape and size of the scour hole and the time scale of scour formation. Later this work was further investigated by Jensen et al. [12]. They investigated experimentally the flow around a pipeline placed initially on a flat, erodible bed at five characteristics stages of the progressive process in 2 The Scientific World Journal currents. The results showed that as the scour develops with time and space, the mean flow field and turbulence around and the forces on a pipeline undergo considerable changes.
Experimental observation [13] revealed that the horseshoe vortices play the basic role in the scour around the bridge piers. The first scouring takes place in the wake of the cylinder. The primary wake vortices and the accelerated side flow are the main cause of this scouring. A CFD study to predict the local scour hole around the cylindrical pier was carried out by [14]. They found that the CFD methods with the powerful flow visualization show the ability of flow representation during local scouring. Authors concluded that the model used is sufficient to predict the complicated flow field around piers that mounted on sandy bed. In a similar numerical study [15], authors used Eulerian twophase model coupled with Euler-Euler governing equations for fluid and solid phases. Investigations into the mechanism of scour reveal that three sediment transport modes (bedload, suspended-load, and laminated-bed) are associated with the scour development.
Roy and Matin [16] experimentally investigated the behavior of scour at floodplain and main channel using three bed materials with three discharges and four lengthwidth ratios. They found that scour behavior at flood plain and main channel for different bed materials, discharges, structure shapes, and length-width ratio is nearly same. They also found that the flood plain is lower in main channel than the flow velocity and as a results, flow velocity is responsible for deeper scour in main channel. The work reported in [17] is a numerical study for simulating a flow field around a circular pier on sandy bed. Reynolds-averaged Navier-Stokes equations (RANS) using standard -model and space averaged Navier-Stokes equations using large eddy simulation (LES) with the standard Smagorinsky subgrid model were employed to simulate the flow field. It was concluded that RANS simulation is an accurate and sufficient model for scouring in engineering applications and requires lesser grids for simulation. There has been substantial amount of work using finite volume method (FVM) for RANS equations in the various flow condition. del Coz Díaz et al. [18] and Nieto et al. [19] in their work used air as the flowing medium and investigate the performance of snow fence and self-weighted metallic roof using FVM.
Kazeminezhad et al. [20] investigated numerically the force components and vortex shedding frequency of a pipe exposed to a steady current in terms of the drag coefficient, lift coefficient, and Strouhal number. It was concluded that the mean force coefficients and the root-mean-square (RMS) lift coefficient are strongly affected by the gap to diameter ratio while the Strouhal number is slightly affected by the gap ratio. Later in [21] they numerically investigated wave-induced tunnel scour beneath marine pipelines. Investigations revealed that the tremendous sediment transport takes place during the tunnel scour stage and the high turbulence intensity.
This study aims to investigate the effect of turbulence models on the flow field behavior at different five scouring phases and study the effect of scouring on flow parameters such as velocity under pipe, pressure coefficient around the pipe, and wall shear stress on the scour's bed. Figure 1 shows the schematic of the two-dimensional (2D) geometrical domain used in the present study along with the corresponding boundary conditions. A logarithmic velocity profile as presented in Figure 3 is created using user-defined functions (UDF) in Fluent based on the following formulation [22]:

Geometrical Structure and Boundary Conditions.
(1) The velocity profile is applied at the inlet. The profile from our CFD model is compared with that from experimental study of Dudley [23] for consistency; see Figure 3. Zero pressure outlet boundary condition is applied at the flow exit. The water surface (top wall) is set as a symmetry boundary condition. No-slip boundary condition is applied on the pipe surface and the scour bed. Gravity acts in the negativedirection.

Governing Equations.
The continuity and the momentum equations for the present case are as given below.
Continuity equation is as follows: -component of the momentum equation is as follows: -component of the momentum equation is as follows:

Turbulence Modeling. Turbulence models of twoequation -, two-equation -models, and five-equation
Reynolds stress models are used in the present investigation and their results are compared with experimental data from literature.

-Models.
Two-equation -models, turbulent kinetic energy and turbulent dissipation , are the simplest and the most widely used models among all turbulence models that aim to study the effect of turbulence in the flow. Two-equation model signifies that it includes two extra transport equations to represent the turbulence properties of the flow. There are three different models that are derived from -model, standard -model, realizable -model, and renormalization group method (RNG). Despite having the two common equations, these turbulence models use   the different ways to calculate the principle form of the eddy viscosity equation. The variants -models approximate the eddy viscosity as The turbulent kinetic energy ( ) and its dissipation rate ( ) for the standard k-model are calculated from The model constants are 1 = 1.44, 2 = 1.92, 3 = −0.33, = 0.09, = 1.0, and = 1.3. The modeled transport equations for ( ) and ( ) in the RNG -model are The model constants are 1 = 1.42, 2 = 1.68, and = 0.084.
The model constants 2 , , and have been established to ensure that the model performs well for certain flows. The model constants are -models, also known as two-equation models, have the same definition for as in -model. However, it differs in the selection of second variable ( ). This model is broadly categorized into two types, the standardmodel and the shear stress transport (SST) model. The turbulent viscosity, , for the -model, is computed by combining ( ) and ( ) as follows: The modeled transport equations for ( ) and ( ) in the standard -model are The modeled transport equations for ( ) and ( ) in the SST -model are

Reynolds Stress Model (RSM).
Unlike the previous two turbulence models, Reynolds stress model (RSM) abandons the calculation of eddy viscosity and solves the transport equation for Reynolds stresses. It provides four additional equations for 2D flow and seven for 3D flow. Having these additional equations, RSM takes comparatively more effort and time to simulate the flow.
The exact transport equation for the Reynolds stresses, ( ), may be written as follows: 2.4. Numerical Methods. The commercial CFD software FLUENT 14.0 [24] which is based on finite volume method (FVM) is used to solve the Reynolds-averaged Navier-Stokes (RANS) equations for the incompressible flow. The transport governing equations are discretized using the second order upwind spatial discretization method. The pressure-implicit with splitting of operators (PISO) scheme was used for the coupling of the pressure and the velocity fields. The underrelaxation factor of all the components, such as velocity components and pressure correction, is kept at 0.3. The scaled residuals of 1 × 10 −6 are set as the convergence criteria for the continuity and momentum equations. Transient model based on implicit scheme with a time step was used in the current numerical investigation. The typical wall treatment function + (= /]) value of the first node in all turbulence models near the bed profile is less than 1.

Simulation Cases.
A total of 11 cases were simulated in the present study; see Table 1. In the first section of this study, the effect of different turbulence models such as -models (standard, RNG, and realizable), -models (SST, standard), and Reynolds stress model (RSM) on horizontal velocity and kinetic energy square root has been investigated. For studying this effect, we have adopted the domain suggested by Jensen et al. [12] to validate the results of Mao [11]. This domain is of 0.  In the second part of this paper, a parametric study has been carried out by using five bed profiles suggested by Mao [11] at 10 min, 30 min, 100 min, 200 min, and 300 min. Five simulation cases (7)(8)(9)(10)(11) were run to obtain the horizontal velocity profile under the pipe, pressure coefficient , wall shear stress , drag coefficient , and lift coefficient around the pipe with different gaps at position / = 0 using standard -turbulence model.

Mesh Independency and Time
Step Test. The domain is meshed using a structured grid; see Figure 2. Mesh independency is carried out using five different types of grids having cell number of 40000, 82500, 130225, 223404, and 315623. A finer grid is used near the bed. The mean velocity was calculated for all selected mesh sizes and plotted against the water depth as shown in Figure 4. The domain with 223,404 cells is selected throughout this study, because it shows reasonable accuracy and the lowest deviations of the velocity profile and turbulent kinetic energy (TKE) for turbulence models compared to Jensen laboratory data [12]. Four different time step sizes (Δ = 0.0002, 0.002, 0.02, and 0.2) were tested and 0.002 is used throughout this current numerical study as there were no deviations observed below this value of time step. The number of iterations for this study is kept between 3000 and 3500.

Results and Discussion
Current numerical study has been carried out in two major sections. First section deals with the effect of various turbulence models on horizontal velocity profile and kinetic energy square root while the second one deals with the effect of scour gap on velocity profile under the pipe, wall shear stress, and others. A large number of simulations have been run to study the effect of turbulence models on unsteady horizontal velocity profile and turbulent kinetic energy square root, and the results are discussed in Section 3.1. The preferred turbulence model identified using above study is used for further investigation of effect of scour gap on velocity profile under the pipe, , wall shear stress, , and which is discussed in Section 3.2. the turbulent kinetic energy (TKE) at some axial length (i.e., / = −3.0, 1.0, 1.5, and 2.0) using different type of turbulence models and the results are presented in Figures 5  and 6, respectively. For comparisons, the experimental results reported in Chang [2] are also presented in the figure. Note that the location of the axial position covers the front and rear part of the pipe (or obstruction). Referring to Figures 5(a) to 5(d), which is presented at / = −3.0, 1.0, 1.5, and 2.0, respectively, the standard -turbulence model prediction is much closer to the experimental data at various water depths; some of them overlap each other. Nearly the same can be said for the realizable -model and RSM model but some deviation seen at some of water depths, for example, at water depths below 0.01 m and those between 0.03 and 0.04 m at / = 2.0 ( Figure 5(d)) for the both models. RNG -model prediction is quite accurate for / = −3.0 and 1.0 (Figures 5(a) and 5(b)), but this is not the case for / = 1.5 and 2.0 (Figures 5(c) and 5(d)); overprediction of the velocity is seen at water depths between 0.01 and 0.03 m. Relatively, standard -and SST -models are the most inaccurate among the turbulence models studied. Significant under prediction of the velocity is seen at water depths from 0 to 0.04 m at / = 1.0, 1.5, and 2.0. Referring to Figure 6, the TKE is well predicted by standard -turbulence with some deviation at / = 2.0 ( Figure 6(d)) in comparison to other turbulence models. RNG -model and RSM model also reasonably predict TKE especially for / = 1.5 and 2.0 and realizable -model has a good prediction only at / = 1.0. On overall, realizable -model, standard -, and SST -models can be regarded as the least accurate among the turbulence models studied in predicting the turbulence kinetic energy.

Effect of Different Turbulence
Among all turbulence models, standard -model gives better predictions for the velocity profile when compared with that of experimental data of Jensen et al. Standardmodel is widely acceptable for scour process modelling because of its better accuracy and reasonable computation time and does not require high computational facilities. However, the -models show a discrepancy with the experimental data and the velocity profile deviates largely near the wake region (highly turbulent region) and less in the farther region in the downstream profile. This deviation is more significant at the bottom (near the wall bed) and around the pipeline and reduces towards the water surface. This deviation may be because -model produces slightly large turbulence in the weak region which is less encountered with the variants -models and the Reynolds stress model (RSM). It may also be because of the absence of the wall treatment function in -models, and they are also very sensitive to the inlet boundary conditions flow. Figure 8 shows the horizontal velocity profile under the pipe (obstruction) at / = 0 for different types of scoured bed profiles as presented in Figure 7. A higher horizontal velocity can be expected for the scoured bed profiles of 10 min and 30 min due to a narrow gap under the pipe to the surface of the bed. The opposite can be said for the scoured bed profiles of 100, 200, and 300 min.

Effect of Scouring Depth.
A reduction of around 18% of horizontal velocity is observed from time phase from 10 min to 300 min. Figure 9 shows the prediction of the pressure coefficient and wall shear stress variation around the pipe at five different circumferential points ( ). The variation of the pressure coefficient with gap depth at different circumferential positions is shown in Figure 9  at that location, whereas the negative value indicates that the pressure drops and the water level decreases. Figure 9(b) shows the variation of the wall shear stress with gap depth. The wall shear stress forms as a sequence of approaching the flow on the pipe parallel to bed surface and it experiences the maximum value max at the flow separation point which causes scouring to begin. At the early stage, max occurs near the bottom in each side of the pipe forming shedding vortex that creates a space on both sides of the pipe which is soon occupied with water and leads the scour formation. Further, the maximum shear stress moves downward approaching the bed and moves the sand particles in the flow direction. It is observed that the higher the shear stress, the deeper the scour hole below the pipe.
When the flow approaches the pipe, it imparts drag and lift forces on the pipe. In the current numerical study, variations in these forces are observed in terms of drag and lift coefficients for five different time phases and presented in Figure 10. Shear stress causes the bed erosion and as soon as it comes into action, the pipe experiences negative lift coefficient ( ) and it decreases with time and consequently increases the scour gap. The negative lift can be attributed to the suction below and behind the pipe caused by the scour gap. For velocity, negative lift can be explained by the position of the stagnation point pipe and the angle of attack of the approaching flow. As a result of lift coefficient elimination, the drag coefficient is reduced with time as well. A reduction of 23.7% in and 51.3% in was observed between the time phases of 10 min and 300 min.

Conclusions
Two-dimensional (2D) CFD analyses were carried to investigate fluid flow over an obstruction under different type of bed profiles using a number of turbulence models. Unsteady horizontal velocity profile and the kinetic energy square root at few axial directions are investigated. The effect of scour depth on the velocity distribution under the pipe, the wall shear on the bed and the pressure coefficient, the drag coefficient, and the lift coefficient of the obstruction body were numerically investigated. The conclusions of the current study are as follows.
(i) The standard -model is able to predict accurate results in comparison to other turbulence models when compared to experimental data for unsteady horizontal velocity and turbulent kinetic energy square root.
(ii) From the velocity behavior under the pipe, results show that the maximum velocity at each phase decreases with increasing time till the scour reaches its equilibrium depth and a reduction of around 18% of horizontal velocity is observed from time phase from 10 min to 300 min. Gap between maximum scoured bed and cylinder (m) : The generation of turbulence kinetic energy due to buoyancy : The generation of turbulence kinetic energy due to the mean velocity gradients : von