Numerical Simulation of an Air-Core Vortex and Its Suppression at an Intake Using OpenFOAM

: A common challenge faced by engineers in the hydraulic industry is the formation of free surface vortices at pump and power intakes. This undesirable phenomenon which sometimes entrains air could result in several operational problems: noise, vibration, cavitation, surging, structural damage to turbines and pumps, energy losses, e ﬃ ciency losses, etc. This paper investigates the numerical simulation of an experimentally observed air-core vortex at an intake using the LTSInterFoam solver in OpenFOAM. The solver uses local time-stepping integration. In simulating the air-core vortex, the standard k − ε , realizable k − ε , renormalization group (RNG) k − ε and the shear stress transport (SST) k − ω models were used. The free surface was modelled using the volume of ﬂuid (VOF) model. The simulation was validated using a set of analytical models and experimental data. The SST k − ω model provided the best results compared to the other turbulence models. The study was extended to simulate the e ﬀ ect of installing an anti-vortex device on the formation of a free surface vortex. The LTSInterFoam solver proved to be a reliable solver for the steady state simulation of a free surface vortex in OpenFOAM. workstation The simulation convergence was assessed using a normalized residual value of order 1 × 10 − 5 . Visualisation of the numerical results was performed in ParaView.


Introduction
A common challenge faced by engineers in the hydraulic industry is the formation of free surface vortices at pump and power intakes, a phenomenon which may entrain air [1][2][3]. This undesirable phenomenon often results in several operational problems: noise, vibration, cavitation, energy losses, efficiency losses, etc. [4,5]. Free surface vortices have also seen beneficial applications in water vortex hydropower plant systems [6][7][8], water treatments [9], vortex drop shafts in sewer systems [10], vortex settling basins [11] and vortex confined chambers [12].
Free surface vortices involving incipient air-entrainment occur at low submergence depth, a depth often referred to as the critical submergence [2]. This implies that ensuring adequate submergence presents an efficient way to avert the occurrence of free surface vortices. Apart from ensuring adequate submergence during the design of intakes, the use of anti-vortex devices also provides another cost-effective means of suppressing the formation of free surface vortices [5,13].
Refs. [14][15][16][17] proposed analytical models to describe the key vortex characteristics. Several validation test cases on these analytical studies also exist in the literature [18][19][20]. Some researchers [21][22][23] have proposed mathematical models and design plots which predict Ref. [14] describes a vortex as being a solid rotating body comprising of an inner core and an outer free vortex. The equations for the tangential velocity of the inner core and outer core vortices are given by Equations (1) and (2), respectively.
where V θ represents the tangential velocity; ω refers to the angular velocity of the vortex centre; r refers to the radius; r m is the radius at maximum tangential velocity; Γ refers to the constant circulation of the outer zone; and the normalised radius, R = r/r m .

Experimental Set-Up
The first set of numerical investigations involving an intake without an anti-vortex device was based on the experimental work by [39], as illustrated in Figure 1. The experimental work entails the formation of an air-core vortex at an intake. To be able to study the free surface vortex occurrence in the set-up shown in Figure 1, water was made to circulate in a closed-loop set-up, such that a constant water level was maintained at the elevation sump. The intake had a diameter (D) of 16 cm and the water level in the sump was adjusted by using the pump to vary the volume of water in circulation. Visualisation of flow within the prototype sump, intake geometry and the tunnel was made possible by the use of clear Perspex in these areas of interest. The set-up was initially made to run for several hours and readings were only taken when steady conditions had been attained at a constant water level. A stable air-core vortex was observed when the relative critical submergence, S/D, and intake Froude number, F r , were set at 1.5 and 0.6, respectively.
where represents the tangential velocity; refers to the angular velocity of the vortex centre; refers to the radius; is the radius at maximum tangential velocity; refers to the constant circulation of the outer zone; and the normalised radius, = ⁄ .

Experimental Set-Up
The first set of numerical investigations involving an intake without an anti-vortex device was based on the experimental work by [39], as illustrated in Figure 1. The experimental work entails the formation of an air-core vortex at an intake. To be able to study the free surface vortex occurrence in the set-up shown in Figure 1, water was made to circulate in a closed-loop set-up, such that a constant water level was maintained at the elevation sump. The intake had a diameter ( ) of 16 cm and the water level in the sump was adjusted by using the pump to vary the volume of water in circulation. Visualisation of flow within the prototype sump, intake geometry and the tunnel was made possible by the use of clear Perspex in these areas of interest. The set-up was initially made to run for several hours and readings were only taken when steady conditions had been attained at a constant water level. A stable air-core vortex was observed when the relative critical submergence, ⁄ , and intake Froude number, , were set at 1.5 and 0.6, respectively. The second set of numerical investigations, which sought to assess the impact of installing an anti-vortex device, was based on the experimental work by [31]. Both studies implemented similar experimental set-ups but with an anti-vortex plate which has a 66% perforated opening installed at the intake of the experimental set-up by [31]. The plate had a thickness of 2 mm and measured 1.5D Figure 1. Experimental set-up by [39].
The second set of numerical investigations, which sought to assess the impact of installing an anti-vortex device, was based on the experimental work by [31]. Both studies implemented similar experimental set-ups but with an anti-vortex plate which has a 66% perforated opening installed at the intake of the experimental set-up by [31]. The plate had a thickness of 2 mm and measured 1.5D by 1D in length and width, respectively. The discharge flow and relative submergence for this experimental case were set at 0.015 m 3 /s and 1.5, respectively.

Governing Equations
The governing Reynolds-Averaged Navier-Stokes (RANS) equations that underline the LTSInterFoam solver are the continuity, transport of phase-fraction and the momentum equations which are given by Equations (6)-(8), respectively. The flow parameters are decomposed using the Reynolds decomposition approach such that u = u + u and p = p + p .
where u refers to the velocity; u is the time-averaged velocity; u refers to the velocity fluctuation; p is the pressure; p represents the time-averaged pressure; p is the pressure fluctuation; α represents the volume fraction which ranges from 0 to 1 (with cells filled with gas given a value of 0 whilst cells filled with liquid are assigned a value of 1); g is the acceleration due to gravity; and µ e f f is the effective viscosity. Equations (9) and (10) are used to compute the density ρ and effective viscosity µ e f f of the fluid, respectively. The subscripts l and g in the equations represent liquid and gas respectively.

OpenFOAM Application and Procedures
The governing differential equations are discretized using the finite volume method (FVM) along with the local time stepping (LTS) scheme for temporal discretization in order to ensure a rapid steady-state solution. The pressure-velocity equations are resolved using the PIMPLE algorithm. The linearUpwind grad (U) was selected for the convection of U. The divergence related to k and ω was computed using the linearUpwind limitedGrad scheme, whilst for that of epsilon, the upwind scheme was selected. For the divergence of alpha, div (phi,alpha), and the compression of the interface, div(phirb,alpha), the van Leer and the linear schemes were applied, respectively. The Gauss linear corrected scheme was used for the Laplacian schemes.
The building and meshing of the geometries were done in CAD and cfMesh, respectively. The cells were mainly hexahedra, comprising of coarse (average cell size 0.08 m), medium (average cell size 0.02 m) and fine (average cell size 0.01 m) cells. The meshed geometry with its associated boundaries is shown in Figure 2.
In CFD analysis, boundary conditions are used to specify the spatial or temporal variable values or behavior required to produce a unique solution. The boundary conditions used for the simulation are provided in Table 1 where U, p_rgh, alpha.water, nut, k and omega refer to the velocity, reduced pressure, phase fraction, turbulent viscosity, turbulence kinetic energy and specific rate of dissipation of turbulence kinetic energy, respectively. The fixedValue represents a Dirichlet Boundary Condition which is specified by the user. The pressureInletOutletVelocity boundary condition applies a zero-gradient for outflow, whilst the inflow velocity is the patch-face normal component of the internal-cell value. The fixedFluxPressure provides an adjustment to the pressure gradient, such that the flux on the boundary is the one specified by the velocity boundary condition. The totalPressure boundary condition is computed as static pressure reference plus the dynamic component due to velocity. The inletOutlet boundary condition provides a zero-gradient outflow condition for a fixed value inflow. The kqRWallFunction is the wall function for the turbulence kinetic energy while the nutkRoughWallFunction is the rough wall function for kinetic eddy viscosity and omegaWallFunction represents the wall function for frequency. The groovyBC boundary condition (contained in the swak4Foam library) was used to impose a constant water level.
Fluids 2020, 5, x FOR PEER REVIEW 5 of 12 the internal-cell value. The fixedFluxPressure provides an adjustment to the pressure gradient, such that the flux on the boundary is the one specified by the velocity boundary condition. The totalPressure boundary condition is computed as static pressure reference plus the dynamic component due to velocity. The inletOutlet boundary condition provides a zero-gradient outflow condition for a fixed value inflow. The kqRWallFunction is the wall function for the turbulence kinetic energy while the nutkRoughWallFunction is the rough wall function for kinetic eddy viscosity and omegaWallFunction represents the wall function for frequency. The groovyBC boundary condition (contained in the swak4Foam library) was used to impose a constant water level. The simulations were run on 28 cores (2 × 14 core Intel Xeon E5-2660 2.00 GHz) of a desktop workstation which has 256 GB RAM. The simulation convergence was assessed using a normalized residual value of order 1 × 10 . Visualisation of the numerical results was performed in ParaView.

Visualisation of the Air-Core Vortex
Using an isosurface of α = 0.96, the stable air-core vortex was depicted as a conical drop extending from the free surface towards the intake, as illustrated in Figure 3. In fluid flow, the presence of a vortex is often depicted by a roughly circular or spiral pattern of streamlines [48,49]. The streamlines, which were drawn with surface line integral convolution (Figure 4), exhibit a spiral pattern, thus confirming the presence of a vortex. The velocity vector plot ( Figure 5) illustrates how the fluid flows spirally into the intake. From the observations in Figures 4 and 5, the fluid flow can  The simulations were run on 28 cores (2 × 14 core Intel Xeon E5-2660 2.00 GHz) of a desktop workstation which has 256 GB RAM. The simulation convergence was assessed using a normalized residual value of order 1 × 10 −5 . Visualisation of the numerical results was performed in ParaView.

Visualisation of the Air-Core Vortex
Using an isosurface of α = 0.96, the stable air-core vortex was depicted as a conical drop extending from the free surface towards the intake, as illustrated in Figure 3. In fluid flow, the presence of a vortex is often depicted by a roughly circular or spiral pattern of streamlines [48,49]. The streamlines, which were drawn with surface line integral convolution (Figure 4), exhibit a spiral pattern, thus confirming the presence of a vortex. The velocity vector plot ( Figure 5) illustrates how the fluid flows spirally into the intake. From the observations in Figures 4 and 5, the fluid flow can be described as moving spirally from the water surface towards the intake and also rotating around the vortex axis. The plot depicted in Figure 6 shows the contour plot of the normalized velocity field expressed as a percentage, (V XZ /V)%, where V XZ refers to the velocity in x-z plane whilst V is the intake flow velocity. The plot is made at a horizontal depth of D (16 cm) above the axis of the intake. The contours of the velocity field indicate the boundaries of the conical flow towards the intake [39].
Fluids 2020, 5, x FOR PEER REVIEW 6 of 12 be described as moving spirally from the water surface towards the intake and also rotating around the vortex axis. The plot depicted in Figure 6 shows the contour plot of the normalized velocity field expressed as a percentage, ( ⁄ )%, where refers to the velocity in x-z plane whilst is the intake flow velocity. The plot is made at a horizontal depth of (16 ) above the axis of the intake.
The contours of the velocity field indicate the boundaries of the conical flow towards the intake [39].   Fluids 2020, 5, x FOR PEER REVIEW 6 of 12 be described as moving spirally from the water surface towards the intake and also rotating around the vortex axis. The plot depicted in Figure 6 shows the contour plot of the normalized velocity field expressed as a percentage, ( ⁄ )%, where refers to the velocity in x-z plane whilst is the intake flow velocity. The plot is made at a horizontal depth of (16 ) above the axis of the intake.
The contours of the velocity field indicate the boundaries of the conical flow towards the intake [39].

Assessment of Different Turbulence Models
During the numerical simulation, the authors observed that the use of the standard − and RNG − models resulted in the development of flow instabilities at the free surface, especially at the upstream end. This phenomenon was particularly pronounced at the high-pressure gradient region, thus causing the numerical simulation to fail abruptly. This situation may be linked to the fact that adverse pressure gradient conditions often cause the standard − and RNG − models to exhibit poor performances [50].
On the other hand, the realizable − , − and the SST − models successfully predicted the occurrence of the air-core vortex. The realizable − model outperformed the other

Assessment of Different Turbulence Models
During the numerical simulation, the authors observed that the use of the standard − and RNG − models resulted in the development of flow instabilities at the free surface, especially at the upstream end. This phenomenon was particularly pronounced at the high-pressure gradient region, thus causing the numerical simulation to fail abruptly. This situation may be linked to the fact that adverse pressure gradient conditions often cause the standard − and RNG − models to exhibit poor performances [50].
On the other hand, the realizable − , − and the SST − models successfully predicted the occurrence of the air-core vortex. The realizable − model outperformed the other

Assessment of Different Turbulence Models
During the numerical simulation, the authors observed that the use of the standard k − ε and RNG k − ε models resulted in the development of flow instabilities at the free surface, especially at the upstream end. This phenomenon was particularly pronounced at the high-pressure gradient region, thus causing the numerical simulation to fail abruptly. This situation may be linked to the fact that adverse pressure gradient conditions often cause the standard k − ε and RNG k − ε models to exhibit poor performances [50].
On the other hand, the realizable k − ε, k − ω and the SST k − ω models successfully predicted the occurrence of the air-core vortex. The realizable k − ε model outperformed the other family of k − ε models (standard k − ε and RNG k − ε) due to its ability to perform well with rotation and separation flow, a situation which is linked to the inclusion of a realizability constraint on the predicted stress tensor of the realizable k − ε model [43].
However, regarding the realizable k − ε model, we observed an amount of unusual free surface deformation along the walls, although the simulation progressed successfully. This observation was generally absent when the k − ω models were used. These observations may be attributed to the fact that the k − ω models often exhibit better numerical prediction at the walls compared to the k − ε models [43,48].
A comparison of results in terms of the normalized tangential velocity distribution from the three (3) different turbulence models (realizable k − ε, k − ω and the SST k − ω models) has been illustrated in Figure 7a, where R refers to the normalized radius. From Figure 7b, the three (3) turbulence models were assessed by comparing the normalized tangential velocity distribution of the turbulence models with some selected analytical vortex models by [14,17,20,38] as well as other related experimental data by [12,20,51]. The plot showed a generally fair agreement with the realizable k − ε model and the k − ω model. However, an amount of over-prediction of the tangential velocity occurred at the outer core of the vortex. Furthermore, at the point of maximum tangential velocity (R = 1), the realizable k − ε model over-predicted the tangential velocity. The ability of the k − ω model to provide a better prediction at the boundary or near-wall layers could be the reason why the model outperformed the realizable k − ε model [43,48]. The best agreement was, however, obtained from the SST k − ω model, which outperformed the k − ω model, and this could be attributed to the fact that even though the k − ω model is suitable for moderate adverse pressure gradients, the model could exhibit some amount of failure with pressure-induced separation [52], which includes the formation of an air-core vortex. Furthermore, the SST k − ω model has been found to be generally robust, with an enhanced performance in terms of flow at walls, adverse pressure gradients as well as flow separation [46,47,[51][52][53][54]. From the illustration in Figure 7b, the ranking of the performance of the three turbulence models beginning from the most suitable model is, thus, SST k − ω model, k − ω model and realizable k − ε model. family of − models (standard − and RNG − ) due to its ability to perform well with rotation and separation flow, a situation which is linked to the inclusion of a realizability constraint on the predicted stress tensor of the realizable − model [43]. However, regarding the realizable − model, we observed an amount of unusual free surface deformation along the walls, although the simulation progressed successfully. This observation was generally absent when the − models were used. These observations may be attributed to the fact that the − models often exhibit better numerical prediction at the walls compared to the − models [43,48].
A comparison of results in terms of the normalized tangential velocity distribution from the three (3) different turbulence models (realizable − , − and the SST − models) has been illustrated in Figure 7a, where refers to the normalized radius. From Figure 7b, the three (3) turbulence models were assessed by comparing the normalized tangential velocity distribution of the turbulence models with some selected analytical vortex models by [14,17,20,38,] as well as other related experimental data by [12,20,51]. The plot showed a generally fair agreement with the realizable − model and the − model. However, an amount of over-prediction of the tangential velocity occurred at the outer core of the vortex. Furthermore, at the point of maximum tangential velocity ( = 1), the realizable − model over-predicted the tangential velocity. The ability of the − model to provide a better prediction at the boundary or near-wall layers could be the reason why the model outperformed the realizable − model [43,48]. The best agreement was, however, obtained from the SST − model, which outperformed the − model, and this could be attributed to the fact that even though the − model is suitable for moderate adverse pressure gradients, the model could exhibit some amount of failure with pressure-induced separation [52], which includes the formation of an air-core vortex. Furthermore, the SST − model has been found to be generally robust, with an enhanced performance in terms of flow at walls, adverse pressure gradients as well as flow separation [46,47,51,52,53,54]. From the illustration in Figure 7b, the ranking of the performance of the three turbulence models beginning from the most suitable model is, thus, SST − model, − model and realizable − model.

Suppression of Air Entrainment Using an Anti-Vortex Device
As reported by [31] and observed in Figure 8, the use of the anti-vortex plate resulted in the suppression of the air-core vortex. The region for the anticipated vortex only showed the formation of a slight depression. Ref. [31] observed a reduction of a Type 6 vortex (a stable air-core vortex) to a Type 1 vortex (coherent surface swirl) when a perforated plate serving as an anti-vortex device was installed. The anti-vortex device suppressed the formation of the air-core vortex by inducing turbulence and friction, which succeeded in cutting the path of the vortex core, thus ensuring a relatively smoother flow towards the intake [31]. The formation of the suppressed vortex is illustrated in Figure 8.

Suppression of Air Entrainment Using an Anti-Vortex Device
As reported by [31] and observed in Figure 8, the use of the anti-vortex plate resulted in the suppression of the air-core vortex. The region for the anticipated vortex only showed the formation of a slight depression. Ref. [31] observed a reduction of a Type 6 vortex (a stable air-core vortex) to a Type 1 vortex (coherent surface swirl) when a perforated plate serving as an anti-vortex device was installed. The anti-vortex device suppressed the formation of the air-core vortex by inducing turbulence and friction, which succeeded in cutting the path of the vortex core, thus ensuring a relatively smoother flow towards the intake [31]. The formation of the suppressed vortex is illustrated in Figure 8.

Suppression of Air Entrainment Using an Anti-Vortex Device
As reported by [31] and observed in Figure 8, the use of the anti-vortex plate resulted in the suppression of the air-core vortex. The region for the anticipated vortex only showed the formation of a slight depression. Ref. [31] observed a reduction of a Type 6 vortex (a stable air-core vortex) to a Type 1 vortex (coherent surface swirl) when a perforated plate serving as an anti-vortex device was installed. The anti-vortex device suppressed the formation of the air-core vortex by inducing turbulence and friction, which succeeded in cutting the path of the vortex core, thus ensuring a relatively smoother flow towards the intake [31]. The formation of the suppressed vortex is illustrated in Figure 8.

Conclusions
Optimizing the intake of a hydropower plant often involves the suppression or prevention of the formation of free surface vortices whose occurrence could jeopardize the operations at the facility. In general, numerical simulations have been found to be relevant to such interventions by providing a key complement to experimental studies. In this regard, the study investigated the numerical simulation of an experimentally observed air-core vortex using different turbulence models as well as the suppression of the air-core vortex at an intake. Based on the results from the study, the following conclusions have been drawn: 1.
The new approach of using the LTSInterFoam solver, as highlighted in this study, proved to be a reliable means for the steady-state simulation of a free surface vortex at an intake. The LTSInterFoam solver, until now, has been associated with hydrodynamic studies involving ships. However, being a steady state solver, it is unable to account for the transient process of evolution of a free surface vortex.

2.
The SST k − ω model provided the best results compared to the other turbulence models. 3.
The approach highlighted in this study can be used to investigate the effectiveness of an anti-vortex device in suppressing the formation of an air-core vortex.
The authors, thus, recommend that future studies explore the use of the interFoam solver in OpenFOAM to simulate the occurrence of free surface vortices in order to account for the transient process associated with the formation of free surface vortices.