3-D Numerical Investigation on Oxygen Transfer in a Horizontal Venturi Flow with Two Holes

In order to investigate the dissolved oxygen increase caused by air suction in a horizontal Venturi flow with two holes, a 3-D computational fluid dynamics model was used to explore the water and bubble mixture flow, coupled with a dissolved oxygen transfer model. A series of experiments were conducted to validate the mathematical model. A relative saturation coefficient correlation was examined factoring in dissolved oxygen concentration at the inlet, water velocity at the inlet, the hole’s diameter, contraction ratio at throat section, and the downstream length of Venturi pipe. It was found that the relative saturation coefficient increases with increasing dissolved oxygen concentration at the inlet and downstream length of Venturi pipe respectively. However, it increases with decreasing water velocity at the inlet and contraction ratio at the throat section to some extent. The hole’s diameter plays a complex role in the relative saturation coefficient. The dimensional analysis method and the least square method were used to deduce a simple formula for the relative saturation coefficient, and this was consistent with related data.


Introduction
Dissolved oxygen (DO) concentration is an important water quality index.A lower value usually causes a serious environmental issue for a body of water, including decreasing the capacity of creatures living within it to grow, reproduce, and survive [1].To improve DO concentration in lakes and rivers and prevent water quality deterioration, aeration is commonly used [2][3][4].As a simple device, the Venturi meter is widely applied in industrial applications to control the flow rate of fluids.The flow characteristics and pressure drop of a gas-coal mixture through the Venturi under high pressure was investigated experimentally [5].The upstream and downstream pressure effects of Venturi flow were studied in terms of the flow rate and critical pressure ratio, both experimentally and numerically.Geometrical parameters such as the throat diameter, throat length, and diffuser angle were also considered [6].The pressure drops between different sections have also been measured and analyzed for a two-phase flow of Venturi, as in the case of a typical nuclear accident [7].A series of Venturi cavitation experiments were conducted with various upstream, throat and downstream lengths, and a model for water vapor void fraction was proposed to investigate the cavitating process [8].He and Bai [9] found that the k-ε model was in agreement with experimental data on wet gas flow in the Venturi meter.Due to the water quality improvement caused by aeration, by using a Venturi flow with holes, the aeration coefficient relationship was examined and found to be well within the geometrical parameters experimentally and numerically [10][11][12].
The aforementioned hydraulic characteristic of a Venturi flow has been investigated extensively.However, a thorough understanding of the DO variation of a Venturi flow with air bubbles is still The aforementioned hydraulic characteristic of a Venturi flow has been investigated extensively.However, a thorough understanding of the DO variation of a Venturi flow with air bubbles is still lacking.The objective of this research was to evaluate DO variation in an air-water two-phase flow through a horizontal Venturi pipe with two holes (HVFTH), experimentally and numerically.In addition, the relative saturation coefficient relationship with the hydraulic parameters and geometrical parameters of the Venturi pipe was investigated.Furthermore, a simple formula for the relative saturation coefficient was developed with the dimensional analysis method and the least square method.

Experimental Setup and Instrumentation
In order to examine the DO variation for HVFTH, a series of experiments were conducted in the Hydraulic Lab at Ocean University of China.Figure 1 shows that the Venturi pipe with two holes consists of an upstream section, a contraction section, a throat section, an expansion section, and a downstream section.At the throat section, two suction holes are drilled through its upper wall and bottom wall, respectively.A PVC circular horizontal pipe with an inner diameter D of 0.021 m and a length of 3.1 m was used to connect two cylindrical water tanks, as shown in Figure 2. The diameter and height for both tanks are 0.58 m and 0.73 m, respectively.Water flowed from the upstream tank to the downstream tank via the A# water pump located 0.5 m downstream of the pipe inlet, and then it was pumped back to the upstream tank by the B# water pump in the downstream tank.The water flow rate in the pipe was controlled by a controlling valve located 0.4 m downstream of the A# water pump.The water velocity inside the pipe was measured using an electromagnetic flowmeter located 0.5 m downstream of the controlling valve.The Venturi pipe with two holes was fixed 0.6 m downstream of the flowmeter, and the center of the holes is defined as zero location in the longitudinal x-direction.The throat diameter of Venturi pipe, Dt, is 0.005 m, and the two circular holes' diameters, d, are both 0.0007 m.Six DO probes numbered 1# to 6#, as shown in Figure 2, were used to measure the DO concentrations and water temperatures, where 1# DO probe was fixed in the upstream tank.Note that it is arduous to insert the other five DO probes into the PVC pipe to investigate the DO concentration variations, due to the undersized diameter of pipe.In order to solve this problem, five thin connection tubes were inserted into the downstream section with a spacing of 0.2 m between each tube.Five glass cups were used to collect the water through the five connection tubes.The 2# to 6# DO probes were submerged in the corresponding cups, and were used to measure the DO concentrations.In particular, the 2# to 6# DO probes approach the connection tube outlets as far as possible to ensure measurement accuracy and avoid the effect of oxygen transfer on the water surface.The corresponding DO concentrations were regarded as the values at the sections of five connection tubes without causing major deviations.

Experimental Procedure and Phenomena
Before the start of each experiment, the two tanks were filled with fresh tap water to 75% of the tanks' height.It was well known that anhydrous sodium sulfite (Na2SO3) is capable of scavenging DO effectively and has long been used in the presence of a catalyst (CoCl2) to scavenge DO in water [13].According to the chemical Equation ( 1), the suitable weight of Na2SO3 was determined by considering the water volume in the system, the present DO concentration, and its desired value.The amount of CoCl2 was suggested to be about 0.1% of Na2SO3 by weight [13,14].Na2SO3 and CoCl2 were put into the two water tanks as uniformly as possible.The water pumps then started to transport the water between the two tanks while mixing the water, Na2SO3 and CoCl2 effectively.The controlling valve was used to control the water flow rate so as to ensure air suction through the two holes at the Venturi pipe, and the water flow rate was measured by the flowmeter.The DO concentrations were measured and recorded once per minute by six DO probes.In order to avoid a possible deviation due to the incomplete reaction between Na2SO3 and DO, DO concentration data were not used until it increased steadily.
The experiment shows that the outside gas moves into the Venturi pipe through the two holes due to the inside and outside pressure difference, and a large number of air bubbles are produced.At the throat section and expansion section, a strong bubbly collision, coalescence and breakup occurs because of the water turbulence and shearing force, and the bubbles are distributed widely.When the water and bubble mixture flows into the downstream section, the buoyancy effect strengthens, the bubbles begin to ascend to some extent, and some bubbles swarm to collide with the upper wall.At the end of the downstream section, the water and bubble mixture is injected into the downstream tank, an inverse flow appears with a counter-clockwise direction vortex, and the bubbles ascend to escape from the water surface in the downstream tank.In terms of the DO concentration variations measured by 2# to 6# DO probes, they increase successively, confirming the ability to enhance the DO concentrations by using a HVFTH.

Governing Equations
Concerning the mathematical model to investigate the turbulence of water and the air bubbles' two-phases, considerable attention has been devoted to discussing their coupling effect.Euler-Euler methods consider the water and the bubbles as a pseudo-continuum and an averaging procedure is applied to both of them [15,16].Euler-Lagrange methods represent only the water motion through the Navier-Stokes equations, while the bubbly motion is obtained by the direct solution of its trajectory using Newton's second law [15].Torti [15] used the Eulerian-Lagrangian one-way

Experimental Procedure and Phenomena
Before the start of each experiment, the two tanks were filled with fresh tap water to 75% of the tanks' height.It was well known that anhydrous sodium sulfite (Na 2 SO 3 ) is capable of scavenging DO effectively and has long been used in the presence of a catalyst (CoCl 2 ) to scavenge DO in water [13].According to the chemical Equation ( 1), the suitable weight of Na 2 SO 3 was determined by considering the water volume in the system, the present DO concentration, and its desired value.The amount of CoCl 2 was suggested to be about 0.1% of Na 2 SO 3 by weight [13,14].Na 2 SO 3 and CoCl 2 were put into the two water tanks as uniformly as possible.The water pumps then started to transport the water between the two tanks while mixing the water, Na 2 SO 3 and CoCl 2 effectively.The controlling valve was used to control the water flow rate so as to ensure air suction through the two holes at the Venturi pipe, and the water flow rate was measured by the flowmeter.The DO concentrations were measured and recorded once per minute by six DO probes.In order to avoid a possible deviation due to the incomplete reaction between Na 2 SO 3 and DO, DO concentration data were not used until it increased steadily.
The experiment shows that the outside gas moves into the Venturi pipe through the two holes due to the inside and outside pressure difference, and a large number of air bubbles are produced.At the throat section and expansion section, a strong bubbly collision, coalescence and breakup occurs because of the water turbulence and shearing force, and the bubbles are distributed widely.When the water and bubble mixture flows into the downstream section, the buoyancy effect strengthens, the bubbles begin to ascend to some extent, and some bubbles swarm to collide with the upper wall.At the end of the downstream section, the water and bubble mixture is injected into the downstream tank, an inverse flow appears with a counter-clockwise direction vortex, and the bubbles ascend to escape from the water surface in the downstream tank.In terms of the DO concentration variations measured by 2# to 6# DO probes, they increase successively, confirming the ability to enhance the DO concentrations by using a HVFTH.

Governing Equations
Concerning the mathematical model to investigate the turbulence of water and the air bubbles' two-phases, considerable attention has been devoted to discussing their coupling effect.Euler-Euler methods consider the water and the bubbles as a pseudo-continuum and an averaging procedure is applied to both of them [15,16].Euler-Lagrange methods represent only the water motion through the Navier-Stokes equations, while the bubbly motion is obtained by the direct solution of its trajectory using Newton's second law [15].Torti [15] used the Eulerian-Lagrangian one-way coupling method to evaluate the water and the bubbly convection separately, where the bubbles were simply Water 2018, 10, 174 4 of 13 considered to be carried by the water and the bubbly effect they had on the water was negligible.However, the algebraic slip mixture model does not assume that there is an interface between these two immiscible phases, and it allows the two phases to interpenetrate as well as move at different velocities [12].Consequently, it was used to model the water and air bubbles' motion in HVFTH, due mainly to its good simplification and computational accuracy.The continuity equation can be expressed as Equation (2).
where t is time, x i is the coordinate in i direction.For i = 1, 2, 3, x i is x, y and z, respectively.The mixture density ρ m and the mixture velocity u mi in i direction can be defined as Equations ( 3) and ( 4).
where α n , ρ n , and u ni are the volume fraction, the density and the i direction velocity of the n-phase, respectively.n = 1 represents the air phase, and n = 2 represents the water phase.
Its momentum equation was expressed as Equation (5).
where µ m is the dynamic viscosity of the water and air bubble mixture, are the dynamic viscosities of air bubbles and water, respectively.g is the gravitational acceleration.p m is the mixture pressure.u dr,n is the diffusion velocity of the n-phase and it can be expressed as Equation (6).The diffusion velocity of a phase is usually caused by density differences, which results in forces on the bubbles differing from those on water: where u qn is the slip velocity between two phases and it can be expressed as Equation ( 7); µ eff is the effective viscosity; and d b is the bubble diameter.The drag function f drag can be defined as Equation (8).
where Re b is the bubble Reynolds number.The air volume fraction equation was used to predict the air volume evolution as Equation (9).
The Renormalization Group (RNG) k-ε model derived by [17] was used to investigate the air-water turbulence due to a satisfactory capacity to model the complex boundary flow considering the buoyancy effect.The turbulent kinetic energy, k, equation was written as Equation (10).
where G k represents the generation of turbulent kinetic energy caused by the average velocity, ε , C µ = 0.0845.The value of C µ is very close to the empirically-determined value of 0.09 used in the standard k-ε model; Pr t is the turbulence Prandtl number, Pr t = 0.07179.β k and β ε are the inverse effective Prandtl numbers for k and ε, respectively.C 1ε = 1.42 and C 2ε = 1.68.Note that although these constants are simple, they have a wide range of applications with commendable accuracy [17,18].
The DO transfer equation was utilized to compute the DO concentration variation, and it was denoted as Equation ( 12) [19,20].
where C is the DO concentration; Γ C is the turbulent diffusion coefficient, Γ C = ρ m ν t /σ t , where ν t is turbulent viscosity, and σ t is turbulent Schmidt number, where σ t = 0.7.S c is the source term representing the oxygen mass transfer across the interface between the water and bubbles, and it can be modeled as Equation ( 13) [21].
where K L,B a B is the bubble transfer coefficient, K L,B is the liquid film mass transfer coefficient, and a B is the specific interfacial area of the bubbles defined as the contact area with water per unit volume of air-water mixture, and it can be expressed as a B = 3α 1 /d b .d b is the air bubble diameter with the assumption of a spherical shape, and it can be expressed as Equation ( 14) [22].
where σ is the surface tension, and U 1 is the superficial air velocity.The expression of d b is relatively simple; however, it was found to work well in a number of studies [23,24].C s is the saturation DO concentration, related to water temperature and pressure [25].In this study, K L,B is calculated as the maximum value between bubbles rising in stagnant liquid [26] and bubbles in turbulent flow [27], and it can be expressed as Equation (15).
where D m is the oxygen molecular diffusion coefficient related to temperature ν 1 and ν 2 are the kinematic viscosity coefficients of mixture, air bubbles and water, respectively.P e is the Peclet number.

Initial and Boundary Conditions
The two sections of I-I and II-II, as shown in Figure 2, were regarded as the inlet and outlet boundaries in the mathematical model, and the computational domain was set from section I-I to section II-II.At the I-I section boundary, the water velocity measured by the flowmeter and DO Water 2018, 10, 174 6 of 13 concentration measured by 1# DO probe were used as the desired values respectively, and the turbulent parameters k and ε were determined approximately using Equations ( 16) and (17).
where U in is the average velocity at I-I cross section, I is the turbulence intensity, and l is the turbulence length scale.In addition, l is a physical quantity related to the size of the large eddies that contains the energy in turbulent flows.In fully-developed pipe flows, l is restricted by the size of the pipe, since the turbulent eddies cannot be larger than the pipe size.An approximate relationship between l and the physical size of the pipe is expressed as Equation (18).
where the factor of 0.07 is based on the maximum value of the mixing length in fully-developed turbulent pipe flow.At the II-II section boundary fully developed turbulence was assumed, and the variation ratios for all the variables were set as 0. At the two holes boundaries, the atmosphere pressure was assumed, and the DO concentration was set as the same value as at the I-I section without high deviation.

Numerical Details
The CFD solver Fluent 16.0 was utilized to solve the governing equations.The Finite Volume Method was used to discretize the governing equations.Pressure, turbulent kinetic energy, turbulent dissipation rate, air volume fraction, and DO concentration are set at the middle point of the cell.The velocity components were set at the cell edges.The convection term was solved with QUICK scheme, the diffusion term was solved with central scheme, and the SIMPLEC algorithm was used to couple pressure and velocity.In FLUENT 16.0 software, the User-Defined Scalar (UDS) method reported by Yin et al. [20] was used to compute DO concentrations by compiling DEFINE_SOURCE and DEFINE_DIFFUSIVITY macros in the C language environment.In addition, the grid acceleration convergence techniques were used to achieve better convergence speed and accuracy.
Solidworks 2012 software (Dassault Systèmes SolidWorks Corporation, Waltham, MA, USA) was used to generate the 3-D geometric body, and the ICEM 14.5 software (ANSYS, Canonsburg, PA, USA) was utilized to partition the computational cells.The refined cells were used to partition the local domain near the two holes.The relatively coarse cells were used to partition the other domain.Grid independence was tested using the grid convergence index (GCI) method, which is an acceptable and recommended method that has been evaluated in over several hundred CFD cases [28][29][30][31].Table 1 shows the discretization error for the water velocities at the centre point of the 2 holes cross section for water flow rate Q w = 2.05 × 10 −4 m 3 /s and DO concentration at the inlet C in = 5 × 10 −4 kg/m 3 .The extrapolated relative error and GCI are computed by Richardson extrapolation.In x, y and z directions, three cell sizes are used in the GCI calculations: 0.00007 m × 0.00015 m × 0.00017 m of Grid 1, 0.00014 m × 0.00031 m × 0.00033 m of Grid 2 and 0.00028 m × 0.00062 m × 0.00066 m of Grid 3, respectively.The total cells number and the refined cells number for Grid 1, Grid 2 and Grid 3 are shown in Table 1.It was found that the discretization error is very small at 1.9% between Grid 2 and Grid 1 with the increasing refined cells number from 172,919 to 1,383,413.Thus, the Grid 2 scheme is relatively independent of the computational results.Therefore, it was selected for the following computations after considering the balance between accuracy and computational speed.

Model Validation
Figure 3 shows the computational DO concentrations and experimental data in the longitudinal direction for Q w = 2.05 × 10 −4 m 3 /s and C in = 5 × 10 −4 kg/m 3 .It was found that the maximum value of the relative deviation is 9.2% at the x = 0.2 m, while the other relative deviations are all smaller than 8%.
relatively independent of the computational results.Therefore, it was selected for the following computations after considering the balance between accuracy and computational speed.

Model Validation
Figure 3 shows the computational DO concentrations and experimental data in the longitudinal direction for Qw = 2.05 × 10 −4 m 3 /s and Cin = 5 × 10 −4 kg/m 3 .It was found that the maximum value of the relative deviation is 9.2% at the x = 0.2 m, while the other relative deviations are all smaller than 8%.The proposed mathematical model was further validated on air flow rate Qa using the experimental data reported by [12].In the experiment by [12], Dt is 0.027 m, D is 0.036 m, d is 0.005 m, and the downstream section length, L, is 0.25 m.Uin was controlled at 1.47 m/s, 1.94 m/s, 2.37 m/s, The proposed mathematical model was further validated on air flow rate Q a using the experimental data reported by [12].In the experiment by [12], D t is 0.027 m, D is 0.036 m, d is 0.005 m, and the downstream section length, L, is 0.25 m.U in was controlled at 1.47 m/s, 1.94 m/s, 2.37 m/s, 3.07 m/s, 4.35 m/s, 4.88 m/s, 5.76 m/s and 6.19 m/s respectively.Figure 4 shows the computational Q a in comparison to experimental data.It was found that Q a increases with increasing U in .The relative deviations between experimental values and computational results are smaller than 10% in general, and are consistent with each other.Figures 3 and 4 show that it is acceptable to predict the hydraulic parameters, such as C and Q a , in a HVFTH using the proposed mathematical model.4 shows the computational Qa in comparison to experimental data.It was found that Qa increases with increasing Uin.The relative deviations between experimental values and computational results are smaller than 10% in general, and are consistent with each other.Figures 3 and 4 show that it is acceptable to predict the hydraulic parameters, such as C and Qa, in a HVFTH using the proposed mathematical model.Water 2018, 10, 174 8 of 13

Results and Discussion
265 scenarios of D = 0.2 m and T = 20 • C were simulated with: d increasing from 0.015 m to 0.020 m, 0.025 m, and 0.030 m; contraction ratio at the throat section (δ = D t /D) increasing from 0.5 to 0.6, 0.7 and 0.8; L increasing from 2 m to 3 m; U in increasing from 2.78 m/s to 6.5 m/s; and C in increasing from 0 kg/m 3 to 0.001 kg/m 3 , 0.002 kg/m 3 , 0.003 kg/m 3 , 0.004 kg/m 3 and 0.005 kg/m 3 .
A relative saturation coefficient, E, is defined for analyzing DO recovery as Equation ( 19) [20].
where C out is the average value of the DO concentration at section II-II, and ∆t is the travel time from the two holes' cross section to section II-II.

Effect of Hydraulic Parameters
Figure 5 shows the E 20 relationship with U in for D = 0.2 m, L = 3 m and C in = 0 kg/m 3 , where E 20 is the E value at 20 • C water temperature.It was found that E 20 decreases with increasing U in .The possible reason is that as ∆t decreases due to the increasing U in , the bubbles flow faster in the Venturi pipe and oxygen has less time to transfer from air bubbles to the water body, resulting in the decrease of E 20 according to Equation (19).Therefore, U in plays a significant role in E 20 .

Results and Discussion
265 scenarios of D = 0.2 m and T = 20 °C were simulated with: d increasing from 0.015 m to 0.020 m, 0.025 m, and 0.030 m; contraction ratio at the throat section (δ = Dt/D) increasing from 0.5 to 0.6, 0.7 and 0.8; L increasing from 2 m to 3 m; Uin increasing from 2.78 m/s to 6.5 m/s; and Cin increasing from 0 kg/m 3 to 0.001 kg/m 3 , 0.002 kg/m 3 , 0.003 kg/m 3 , 0.004 kg/m 3 and 0.005 kg/m 3 .
A relative saturation coefficient, E, is defined for analyzing DO recovery as Equation ( 19) [20].
where Cout is the average value of the DO concentration at section II-II, and Δt is the travel time from the two holes' cross section to section II-II.

Effect of Hydraulic Parameters
Figure 5 shows the E20 relationship with Uin for D = 0.2 m, L = 3 m and Cin = 0 kg/m 3 , where E20 is the E value at 20 °C water temperature.It was found that E20 decreases with increasing Uin.The possible reason is that as Δt decreases due to the increasing Uin, the bubbles flow faster in the Venturi pipe and oxygen has less time to transfer from air bubbles to the water body, resulting in the decrease of E20 according to Equation (19).Therefore, Uin plays a significant role in E20., and it increases as a result.Therefore, Cin also plays an important role in E20.In addition, it was interesting to note that E20 sensitivity decreases with increasing Cin and δ.   19) can be rewritten as E = 1 + C out −C s C s −C in , and it increases as a result.Therefore, C in also plays an important role in E 20 .In addition, it was interesting to note that E 20 sensitivity decreases with increasing C in and δ.  19) can be rewritten as , and it increases as a result.Therefore, Cin also plays an important role in E20.In addition, it was interesting to note that E20 sensitivity decreases with increasing Cin and δ.

Effect of Geometric Parameters
Figure 7 shows that E20 increases with increasing d to some extent, especially for small δ scenarios, and an optimum d exists corresponding to the maximum E20.For large δ and d scenarios, E20 does not increase much with the continual increase of d.Therefore, the diameter of the hole plays a complex role in E20 and a simple comment on its effect is warranted.With a given water velocity at the inlet, the negative pressure usually occurs near the throat section of a conventional Venturi flow without hole(s).With the increase of the diameter of the holes from 0, the air flow rate through the holes increases, and the air volume fraction also increases, contributing to the increasing E20.With the continual increase of the holes' diameter, the increasing rate of Qa slows due mainly to the negative pressure attenuation.In the case that the local pressure near the holes reaches the outside atmosphere pressure, the air suction phenomenon vanishes.With a large enough hole diameter, the local pressure becomes positive and the inside water will jet into the outside atmosphere through the holes, instead of air suction.Therefore, an optimum value of the diameter of the holes can be investigated for a given scenario.In addition, E20 increases with decreasing δ as shown in Figure 7.A possible reason is that the negative pressure near the holes decreases due to the decreasing δ, and the air flow rate increases as a result, contributing to the increasing E20.Furthermore, E20 sensitivity attenuates in δ with the increase of the diameter of the holes.

Effect of Geometric Parameters
Figure 7 shows that E 20 increases with increasing d to some extent, especially for small δ scenarios, and an optimum d exists corresponding to the maximum E 20 .For large δ and d scenarios, E 20 does not increase much with the continual increase of d.Therefore, the diameter of the hole plays a complex role in E 20 and a simple comment on its effect is warranted.With a given water velocity at the inlet, the negative pressure usually occurs near the throat section of a conventional Venturi flow without hole(s).With the increase of the diameter of the holes from 0, the air flow rate through the holes increases, and the air volume fraction also increases, contributing to the increasing E 20 .With the continual increase of the holes' diameter, the increasing rate of Q a slows due mainly to the negative pressure attenuation.In the case that the local pressure near the holes reaches the outside atmosphere pressure, the air suction phenomenon vanishes.With a large enough hole diameter, the local pressure becomes positive and the inside water will jet into the outside atmosphere through the holes, instead of air suction.Therefore, an optimum value of the diameter of the holes can be investigated for a given scenario.In addition, E 20 increases with decreasing δ as shown in Figure 7.A possible reason is that the negative pressure near the holes decreases due to the decreasing δ, and the air flow rate increases as a result, contributing to the increasing E 20 .Furthermore, E 20 sensitivity attenuates in δ with the increase of the diameter of the holes.Figure 8 shows that E20 increases with increasing L for the scenarios of D = 0.2 m and Cin = 0 kg/m 3 .The reason is that with increasing L, Δt increases as a result, and it has more time to transfer oxygen from air bubbles to the water body, resulting in the increase of E20 according to Equation ( 19).In addition, E20 increases with increasing d and decreasing δ to some extent.Figure 8 shows that E 20 increases with increasing L for the scenarios of D = 0.2 m and C in = 0 kg/m 3 .The reason is that with increasing L, ∆t increases as a result, and it has more time to transfer oxygen from air bubbles to the water body, resulting in the increase of E 20 according to Equation (19).In addition, E 20 increases with increasing d and decreasing δ to some extent.
(c) (d) Figure 8 shows that E20 increases with increasing L for the scenarios of D = 0.2 m and Cin = 0 kg/m 3 .The reason is that with increasing L, Δt increases as a result, and it has more time to transfer oxygen from air bubbles to the water body, resulting in the increase of E20 according to Equation (19).In addition, E20 increases with increasing d and decreasing δ to some extent.

A Correlation for the Relative Saturation Coefficient
It is not easy to calculate the relative saturation coefficient by Equation ( 19), due to the difficulty determining k L,B a B t. Consequently, it is necessary to deduce a simple formula for the relative saturation coefficient considering the effects of D, D t , d, L, U in , C in , C s and ν 2 , and E 20 can be written as Equation (20).
where Re = U in D/ν 2 .The 249 groups of numerical data from the aforementioned 265 scenarios were used randomly to deduce the following Equation ( 22) with the least squares method, The correlation coefficient is 0.88.In order to validate Equation ( 22), the other 16 scenarios were used as shown in Table 2. Figure 9 shows the E 20 validation between Equation ( 22) and the corresponding results yielded from Equation (19).It was found that they agree satisfactorily with each other; all the relative deviations between them are smaller than 20%, and most of them are smaller than 10%.Therefore, it is reliable to compute the relative saturation coefficient using Equation ( 22) without causing large deviations.Note that Equation (22) has good applicability to predict our results, while it still needs further validation for practical applications.

Conclusions
A 3-D two-phase flow mathematical model, coupled with a DO transfer model and a series of experiments, were conducted to investigate the DO recovery at a HVFTH.The computational DO concentration and air volume fraction were validated with experimental data, and they agree with each other well.In addition, the validated mathematical model was used to compute the hydraulic parameters, DO concentration and relative saturation coefficient within different scenarios.It was found that with increasing DO concentration at the inlet and downstream section length, the relative saturation coefficient increases.However, it decreases with the increase of water velocity at the inlet and throat section contraction ratio, respectively.The diameter of the holes plays a complex role in

Conclusions
A 3-D two-phase flow mathematical model, coupled with a DO transfer model and a series of experiments, were conducted to investigate the DO recovery at a HVFTH.The computational DO concentration and air volume fraction were validated with experimental data, and they agree with each other well.In addition, the validated mathematical model was used to compute the hydraulic parameters, DO concentration and relative saturation coefficient within different scenarios.It was found that with increasing DO concentration at the inlet and downstream section length, the relative saturation coefficient increases.However, it decreases with the increase of water velocity at the inlet and throat section contraction ratio, respectively.The diameter of the holes plays a complex role in the relative saturation coefficient, and an optimum value may exist when the relative saturation coefficient reaches the maximum.The dimensional analysis method and the least square method were used to deduce a simple formula for the relative saturation coefficient considering the effect of geometric sizes and the hydraulic parameters, and it was validated well with the results of Equation (19).

Figure 3 .
Figure 3. Measured DO concentrations and predicted data in the longitudinal direction.

Figure 3 .
Figure 3. Measured DO concentrations and predicted data in the longitudinal direction.

Figure 6
Figure 6 shows the E20 relationship with Cin for D = 0.2 m, d = 0.02 m and L = 3 m.It was found that E20 increases with increasing Cin.For a given Cs in hypoxic water, with the increasing Cin, Cout increases,

Figure 6
Figure 6 shows the E 20 relationship with C in for D = 0.2 m, d = 0.02 m and L = 3 m.It was found that E 20 increases with increasing C in .For a given C s in hypoxic water, with the increasing C in , C out increases, C s − C in decreases, and C out − C s increases.Equation (19) can be rewritten as E = 1 + C out −C s C s −C in , and it increases as a result.Therefore, C in also plays an important role in E 20 .In addition, it was interesting to note that E 20 sensitivity decreases with increasing C in and δ.

that
E20 increases with increasing Cin.For a given Cs in hypoxic water, with the increasing Cin, Cout increases,

Water 2018 ,
10, x FOR PEER REVIEW 12 of 14

Table 1 .
Analysis of discretization error.

Table 1 .
Analysis of discretization error.