Numerical Modeling of Venturi Flume

: In order to measure ﬂow rate in open channels, including irrigation channels, hydraulic structures are used with a relatively high degree of reliance. Venturi ﬂumes are among the most common and efﬁcient type, and they can measure discharge using only the water level at a speciﬁc point within the converging section and an empirical discharge relationship. There have been a limited number of attempts to simulate a venturi ﬂume using computational ﬂuid dynamics (CFD) tools to improve the accuracy of the readings and empirical formula. In this study, simulations on different ﬂumes were carried out using a total of seven different models, including the standard k– ε , RNG k– ε , realizable k– ε , k– ω , and k– ω SST models. Furthermore, large-eddy simulation (LES) and detached eddy simulation (DES) were performed. Comparison of the simulated results with physical test data shows that among the turbulence models, the k– ε model provides the most accurate results, followed by the dynamic k LES model when compared to the physical experimental data. The overall margin of error was around 2–3%, meaning that the simulation model can be reliably used to estimate the discharge in the channel. In different cross-sections within the ﬂume, the k– ε model provides the lowest percentage of error, i.e., 1.93%. This shows that the water surface data are well calculated by the model, as the water surface proﬁles also follow the same vertical curvilinear path as the experimental data.


Introduction
The Parshall flume is a simple static measuring device with no moving parts that is used to determine the flow rate in an open channel where a constant recording of discharge is required. The initial idea by Ralph Parshall in designing the Parshall flume was to make it easier for water users, like farmers, who do not have access to sophisticated types of equipment, to be able to determine how much water is delivered to them with an acceptable level of accuracy [1]. Currently, Parshall flumes are mainly used in irrigation and sewer systems to measure the flowrate [2]. In general, they are designed to generate a critical flow within the throat section, which affects the water level along the converging section upstream, implementing an empirical relationship between water surface elevation and discharge results in finding the discharge value at a specific time from the water surface elevation.
Variations of Parshall flumes are typically restricted to the dimensions of the geometries proposed by Ralph Parshall in 1936 for limited size numbers, i.e., 16 sizes, which vary in the opening from 1" to 144". The arbitrary dimensions of the flumes used in various open channels in recent years have not been comprehensively studied. Since then, technology has advanced significantly, especially in the development of computational tools. Numerical simulations have introduced a new revolutionary chapter to the design of hydraulic structures, allowing engineers to extend their full potential by designing a variety of different hydraulic structures with new arrangements in their dimensions and shapes. They allow designers to optimize the specifications of their design to suit the need of the project's objective. Computers and modeling software offer possibilities to perform extensive calculations, and they allow new designs to be tested with relatively low cost. For instance, a proposed hydraulic structure can be simulated to study its behavior under working conditions. This provides engineers with extensive variations in geometries and dimensions when designing a structure. Using computer models makes it possible to implement Parshall flumes with complex shapes and dimensional limitations where needed.
The authors of [3] used a numerical model to develop an alternative rating equation to be implemented at low discharge for different flume sizes. In a study on a submerged Montana flume, to prove the reliability of the CFD programs' accuracy, the flow rate was calculated with FLOW-3D. It is shown that the numerical results closely matched the experimental results. Furthermore, it is revealed that a free-flow equation for a Parshall flume would also be a good fit for a Montana flume [4]. To determine the accuracy of a fieldscaled Parshall flume at a wastewater system in Minneapolis, Minnesota, [5] implemented the large-eddy simulation (LES) and level-set method to compute the turbulent flow under two-phase flow conditions. A three-dimensional finite difference code, SOLA-FLUMP, is presented by [6] to assess the various effects due to the sloped channel, the upstream velocity profile distortion, and the geometry of the flume.
A parametric study was performed in terms of stability and accuracy on Reynoldsaverage Navier-Stokes (RANS) and hybrid RANS/LES turbulence models and numerical schemes offered in openFOAM, an open-source software, by [7]. From the results, the second-order upwind scheme and limiters were found to be the most stable with the lowest computational cost, thereby providing the highest level of accuracy for RANS models.
It is important to determine the degree of reliability of the model that is going to be used for a specific simulation. One of the most dependable approaches is to evaluate the simulated results with actual experimental data. In this paper, different simulated datasets generated from seven dissimilar methods, i.e., standard k-ε, RNG k-ε, realizable k-ε, and k-ω SST; large-eddy simulation (LES); and detached eddy simulation (DES) are examined versus four physical datasets obtained from physical experiments on different flumes with different discharges.
A high-order partial differential equation, like the momentum transport equation, that includes nonlinear terms cannot be solved analytically to obtain a general solution. Numerical solutions require discretization techniques to reshape the continuous partial differential equation into a discrete equation. Computational fluid dynamics (CFD) using the finite volume method (FVM) [8] can be very useful for this purpose. With the rapid increase in the computational capability of computers, large-eddy simulation (LES) is increasingly embedded in the CFD models used by researchers and engineers to solve turbulent flows [9]. LES is a compromise between the efficiency of Reynolds-average Navier-Stokes (RANS) and the prohibitive computational cost of direct numerical simulation (DNS). Approaches of LES or the variant hybrid family like RANS/LES (DES) are progressively taking over the computationally expensive DNS approach to solve problems with compound geometries and flow properties [10].
The authors of [11] tried to introduce a correction coefficient to a 24" Parshall flume where the positions of the staff gauges were mislocated and the condition of the flume entrance was set up differently from the one introduced by [1]. They used numerical modeling to implement the correction factor for other sizes of Parshall flumes. Based on this study, a part of the small deviations, the physical model, and numerical simulation were aligned with one another.
In order to decrease the head loss in a curved flume, three flumes were studied on the basis of critical flow by [12]. The study was conducted on laboratory experimental results versus numerical simulation data. The hydraulic parameters such as velocity of the free surface and the depth of water were analyzed and compared. A maximum error of 4.7% in the water depth was obtained. It was shown that a good consistency was achieved between numerical simulation and experimental data.
The objective of this paper is to solve the Navier-Stokes equations using OpenFOAM to study the behavior of various numerical models in simulating Parshall flumes similar in size and flowrate to the flumes used in physical experiments by Dursun (2016).
Based on the literature review at the beginning, there is no comprehensive numerical solution that has been conducted on the models simulating the hydraulic structures, i.e., Parshall flume. The reliability of RANS LES and DES turbulence models in simulations has been studied in this research by evaluating the performance of seven turbulence models against the experimental data of different Parshall flume structures. Subsequently, the consistency of the simulations was determined in various scenarios in relation to the experimental data. Based on what was discussed in the literature review, seven turbulence models, i.e., RAS models (including standard k-ε, RNG k-ε, realizable k-ε, and k-ω SST) and hybrid RANS/LES models (such as k-ω, SST-DES, and an LES models, namely, the Smagorinsky method and dynamic K LES method), respectively, were selected due to their wide usage and advantages compared to other turbulence models.
This paper is organized as follows: the methodology is explained in the next section, including governing equations, turbulence models (including RANS, LES, and DES models), numerical setup (including initial and boundary conditions), mesh analysis, and data; the results and discussion are then presented, in which the performance of numerical models is discussed. Some concluding remarks complete the study.

Methodology
By the increasing power of processors, computational fluid dynamics (CFD) has become the most convenient tool to simulate fluid motion. In CFD, the fluid's flow establishment follows physical parameters, including pressure, viscosity, velocity, and temperature. In order to simulate a physical case related to the fluid flow, the physical properties should be taken into account accurately.
CFD approaches are used in solving the fluid flow equations as well as fluid interaction with solid bodies. The Euler equation for inviscid fluid and the Navier-Stokes equation for viscous fluid can be derived in their integral arrangement with respect to the conservation of energy, mass, and momentum [13].
OpenFOAM is one of the open-source solvers for CFD that is widely used for simulations. It is a platform including numerous C++ libraries and applications that are able to solve numerically the continuum mechanics problems [14]. It uses a tensorial method that implements an object-oriented programing approach and employs the finite volume Method (FVM).

Governing Equations
A viscous incompressible fluid flow is governed by a general three-dimensional system of equations, called the Navier-Stokes system, that consists of momentum and continuity equations. The system is described as follows [15,16]: where ρ denotes density; p represents total pressure; u, v, and w represent the velocity in three different directions, i.e., x, y, and z; t is used for time; and gravitational acceleration is denoted by g. ρ is obtained using the following equation: Here, ρ 1 and ρ 2 represent the air and water densities, the two phases of the involved fluid. The α value varies from 1 to 0 depending on the location, where 1 denotes the presence of water and 0 shows the presence of air. Any number between these two values represent the interface. Finally,

Equation of the Free Surface
With respect to the zero pressure at the surface, the free surface was analyzed with the volume-of-fluid (VoF) method. The following equation is used by VoF: As stated earlier in this paper, a powerful open-source CFD, i.e., OpenFOAM application, was used for the purpose of the numerical simulations.
The flow motion in a Parshall flume was simulated in this study with the help of seven different turbulence models: standard, realizable, and RNG k-ε models; k-ω SST models; and detached eddy simulation (DES) models, such as k-ω SST-DES; as well as LES methods, including the Smagorinsky LES model and dynamic K LES model. In the following, a brief description of each selected model is provided.

Reynolds-Average Navier-Stokes (RANS) Approach
The Reynolds-average Navier-Stokes Model is currently the most popular approach for the simulation of fluid flow. This approach essentially uses a viscosity term to approximate the turbulence equations. K is a term in these models that represents the fluctuations of the turbulence kinetic energy per unit mass.

Standard k-ε Model
This model requires two additional transport equations: one for turbulent kinetic energy (k) and another for energy dissipation (ε). Apart from its poor performance in large adverse pressure gradient cases, it is known as one of the most popular turbulence models [17]. This model comes from the Reynolds-averaged Navier-Stokes (RANS) category where modeling is applied to all properties of fluid motion.
The equation for turbulent kinetic energy k and dissipation ε is shown below: ε is the component that controls the turbulence scale where k represents the turbulence kinetic energy. The reader is referred to [15] for further details and values of the coefficients.

Realizable k-ε
The latest improved form of the three k-epsilon models is the realizable k-epsilon model [18,19]. There are two significant differences when this model is compared to the standard k-epsilon model. Firstly, formulation of the turbulence viscosity has been revised. Secondly, the dissipation rate transport equation is explained based on the equation of transport of the mean-square vorticity [19,20].
The reader is referred to [21] for further details and value of the coefficients.

LES Approach Smagorinsky LES Model
This model was originally developed within the metrological community to simulate atmospheric air currents [22]. As a well-known subgrid-scale model according to [23], the Smagorinsky model estimates the shear as where The reader is referred to [22] for further details and values of the coefficients. One of the main disadvantages of this model is the lack of ability to predict the energy transfer from subgrid-scale structures to the greater resolved scales; thus, the model is totally dissipative. Another problem with the Smagorinsky model is that its coefficient has to be adjusted for every flow field. However, it is still one of the well-known models in the field of CFD. The Smagorinsky model is easy to use in numerical simulations. If the Smagorinsky coefficient is adjusted based on the local characteristics of the fluid motion, it can generate more accurate results [23].

The Dynamic One-Equation Model
The SGS stresses determine how successful an LES model can be. As a simple model, in the Smagorinsky model, the factor of the proportionality is a fixed value that has to be determined before running the model. In reality, the factor is a flow-dependent value and is not defined as a single universal constant. The weak point of the model comes from this section. There have been attempts to improve this model [24,25]. Moreover, this model is completely dissipative, and there is always a transformation of large-to-small scale for the energy.
On the other hand, dynamic models are the best choice to substitute the Smagorinsky model. In this model, the C value of the subgrid eddy viscosity is determined while the simulation is computed [26]. In recent decades, a one-equation dynamic model has been presented [27]. The equation of the dynamic model is presented below.
where the C e and S k coefficients are derived from local flow properties [21]. The reader is referred to [28] for further details.

DES Approach k-ω SST-DES
While the RANS model is derived through Reynolds temporal averaging, the LES model is the result of spatial filtering. Between these two methods, the difference is the magnitude of the generated eddy viscosity. This is one of the main reasons for the development of the DES model, which has the ability to cover the weaknesses of LES when it comes to treating wall regions with very fine mesh [29].
The RANS approach is used in the near-wall region by the DES method where, at the same time, the LES model is applied to the rest of the region excluding the wall region. The region associated with the LES model is usually the core turbulent area, where large-scale turbulences play a major role. Within this area, the DES approach uses a LES subgrid-scale model, while for the near-wall region, it uses the RANS model [30].
A DES-improved form of the k-ω SST method is the k-ω SST DES approach [31,32]. Recently, in the aerodynamic field, DES has been widely implemented due to its computational speed and quality of results. It is proven to be less computationally expensive, and it generates better results than steady RANS [33,34].
The turbulence specific dissipation rate equation is given by (16) and the turbulence kinetic energy is calculated as the length scale, d, is given by and the turbulence viscosity is obtained using The reader is referred to [35] for further details and value of the coefficients.

Numerical Solution Details
The combination of the finite volume method with the VoF method was used in the model. To solve the governing equation of motion, the "interFoam" solver was implemented in OpenFOAM. The temporal term was discretized with the help of a Eulerian scheme, besides the Gauss linear method, which is used for the gradient term. For the Laplacian, the corrected Gauss linear method was applied. The divergence terms were discretized using a Gauss vanLeer plus Gauss linear scheme. The linear scheme was used to discretize the interpolation term.
The tolerance level was defined for each individual variable, while the desired convergence was expected to be achieved where an iterative solver was used in the processing. The Gauss-Seidel technique was applied with a level of accuracy of 10 −5 for the fraction of liquid (α), 10 −8 for pressure, and 10 −8 for the velocity.

Initial Conditions
The flume had a fixed inlet flow velocity for each case, i.e., a discharge of 20 L/s or 30 L/s as the initial conditions. There was no initial acceleration or dissipation defined in the model. There was no flow across the walls, i.e., there was no liquid coming in or getting out through the wall boundaries at the defined locations. Initial water surface was set to be constant everywhere.

Boundary Conditions
There were different boundary conditions employed in this simulation as illustrated in Figure 1. Hydraulically smooth walls were considered in this study, and standard wall functions were employed. Water discharge was specified at the inlet. Zero gradient condition was considered at the outlet. The free surface was tracked based on the volumeof-fluid method based on a zero pressure condition at the interface of air and water.

Mesh Analysis
In this study, in order to determine the optimum size of the mesh, mesh sensitivity analysis was conducted. The mesh grid used in this paper is a structured mesh. The purpose of this analysis is to find the finest mesh size for which the results will not be affected further. Figure 2 shows the mesh sensitivity analysis performed on the Parshall flume. The maximum cell number was 263,700 cells for this structure. The optimal case in terms of computational cost and changes in results is that illustrated in Figure 2b, with a total of 74,496 cells. Greater reduction applied on the cell size did not significantly change the numerical results. In Figure 3, the results of the mesh sensitivity analysis are provided for the crosssections 4-6. As the size of the original coarse mesh became finer, the simulated water levels approached the values of the experimental data. In this study, the initial number of cells was 52,200, and it was gradually increased to 74,496 and further to 263,700 cells. The results from the last two finer meshes show negligible difference; hence, the mesh containing 74,496 cells was selected and further used to implement the remaining turbulence models.

Data
The experimental tests were conducted in the hydraulic laboratory of Firat University, Elazig, Turkey (Dursun, 2016). The tests were performed in a rectangular channel with the dimensions of 0.4 by 5 by 0.6 m in width, length, and depth, respectively. The purpose of these experiments was to determine the changes in the quantity of the dissolved oxygen of the stream before and after a Parshall flume structure was introduced. The dimensions of the Parshall flumes that was used in this research were the same as the 3-inch (7.62 cm) Parshall flume with a 45 degree wing wall. The discharge in the experiment was measured using an electromagnetic flow meter.
The model was run for 300 time-steps where, after 120 steps, water levels became steady; water surface fluctuations, for instance, were limited to +/−0.003 m. The water level is shown using a post-data analysis software called ParaView. At the point of interest, two vertical planes (perpendicular to each other) were introduced, where the intersection line between the two planes was set to pass the desired point. The Y-coordinate (water level) of the line was extracted using the calculator filter; for example, the water level at each desired time step is shown afterwards in the tabulated mode in a separate window in ParaView.

Results
The numerical simulation for the experimental case conducted by [2] is performed in this study. The Parshall flume used in this study is a 3-inch Parshall flume modified to meet the experimental criteria. OpenFOAM was used as an open source CFD tool to carry out the numerical simulation of the Parshall flumes.
Switching from a coarse size mesh grid to a finer grid size in the mesh sensitivity analysis led to some change seen in the numerical simulation model results. As the flow rate decreased from 30 to 10 L/s, the model tended to produce better results. This pattern is observed in the k-ε model as well as the other RANS models.
Water surface elevations in seven cross-sections were compared with the experimental results with three different flowrate values, i.e., 10, 20, and 30 L/s. Figure 4 shows water surface elevation for the flowrate of 10 L/s using seven various models. As shown in Figure 4, all models follow the same pattern as the experimental data. Figures 5 and 6 are introduced to show more clarified view of the water level data. Based on the results of various turbulence models used in simulations, the k-ε turbulence model provides the most accurate simulation compared to the other models used in this study.   The locations of the cross-sections in the Parshall flume are denoted in Figure 7 where at all critical locations, a cross-section is introduced. The reason for choosing the crosssection locations, as shown in Figure 7, is that in the laboratory experiment conducted by Dursun in 2016, the same locations were chosen; hence, water level data were also available for them. In Figure 8, the fluid flow in the Parshall flume with a flowrate of 10 L/s is illustrated. As the stream lines indicate, the velocity of the fluid prior to the flume throat displays a maximum of 1 m/s, and by the time it reaches the narrow section with a declining floor, the speed rapidly increases by 50%. This is the section where the fluid experiences supercritical flow. It continues to increase by the end of the divergence section where maximum velocity is reached (red arrows), i.e., 2 to 2.5 m/s. Once the fluid reaches the inclined slope in the divergence section downstream, the flume forces it to develop a hydraulic jump. The velocity profiles of all seven turbulence models at cross-section 5, where the numerical models exhibited the maximum velocity, are shown in Figure 9. The comparison between different velocity profile shows that in the seven different turbulence models, the values are almost identical, with a difference margin of 0.04 m/s at the maximum velocity points. The maximum recorded velocity in the flume for different models is within the range of 1.26-1.28 m/s in the different models. The "y" axis is the local coordinate of the flume throat defined in the OpenFOAM computational domain. Figure 9. Comparison of the velocity profile at cross-section 5 and the results obtained using the investigated turbulence models.
As illustrated in Figure 10, the velocity is constant upstream of the midsection of the converging area. At the end of this section, the flow is forced to increase its velocity due to the narrow design of the throat combined with a sharply slopped bed. The velocity streamlines after the throat section of the flume show that the velocity magnitude is greater than that observed throughout the rest of the structure. Once the flow immediately exits the divergence section, it attains maximum velocity. The flume's pressure field is represented in Figure 11. In the section between crosssections 5 and 6, the pressure is negative. As illustrated in this figure, due to the throat section, the pressure is built up in the upstream, while when passing cross-section 5 toward downstream, the value of the pressure drops rapidly.

Discussion
The percentage difference of water level values of the numerical models and the experimental data is specified in Table 1 below. The calculation is performed based on the following relationships [15]: Here, the root mean square error (RMSE), centered root mean square error (CRMSE), and correlation coefficient (R) are used in order to compare the simulation results with experimental data. The results from the above formulas are presented in Table 1, where different types of error analysis are applied in order to determine the reliability of the numerical simulations. Among the seven different turbulence models applied in this study, the Smagorinsky and dynamic K equation models provide the least average error values. This model provides the minimum error percentages not only in the average value of all of the crosssections but also within the first three, i.e., 1, 2, and 3.
By calculating the root mean square error for all the turbulence models, the results show that two models, k-ε and dynamic K LES, provide the lowest values, 1.93% and 2.08%, respectively, while the k-ω SST model has a greater value compared to the others. However, for comparison of the R 2 values of this turbulence model, the Smagorinsky and dynamic K LES models have the closest value to 1, i.e., 0.998.
The ranges of the average error percentage for all of the turbulence models are within a narrow domain a minimum value of 1.93% and a maximum value is 2.53%. While most of the turbulence models used in this study provided reliable error percentages, the k-ω SST model was the only one that generated results with a high average error of 2.53% and RMSE of 0.64%, and the lowest R 2 was also recorded for this model. Table 1 and Figure 12 tabulate the results of the different numerical simulations at various locations in the flume, i.e., at the seven cross-sections. The average values and the statistical analysis that were discussed earlier in this section are also included here. At cross-section 1, the error percentages for all the models were larger than those calculated at the next cross-section (cross-section 2). Since the location of the first crosssection was chosen very close to the defined inlet of the flume, the unsteady water surface affected the results at this point.

RANS Models
The best average performance was obtained using the standard k-ε turbulence model. The other turbulence models that fall into this family showed less accuracy than the LES models. However, the difference percentages at Sections 2-4 are better than those obtained using the LES and DES models.

DES Model
The only DES model used in this study, k-ω SST DES , provided almost the same average result as those obtained using the k-ω SST , which has the greatest average error percentage among the models. The second highest average error percentage was observed for the DES family of turbulence models.

LES Model
Two LES models were used in this study. The water surface elevation results from the Smagorinsky model and dynamic K LES are similar to those of the RANS model family.
The average error value, as mentioned earlier, is the lowest for dynamic K LES following the standard k-ε model.
The pressure and velocity fields for all the models follow same pattern, where the low velocity of the flow is converted to supercritical flow just after entering the throat section of the flume. The pressure field exhibits opposite behavior as the initial high-pressure flow enters the throat of the flume: the pressure drops and reaches a negative value in the throat over a short distance.
At the sections where the velocity is locally maximum, the pressure value drops dramatically. This is due the higher viscous losses occurring with the high velocity flow.

Conclusions
The main objective of this study was the numerical examination of fluid motion in a Parshall flume. Assessing the results of different turbulence models, such as standard k-ε, realizable k-ε, RNG k-ε, k-ω SST, k-ω SST DES, Smagorinsky, and the dynamic K LES. reveals that despite some poor performance of k-ω SST within two cross-sections, the same as the rest of the turbulence models except DES family models, this method estimates the water level accurately enough overall.
This study shows that the results from a turbulence model, in general, a CFD model like OpenFOAM, can provide reliable solutions to Parshall flume design problems. CFD models are able to simulate Parshall flumes to find the optimum design specifications for different scenarios like design changes to increase the efficiency of Parshall flumes. The results are provided with the lowest cost compared to experimental methods.
At the first location chosen to collect data (cross-section 1), the quality of simulated results was not adequate, as discussed in the previous section. To eliminate this issue, the authors recommend conducting further studies on this type of structure.
This study can be further expanded by proposing some slight changes in the design of the Parshall flume, for example, increasing the reliability in terms of the gauge reading of the stilling well at the converging section and also providing more accurate empirical relationships with respect to these changes.
It is recommended as a continuation of this study to conduct a detailed investigation into finding a proper correction factor for the numerical simulated models when different higher flow rates are used. In this study, the simulation data with a 10 L/s flowrate achieved the lowest error value compared to the experimental results, while those with larger flowrates, such as 20 or 30 L/s, where experimental data were available, demonstrated higher error values.
In order to conduct more accurate assessments of numerical model performance, it is recommended that future studies use more sources of experimental data to reduce the impact of experimental data error on the assessment of simulated results.