Numerical simulation of an air-core vortex at a hydraulic intake using OpenFOAM

A vortex is a ubiquitous everyday phenomenon that is observed in nature and it is formed due to the rotational motion of fluid around an axis perpendicular to the free surface. Free surface vortices are a common unwanted occurrence at hydraulic intakes which can cause serious detrimental impacts on mechanical devices such as turbines and pumps. In this paper, an experimentally observed air-core vortex is numerically simulated using the OpenFOAM LTSInterFoam solver. The LTSInterFoam solver has hitherto been mainly used for hydrodynamic studies relating to ship manoeuvrability by researchers. This solver uses a local time stepping approach to speed up convergence towards steady state conditions thus overcoming some of the challenges associated with the use of the conventional interFoam solver for the simulation of free surface vortices. The Shear Stress Transport (SST) k − ω Model was used for the simulation. There was generally good agreement when results from the study were compared with other vortex-related analytical models and experimental data. Overall, the study concludes that the OpenFOAM LTSInterFoam solver is capable of simulating free surface vortices at hydraulic intakes. However, being a steady state solver, the solver cannot account for the transient process involved in the evolution of free surface vortices. Crown Copyright © 2020 Published by Elsevier B.V. on behalf of African Institute of Mathematical Sciences / Next Einstein Initiative. This is an open access article under the CC BY license. ( http://creativecommons.org/licenses/by/4.0/ )


INTRODUCTION
A vortex is described as a rotating motion of fluid about a common centre [19] . Free surface vortices are a common unwanted occurrence at hydraulic intakes [ 16 , 36 , 37 ]. Hydraulic intakes are designed bearing in mind as much as possible the need to permit a smooth transition from open channel flow into pressurised flow thus ensuring a steady flow to mechanical devices such as turbines and pumps [35] .
If this is not properly achieved, free surface vortices may develop and thus the efficiency of devices such as pumps and turbines could be negatively impacted [3] . It is reported that a 1% volume air entrainment as a result of free surface vortices can be responsible for as much as 15% efficiency reduction in centrifugal pumps [6] . Aside from the loss in efficiency of hydraulic devices, other detrimental impacts such as debris entrainment, local corrosion damage, cavitation, reduced discharge into intakes, surging, and associated reduction in life-span of turbines and pumps could also be the result of the occurrence of free surface vortices [ 16 , 35 ]. The major causes of vortices include high velocity gradients, rotational wakes and non-uniform approach flow [7] . In hydropower plants, especially low-head plants, low water levels could initiate the development of vortices [1] .
Vortex dynamics is generally a complex phenomenon. The complexity of the vortex dynamics stems from the fact that its occurrence is affected by a number of factors such as density variation between air and water, surface tension forces, gravity and buoyancy [1] . In response to this and as a result of the intriguing nature of its occurrence and dynamics, the study of vortices at hydraulic intakes has become a popular subject among engineers and researchers with several theoretical and experimental studies by various authors [ 23 , 31 , 40 ]. The formation of free surface vortices is dependent on a number of factors. Prominent among these factors are the Froude number ( F r ) , Critical submergence ( S c ) , Circulation number ( ) and Reynolds number ( R e ) . Critical submergence is the submergence below which incipient vortices are formed at intakes [17] .
In recent times and with the continuous increase in computational power, numerical simulation of vortices has gained popularity among researchers. Numerical simulations have an added advantage of overcoming many of the limitations such as cost and time constraints associated with experimental studies [43] . Whereas most researchers have utilised commercial Computational Fluid Dynamics (CFD) software tools such as FLOW-3D, STAR CCM + and ANSYS (CFX and FLUENT) in simulating free surface vortices, it appears from our literature survey that very few researchers have utilised open source CFD codes such as OpenFOAM in modelling free surface vortices at hydraulic intakes. OpenFOAM is an open source CFD tool written in C ++ library and it is based on the Finite Volume Method (FVM). The OpenFOAM CFD package comprises meshing, pre-processing and post-processing tools for solving fluid flow problems [42] .
Previous authors [ 11 , 43 ] who simulated air core-vortices in pumps using OpenFOAM, used the transient two-phase interFoam solver. The researchers reported two major challenges: development of waves at the inlet as well as long computational time even with the use of a computer cluster. According to Li et al. [18] , it is possible to simulate air-core vortices using steady state solvers in CFD. The LTSInterFoam solver, which is one of the numerous solvers in OpenFOAM, is a steady state Volume of Fluid (VOF) solver that has been used by previous authors [ 5 , 12 , 30 ] for studies on the manoeuvrability of ships and boats. The solver uses the PIMPLE algorithm alongside a local time stepping integration. The purpose of the local time stepping technique is to speed up convergence towards steady state conditions. The technique manipulates the time-step for each cell rendering it as high as possible to enhance steady state convergence [ 9 , 28 ]. For further details on the implementation of the local time stepping approach readers are referred to Espinoza et al. [9] .
The Shear Stress Transport (SST) k − ω model has been found to be suitable for numerical modelling of free surface vortices at hydraulic intakes [ 1 , 25 , 29 ], hence, its selection for the study. The model's suitability is as a result of the fact that the SST k − ω model is robust and also exhibits better performance regarding the prediction of flows at walls and adverse pressure gradient in comparison to other eddy-viscosity models ( [20][21][22] . In this regard, the study explored the use of the steady state LTSInterFoam solver in OpenFOAM along with the SST k − ω model.
In this paper, an air-core vortex observed in an experimental set-up by Sarkardeh et al. [34] is numerically simulated and validated against appropriate vortex-related analytical models and experimental data. Also, a means of applying the proposed model to a typical hydropower facility facing challenges with generating adequate power during low water levels has been highlighted. The paper is structured as follows: firstly, a discussion on some related vortex models is presented, followed by an overview of the experimental set-up and the numerical approach, a discussion of results and eventually a conclusion based on the findings of the study.

Vortex models
In order to provide validation for the CFD results, comparison with analytical and empirical models for vortex flow can be used. Accordingly, in this section, appropriate mathematical vortex models will be reviewed.
Rankine [31] proposed a model for a vortex as having a solid rotating body made up of an inner core and an outer irrotational field. The relations for the tangential velocity of the inner core and outer core vortices are presented in Eqns Vortex models proposed by other authors, Mih [23] , Vatistas et al. [39] , and Sun and Liu [37] are given in Eqns. 3 , 4 & 5 respectively.
Regarding the water surface profile associated with the formation of an air-core vortex, Hite and Mih [13] obtained a relation ( Eqn. 6 ) based on a hydrostatic pressure assumption in the vertical direction. Vatistas et al. [39] also proposed an equation ( Eqn. 7 ) for the free surface profile with the same assumption of hydrostatic pressure distribution but using a different approach whereas Sun and Liu [37] proposed Eqn. 8 .
where V θ is the tangential velocity, ω is the angular velocity of the vortex centre, r is the radius, r m is the radius at maximum tangential velocity, is constant circulation of the outer zone, H r refers to the water surface elevation at r, H o is the water surface elevation at the centre, h refers to the depth of the air-core as R approaches infinity, and the normalised radius, R = r/ r m .

Experimental set-up
The experimental set-up under investigation is reported in Sarkardeh et al. [34] and it is illustrated in Fig. 1 . The closedcircuit experimental set-up comprised of a constant elevation sump with a 16 cm diameter (D ) intake which is connected to a pump. Apart from providing the intake flow rate, the pump also enables adjustment of the water level in the tank by varying the volume of water in circulation. The walls of the prototype, intake geometry and tunnel were made of clear Perspex to aid visualization of the flow phenomenon. The set-up was made to run for long hours until steady conditions were obtained at a constant water level. The results presented here are for flow conditions at the relative Critical submergence, S/D = 1 . 5 and intake Froude number, F r = 0.6 where an air-core vortex was observed upstream of the intake.

Governing equations
The governing Reynolds-averaged Navier-Stokes (RANS) equations underlining the LTSInterFoam solver are the Continuity ( Eqn. 9 ), Transport of phase-fraction ( Eqn. 10 ), and Momentum equations ( Eqn. 11 ). The RANS equations implement the Reynolds decomposition of flow parameters such that u = ū + u and p = p + p .
where u is the velocity, ū is the time-averaged velocity, u is the velocity fluctuation, p is the pressure, p is the timeaveraged pressure, p is the pressure fluctuation, α is the volume fraction which ranges from 0 to 1 (with cells fully filled with gas assigned 0 whilst cells fully filled with liquid are assigned a value of 1; the location of the interface is captured by cells with volume fraction 0 <α< 1), g is the acceleration due to gravity, μ e f f is the effective viscosity, F s is the surface tension force.
The density ρ and viscosity μ of the fluid are computed from the relations in Eqns. 12 and 13 . The subscripts l and g in the equations refer to liquid and gas respectively. μ In the SST k − ω model, the viscosities are calculated from Eqns. 14 , 15 & 16 .
μ e f f = μ + μ turb (14) μ turb = where μ refers to the viscosity of the fluid, μ turb is the turbulent viscosity, k is the turbulent kinetic energy, ω is the specific rate of dissipation of turbulence kinetic energy, a 1 is the experimental constant which equals 0.31, S refers to the strain rate magnitude, F 2 represents the blending function which is computed using Eqn. 16 and β * = 0 . 09 [ 22 , 24 ]. Additional transport equations for k and ω can be obtained from Menter et al. [22] . In order to address challenges associated with numerical diffusion at the interphase as well as attain a sharp interphase, OpenFOAM implements an artificial velocity u r which modifies Eqn. 10 to obtain Eqn. 17 [15] .
Numerical solution procedure The finite volume method is used to discretise the governing differential equations. As already indicated, the solver uses the Local Time Stepping (LTS) scheme for temporal discretisation in order to attain a rapid steady state solution along with the PIMPLE algorithm. Regarding the convection schemes, the linearUpwind grad(U) was used for U whilst the linearUpwind limitedGrad scheme was selected for the divergence related to k and ω. For the divergence of alpha, div(phi,alpha), and for the compression of the interface, div(phirb,alpha), the vanLeer and the linear schemes were respectively used. For the Laplacian schemes, the Gauss linear corrected scheme was implemented.

CFD Modelling
The geometry of the model was built in CAD and meshed using the open source meshing tool cfMesh. In view of the computational cost and the results of the output, a meshed geometry of approximately 9.6 million cells which are mainly hexahedra cells was used for the simulation. The region without water flow was assigned a coarse (average cell size 0.08m) mesh, while medium (average cell size 0.02m) and fine (average cell size 0.01m) meshes were provided for regions with flow of water and the air-water interface respectively. The boundary conditions used for the simulation are provided in Table 1 .
In order to impose a constant water level in the numerical study, the groovyBC boundary condition was implemented at the inlet. This functionality is contained in the swak4Foam library and it enables the setting-up of arbitrary boundary condition based on expressions. The LTSInterFoam solver with the SST k − ω model was used for the simulation. The simulation was run on 28 cores (2 × 14 core Intel Xeon E5-2660 2.00GHz) of a desktop workstation with 256GB RAM. The simulation time was approximately 18 hours. A normalised residual value of order 1 × 10 −5 was used to assess simulation convergence. The numerical results were post-processed in ParaView for visualisation after the simulation converged after 50 0 0 iterations. After the end of the simulation, an air-core vortex was formed in front of the intake reproducing the experimental observation.

The case of the Akosombo Hydropower Plant
The major problem facing the Akosombo Hydropower Plant is the inability to generate adequate power during low water levels, a situation which is caused by low inflows. This recurrent problem often compels operators of the dam facility to operate less than the installed number of turbine units resulting in national electricity rationing in certain years including 1983, 1986, 1998, 1999, 20 06, 20 07, 2015 and 2016 [ 2 , 8 ].
According to Obeng and Fiagbe [27] , the estimated average plant efficiency below the minimum operating level of the Akosombo Hydropower Plant is 84.5%. Meanwhile, the turbines at the facility are expected to operate at 95% efficiency even at low water levels [10] . This implies that the possible formation of free surface vortices at low water levels could be responsible for the reduced plant efficiency. Due to the peculiar feature of this hydropower plant, a means of applying the proposed model to the case of the Akosombo Hydropower Plant has been discussed at the results and discussions section.

RESULTS AND DISCUSSIONS
The results and discussions presented below relate to the vortex identification and visualisation as well as the tangential velocity and free surface profile plots. Fig. 2 shows an air-core vortex identified using an isosurface of α = 0.96. The threshold value of α = 0.96 is in agreement with the threshold value that was used by Yamasaki et al. [43] to detect free surface vortices in a numerical simulation involving a pump sump. Inferring from Fig. 2 , the vortex has a conical shape that extends from the free surface towards the intake as observed in the experimental work by Sarkardeh et al. [34] . The streamlines associated with the vortex formation are illustrated in Fig. 3 and they were drawn using the Line Integral Convolution technique [4] . The streamlines descend and then start to converge spirally at the centre of the vortex and then towards the intake in a similar manner to that observed by Suerich-Gulick et al. [35] .

Tangential velocity profile
A comparison between the numerical results and other proposed vortex models [ 23 , 31 , 37 , 39 ] as well as experimental data [ 32 , 37 , 40 ] in terms of the normalised tangential velocity is illustrated in Fig. 4 . From Fig. 4 , there is an observed general increase in tangential velocity up to the inner core radius in order to conserve angular momentum. It is a generally known fact that the Rankine model over-estimates the tangential velocity at R = R m and it is also non-differentiable at this point [41] . The model proposed by Mih [23] exhibits a slight deviation around the vortex core. This is as a result of the fact that the maximum tangential velocity for the Mih's model occurs at R = 0 . 71 instead of R = 1 [37] . The rest of the models, as well as results from the numerical simulation, show good agreement with the experimental data.

Water surface profile
As reported by Hite and Mih [13] and demonstrated by Fig. 5 , the water surface profile associated with the free surface vortex exhibits a hyperbolic drop from the far-field to the half-depth of the water surface depression. Beyond that, there is a parabolic decrease up to the vortex centre. From Fig. 5 , the numerical prediction performed well with free surface vortex models reported by Hite and Mih [13] , Sun and Liu [37] and Vatistas et al. [39] as well as experimental data by Julien [14] , Newman [26] , Roschke [33] and Vatistas [38] available in the literature.

A proposed application of the model
The use of the OpenFOAM LTSInterFoam solver to predict the occurrence of free surface vortices at hydropower intake systems as highlighted in this paper will be crucial for addressing flow problems at dam intakes. The model could be used for the assessment of the impacts of the formation of free surface vortices on energy generation from hydropower facilities facing challenges with perennial low water levels such as the Akosombo Hydropower Plant in Ghana. Such a numerical simulation of a hydropower facility like the Akosombo Hydropower Plant will require the construction and experimental study of a reduced-scale physical hydraulic model of the dam intake structure. A numerical model of the physical hydraulic model will then be set-up and simulated as discussed in this paper. The results of the simulation will then be validated against data obtained from the experimental scale model. Information from the validated numerical model can then be used to optimise relevant dam structures in order to attenuate the impact of the occurrence of vortices at low water levels if  required. Such remedial interventions could include the recommendation and installation of appropriate anti-vortex devices to attenuate the development of free surface vortices.

CONCLUSIONS
In this paper, an experimentally observed air-core vortex is numerically simulated using the OpenFOAM LTSInterFoam solver. This solver uses a local time stepping approach to speed up convergence towards steady state conditions thus overcoming some of the challenges associated with the use of the conventional interFoam solver for the simulation of free surface vortices. The SST k − ω model was used for the simulation. There was generally a good agreement when results from the study were compared with other vortex-related analytical models and experimental data available in literature. Overall, the study concludes that the OpenFOAM LTSInterFoam solver is capable of simulating free surface vortices at hydraulic intakes. However, being a steady state solver, the solver cannot account for the transient process involved in the evolution of free surface vortices.
It is recommended that a transient solver such as the interFoam solver be utilised in future studies so that the transient process involved in the evolution of free surface vortices could be accounted for. Also, future studies could look into the possibility of applying the model to a typical hydropower facility such as the Akosombo Hydropower Plant to assess the impact of the formation of free surface vortices at low water levels. Information from such a study could be used to optimise energy generation from the hydropower plant.

Declaration of Competing Interest
None.