Turbulent free surface flow over semicircular and circular obstacle in a duct

Turbulent free surface flows subject to complex boundaries are frequently encountered in nature and engineering, such as channel flows over uneven bed and submerged obstacles, water wave propagations in coastal regions and interactions with structures. Focusing on the difficulties from the combination of free surface, turbulence and complex boundary configuration, researchers made lots of efforts in the last several decades. Many techniques have been developed and used successfully over the years. These methods can be classified into two main types: interface tracking and interface capturing. The former computes the liquid flow only, using a numerical grid that adapts itself to the shape and position of the free surface. The free surface is represented and tracked explicitly either by marking it with special marker points, or by attaching it to a mesh surface. Various surface fitting methods for attaching the interface to a mesh surface were developed during the past decades [1-4]. In the interface-tracking methods, the free surface is treated as a boundary of the computational domain, where the kinematic and dynamic boundary conditions are applied. These methods cannot be used if the interface topology changes significantly [4]. The second possible approach is given by the so-called interface-capturing methods [5-10]. One may also divide these methods into three groups: the Lagrangian method, the Eulerian method and the arbitrary Lagrangian Eulerian method. Generally speaking, the capturing method has advantages over the tracking methods from the viewpoint of generality although they are not as accurate as the latter one. The Volume of Fluid (VOF) method was proposed and widely used for numerous problems due to its small requirement on computer storage and inherent abilities to dealing with complex free surface process, such as large deformation, breaking, overlapping and splashing [11-13]. One of the main concerns with the interface capturing approaches has been sustaining global mass conservation in long-time integrations [14]. In order to make up the limitation of VOF mentioned above, some researchers have made various efforts [14-17]. Turbulence is an important characteristic of fluid flows and remains still a great challenge for computational fluid dynamics [18]. Nowadays, most of the existing numerical models for free surface flows employ the Reynolds-Averaged Navier-Stokes (RANS) procedure to account for the turbulence. The increasing power of computers allows a more advanced method to be integrated in the VOF based on numerical models. As an alternative of Direct Numerical Simulation (DNS) [6], the Large Eddy Simulation (LES) technique is a promising approach [1921]. Lu et al. [22] proposed a two-dimensional hybrid numerical model, FEM-LES-VOF, for free surface flows, which is a combination of three-step Taylor-Galerkin finite element method, large eddy simulation with the Smagorinsky sub-grid model and Computational LagrangianEulerian Advection Remap Volume of Fluid (CLEARVOF) method. Löhner et al. [4] developed a twodimensional computer model to simulate free surface flow interaction with a moving body. VOF technique in coupling with an incompressible Euler/Navier–Stokes solver operating on adaptive, unstructured grids has been used to simulate the interactions of extreme waves and threedimensional structures [4]. To the best of our knowledge no comprehensive model has been published yet to study the free surface flow over semicircular and circular cylinder obstruction. In this article, the governing equations are numerically solved in the domain by the control volume approach based on the SIMPLE technique. In this study, effects of semicircle and circle obstruction on free surface flow are investigated and the final profile of free surface in the vicinity of obstruction is comprised with the final profile in a duct without obstruction. The remainder of the paper is organized as follows: Section 2 summarizes the Problem formulation and numerical modelling; Sections 3 describes the Numerical method; code validation is shown in Section 4; Sections 5 discusses the results and finally, some conclusions are given in Section 6.


Introduction
Turbulent free surface flows subject to complex boundaries are frequently encountered in nature and engineering, such as channel flows over uneven bed and submerged obstacles, water wave propagations in coastal regions and interactions with structures.Focusing on the difficulties from the combination of free surface, turbulence and complex boundary configuration, researchers made lots of efforts in the last several decades.
Many techniques have been developed and used successfully over the years.These methods can be classified into two main types: interface tracking and interface capturing.The former computes the liquid flow only, using a numerical grid that adapts itself to the shape and position of the free surface.The free surface is represented and tracked explicitly either by marking it with special marker points, or by attaching it to a mesh surface.Various surface fitting methods for attaching the interface to a mesh surface were developed during the past decades [1][2][3][4].In the interface-tracking methods, the free surface is treated as a boundary of the computational domain, where the kinematic and dynamic boundary conditions are applied.These methods cannot be used if the interface topology changes significantly [4].The second possible approach is given by the so-called interface-capturing methods [5][6][7][8][9][10].One may also divide these methods into three groups: the Lagrangian method, the Eulerian method and the arbitrary Lagrangian -Eulerian method.Generally speaking, the capturing method has advantages over the tracking methods from the viewpoint of generality although they are not as accurate as the latter one.
The Volume of Fluid (VOF) method was proposed and widely used for numerous problems due to its small requirement on computer storage and inherent abilities to dealing with complex free surface process, such as large deformation, breaking, overlapping and splashing [11][12][13].One of the main concerns with the interface capturing approaches has been sustaining global mass conservation in long-time integrations [14].In order to make up the limitation of VOF mentioned above, some researchers have made various efforts [14][15][16][17].
Turbulence is an important characteristic of fluid flows and remains still a great challenge for computational fluid dynamics [18].Nowadays, most of the existing numerical models for free surface flows employ the Reynolds-Averaged Navier-Stokes (RANS) procedure to account for the turbulence.The increasing power of computers allows a more advanced method to be integrated in the VOF based on numerical models.As an alternative of Di-rect Numerical Simulation (DNS) [6], the Large Eddy Simulation (LES) technique is a promising approach [19][20][21].
Lu et al. [22] proposed a two-dimensional hybrid numerical model, FEM-LES-VOF, for free surface flows, which is a combination of three-step Taylor-Galerkin finite element method, large eddy simulation with the Smagorinsky sub-grid model and Computational Lagrangian-Eulerian Advection Remap Volume of Fluid (CLEAR-VOF) method.Löhner et al. [4] developed a twodimensional computer model to simulate free surface flow interaction with a moving body.VOF technique in coupling with an incompressible Euler/Navier-Stokes solver operating on adaptive, unstructured grids has been used to simulate the interactions of extreme waves and threedimensional structures [4].
To the best of our knowledge no comprehensive model has been published yet to study the free surface flow over semicircular and circular cylinder obstruction.In this article, the governing equations are numerically solved in the domain by the control volume approach based on the SIMPLE technique.In this study, effects of semicircle and circle obstruction on free surface flow are investigated and the final profile of free surface in the vicinity of obstruction is comprised with the final profile in a duct without obstruction.
The remainder of the paper is organized as follows: Section 2 summarizes the Problem formulation and numerical modelling; Sections 3 describes the Numerical method; code validation is shown in Section 4; Sections 5 discusses the results and finally, some conclusions are given in Section 6.

Computational domain
Consider a two-dimensional computational domain with the opening length L = 0.6 m and, height H = 0.125 m through which water flows in to the domain.Inlet velocity of the water is considered U in = 0.075 m/s.In the middle of the domain, an obstruction is situated and the effect of this obstruction is seen on free surface flow.As it could be seen later, two types of obstruction are investigated in this study a semi-circular and a fully circular cylinder.Diameter of the circle is considered to be equal to 0.06 m.In case of fully circular cylinder, the height of computational domain is considered to be 0.15 m.
For the computational domain the grid independency test was performed with 4 grid levels 30000, 40000, 50000 and 60000.In all cases a C type mesh is used for simulating as shown in Fig. 1.Also a unified grid clustering factor near to the circle has been used.It is observed that a 40000 grid, adequately resolves the physics of the problem and further refinement is unnecessary.

Governing equations
The RANS equations can be written for incompressible flow in tensor notation as follows Upper case denotes ensemble-mean quantities and lower case fluctuating or turbulence quantities; P is pressure and ρ is density.An over bar is used to denote Reynolds averaging.(5) where Re t and y * are turbulent Reynolds number and dimensionless distance from the wall, respectively.k is the turbulent kinetic energy and  is the dissipation rate.

Transport equations for
k and  are where the turbulent kinetic energy production rate is Application of the Boussinesq assumption, Eq. ( 3), yields and mean strain tensor is defined by 1 2

Free surface modelling
f which is defined as fluid volume indicator across the interface is a function of space and time, which is advected in the flow field according to The fluid volume indicator function f is set to 1 in the liquid region and 0 in the gaseous region.For the interface between liquid and gas, f is between 0 and 1.

Boundary conditions
At the inflow boundary, with a uniform velocity U as an input, the turbulence kinetic energy is set to a low value, At an outflow boundary the convective boundary condition may be applied for each variable  which is the Neumann boundary condition.
Enforcement of no-slip condition on the embedded boundaries is the main difficulty of the method.The values of the flow field variables have been arranged to be the solution at the embedded boundaries grid points by Great-Source-Term.

Numerical method
Our code uses Finite Volume method and the SIMPLE algorithm for discretizing the governing equa-tions of flow and resolving the pressure-velocity coupling system.In addition, all the variables are stored in the same nodes by using a collocated grid.Collocated grid has various advantages over the staggered grid, e.g. the control volumes for all variables coincide with the boundaries of the solution domain, facilitating the specification of boundary conditions, also geometrical data need to be calculated only for one set of the control volumes.The diffusion term of the Eqs.6 and 7 is discretized using a central difference algorithm.

Validation
In order to assess the accuracy of our numerical procedure for free surface flow, we have tested our algorithm with dam breaking problem.The results obtained for the horizontal location of the free surface along the bottom wall are compared to the experimental values of Martin and Moyce [24] and the numerical results obtained by Hansbo [25], Walhorn [26] and Koelke [27] in Fig. 2. As shown in this figure, the agreement is very good.Fig. 3 shows the comparison of the free surface profile in the vicinity of the semicircular obstruction at the channel from present work and those of the Ashgriz et al. [28] and also shows a very good agreement.

Results and discussion
At the first step, flow in a wide channel is simulated.Because of the width of the channel, we can assume that the flow is 2D.Fluid enters the computational domain with a constant velocity U in and flows over the bottom surface of the channel to reach the exit boundary at x = 0.6 m.Although the channel height is assumed to be h = 0.125 m, but the width of the input region is limited to h = 0.075 m.Gradually fluid passes the channel length and reaches the exit boundary approximately at the t = 0.45 s.Afterwards, the free surface moves up at the channel to reach its final height which can be reached nearly at t = 2.0 s and then flow is steady.Other than input region, the height of the fluid is uniform all over the channel and the free surface is very smooth, at steady state condition.Because of the absence of any obstacle in the flow path, all stream lines would be smooth.Fig. 4 shows the detail flow characteristics of this problem at 4 different time steps.
As a new case, we put an obstacle in the middle of the channel to study the flow pattern over it.A semicircle is chosen as the obstruction.In the present case, the upstream Froude number is 0.373 which is less than 1.0 and leads to a sub-critical flow, which means the reflected waves due to the existence of obstruction will be generated and travel upstream.On the other hand, a supercritical flow takes place downstream and converges to a steady state condition quickly with respect to the sub-critical flows upstream because the disturbance may be washed out quickly by the higher-velocity flow downstream.The present problem was also examined by Lu et al. [22] numerically and by Forbes et al. [17] with an analytical method based on the potential theory which neither the evolution of free surface nor the turbulence was presented.At the early time steps and before hitting the obstacle, the situation is similar to the dam breaking and previous problem (Fig. 4).Water strikes the semicircle at about t = 0.25 s, passes over it and flows to the end of the channel which is reached at t ≈ 0.5 s.When the flow passes the top of the semicircle, a jet or separation of the fluid can be seen which causes some air to be trapped at the tail of the semicircle.This trapped air moves slowly and leaves the domain at t = 1.8 s after which, the free surface at the downstream is smoothened.A new vortex appears at the upstream and just at the corner of the semicircle.This vortex develops further and enlarges significantly as shown in Fig. 5.
conditions are similar to that of the channel without any obstacle (Fig. 4).As the flow reaches the circle, bubbles of air are trapped by the fluid at the lower corner of the obstacle.Fluid jumps over the circle like a cascade and passes it (see Fig. 6 at time t = 0.3, 0.45 s) and hence a large amount of air is trapped by this cascade.
Simultaneously the trapped air in the upstream makes a wake that moves to the entrance area and after reaching this, wavy free surface is damped and gradually the trapped bubbles of air vanish.At t = 1.9 s the conditions in the upstream, reach the steady state and the free surface upstream the circle obtains its final smooth profile.
At the downstream, the fluid level moves up and the captured air moves forward and exit from the flow domain, like a vortex.Flow field reached its steady state after t = 4 s as it can be seen in the Fig. 6.Despite of the previous problem, the captured air does not leave the right side of the circle completely which may be a result of the size of the captured bubble, as you can see at t = 20 s (Fig. 6) which is greater than the steady state time.Stream lines of the last three problems are compared in the Fig. 7.The three Stream line patterns are similar, but a small difference could be seen.In the second scheme (flow over a semicircle) a small vortex exists at the right side of the obstacle, but at the third scheme, in addition to the downstream vortex which is greater than that of the previous problem, a new vortex is seen at the upstream of the circle which does not vanish.
Fig. 8 shows the height of the fluid along the channel at three different previous problems.As the obstacle size grows, the level of the free surface at the upstream grows, and the height of the free surface at the downstream, decreases.Since the entered mass flow is constant in three cases, flow average velocity and Froud number at the exit boundary increases with increasing the size of the obstruction.

Conclusions
A numerical model for simulating turbulent free surface flow is developed in this article, which is based on the finite volume on collocated grid using the scheme suggested by Rhie and Chow.Furthermore, Reynolds-Averaged Navier-Stokes (RANS) equations and Volume of Fluid (VOF) method is used for modelling the turbulence effects and free surface capturing, respectively.
The comparisons and validations show that the proposed FVM-RANS-VOF model produces reliable and reasonable predictions for the fluid flows involving complicated interface and boundary configurations.
Present results show that as the obstacle size grows, the level of the free surface at the upstream grows, and the height of the free surface at the downstream, decreases.Since the entered mass flow is constant in three cases, the flow average velocity at the exit plane increases with increasing the size of the obstruction.This means that Froud number in downstream increases with increasing the obstacle size.

Fig. 2 Fig. 3
Fig. 2 Comparison of the horizontal location of the free surface along the bottom wall obtained from present study and the literature

Fig. 4
Fig. 4 Water free surface profile at different time steps for flow in a channel without any obstacle

7 Fig. 8
Fig. 8 Final free surface profile in the vicinity of obstruction The momentum equations are not closed.
ij uu  .i  is the body force per unit volume on each direction.