Three-Dimensional Flow of a Vortex Drop Shaft Spillway with an Elliptical Tangential Inlet

Vortex drop shaft (VDS) spillways are eco-friendly hydraulic structures used for safely releasing flood. However, due to the complexity of the three-dimensional rotational flow and the lack of suitable measuring devices, current experimental work cannot interpret the flow behavior reliably inside the VDS spillway, consequently experimental and CFD study on a VDS spillway with an elliptical tangential inlet was conducted to further discern the interior three-dimensional flow behavior. Hydraulic characteristics such as wall pressure, swirl angle, annular hydraulic height and Froude number of the tapering section are experimentally obtained and acceptably agreed with the numerical prediction. Results indicated that the relative dimensionless maximum height of the standing wave falls off nearly linearly with the increasing Froude number. Nonlinear regression was established to give an estimation of the minimum air-core rate. The normalized height of the hydraulic jump depends on the flow phenomena of pressure slope. Simulated results sufficiently reveal the three-dimensional velocity field (resultant velocity, axial velocity, tangential velocity and radial velocity) with obvious regional and cross-sectional variations inside the vortex drop shaft. It is found that cross-sectional tangential velocity varies, resembling the near-cavity forced vortex and near-wall free vortex behavior. Analytic calculations for the cross-sectional pressure were developed and correlated well with simulated results.


Introduction
Vortex drop shaft (VDS) spillways are widely used discharge structures in China owing to various excellent performances, such as high-energy dissipation, steady flow regime and good terrain adaptability [1]. Additionally, compared with traditional ski-jump energy dissipaters, VDS spillways can transfer the flood discharge and energy dissipation task into the tunnel, which effectively avoids outlet atomization formed by the interaction between air boundary water and reduces the impact on the ecological environment [2,3]. Examples of existing VDS spillway include Shapai, China [4], Xiaowan, China [5], Jiayan, China [6], etc. While for other countries, such structures are now commonly applied in urban drainage systems to convey sewage or surface runoff from different levels [7,8]. Regardless of different application purpose, the investigation focused on the general structure, and flow characteristics have many similarities.
Typical design components for the VDS include an inlet structure, a vertical structure and an outlet structure. Inlet structures are mainly classified into the three forms: (a) tangential inlet [9], (b) vortex slot inlet [10], and (c) spiral inlet [11]. The intake structure aims to directly transform horizontal or sloping flow into a stable helicoidal flow and to maintain a stable air core, leaving enough space for the air to escape [12,13]. Extensive experimental and theoretical studies have been conducted to determine the flow properties of the inlet structures. Early Slisskii and Kuznetsov [14] reported a simple hydraulic model able to reproduce the experimental flow phenomena while straightforward obtaining the detailed hydraulic parameter inside the VDS. Recently, Liu et al. [2] have reported the successful use of the RNG k-ε turbulence model on a newly developed VDS spillway with the spiral-flow-generating piers. However, there are also no detailed three-dimensional flow characteristics.
Hence, in the present study, physical model experiments with corresponding numerical simulation were conducted to further discern the three-dimensional flow behavior of a VDS spillway, comprising a tangential vortex intake with 1/4 ellipse guiding wall, and the relevant computational results are verified by experimental data. Meanwhile, experimental and simulated results such as the performance of the standing wave, air-core ratio, cross-sectional pressure and three-dimensional velocity were comprehensively analyzed. Furthermore, the combined vortex rule was introduced to predict the cross-sectional pressure distribution of the swirling flow in the VDS. The outcomes of the current research can provide comprehensive insights for the research and application of similar engineering.

Physical Model and Measurement
A 1:40 Froude model fabricated from transparent plexiglass was made in the College of Water Resources and Architectural Engineering, Northwest A&F University, to provide a useful reference for future prototype VDS spillways serving in China. Generally, the scale effect can be ignored for some time-averaged hydraulic characteristics [36]. The physical model consists of an elliptical tangential inlet (include a tapering section and a vortex chamber), a gradient (contracted) section, a vertical shaft, an energy dissipation well, a pressure slope and an outlet tunnel, as shown in Figure 1, in which the vortex chamber, gradient section, vertical shaft, dissipation well can be collectively referred to as the VDS. 021, 13, x FOR PEER REVIEW 3 of 24 lack of suitable measuring devices, detailed information such as the three-dimensional velocity field inside the VDS spillway is hitherto lacking. Fortunately, numerical simulation is able to reproduce the experimental flow phenomena while straightforward obtaining the detailed hydraulic parameter inside the VDS. Recently, Liu et al. [2] have reported the successful use of the RNG k-ε turbulence model on a newly developed VDS spillway with the spiral-flow-generating piers. However, there are also no detailed three-dimensional flow characteristics. Hence, in the present study, physical model experiments with corresponding numerical simulation were conducted to further discern the three-dimensional flow behavior of a VDS spillway, comprising a tangential vortex intake with 1/4 ellipse guiding wall, and the relevant computational results are verified by experimental data. Meanwhile, experimental and simulated results such as the performance of the standing wave, air-core ratio, cross-sectional pressure and three-dimensional velocity were comprehensively analyzed. Furthermore, the combined vortex rule was introduced to predict the cross-sectional pressure distribution of the swirling flow in the VDS. The outcomes of the current research can provide comprehensive insights for the research and application of similar engineering.

Physical Model and Measurement
A 1:40 Froude model fabricated from transparent plexiglass was made in the College of Water Resources and Architectural Engineering, Northwest A&F University, to provide a useful reference for future prototype VDS spillways serving in China. Generally, the scale effect can be ignored for some time-averaged hydraulic characteristics [36]. The physical model consists of an elliptical tangential inlet (include a tapering section and a vortex chamber), a gradient (contracted) section, a vertical shaft, an energy dissipation well, a pressure slope and an outlet tunnel, as shown in Figure 1, in which the vortex chamber, gradient section, vertical shaft, dissipation well can be collectively referred to as the VDS. Unlike the previous tangential slot inlet where the tapering section and the vortex chamber are directly tangential connected in a straight line is widely applied in urban drainage systems. In this study, an elliptical tangential inlet equipped with a 1/4 elliptical guiding wall (  Unlike the previous tangential slot inlet where the tapering section and the vortex chamber are directly tangential connected in a straight line is widely applied in urban drainage systems. In this study, an elliptical tangential inlet equipped with a 1/4 elliptical guiding wall (x 2 /10.5 2 + y 2 /18.5 2 = 1) was applied, as shown in Figure 1. Additionally, a deflector was set in the vortex chamber to prevent the backwater phenomenon. The inclination angle (β) of the tapering section is 5.71 degrees, and the entrance width (E) of the vortex chamber is 0.162 m. The diameter of the vortex chamber (D 1 ) and vertical shaft (D 2 ) is 0.275 m and 0.2 m separately. The total height of the VDS is 1.91 m. A more detailed dimension is shown in Figure 1. In order to better understand the hydraulic behavior, a total of five experimental flow rates were considered, namely 45.46 L/s, 40.52 L/s, 24.11 L/s, 20.26 L/s and 17.59 L/s, corresponding to flood frequency P = 0.05%, P = 0.1%, P = 1%, Water 2021, 13, 504 4 of 23 P = 2% and P = 3.3%, respectively, and the first 3 test-runs (P = 0.05%, 0.1%, 1%, respectively) were simulated.
In this experiment, a steady discharge generated by a recirculating flow system was provided to guarantee the rigor and repeatability of the experiments, as simply shown in Figure 1. The discharge was measured by a rectangular weir. Piezometers were used for acquiring time-average pressure values of the wall. Due to the pulsation of water flow at dissipation well, a possibly largest pressure error of ±10% was estimated, and the error range of the remaining points was estimated to be within plus or minus 3%. The location of measuring points is as shown in Figure 1. The velocity at the tapering section was detected using a propeller-type current meter (relative error ≤ ±5%). The measuring points are, respectively, taken at a distance of about 1.5 cm from the water surface and the bottom and taken in the middle. The average velocity can be obtained through the above three sets of data. The near-wall swirl angle was measured by the tracer method with a protractor, and the accuracy is 1 degree. Given the available measuring devices and laboratory limitations, three sets of measurements were taken, and the experimental data were averaged.
Inside the VDS, the hydraulic behavior of the numerical model is similar. However, in that region, the laboratory measurements such as the pitot tube may be affected by the turbulently uneven water thickness and misaligned direction, which might result in an error up to 20-30% [31]. Hence, the flow velocity in this research was supplemented by numerical simulation.

Governing Equations
Although the k-ε turbulent models have been widely applied to flow simulations because of the fine convergence and computational efficiency [37][38][39], the RNG k-ε turbulence model proposed by Yakhot and Orzag, by providing an additional coefficient correcting the turbulent viscosity and taking into account the average rotational motion, can effectively process the strong swirling flow or curved wall jet under high Reynolds number [40]. Furthermore, regarding to previous research [2,24,41] were found that the RNG k-ε turbulence model can give the best prediction of the jet flow in such shaft spillway. Accordingly, this turbulence model was utilized to model the intricate flow in the VDS spillways via ANSYS 16.0 software, Fluent. The basic control equations are as follows: Continuity Equation: Momentum Equation: k Equation: ε Equation: ∂(ρε) ∂t where u i is the velocity component in x i direction; t is the time; ρ and v are the density and the kinematic viscosity coefficient of water, respectively. P is the corrected pressure; f i is the external volume force component; ∂k and ∂ε are the turbulent Prandtl constant; G k is the turbulent kinetic energy generated by the velocity gradient; µ e f f is the corrected effective Water 2021, 13, 504 [42].
The VOF method presented by Hirt and Nichols [43] was used to catch the free surface of the high-speed flow, which is achieved by resolving the following Equation: where a w is the volume fraction of water.

Mesh and Boundary Conditions
The simulated domain with the same size as the physical model is shown in Figure 2. The origin of the coordinate system is set at the center of the bottom dissipation well. In the present simulation, a combination of structured and unstructured grids was utilized to mesh different parts of the VDS spillways. The well adaptable unstructured grids were applied for the VDS, and the grids were properly refined to ensure accuracy. The structural mesh was used for other parts. The partial meshing of the calculated domain is shown in Figure 2.   Mesh convergence is tested under four grid densities at the flood frequency of 0.05%; the results are shown in Table 1. It can be seen that the predicted outlet average velocity between mesh Scheme 3 and Scheme 4 shows less than a difference of 7.5% to the laboratory data and a relative error of 2.3%, suggesting that both grid densities are practically available. Scheme 3 was eventually selected to save calculation time. Table 2 shows the comparison of experimental and simulated outlet average velocity based on the grid Scheme 3 under different flood conditions. The small deviations that are also visible once again prove that grid Scheme 3 is reliable. The boundary conditions are set as shown in Figure 2. For the water inlet, the flow depth and velocity were selected while the pressure-outlet prescribed with 0 gauge pressure was defined. Additionally, the pressure-inlet was adopted for the two air vents. The wall boundary was based on no-slip conditions, and a standard wall function was used for processing the area near the wall.
The governing equations are solved using the finite volume method; the pressurevelocity coupling adopts the well-established PISO algorithm; a quick scheme is used for momentum, a geo-reconstruct scheme for the volume fraction, and a second-order upwind scheme for the turbulent kinetic energy. The density and the kinematic viscosity coefficient of water are 998.2 kg/m 3 and 1.003 10 −3 N·s/m 2 , respectively. A transient simulation using a variable time step was performed until the residual error was within 10 −4 . A group of successful simulation required about 50-60 s (approximately runs for 12 days), and three sets of stable computational results were averaged for different parameters.  Figure 4, which reduces with increasing discharge. It can be noted that the experimental results are in acceptable agreement with the numerical prediction (percentage error less than 9.12%). For all flow cases, the approach flow is supercritical within the entire tapering section, where Fr is larger than 1. tangential inlet are shown in Figure 3. The figure reveals that the flow at long tapering section can be viewed as the rectangular open channel flow. The Froude number (Fr) calculated by F / r v gh = is also plotted in Figure 4, which reduces with increasing discharge. It can be noted that the experimental results are in acceptable agreement with the numerical prediction (percentage error less than 9.12%). For all flow cases, the approach flow is supercritical within the entire tapering section, where Fr is larger than 1.    After the flow entered the vortex chamber with a 1/4 elliptical guiding wall, the water surface elevation had a rising process due to flow deflection, generating one visible standing shock wave along the guiding wall, as observed in Figure 3. According to the experimental data in previous research [11], the shock wave was characterized by the maximum shock wave elevation (hm) for the spiral inlet, whereas few experimental data are available for the elliptical tangential intake. In the present investigation, the After the flow entered the vortex chamber with a 1/4 elliptical guiding wall, the water surface elevation had a rising process due to flow deflection, generating one visible standing shock wave along the guiding wall, as observed in Figure 3. According to the experimental data in previous research [11], the shock wave was characterized by the maximum shock wave elevation (h m ) for the spiral inlet, whereas few experimental data are available for the elliptical tangential intake. In the present investigation, the performance of the maximum height of the standing wave is expressed as a dimensionless parameter η m (relative maximum height η m = (h m − h 0 )/D 1 , where h m , h 0 and D 1 is the maximum standing wave elevation, the water depth at the beginning of the elliptical guiding wall and the diameter of the vortex chamber, respectively, as shown in Figure 5. After the flow entered the vortex chamber with a 1/4 elliptical guiding wall, the water surface elevation had a rising process due to flow deflection, generating one visible standing shock wave along the guiding wall, as observed in Figure 3. According to the experimental data in previous research [11], the shock wave was characterized by the maximum shock wave elevation (hm) for the spiral inlet, whereas few experimental data are available for the elliptical tangential intake. In the present investigation, the performance of the maximum height of the standing wave is expressed as a dimensionless parameter ηm (relative maximum height where hm, h0 and D1 is the maximum standing wave elevation, the water depth at the beginning of the elliptical guiding wall and the diameter of the vortex chamber, respectively, as shown in Figure 5. The dimensionless parameter ηm at the elliptical tangential inlet is plotted as a function of Fr0 at the beginning of the elliptical guiding wall in Figure 6. It is shown from this figure that with increasing the Fr0 (approximately Fr0 = 1-4) or equivalently the decrease of discharge, the ηm tends to linearly reduce. The linear regression is presented in Equation (6) with R 2 = 0.987. For spiral inlet with the supercritical flow, Crispino [17] stated that the shock wave maxima also depends on the Fr0 with a linear trend when not considering other geometric parameters, expressed as  The dimensionless parameter η m at the elliptical tangential inlet is plotted as a function of Fr 0 at the beginning of the elliptical guiding wall in Figure 6. It is shown from this figure that with increasing the Fr 0 (approximately Fr 0 = 1-4) or equivalently the decrease of discharge, the η m tends to linearly reduce. The linear regression is presented in Equation (6) with R 2 = 0.987. For spiral inlet with the supercritical flow, Crispino [17] stated that the shock wave maxima also depends on the Equation (6) is based on the present geometric model and data, and it can be used to provide a reference for the tangential inlet with an elliptical guiding wall. Meanwhile, more experimental work by varying discharge, bottom slope and elliptical form is necessary to further discern this behavior.

Air Core
The vortex chamber forces the approach flow to the spiral down and adheres to the wall, forming a stable cavity along the vertical axis. Figure 7 shows the experimental and predicted air core. The numerical prediction is able to overcome the test limitation inside the VDS and clearly capture the vertical variation, as shown in Figure 7b. The air core was hardly to be symmetrical because of the non-axisymmetric tangential inlet, resulting in uneven distribution of water layers in different cross-sections. Accordingly, the rotary flow can be divided into the main flow with thicker water layers and the non-main flow Equation (6) is based on the present geometric model and data, and it can be used to provide a reference for the tangential inlet with an elliptical guiding wall. Meanwhile, more experimental work by varying discharge, bottom slope and elliptical form is necessary to further discern this behavior.

Air Core
The vortex chamber forces the approach flow to the spiral down and adheres to the wall, forming a stable cavity along the vertical axis. Figure 7 shows the experimental and predicted air core. The numerical prediction is able to overcome the test limitation inside the VDS and clearly capture the vertical variation, as shown in Figure 7b. The air core was hardly to be symmetrical because of the non-axisymmetric tangential inlet, resulting in uneven distribution of water layers in different cross-sections. Accordingly, the rotary flow can be divided into the main flow with thicker water layers and the non-main flow with the relatively thin water layers. According to the location change of the main flow in Figure 7b, it can be estimated that water rotated close to about 1 circle in the VDS. Overall, the simulation was able to reproduce the experimental flow phenomena while straightforward showing the flow regime of the rotary flow inside the VDS.

Air Core
The vortex chamber forces the approach flow to the spiral down and adheres to the wall, forming a stable cavity along the vertical axis. Figure 7 shows the experimental and predicted air core. The numerical prediction is able to overcome the test limitation inside the VDS and clearly capture the vertical variation, as shown in Figure 7b. The air core was hardly to be symmetrical because of the non-axisymmetric tangential inlet, resulting in uneven distribution of water layers in different cross-sections. Accordingly, the rotary flow can be divided into the main flow with thicker water layers and the non-main flow with the relatively thin water layers. According to the location change of the main flow in Figure 7b, it can be estimated that water rotated close to about 1 circle in the VDS. Overall, the simulation was able to reproduce the experimental flow phenomena while straightforward showing the flow regime of the rotary flow inside the VDS.   Figure 8 shows the predicted air-core ratio λ (ratio of the air-core area to cross-sectional area) in the contracted section and vertical shaft, showing an increasing trend of λ with decreasing discharge. It can be seen that λ was the smallest around the junction with a value larger than the minimum engineering requirement of 0.25, indicating that the structure can meet the maximum flood discharge requirements and ensure sufficient aeration conditions [7,26]. With the rapid rotation and diffusion of flow in the vertical shaft, the water thickness decreases and equivalently λ increases. However, the increasing effect seems to flatten off remarkably as decreasing level z.
The minimum λ (λm), a vital design parameter, can be plotted against a dimensionless discharge parameter F a (F a = 4(Q 2 e/(gπ 3 D 1 6 cos β)) 1/3 /(1 − e/D 1 )) on previous experimental investigations [15] for the tangential slot inlet. Figure 9 shows the numerical predictions of λm against Fa (according to [1], a value of e should be given as 0.9E for the elliptical tangential intake); the data of Yu and Lee [15] and Jain [7] are also extracted. It is seen that the results nearly match on one line. To estimate the minimum air-core ratio for the tangential intake, nonlinear regression was founded to fit collective data and proposed as follows: cross-sectional area) in the contracted section and vertical shaft, showing an increasing trend of λ with decreasing discharge. It can be seen that λ was the smallest around the junction with a value larger than the minimum engineering requirement of 0.25, indicating that the structure can meet the maximum flood discharge requirements and ensure sufficient aeration conditions [7,26]. With the rapid rotation and diffusion of flow in the vertical shaft, the water thickness decreases and equivalently λ increases. However, the increasing effect seems to flatten off remarkably as decreasing level z. The minimum λ (λm), a vital design parameter, can be plotted against a dimensionless discharge parameter a F ( ) on previous experimental investigations [15] for the tangential slot inlet. Figure 9 shows the numerical predictions of λm against Fa (according to [1], a value of e should be given as 0.9E for the elliptical tangential intake); the data of Yu and Lee [15] and Jain [7] are also extracted. It is seen that the results nearly match on one line. To estimate the minimum air-core ratio for the tangential intake, nonlinear regression was founded to fit collective data and proposed as follows: Equation (7) can be expected to give a satisfactory estimation of the minimum air-core rate for the VDS spillway with the tangential inlet.

Annular Hydraulic Jump
Finally, the helicoidal jet dropped into the milky hydraulic jump and formed a water cushion, consuming much energy, as shown in Figure 10. In the present experimental discharge, the pressure slope was not always submerged by the air-water mixture flow, showing a changing performance in relation to the discharge. The experimental observation shows that the condition of flow in pressure slope altering from weir flow (P = 3.33%, Q = 17.59 L/s) to orifice flow phenomenon (P = 1%, Q = 24.11 L/s). With a further increase in discharge, the free surface profile of the hydraulic jump is on the rise. Equation (7) can be expected to give a satisfactory estimation of the minimum air-core rate for the VDS spillway with the tangential inlet.

Annular Hydraulic Jump
Finally, the helicoidal jet dropped into the milky hydraulic jump and formed a water cushion, consuming much energy, as shown in Figure 10. In the present experimental discharge, the pressure slope was not always submerged by the air-water mixture flow, showing a changing performance in relation to the discharge. The experimental observation shows that the condition of flow in pressure slope altering from weir flow (P = 3.33%, Q = 17.59 L/s) to orifice flow phenomenon (P = 1%, Q = 24.11 L/s). With a further increase in discharge, the free surface profile of the hydraulic jump is on the rise. ter cushion, consuming much energy, as shown in Figure 10. In the present experimental discharge, the pressure slope was not always submerged by the air-water mixture flow, showing a changing performance in relation to the discharge. The experimental observation shows that the condition of flow in pressure slope altering from weir flow (P = 3.33%, Q = 17.59 L/s) to orifice flow phenomenon (P = 1%, Q = 24.11 L/s). With a further increase in discharge, the free surface profile of the hydraulic jump is on the rise.  To examine the flow behavior of the annular hydraulic jump in more detail, the normalized height of the hydraulic jump h j /D 2 (where h j is the elevation of the hydraulic jump, D 2 is the diameter of the vertical shaft.) is plotted against Q/ gD 2 5 in Figure 11. In the numerical test, the annular hydraulic jump occurred at h j /D 2 = 3.60 (P = 0.05%), 3.01 (P = 0.1%) and 2.4 (P = 1%), which agreed well with experimental data under a negligible percentage error (2.1-3.6%). For different discharges, h j /D 2 increases as the Q/ gD 2 5 increases. However, h j /D 2 does not shows a monotonic increasing tendency with respect to Q/ gD 2 5 . It seems that h j /D 2 initially increases slowly resembling with a state of weir flow when the discharge is lower than the critical value (near Q/ gD 2 5 = 0.42, Q = 24.11 L/s, P = 1%), whereas this increase becomes quickly with a trait of parabolic shape for discharge larger than these ranges. This means that the face of orifice flow is generated. The results are exactly as shown in Figure 10. , Q = 24.11 L/s, P = 1%), whereas this increase becomes quickly with a trait of parabolic shape for discharge larger than these ranges. This means that the face of orifice flow is generated. The results are exactly as shown in Figure 10.

Wall Pressure Distribution
The pressure corresponding to the shaft wall is one of the vital parameters to evaluate the operation safety of the VDS spillways. In order to prevent cavitation occurrence and cavitation erosion, the negative pressure area should be avoided or reduced as much as possible [24]. Figure 12 shows the measured data with corresponding simulated data of wall pressure in the VDS. It can be clearly seen that the experimental and predicted pressure show similar behavior with an acceptable deviation (mean error is equal to 8.07%), which further verifies the accuracy of the numerical prediction. However, there is a significant deviation below the elevation of the annular hydraulic jump (a percentage error of 5.89-33.69%). A reasonable explanation for this phenomenon was the  The pressure corresponding to the shaft wall is one of the vital parameters to evaluate the operation safety of the VDS spillways. In order to prevent cavitation occurrence and cavitation erosion, the negative pressure area should be avoided or reduced as much as possible [24]. Figure 12 shows the measured data with corresponding simulated data of wall pressure in the VDS. It can be clearly seen that the experimental and predicted pressure show similar behavior with an acceptable deviation (mean error is equal to 8.07%), which further verifies the accuracy of the numerical prediction. However, there is a significant deviation below the elevation of the annular hydraulic jump (a percentage error of 5.89-33.69%). A reasonable explanation for this phenomenon was the measurement error caused by the rapidly up and down fluctuation of the hydraulic jump. Another explanation lied in the limitations of the VOF method for simulating vigorous air-water mixed-flow [2]. There was no negative pressure on the wall; the pressure distribution was asymmetrical because of uneven water thickness. The pressure initially increased from the vortex chamber to the transition section; this is due to the reason that the tapered gradient section shrinks, resulting in the increase of water thickness and tangential velocity. After this, the wall pressure above the annular hydraulic jump reduced, and the overall value was small. Finally, owing to the thickness and impingement of the water cushion, the pressure below the hydraulic jump sharply increased until the pressure reached its maximum value in the bottom of the well. The maximum value in the test was 10.76 kPa on the left wall (P = 0.05%) and 9.82 kPa on the right wall (P = 0.05%), corresponding to the simulated values of 12.34 kPa and 9.23 kPa. After being converted to prototype, it was a relatively large value and should be paid attention to. This tendency recorded was similar to the result obtained by He et al. [27] who simulated a large-scale vortex shaft spillway with superhigh water head and large flood discharge (Q > 1400 m 3 /s).

Cross-Sectional Pressure Distribution
Pressure results for the different parts were obtained based on the three representative cross-sections, which are located in the vortex chamber (z = 1.55 m), the gradient section (z = 1.40 m), and the vertical shaft (z = 1.10 m), respectively. The general radial variation of the pressure cannot be obtained by experiment work but successfully captured by the numerical simulation. Figure 13 shows the cloud diagram of cross-sectional pressure, and Figure 14 drawn from the extracted node data, is the pressure distribution on the right side of the shaft axis. To avoid redundancy, the left side result that has the same variation as the right side is not considered. It can be seen from Figure 14 that the variation of pressure distribution is represented as being the smallest at the air core, approximately 0, then gradually increasing with the increasing radius for the rotary water, and finally reaching the maximum value on the wall. Additionally, the pressure of the same cross-section with the same radius is not symmetrical because of the uneven distribution of water thickness, as seen in Figure 13. Interestingly, the pressure isoline density near the cavity is slightly larger than that near the shaft wall. This rather complicated behavior can be explained in terms of the adequately turbulent mixed flow at the airwater interface. There was no negative pressure on the wall; the pressure distribution was asymmetrical because of uneven water thickness. The pressure initially increased from the vortex chamber to the transition section; this is due to the reason that the tapered gradient section shrinks, resulting in the increase of water thickness and tangential velocity. After this, the wall pressure above the annular hydraulic jump reduced, and the overall value was small. Finally, owing to the thickness and impingement of the water cushion, the pressure below the hydraulic jump sharply increased until the pressure reached its maximum value in the bottom of the well. The maximum value in the test was 10.76 kPa on the left wall (P = 0.05%) and 9.82 kPa on the right wall (P = 0.05%), corresponding to the simulated values of 12.34 kPa and 9.23 kPa. After being converted to prototype, it was a relatively large value and should be paid attention to. This tendency recorded was similar to the result obtained by He et al. [27] who simulated a large-scale vortex shaft spillway with superhigh water head and large flood discharge (Q > 1400 m 3 /s).

Cross-Sectional Pressure Distribution
Pressure results for the different parts were obtained based on the three representative cross-sections, which are located in the vortex chamber (z = 1.55 m), the gradient section (z = 1.40 m), and the vertical shaft (z = 1.10 m), respectively. The general radial variation of the pressure cannot be obtained by experiment work but successfully captured by the numerical simulation. Figure 13 shows the cloud diagram of cross-sectional pressure, and Figure 14 drawn from the extracted node data, is the pressure distribution on the right side of the shaft axis. To avoid redundancy, the left side result that has the same variation as the right side is not considered. It can be seen from Figure 14 that the variation of pressure distribution is represented as being the smallest at the air core, approximately 0, then gradually increasing with the increasing radius for the rotary water, and finally reaching the maximum value on the wall. Additionally, the pressure of the same cross-section with the same radius is not symmetrical because of the uneven distribution of water thickness, as seen in Figure 13. Interestingly, the pressure isoline density near the cavity is slightly larger than that near the shaft wall. This rather complicated behavior can be explained in terms of the adequately turbulent mixed flow at the air-water interface.

Three-Dimensional Velocity Field in the Vortex Drop Shaft
The flow velocity is an important parameter when analyzing the energy dissipation effect and calculating the cavitation number. The laboratory measurement such as pitot tube has great limitations on measuring the resultant velocity in the VDS. However, Nan et al. [31] has validated the reliability for simulated velocity field When the numerical results are able to well predicts the flow pattern, the pressure, and the swirl angle (mentioned below). For this reason, the following analysis is based on numerical results.

Velocity Field Near the Wall
Generally, the rotational motion within the VDS not only depends on geometric parameters such as the shape of the inlet structure, the shrinkage degree of the gradient section and the size of the vertical shaft, but also is closely linked to the turbulence and shear of flow, the mixing of the air and water, and the centrifugal force. Previous studies [2,24,27] focused more on the resultant velocity rather than each velocity component. As a typical three-dimensional motion, the rotary flow motion can be decomposed into axial motion, tangential motion and radial motion, satisfying the expression as follows: where v denotes the resultant velocity, z v , t v and r v are the axial velocity, the tangential velocity and the radial velocity, respectively.

Three-Dimensional Velocity Field in the Vortex Drop Shaft
The flow velocity is an important parameter when analyzing the energy dissipation effect and calculating the cavitation number. The laboratory measurement such as pitot tube has great limitations on measuring the resultant velocity in the VDS. However, Nan et al. [31] has validated the reliability for simulated velocity field When the numerical results are able to well predicts the flow pattern, the pressure, and the swirl angle (mentioned below). For this reason, the following analysis is based on numerical results.

Velocity Field Near the Wall
Generally, the rotational motion within the VDS not only depends on geometric parameters such as the shape of the inlet structure, the shrinkage degree of the gradient section and the size of the vertical shaft, but also is closely linked to the turbulence and shear of flow, the mixing of the air and water, and the centrifugal force. Previous studies [2,24,27] focused more on the resultant velocity rather than each velocity component. As a typical three-dimensional motion, the rotary flow motion can be decomposed into axial motion, tangential motion and radial motion, satisfying the expression as follows: where v denotes the resultant velocity, z v , t v and r v are the axial velocity, the tangential velocity and the radial velocity, respectively.

Three-Dimensional Velocity Field in the Vortex Drop Shaft
The flow velocity is an important parameter when analyzing the energy dissipation effect and calculating the cavitation number. The laboratory measurement such as pitot tube has great limitations on measuring the resultant velocity in the VDS. However, Nan et al. [31] has validated the reliability for simulated velocity field When the numerical results are able to well predicts the flow pattern, the pressure, and the swirl angle (mentioned below). For this reason, the following analysis is based on numerical results.

Velocity Field near the Wall
Generally, the rotational motion within the VDS not only depends on geometric parameters such as the shape of the inlet structure, the shrinkage degree of the gradient section and the size of the vertical shaft, but also is closely linked to the turbulence and shear of flow, the mixing of the air and water, and the centrifugal force. Previous studies [2,24,27] focused more on the resultant velocity rather than each velocity component. As a typical three-dimensional motion, the rotary flow motion can be decomposed into axial motion, tangential motion and radial motion, satisfying the expression as follows: where v denotes the resultant velocity, v z , v t and v r are the axial velocity, the tangential velocity and the radial velocity, respectively. Figure 15 shows the resultant velocity with corresponding velocity components of the swirling flow cling to the wall (the velocity extracting points are about 5-8 mm from the wall, and elevation corresponds to the pressure measurement point). With the exception of local fluctuation caused by shrinkage effects in the gradient section, the changing trend on both sides is basically similar under the tested discharge. The unequal values on both sides at a fixed elevation are probably due to the flow thickness and the wall friction. It is well known that with decreasing shaft elevation, gravity potential energy is gradually converted into kinetic energy, resulting in the increase in v and v z . However, once the helicoidal jet is distant from the top of the vertical shaft, the increasing rate of v and v z is no longer obvious. Following Hager [19], the rotational flow starts maintaining quasi-uniform flow with a constant maximum velocity. This behavior seems to occur for the small discharge (P = 1%), in which a stable velocity value is around 3.8 m/s, as shown in Figure 15. Dong [1] provided an approximation of maximum velocity v m , as v m = 8gQ nπD 2 1/3 (9) Water 2021, 13, x FOR PEER REVIEW 15 of 24 Figure 15 shows the resultant velocity with corresponding velocity components of the swirling flow cling to the wall (the velocity extracting points are about 5-8 mm from the wall, and elevation corresponds to the pressure measurement point). With the exception of local fluctuation caused by shrinkage effects in the gradient section, the changing trend on both sides is basically similar under the tested discharge. The unequal values on both sides at a fixed elevation are probably due to the flow thickness and the wall friction. It is well known that with decreasing shaft elevation, gravity potential energy is gradually converted into kinetic energy, resulting in the increase in v and z v . However, once the helicoidal jet is distant from the top of the vertical shaft, the increasing rate of v and z v is no longer obvious. Following Hager [19], the rotational flow starts maintaining quasi-uniform flow with a constant maximum velocity. This behavior seems to occur for the small discharge (P = 1%), in which a stable velocity value is around 3.8 m/s, as shown in Figure 15. Dong [1] provided an approximation of maximum velocity m v , as Herein n is the roughness coefficient; it is generally considered to be greater than the value calculated by the Manning formula [26]. This coefficient was also experimentally summarized by Dong [1] as: With respect to the tangential velocity, its value roughly increases at the gradient section but decreases at the vertical shaft. The reason for this increased behavior was attributed to the blocking and shrinkage effect of the gradient section, leading to the centrifugal force increases to maintain the conservation of angular momentum. After this, the tangential velocity behaves to decrease because of wall friction. In addition, it can be noted from Figure 15 that the radial velocity is negligible. Alternatively speaking, the rotational flow field in the VDS is dominated by axial and tangential movement. Herein n is the roughness coefficient; it is generally considered to be greater than the value calculated by the Manning formula [26]. This coefficient was also experimentally summarized by Dong [1] as: With respect to the tangential velocity, its value roughly increases at the gradient section but decreases at the vertical shaft. The reason for this increased behavior was attributed to the blocking and shrinkage effect of the gradient section, leading to the centrifugal force increases to maintain the conservation of angular momentum. After this, the tangential velocity behaves to decrease because of wall friction. In addition, it can be noted from Figure 15 that the radial velocity is negligible. Alternatively speaking, the rotational flow field in the VDS is dominated by axial and tangential movement.

Swirl Angle
So far, few analytical approaches for rotational flow angle (α) are available up to the authors' knowledge [2,17]. The swirl angle can be estimated as α = arctan(v t /v z ) for the reason that v r the aforementioned is negligible. According to the results of swirl angle, as shown in Figure 16, the simulated values are very close to those actually obtained in the laboratory, with an acceptable error of 3.3-11.5%. The swirl angle decreases with increasing the incoming discharge at a fixed elevation and varied approximately in the range of 15-80 • . Specifically, the patterns of swirl angle distribution are as follows: (1) The large swirl angle appeared at the vortex chamber, where the tangential velocity was much larger than the axial flow velocity, according to Figure 15. As the initial tangential dominant movement was transformed into a compound movement of the tangential and axial direction, the swirl angle gradually decreased along the wall; (2) The swirl angle in the gradient section tended to fluctuate due to the complex and unstable flow field herein. This phenomenon corresponded to the noticeable fluctuation of v t and v z in Figure 15, which was beneficial to increase the local water head loss to improve the energy dissipation rate; (3) In the vertical shaft, the tangential velocity along the shaft was gradually consumed due to friction effects, whereas the axial velocity increased because of gravity, which together reduced the swirl angle. So far, few analytical approaches for rotational flow angle (α) are available up to the authors' knowledge [2,17]. The swirl angle can be estimated as arctan( / ) for the reason that r v the aforementioned is negligible. According to the results of swirl angle, as shown in Figure 16, the simulated values are very close to those actually obtained in the laboratory, with an acceptable error of 3.3-11.5%. The swirl angle decreases with increasing the incoming discharge at a fixed elevation and varied approximately in the range of 15-80º. Specifically, the patterns of swirl angle distribution are as follows: (1) The large swirl angle appeared at the vortex chamber, where the tangential velocity was much larger than the axial flow velocity, according to Figure 15. As the initial tangential dominant movement was transformed into a compound movement of the tangential and axial direction, the swirl angle gradually decreased along the wall; (2) The swirl angle in the gradient section tended to fluctuate due to the complex and unstable flow field herein. This phenomenon corresponded to the noticeable fluctuation of t v and z v in Figure 15, which was beneficial to increase the local water head loss to improve the energy dissipation rate; (3) In the vertical shaft, the tangential velocity along the shaft was gradually consumed due to friction effects, whereas the axial velocity increased because of gravity, which together reduced the swirl angle.

Cross-Sectional Velocity Field
(1) Resultant velocity distribution The example representing the resultant velocity distribution along the radial direction is shown in Figure 17 at the flood frequency of 0.5%. Figure 18 shows the velocity distribu-tion of the rotational flow (fraction of water accounting for 80% is taken as the air-water interface). Due to the similar variation of left side, only the right results are shown. It can be seen that the velocity behavior in the swirling flow area and the air core area exhibit obvious spatial variations.
from the air-water interface due to viscous dissipation, subsequently plummets to 0 rapidly near the wall affected by the boundary condition, as shown in Figure 18. On the other hand, the distribution of resultant velocity at the vertical shaft where the water layer is relatively thin is distinctly different from that aforementioned. The resultant velocity value no longer declines at first. Instead, its value first increases and then decreases with an insignificant variation, which persuasively proves why it is reasonable to use the near-wall velocity to replace average velocity at the vertical shaft regarding the previous studies [2,24].
With respect to the air, its movement is mainly caused by the entrainment of swirling water [33]. The density of velocity isoline is the largest near the air-water interface, probably because of the highly turbulent and violently mixed air and water here. Additionally, the air velocity in the cavity approximates the maximum at the center, as seen in   from the air-water interface due to viscous dissipation, subsequently plummets to 0 rapidly near the wall affected by the boundary condition, as shown in Figure 18. On the other hand, the distribution of resultant velocity at the vertical shaft where the water layer is relatively thin is distinctly different from that aforementioned. The resultant velocity value no longer declines at first. Instead, its value first increases and then decreases with an insignificant variation, which persuasively proves why it is reasonable to use the near-wall velocity to replace average velocity at the vertical shaft regarding the previous studies [2,24]. With respect to the air, its movement is mainly caused by the entrainment of swirling water [33]. The density of velocity isoline is the largest near the air-water interface, probably because of the highly turbulent and violently mixed air and water here. Additionally, the air velocity in the cavity approximates the maximum at the center, as seen in Figure 17.   For the rotational flow area, the resultant velocity gradient is much smaller than that of the air in Figure 17. On one hand, for the vortex chamber and gradient section with a thick water layer, the velocity distribution is represented as that it linearly decreases from the air-water interface due to viscous dissipation, subsequently plummets to 0 rapidly near the wall affected by the boundary condition, as shown in Figure 18. On the other hand, the distribution of resultant velocity at the vertical shaft where the water layer is relatively thin is distinctly different from that aforementioned. The resultant velocity value no longer declines at first. Instead, its value first increases and then decreases with an insignificant variation, which persuasively proves why it is reasonable to use the near-wall velocity to replace average velocity at the vertical shaft regarding the previous studies [2,24].
With respect to the air, its movement is mainly caused by the entrainment of swirling water [33]. The density of velocity isoline is the largest near the air-water interface, probably because of the highly turbulent and violently mixed air and water here. Additionally, the air velocity in the cavity approximates the maximum at the center, as seen in Figure 17.
(2) Tangential velocity distribution The tangential velocity distribution on a level rotary flow discharge tunnel is similar to that of the hydrocyclone, which is basically in line with the combined vortex distribution that inside and outside swirling flow are corresponding to the forced vortex distribution and free vortex distribution, respectively [30,31,33]. Whereas such phenomena for the VDS spillway is not clear. Therefore, the author tried to analyze and discuss the tangential velocity distribution of the rotational flow in the VDS. According to the test of Bradley and Pulling [44], the tangential velocity distribution of the forced vortex accords with the following law: v t r k = C The distribution of free vortex conforms to the following expression: where the k is an index reflecting the degree to which the swirling flow conforms to the absolute forced vortex or free vortex behavior, 0 < k < 1, the closer k is to 1, the more consistent it is. C is a constant on the same cross-section, but different-in-different crosssections. Figure 19 visualizes the distribution of the tangential velocity of the swirling flow. For the vortex chamber and the contracted section, the tangential velocity increases first, with an increasing radius near the cavity, which is suggestive of a forced vortex behavior, then tends to fall off resembling a free vortex distribution near the wall, ultimately plummets to 0 affected by boundary condition. Moreover, it can be remarkably noted that the range of increase in the tangential velocity is smaller than that of decrease of the tangential velocity. This means that the tangential velocity mainly presents the free vortex distribution here. For instance, in the case of P = 0.05%, the ratio of increasing range to decreasing range is 1.75:1 (z = 1.55 m) and 1.83:1 (z = 1.40 m). However, for the vertical shaft, this trend remains opposite; the forced vortex area is wide because of the relatively large range of increase in the tangential velocity. For instance, in the case of P = 0.05%, the ratio of increasing range to decreasing range is 0.54:1 (z = 1.10 m). This finding is contrary to the previous assumption established by Jain [20], who considered the tangential velocity distribution at the vertical shaft as a free vortex region. (2) Tangential velocity distribution The tangential velocity distribution on a level rotary flow discharge tunnel is similar to that of the hydrocyclone, which is basically in line with the combined vortex distribution that inside and outside swirling flow are corresponding to the forced vortex distribution and free vortex distribution, respectively [30,31,33]. Whereas such phenomena for the VDS spillway is not clear. Therefore, the author tried to analyze and discuss the tangential velocity distribution of the rotational flow in the VDS. According to the test of Bradley and Pulling [44], the tangential velocity distribution of the forced vortex accords with the following law: The distribution of free vortex conforms to the following expression: where the k is an index reflecting the degree to which the swirling flow conforms to the absolute forced vortex or free vortex behavior, 0< k <1, the closer k is to 1, the more consistent it is. C is a constant on the same cross-section, but different-in-different cross-sections. Figure 19 visualizes the distribution of the tangential velocity of the swirling flow. For the vortex chamber and the contracted section, the tangential velocity increases first, with an increasing radius near the cavity, which is suggestive of a forced vortex behavior, then tends to fall off resembling a free vortex distribution near the wall, ultimately plummets to 0 affected by boundary condition. Moreover, it can be remarkably noted that the range of increase in the tangential velocity is smaller than that of decrease of the tangential velocity. This means that the tangential velocity mainly presents the free vortex distribution here. For instance, in the case of P = 0.05%, the ratio of increasing range to decreasing range is 1.75:1 (z = 1.55 m) and 1.83:1 (z = 1.40 m). However, for the vertical shaft, this trend remains opposite; the forced vortex area is wide because of the relatively large range of increase in the tangential velocity. For instance, in the case of P = 0.05%, the ratio of increasing range to decreasing range is 0.54:1 (z = 1.10 m). This finding is contrary to the previous assumption established by Jain [20], who considered the tangential velocity distribution at the vertical shaft as a free vortex region. It is significant to say that the behavior of the tangential velocity changing from quasi-forced vortex to quasi-free vortex is probably due to the reason that the violent It is significant to say that the behavior of the tangential velocity changing from quasi-forced vortex to quasi-free vortex is probably due to the reason that the violent turbulence and full mixing of the air and water at the air-water interface, making the nearcavity rotational flow more affected by the air, consequently presenting a forced vortex distribution, especially obvious when the water layer is thinner, and the cavity is larger. Hence, the tangential velocity at the vertical shaft is mainly presented as the forced vortex distribution. Whereas the wall friction and the water layer shear will cause the tangential velocity to be consumed, thus the near-wall tangential velocity decays freely and presents a free vortex distribution.

Theoretical Analysis between Cross-Sectional Pressure and Tangential Velocity
Herein, a predicted approach for the relationship between tangential velocity and pressure is proposed. The following calculations and results refer to the present geometry of the VDS spillway, with an elliptical tangential intake at the drop shaft head. Now using a swirling micro-surface, as shown in Figure 20. The motion Equation along the radial direction is established as follows: where r is the radius of micro rotary flow, dr is the thickness of micro rotary flow, ∂p ∂r is the pressure change, dv r /d t is the radial acceleration, and ρ is the density of water. Since the radial velocity is extremely small and can be ignored, Equation (13) can be simplified to: Water The k value of different positions was calculated by the simulated tangential velocity. Then combined with Equations (17) and (18), the theoretical pressure at each position can be obtained. The comparison of pressure between the theoretical value and the simulated value is shown in Figure 21, where the dotted line indicates the boundary between the forced vortex and the free vortex region. It can be seen that there is a significant difference between the calculated value and the simulated value near the boundary and that the deviation within the forced vortex region is greater than that within the free vortex region. The possible reason lies in that there is a rapid adjustment of kinetic energy and pressure energy in the forced vortex region [32]. For the near-cavity swirling flow resembling forced vortex behavior, C is a constant of the same section, the Formula (11) can be expressed as: where v tc is the tangential velocity of the swirling flow at the air-water interface, and r c is the corresponding radius. Substituting Formula (15) into Formula (14) and performing integration yields: where C 1 is the integral constant. When r = r c , pressure is equal to 0, Calculate C 1 = −ρv tc 2 /2k, and substitute it into Equation (16) to obtain the pressure expression of the forced vortex distribution area as follows: For the near-wall rotational flow resembling free vortex behavior, only Equation (12) needs to be expressed as v t = v tR (R/r) k . Meanwhile, the boundary conditions that r = R, p = p R are adopted, then the pressure expression can be obtained as: where p R is the wall pressure, v tR is the corresponding tangential velocity, and R is the radius of the VDS. The k value of different positions was calculated by the simulated tangential velocity. Then combined with Equations (17) and (18), the theoretical pressure at each position can be obtained. The comparison of pressure between the theoretical value and the simulated value is shown in Figure 21, where the dotted line indicates the boundary between the forced vortex and the free vortex region. It can be seen that there is a significant difference between the calculated value and the simulated value near the boundary and that the deviation within the forced vortex region is greater than that within the free vortex region. The possible reason lies in that there is a rapid adjustment of kinetic energy and pressure energy in the forced vortex region [32]. The k value of different positions was calculated by the simulated tangential velocity. Then combined with Equations (17) and (18), the theoretical pressure at each position can be obtained. The comparison of pressure between the theoretical value and the simulated value is shown in Figure 21, where the dotted line indicates the boundary between the forced vortex and the free vortex region. It can be seen that there is a significant difference between the calculated value and the simulated value near the boundary and that the deviation within the forced vortex region is greater than that within the free vortex region. The possible reason lies in that there is a rapid adjustment of kinetic energy and pressure energy in the forced vortex region [32].  The simulated and theoretical results exhibit a very similar pattern with an acceptable overall deviation at the vortex chamber (z = 1.55 m) and the gradient section (z = 1.40 m). The mean error between the calculated and simulated values was 6.96% and 7.92%, respectively, while for the vertical shaft (z = 1.10 m), The mean error is 15.11%, suggesting that the theoretical formula based on the combined vortex rule is feasible for the cross-sectional pressure prediction in the VDS spillway, especially in the free vortex region.

Conclusions
The complex flow field of a VDS spillway with an elliptical tangential inlet is tested experimentally and simulated. The results of the model experiment with corresponding numerical prediction were similar in terms of the flow pattern. The wall pressure and the swirl angle flow regime of the VDS spillway is characterized by the standing shock wave, the air core and the annular hydraulic jump. One standing wave was detected under supercritical flow conditions with a moderate Froude number. The dimensionless maximum wave height tends to linearly grow with increasing discharge, but the increasing effect is insignificantly improved. A high-axisymmetric air core along the shaft can be well captured by numerical prediction. On the basis of previous studies, a quadratic equation was fitted to give an estimation of the minimum air-core rate for the VDS spillway with the tangential inlet. As regards the annular hydraulic jump, its dimensionless height, showing variation performance of first increasing slowly and then rapidly increasing with increasing discharge, was related to the pressure slope where the flow condition varies from weir flow to orifice flow simultaneously.
The swirling flow velocity is decomposed into the three-dimensional components for detailed research. The end rotational flow seems to start maintaining quasi-uniform flow with the constant maximum resultant velocity for small discharge. The cross-sectional resultant velocity value linearly decreases with an increasing radius for the vortex chamber and gradient section while initially increases and then decreases with an insignificant variation for the vertical shaft. Moreover, the rotational flow field in the VDS is dominated by the axial and tangential movement because radial velocity with an extremely small value is negligible. The cross-sectional tangential velocity distribution accords with the combined vortex behavior. The free vortex region is wider in the vortex chamber and gradient section, but the vertical is opposite. The theoretical calculation for the cross-sectional pressure based on the combined vortex rule was developed, especially feasible in the free vortex region. Considering that the eco-friendly VDS spillway has become an important research hotspot in the field of flood discharge and energy dissipation, the current work can provide useful insights for future research and application of similar engineering. Furthermore, the relevant prototype measurements work is necessary to conduct in the future.
Author Contributions: Physical experiment, Z.Y., G.X. and H.Y.; computer and numerical simulation, Z.Y., J.Y. and Y.L.; methodology and Chart, Z.Y., J.Y. and G.X.; writing-review and editing, Z.Y., J.Y. and Z.L.; funding support, Z.L. All authors contributed to the overall framing and revision of the manuscript at multiple stages. All authors have read and agreed to the published version of the manuscript.
Funding: This research was funded by the National Natural Science Foundation of China, grant number (No. 51509212).

Institutional Review Board Statement: Not applicable.
Informed Consent Statement: Informed consent was obtained from all subjects involved in the study.
Data Availability Statement: All data, models, or code generated or used during the study are available from the corresponding author by request.

Conflicts of Interest:
The authors declare no conflict of interest.