Impact Statement
Insights about vertical weir and sill overflow are laid out to inform the improved design and operation of weirs to enhance flow control and measurement. The derivation of the weir discharge coefficient is given a new treatment and practically applied as a type of Froude number, giving future adopters of this definition a more intuitive understanding of what contributes to the discharge capacity of an overflow structure. Thresholds for different regimes of weir flow measurement are set forward and justified to provide closure to what constitutes ideal operating conditions. Application of this research to other types of structures beyond vertical thin weirs is also possible, especially in consideration of flow over a tilting weir – which represents a translational structure that operates dynamically between the bounding cases discussed in this study.
1. Introduction
1.1 Curvilinear weir flow
Free-surface flow over a vertical weir has been a problem of considerable practical as well as theoretical interest for many centuries (Eckert, Reference Eckert2024). It has proven useful for determining the efflux of water in hydraulic engineering applications, while also representing an amalgam of several canonical flows that are of interest to hydrodynamicists. Upstream of the weir, essentially one-dimensional hydrostatic open-channel flow is present, which transitions in the approach to the weir to convergent flow through a sharp-edged orifice (Perry, Reference Perry1949) and finally to a gravity-driven free jet (Tuck, Reference Tuck1976). Additionally, a turbulent boundary layer subjected to an adverse pressure gradient (Na & Moin, Reference Na and Moin1998; Kitsios et al., Reference Kitsios, Sekimoto, Atkinson, Sillero, Borrell, Gungor, Jiménez and Soria2017), created by the stagnation along the vertical weir face, undergoes separation and results in a zone of recirculating flow at the base of the weir. The combination of phenomena in this flow field rewards it with the designation of being ‘an exceedingly complex problem at best’ (Street et al., Reference Street, Watters and Vennard1996) – making it a suitable candidate for fundamental research and empirical simplifications when prudent.
The viability of the weir as a reliable device for measuring discharge is dependent on the following two main conditions: (1) the flow field can be well approximated as two-dimensional in the streamwise and vertical directions; and (2) the falling nappe is sufficiently ventilated so that atmospheric pressure is present along both the upper and lower free surfaces.
The dynamic parameter describing a particular weir flow condition is given by the ratio
$h/p$
. This inertial condition with respect to the size of the weir obstacle is defined by the ratio of the pressure head (
$h$
) measured at the gauging section, which drives the flow over the weir, to the height of the weir crest relative to the channel bottom (
$p$
) (see figure 1).
Figure 1. Schematic of curvilinear flow over a thin weir. It is assumed that the flow field can be approximated as two-dimensional and that the overflowing nappe is fully supported by atmospheric pressure. All dimensions plotted are scaled from 3-D numerical results for a specific flow case where
$h/p=0.68$
(case N0.68 in table 2). A separation streamline is plotted to indicate the locations where the boundary layer becomes detached due to the adverse pressure gradient created by the weir and where it then eventually reattaches. The shape of the separation zone will differ under varying flow conditions. EGL is the energy grade line, given by the dashed line.
Weir flow can be considered to have two bounding cases as defined by
$h/p$
. The wall overflow condition occurs in the limit as
$h/p \rightarrow$
0, occasionally called an ‘infinitely high’ weir (Lauck, Reference Lauck1925; Rouse, Reference Rouse1932), and is applicable to when the structure is placed on the crest of a dam spillway. As
$h/p \rightarrow \infty$
, the obstructive effect of the weir vanishes and the flow becomes a free overfall (Rouse, Reference Rouse1936; Rajaratnam & Muralidhar, Reference Rajaratnam and Muralidhar1968). The obstacle is often described as a sill rather than a weir when approaching this upper limit (Ramamurthy et al., Reference Ramamurthy, Tim and Rao1987). Typical weir channel flows observed in hydraulic engineering practice have
$h/p =\mathcal{O}(10^{0})$
(Ackers et al., Reference Ackers, White, Perkins and Harrison1978). In the gauging section, the volumetric discharge per unit channel width is defined as
$q$
, with the characteristic velocity scale in this section being
$U_{0} = q/D$
.
The critical section, or vena contracta, is the location of minimum specific energy head (
$E$
), identified as the longitudinal position where the maximum elevation of the lower nappe surface occurs (Montes, Reference Montes1998). This definition of the critical section necessarily differs somewhat from that used for hydrostatic open-channel flows. For curvilinear weir flow, the pressure is no longer hydrostatic in the transition from subcritical to supercritical and the velocity is no longer uniform. Therefore, the critical section cannot be simply defined here as where the mean flow speed equals the small wave celerity. Instead, it is defined more generally as the location where the derivative of specific energy head with respect to flow depth equals zero. For weir flows, this has been shown to be where the elevation of the bottom boundary reaches a maximum (Castro-Orgaz & Hager, Reference Castro-Orgaz and Hager2017).
The thickness of the nappe in this section is noted as
$t$
and the average velocity in the nappe at this section is given by
$\overline {u} = q/t$
. The relative effects of draw-down on the upper nappe surface, and deflection of the lower nappe surface above the crest elevation, are given by
$\delta$
and
$\varepsilon$
, respectively, so that
$h = \varepsilon + t + \delta$
. Here,
$\eta$
is chosen as an arbitrary parameter to represent the vertical distance of a fluid layer below the free surface at the gauging section.
1.2 Governing equations and review of previous work
1.2.1 Hydrodynamics perspective
The physics governing the motion of fluid over a weir is described by the steady, incompressible Navier–Stokes equations. For purposes of determining steady discharge, an assumption of two-dimensional (2-D) flow is made, meaning that the three-dimensional (3-D) flow field is taken to be uniform in the transverse (
$y$
) direction (the dimension over which the crest spans the channel). This assumption allows boundary drag along the channel side walls to be neglected, along with lateral contraction and detachment of the nappe. This assumption is taken to be valid as long as the approach channel width (
$B$
) is at least eight times
$h$
(Bureau of Reclamation, 2001). With these simplifications, the coupled pressure-velocity field is given by (1.1), with the continuity condition of (1.2):
\begin{equation} \boldsymbol {u} \cdot \nabla\boldsymbol {u} = -\frac {1}{\rho }\nabla P - g\hat {k} + \nu \nabla ^{2} \boldsymbol {u}\:, \end{equation}
\begin{equation} \nabla \cdot \boldsymbol {u} = 0\:, \end{equation}
where
$\boldsymbol {u}$
is the two-dimensional velocity field (
$\boldsymbol {u} = u \hat {i} + w\hat {k}$
);
$P$
is pressure;
$\rho$
and
$\nu$
are the fluid density and kinematic viscosity, respectively; and
$g$
is the gravitational constant (equal in this study to 9.81 m s
$^{-2}$
).
It was recognised early that the wall overflow case is a suitable candidate for analyses invoking the inviscid fluid assumption regarding the convergent nature of the flow field and assumed negligible effect of boundary friction. Approximate analytical solutions can be found via complex analysis of the potential velocity, which uses numerical integration methods to solve for the two-dimensional fluid velocity field described by the stream function (
$\psi$
) and velocity potential (
$\phi$
). For steady, inviscid, incompressible and irrotational flow, boundary conditions for this approach are given by the Bernoulli equation:
\begin{equation} \frac {P}{\rho } + \frac {1}{2}(u^{2} + w^{2}) + gz = \text{constant}. \end{equation}
See Dias & Tuck (Reference Dias and Tuck1991) for a good formulation of the complex-potential flow analysis problem. Early notable advancements in the field taking this approach were those of Lauck (Reference Lauck1925), who approximated the shape of the free-surface and estimated the contraction coefficient (
$C_{c}=t/h$
) for the wall overflow. For this case, Lauck (Reference Lauck1925) found
$C_{c}$
to be the same as Kirchhoff’s coefficient for a planar free jet issuing from a sharp-edged orifice
$[(\pi/(\pi + 2) \approx 0.611]$
. Notably, Strelkoff (Reference Strelkoff1964) expanded the complex-potential flow analysis for the weir problem by putting forward a formulation that proposed a solution for the free-surface profile and discharge characteristics for any value of
$h/p$
. Similar approaches have also been fruitfully applied to approximate the profiles of the free overfall (Tuck, Reference Tuck1976; Dias & Tuck, Reference Dias and Tuck1991), as well as subcritical (Vanden-Broeck & Keller, Reference Vanden-Broeck and Keller1987) and supercritical (McLean et al., Reference McLean, Bowles, Scheichl and Vanden-Broeck2022) weir flows.
As elegant as these analytical solutions of potential flow are, it is not often that these results, and the assumptions of inviscid and irrotational flow implicit within them, are validated against experimental observations. The flow physics explored in these analyses can also sometimes lack practical value, such as imposing a supercritical upstream boundary condition for flow over a weir (McLean et al., Reference McLean, Bowles, Scheichl and Vanden-Broeck2022). It was noted by Rouse (Reference Rouse1964) that the inviscid results of Strelkoff (Reference Strelkoff1964) generally agreed well with experimental results for small values of
$h/p$
, but began to deviate as
$h/p$
increased towards the free overfall limit. This occurs because in this limit, the size of the turbulent boundary layer in relation to the overfall flow field renders invalid the assumptions of inviscid and irrotational flow. Here, we are met again with the infamous injunction to resort to empirical investigations when dealing with the complexities of turbulent flows. In addition to this, the analytical methods previously described proved overly tedious for the hydraulic engineer, whose primary concern was a simple and reliable estimation of discharge. For these reasons, we will now trace the development of the weir and sill flow studies through the practice of hydraulics.
1.2.2 Hydraulics perspective
In the early days, the Italian marquis Giovanni Poleni formulated what is now known as the weir-discharge or Poleni equation by conceptualising the flow over a thin vertical weir as an orifice flow capped by the free surface and driven by an upstream reservoir of hydrostatic stagnation pressure head (
$h$
). Using the Torricelli principle (see (
1.3)) to estimate the velocity of fluid layers within the overflowing jet and integrating from the crest elevation,
$p$
, to
$h+p$
, the following equation (sometimes referred to as the Poleni equation) was eventually put forward as a means for estimating volumetric discharge per unit channel width (
$q$
) passing over the crest of a vertical weir:
\begin{equation} q \approx \frac {2}{3}\sqrt {2g}h^{3/2} \rightarrow q = C_{d}\frac {2}{3}\sqrt {2g}h^{3/2}\:. \end{equation}
Like (
1.3), (
1.4) is derived from assumptions of inviscid and irrotational flow. It also neglects any contraction within the nappe, as well as the velocity head at the gauging section,
${U_{0}}^{2}/2g$
. It is well known in practice that these assumptions are not strictly true. For this reason, the introduction of the discharge coefficient (
$C_{d}$
) is found in (
1.4) to account for these violated assumptions. A reliable empirical equation for
$C_{d}$
as a function of practically measurable parameters,
$h$
and
$p$
, was the goal of much of the previous hydraulics studies. Foundational early research can be credited to Bazin (Reference Bazin1898) and Schoder & Turner (Reference Schoder and Turner1929) for weir flows – the latter of which also examined the influence of the velocity distribution in the approach channel on the discharge characteristics. A robust empirical equation relating the magnitude of
$C_{d}$
to
$h/p$
for weir flows of
$h/p =\mathcal{O}(10^{0})$
was put forward by Theodor Rehbock at the Karlsruhe Institute of Technology, described in English in his discussion of the results of Schoder & Turner (Reference Schoder and Turner1929) in King et al. (Reference King, Nagler Floyd, Streiff, Parshall, Pardoe, Ballester, Williams Gardner, Th., Lindquist Erik and Clemens1929). Rouse (Reference Rouse1936) extended the discharge-measuring capacity of these types of structures further by applying a similar analysis to the free overfall. Kindsvater & Carter (Reference Kindsvater and Carter1957) also notably included a discussion of scale effects due to viscosity and surface tension at low values of
$h$
, suggesting that these can result in a greater than 2 % variation in the prediction of discharge when
$h\lt 5$
cm (assuming the fluid is water).
A link between the results of Schoder & Turner (Reference Schoder and Turner1929) for weir flows and Rouse (Reference Rouse1936) for the free overfall was first provided by Kandaswamy & Rouse (Reference Kandaswamy and Rouse1957), who examined the discharge characteristics of weir and sill flows spanning the full range of
$h/p$
. A commonly used weir-discharge formulation found today is that of Kandaswamy & Rouse (Reference Kandaswamy and Rouse1957), where for
$h/p \leq 6$
:
$C_{d}=0.611 + 0.075h/p$
and for
$ h/p \geq 16.67$
:
$C_{d} = 1.06(1+p/h)^{3/2}$
. This formulation features a non-monotonic relationship between
$h/p$
and
$C_{d}$
, which Kandaswamy & Rouse (Reference Kandaswamy and Rouse1957) hypothesised was due to a maximum in the relative nappe thickness as the flow changes from a weir to sill regime. The reason for this phenomenon has so far not been elaborated and represents one of the main goals of this paper.
1.2.3 Bridging the gap
Although the dependence of
$C_{d}$
on
$h/p$
has received much attention in the literature, a physical description of
$C_{d}$
has remained elusive. It is not immediately apparent which of the aforementioned assumptions regarding irrotationality, contraction and upstream velocity head are a fair approximation in the majority of cases of weir flow, the conditions under which these approximations might break down and how this might inform a physical description of
$C_{d}$
. A survey of common textbooks used in undergraduate fluid mechanics courses in the United States reveals that
$C_{d}$
is described as a ‘lumped parameter’, accounting for friction loss between the gauging section and crest section, contraction of the overflowing nappe, contribution of the upstream velocity head, and scale effects due to surface tension and viscosity (Sturm, Reference Sturm2001; Çengel & Cimbala, Reference Çengel and Cimbala2006). However, we believe what has been lacking in the literature until this point is a cogent description of the contributing factors underlying the physical meaning of
$C_{d}$
through its decomposition from experimental and numerical results to show which of the aforementioned effects are most dominant. For this reason, we choose in this paper to develop a new semi-empirical expression for weir discharge using dimensional analysis in the hopes of demystifying the meaning of
$C_{d}$
.
Elucidation of the nature of
$C_{d}$
and its variability can be sought by examining the flow field dynamics – namely the pressure and velocity profiles at important sections. Earlier work on this was completed by Scimemi (Reference Scimemi1930) and Rouse (Reference Rouse1932), who measured the free surface profiles for several cases of
$h/p$
, pressure profiles at various points using manometers and velocity using Pitot tubes. These results have been critical for validating the output of numerical simulations and analytical approximations. This work was furthered by the results of Rajaratnam & Muralidhar (Reference Rajaratnam and Muralidhar1968, Reference Rajaratnam and Muralidhar1971) and Ramamurthy et al. (Reference Ramamurthy, Tim and Rao1987). However, the flow field dynamics for the full
$h/p$
range from the wall overflow to free overfall have yet to be fully characterised and rightly contextualised for the practical weir flow problem. This represents the main goal of this paper – a revised definition and understanding of
$C_{d}$
gained from new data that reconciles traditional insights from both hydrodynamics and hydraulics.
The need to rejoin these traditionally divergent fields of study was recently highlighted by Michael Eckert in his historical review of the more general efflux problem: ‘From an epistemic and historiographic perspective, this study calls for inquiries that combine the history of science with that of engineering. The history of fluid mechanics may not be approached from one or the other side alone, and specific cases such as the efflux problem provide suitable probes in this quest’ (Eckert, Reference Eckert2024). In this paper, we endeavour to do just that, through the following aims.
-
1. Revise the formulation of the weir-discharge equation to clarify the meaning of
$C_{d}$
. -
2. Provide additional data on the discharge characteristics and flow field dynamics of viscous flow over weirs and sills, especially towards the upper limit of the free overfall.
-
3. Elucidate the reason for the non-monotonic trend between
$h/p$
and
$C_{d}$
observed by Kandaswamy & Rouse (Reference Kandaswamy and Rouse1957) and Ramamurthy et al. (Reference Ramamurthy, Tim and Rao1987). -
4. Clarify the delineation between the respective weir and sill flows, along with a description of other regimes to inform reliable flow measurement.
2. Dimensional analysis for the weir flow problem
Since the development of (1.4), the field of hydraulics has advanced through dimensional analysis using the Buckingham
$\Pi$
theorem. This methodology has been applied to the weir flow problem in the past (Kindsvater & Carter, Reference Kindsvater and Carter1957; Ackers et al., Reference Ackers, White, Perkins and Harrison1978), with insights informing the physical parameters influencing the non-dimensional discharge. For completeness, we show the dimensional analysis for fully ventilated weir flow here. The relevant parameters for determining
$q$
are
\begin{equation} q = f(h,p,g,\rho ,\mu ,\sigma ), \end{equation}
where previously undefined parameters are the fluid properties of dynamic viscosity,
$\mu$
; and surface tension,
$\sigma$
. In the classical dimensional analysis of weir flow, as found from Kindsvater & Carter (Reference Kindsvater and Carter1957),
$h$
is chosen as the relevant dynamic flow parameter and length scale when predicting the discharge. This is because
$h$
is often readily measured from the hydrostatic pressure relationship, and for subcritical flows over a weir,
$U_{0} = f(h)$
(Vanden-Broeck & Keller, Reference Vanden-Broeck and Keller1987). The value of
$h$
is typically measured at a longitudinal distance of approximately 3–4 times
$h$
upstream of the weir crest (Bos, Reference Bos1976). The other important geometric parameter is
$p$
, representing the height of the obstruction. Physical fluid properties must be included to account for scale effects. The basic dimensions can be represented using appropriate scaling variables such that
$M \equiv \rho h^{3}$
,
$L \equiv h$
,
$T \equiv h/\sqrt {gh}$
and
$V \equiv \sqrt {gh}$
. This scaling will result in the following
$\Pi$
terms:
$\Pi _{1} = { Fr}_{h} = q/\left (\sqrt {g}h^{3/2}\right )$
,
$\Pi _{2} = h/p$
,
$\Pi _{3} = { Re}_{h} = \rho \left (\sqrt {g}h^{3/2}\right )/\mu$
and
$\Pi _{4} = { We}_{h} = \left (\rho g h^{2}\right )/\sigma$
.
Three common dimensionless numbers appear: the Froude number,
${ Fr}_{h}$
, the Reynolds number,
${ Re}_{h}$
, and the Weber number,
${ We}_{h}$
. These are given the subscript
$h$
to differentiate them from more common definitions for free-surface flow using other length and velocity scales. Assuming that fluid properties remain constant, it is seen that the only dynamic parameter in
$\Pi _{3}$
and
$\Pi _{4}$
is
$h$
. Since
${ We}_{h} \propto { Re}_{h}^{4/3}$
, only
${ Re}_{h}$
is necessary to account for scale effects at low values of
$h$
. If it is assumed these scale effects are negligible above a certain threshold value of
${ Re}_{h}$
, the simplified dimensionless equation for weir discharge becomes
\begin{equation} { Fr}_{h} = \phi (h/p)\quad \textrm {for}\:\: { Re}_{h} \gg 1. \end{equation}
Thus, a simplified form of the weir discharge equation is
\begin{equation} q = \phi (h/p)\sqrt {g}h^{3/2}\:, \end{equation}
where the function
$\phi (h/p)$
approximates
${ Fr}_{h}$
within a normal operating range of
$h/p$
in which scale effects are negligible. Previous studies have indicated that a suitable value for this threshold is
${ Re}_{h}=3.5\times 10^{4}$
(Kindsvater & Carter, Reference Kindsvater and Carter1957; Pugh et al., Reference Pugh, Venayagamoorthy, Gates, Berni and Rastello2024), equivalent to a minimum
$h$
of
$\approx$
5 cm to drive the water flow. By comparing (
1.4) and (
2.3), it can be seen that
${ Fr}_{h}$
is simply
$(2\sqrt {2}/3) C_{d}$
. We choose here to refer to
${ Fr}_{h}$
as a weir Froude number, intended to serve as an alternative to the classical weir discharge coefficient
$C_{d}$
. The result of this approach is that a physical understanding of the dimensionless coefficient on the front-end of the weir discharge equation is obtained, in that it is best understood as a ratio of inertial to gravitational effects. To our knowledge, the only other work to propose an understanding of the weir discharge coefficient as a type of Froude number is that of Street
et al.
(Reference Street, Watters and Vennard1996). With a fundamental definition of
$C_{d}$
as a ratio of inertial to gravitational effects now established, we will proceed to its derivation from (1.3) to reveal its constituent parts.
3. Revisiting the discharge coefficient
As discussed earlier, Poleni first derived his equation for weir discharge by envisioning the flow over a weir as an integrated sum of Torricelli’s problem of flow through a sharp-edged orifice at the bottom of a reservoir, where the jet freely enters into the atmosphere and the flow is driven by the static pressure head (
$h^{\prime }$
) above the jet opening of width
$d$
. The unit volumetric discharge of the jet,
$q^{\prime }$
, is given by the general orifice equation (Streeter, Reference Streeter1985)
\begin{equation} q^{\prime } = {\overline {u}}^{\prime }t^{\prime } = ({C_{v}}^{\prime }\sqrt {2gh^{\prime }})({C_{c}}^{\prime }d)\:, \end{equation}
where
${\overline {u}}^{\prime }$
is the average velocity in the jet over thickness
$t^{\prime }$
. The term
$\sqrt {2gh^{\prime }}$
is the theoretical velocity given by (
1.3) from assumptions of steady, inviscid and irrotational flow. Then,
${C_{v}}^{\prime }$
is a correction coefficient for any local energy losses that occur due to separation as the flow navigates the orifice opening. As the jet exits the orifice and the streamlines become horizontal, the jet diameter contracts. Thus, the jet thickness in the vena contracta,
$t^{\prime }$
, divided by the height of the orifice opening,
$d$
, is a contraction coefficient
${C_{c}}^{\prime }$
.
Returning now to the derivation of (
1.4), it can be seen that if the Bernoulli equation (
1.3) is applied between the free surface at the gauging section and the crest section (see figure 1), the velocity (
$u(\eta )$
) within the overflowing jet of a fluid layer at a distance
$\eta$
below the free surface can be estimated as
\begin{equation} h+\frac {{U_{0}}^{2}}{2g} = h - \eta + \frac {{u(\eta )}^{2}}{2g} \Rightarrow u(\eta ) = \sqrt {2g}\left ( \eta + \frac {{U_{0}}^{2}}{2g} \right )^{1/2}\:. \end{equation}
Here,
$q$
is estimated by integration of (
3.2) from
$\eta =0$
to
$\eta =h$
:
\begin{equation} q \approx \sqrt {2g}\int _{0}^{h} \left (\eta + \frac {{U_{0}}^{2}}{2g}\right )^{1/2}{d}\eta = \frac {2}{3}\sqrt {2g} \left [ \left (h+\frac {{U_{0}}^{2}}{2g}\right )^{3/2} - \left (\frac {{U_{0}}^{2}}{2g}\right )^{3/2}\right ]. \end{equation}
The theoretical average velocity in the jet,
${\overline {u}}_{t}$
, predicted by (
3.3) is then
\begin{equation} {\overline {u}}_{t} = \frac {1}{h}\frac {2}{3}\sqrt {2g} \left [ \left (h+\frac {{U_{0}}^{2}}{2g}\right )^{3/2} - \left (\frac {{U_{0}}^{2}}{2g}\right )^{3/2}\right ], \end{equation}
and the orifice flow equation ( 3.1) can be written for the weir flow problem as
\begin{equation} q = C_{v}{\overline {u}}_{t}(C_{c}h), \end{equation}
where
\begin{equation} C_{v} = \frac {\overline {u}}{{\overline {u}}_{t}}\quad \text{and}\quad C_{c} = \frac {t}{h}. \end{equation}
Here,
${\overline {u}}_{t}$
is given by (
3.4). As mentioned earlier,
$\overline {u} = q/t$
, and the correction coefficients in (
3.5) are commonly lumped into a single coefficient so that the need for a priori knowledge of
$U_{0}$
in (
3.4) is eliminated. Thus, the relation between (
3.5) and the commonly used Poleni equation (1.4) is given by
\begin{equation} q = \left (C_{v}C_{c}\right ) \frac {2}{3}\sqrt {2g} \left [ \left (h+\frac {{U_{0}}^{2}}{2g}\right )^{3/2} - \left (\frac {{U_{0}}^{2}}{2g}\right )^{3/2}\right ] = C_{d}\frac {2}{3}\sqrt {2g}h^{3/2}\:. \end{equation}
Examining (
3.7), it can be seen that
$C_{d}$
is then the product of three terms:
\begin{equation} C_{d} = C_{v}KC_{c}, \end{equation}
where
$K$
is a correction factor incorporated to account for the missing gauging section velocity head (
${U_{0}}^{2}/2g$
) in (
1.4) (Rouse, Reference Rouse1946):
\begin{equation} K = \frac {1}{h^{3/2}} \left [\left (h + \frac {{U_{0}}^{2}}{2g}\right )^{3/2} -\left (\frac {{U_{0}}^{2}}{2g}\right )^{3/2}\right ] = \left [\left (1 + \frac {{U_{0}}^{2}}{2gh}\right )^{3/2} -\left (\frac {{U_{0}}^{2}}{2gh}\right )^{3/2}\right ]\:. \end{equation}
Thus, for normal operating regimes where scale effects may be neglected, the discharge coefficient,
$C_{d}$
, can be understood as a correction factor that accounts for viscous friction losses, upstream velocity head and contraction. Returning to the dimensional analysis in § 2, the simpler characterisation of dimensionless discharge,
${ Fr}_{h}$
, which retains the intuitive definition of being a ratio between inertial and gravitational effects, can also be easily decomposed as
\begin{equation} { Fr}_{h} = \left (\frac {2\sqrt {2}}{3}\right ) C_{v}KC_{c}\:. \end{equation}
We will now examine the contribution of the three component terms on the right-hand side of (
3.10) towards the overall behaviour of
${ Fr}_{h}$
as
$h/p$
varies in pursuit of the proper form of the function
$\phi (h/p)$
in (
2.3).
4. Research methodology
This analysis synthesises previously published data on sharp-crested weirs, along with new experimental data collected by the authors. A summary of these data is provided in table 1. In general, the experimental data sets of table 1 were completed within laboratory flumes where the weir crest was placed on the downstream end of the flume to ensure full ventilation of the nappe. Flow depth measurements were typically made using a point gauge, hook gauge or similar device. Discharge measurements were made using a previously calibrated weir, a flow meter in the inlet pipe to the flume or a weighing tank mechanism. New experimental and numerical results of the authors include supplemental results of the steady two-dimensional flow field, summarised in table 2, produced using particle image velocimetry (PIV) for one case where
$h/p=0.68$
(E0.68), along with eleven numerical simulations using FLOW-3D (Flow Science, Inc., 2023).
Table 1. Description of compiled historical discharge data for thin weir and sill flows, along with that of the current study.
$n$
is the number of observations in the dataset. The total number of historical observations is 436
Table 2. Details for numerical and experimental flow cases of the current study where the two-dimensional velocity and pressure field were analysed. The numeric value in the case name corresponds to the
$h/p$
value for that case
4.1 Experimental measurements
Initially, observations of head (
$h$
) and steady discharge (
$Q$
) were made over a physical model of a sharp-crested weir at the Environmental Fluid Mechanics Laboratory (EFML) at Colorado State University (see figure 2a). The height of the model was
$p=150$
mm. The length of the flume was 5 m in total, with a width of 0.3 m and a horizontal slope. The weir was placed 3 m downstream of the inlet and PVC tubes were installed beneath the weir crest on the downstream side to provide proper ventilation of the nappe per specifications of Bos (Reference Bos1976). Flow straighteners and a boundary layer tripping tool at the flume inlet were also installed to achieve as-close-as-possible a fully developed turbulent flow profile. Flow depths were measured using a point-gauge placed 1 m upstream of the model and steady discharges were measured using an electromagnetic flow meter installed in the inlet pipe of the flume.
Additional head-discharge measurements were completed at the Hydraulics and Hydromorphology Laboratory (HHLab) at INRAE Lyon, France by the first author. This experimental set-up featured an adjustable-crest weir at the end of an 18 m flume with a width of 1 m. A sharp-edged crest piece was affixed to the top of the weir to ensure springing flow (figure 2b). The range of crest heights examined was 50–150 mm. Flow depth measurements were made using ultrasonic water-level sensors and steady discharge was measured using an electromagnetic flow meter.
Two-dimensional planar PIV observations also were perfdormed at the EFML using a dual head Nd:YLF 527 nm green laser (Photonics model no. DM-527-30), and a high-speed CMOS camera with a 2560 × 1600 pixel resolution, 12-bit depth and maximum capture rate of 800 frames per second (Phantom model no. VEO-E 340L). A LaVision programmable timing unit (PTU 10) was used as a synchroniser, and the LaVision Davis 10 software was used for PIV velocity field calculations on a 16-core Supermicro machine with an Intel Xeon E5-2620 v4 processor and 64 GB of RAM. Pressure-from-PIV calculations also were performed using an in-house code that solved the Poisson equation implicitly, although these results were limited in accuracy due to the difficulty of correctly defining the boundary conditions.
Figure 2. (a) Experimental facility at the EFML. (b) Weir model used at the HHLab facility. (c) Numerical domain for flow case N0.68. The upstream boundary of the domain in the
$x$
-direction has been truncated from −10 to −0.5 m to show details of the near-crest region. Red dots represent pressure probes at the gauging section and just below the weir crest. Text boxes with arrows indicate boundary conditions.
4.2 Numerical simulations
Numerical simulations of fluid flow were completed using FLOW-3D version 22.2 on a 32-core machine with an Intel Xeon E5-2680 v3 processor and 128 GB of RAM. A Reynolds-averaged Navier–Stokes (RANS) modelling approach with the
$k{-}\omega$
turbulence model (Wilcox, Reference Wilcox1988) was chosen due to the desire to model viscous flow phenomena subjected to an adverse pressure gradient and its cheaper computational cost compared with large eddy simulation (LES). Flow development was achieved by allowing the simulation to run for at least 45 s of flow time, along a distance of 10 m in the approach to the crest. The width of the numerical domain was 0.6 m, with a rigid wall of width 0.1 m placed on one side of the domain to allow for a free-falling nappe (see figure 2c). This resulted in an effective approach channel width of 1 m, due to a symmetry boundary condition being placed on the mesh face opposite the solid side wall. The longitudinal boundary conditions were an upstream stagnation pressure given by a fluid elevation that represented approximately 10.1 cm of head above the crest elevation and a pressure outflow boundary condition with a fluid elevation of 0.15 m. For all simulations, the orthogonal grid size was 5 mm in the longitudinal (
$x$
) and vertical (
$z$
) directions, and 10 mm in the transversal (
$y$
) direction. An additional mesh with 1 mm resolution in the
$x$
- and
$z$
-directions was included to enclose the weir geometry so that its sharp interfaces were well resolved. Three-dimensional simulations were completed to ensure the numerical channel width was wide enough to neglect side-wall effects and to allow for sufficient ventilation of the falling nappe. Analysis was completed on a two-dimensional (
$x{-}z$
) planar slice taken from the boundary opposite the channel side-wall. A grid sensitivity analysis for the resolution of the
$x{-}z$
plane was completed for grid sizes of 20, 10, 5 and 2 mm; with sufficient convergence being seen at the 5 mm resolution.
The gauging section was located so that the magnitude of the longitudinal distance,
$x$
, from the crest was equal to
$4(D(x)-p)$
, in accordance with measurement guidelines (Bos, Reference Bos1976; Bureau of Reclamation, 2001). Measurements of
$D$
, and thus
$h$
, were made by evaluating the hydrostatic pressure head at the gauging section.
Numerical results were validated against the experimental results of Scimemi (Reference Scimemi1930), Rouse (Reference Rouse1932), Rajaratnam & Muralidhar (Reference Rajaratnam and Muralidhar1968) and case E0.68. Qualitative comparisons of the velocity field and nappe shapes between E0.68 and N0.68 showed good agreement (see figure 3a,b). The values of
${ Fr}_{h}$
, along with the pressure (figure 3c) and velocity (figure 3d) profiles in the critical section produced from numerical simulations were found to be in good agreement with the available experimental data for similar cases of
$h/p$
. The normalised root mean squared error (NRMSE) values in the velocity profiles shown in figure 3c in comparison to the experimental profiles were 9.57 % for case N0.27, and 8.17 % for case NFO. In comparing the pressure profiles in figure 3d, the NRMSE values were found to be 11.47 % for case N0.27, and 9.08 % for case NFO. It should also be noted that some deviation from the experimental velocity profile was seen in the numerical results near the lower nappe boundary. The location of maximum velocity for the numerical results tended to be slightly above the lower nappe boundary, whereas experimental and theoretical results have suggested that the maximum velocity in the profile should occur at this lower boundary. This discrepancy is likely due to the difficulty in reproducing numerically, and measuring experimentally, the shear layer that forms along the air–water interface of the lower nappe. However, these differences were confined to a thin region of the overall nappe. Overall, our numerical simulations did not produce significant deviations in the trend of results compared with the available experimental data.
Figure 3. Comparison of dimensionless velocity magnitude field (
$V = \sqrt {u^{2}+w^{2}}$
) from (a) PIV experiments (case E0.68) and (b) numerical simulation (case N0.68). (c) Comparison of streamwise velocity profiles in the critical section. (d) Comparison of pressure profiles in the critical section between experimental of previous authors and numerical results of the current study. For the data of Scimemi (Reference Scimemi1930),
$h/p$
= 0.25, and the data of Rouse (Reference Rouse1932) and Rajaratnam & Muralidhar (Reference Rajaratnam and Muralidhar1968) are both for the free overfall case.
5. Results
5.1 Revisiting the weir-discharge equation
Figure 4. (a)
${ Fr}_{h}$
versus
$h/p$
for the entire weir–sill range, with the historical data (table 1) plotted. Classical empirical equations relating
$h/p$
and
${ Fr}_{h}$
are also shown by a solid line, after Kandaswamy and Rouse (Reference Kandaswamy and Rouse1957), with dash-dotted lines representing
$\pm$
5 % variation. For all data,
${ Re}_{h}\gt 3.5 \times 10^{4}$
. For cases shown in table 2: (b)
$C_{v}$
versus
$h/p$
; (c)
$K$
versus
$h/p$
; and (d)
$C_{c}$
versus
$h/p$
. Approximate values of
$C_{v}$
,
$K$
and
$C_{c}$
for the upper and lower limits of
$h/p$
, representing the respective cases of the free overfall and wall overflow, are given by blue horizontal lines.
With a revised presentation of
$C_{d}$
as the weir Froude number,
${ Fr}_{h}$
, and the product of three component terms (
${ Fr}_{h} = (2\sqrt {2}/3)C_{v}KC_{c}$
) we now turn to the form of the function
$\phi (h/p)$
in (
2.3). Along with the current results from numerical simulation, we combine new experimental data with historical data found in the literature to examine the discharge-rating equation for thin weirs and sills. Figure 4 shows, quite remarkably, that a simple linear form of the weir-discharge equation provides an excellent fit to the data within the weir regime. This fit takes the form of the classical Rehbock equation without consideration of scale effects (Kandaswamy & Rouse, Reference Kandaswamy and Rouse1957),
$C_{d} = 0.611 + 0.075h/p$
, multiplied by
$(2\sqrt {2}/3)$
so that
\begin{equation} { Fr}_{h} = \phi (h/p) \simeq 0.576 + 0.071h/p\:, \end{equation}
which is applicable in the range
${ Re}_{h}\gt 3.5\times 10^{4}$
. A non-monotonic trend in
${ Fr}_{h}$
with
$h/p$
can be observed in figure 4, with a maximum value of
${ Fr}_{h}$
occurring at
$h/p =\mathcal{O}(10^{1})$
. The flow cases beyond this peak value and approaching the limit of the free overfall are considered to be in the sill regime. Kandaswamy and Rouse (Reference Kandaswamy and Rouse1957) suggested a fit to the data in the sill regime, where
${ Fr}_{h} \propto (h/p)^{-3/2}$
. However, uncertainty in measurements is much greater in this regime due to variability in the water surface caused by friction losses and surface waves.
The reason for the non-monotonicity in
${ Fr}_{h}$
with
$h/p$
is observed in figure 4d to be due to the behaviour of
$C_{c}$
, as was hypothesised earlier by Kandaswamy and Rouse (Reference Kandaswamy and Rouse1957). It is also seen that the coefficient
$K$
has the largest relative effect on the magnitude of
${ Fr}_{h}$
, and
$C_{v}$
the smallest. For planar orifice flows,
$C_{v}$
was found to range between 0.95 and 0.99, and accounts for local losses that occur in the shear layer formed near the sharp crest (Streeter, Reference Streeter1985). Thus, the value of
$C_{v}$
in the wall overflow limit shown in figure 4b is taken to be 0.95. Additionally,
$C_{v}$
in the free overfall limit is taken to be 1.0, which is akin to stating that local losses are negligible in this case. The trend in values shown in figure 4b and table 2 supports this. Rouse (Reference Rouse1932) stated that the behaviour of
$K$
in the wall overflow and free overfall limits should be 1 and 1.485, respectively, which can be shown to be true when it is recognised that the free overfall case is the critical flow condition where the channel Froude number,
${ Fr} = q/\left (\sqrt {g}D^{3/2}\right )$
, is unity.
From an analysis of orifice flows, it can be seen that for viscous fluids where
$h/p \ll 1$
, the value of
$C_{c}$
in the wall overflow limit should approach 0.64, so that
$C_{d}=0.611$
with
$C_{v}=0.95$
and
$K=1$
. For the free overfall, Rouse (Reference Rouse1936) showed that the value of this coefficient should be approximately 0.715, which again is derived from the assumption of critical flow conditions. Overall, the behaviour of these three constituent coefficients from the results of the current study shown in figure 4 agrees well with the current understanding in the literature.
5.2 Examining flow dynamics
To elucidate the data trends observed in figure 4, pressure and velocity profiles within the critical section are shown in figure 5 for flow cases given in table 2. It can be seen in figure 5a that a monotonic trend in streamwise inertia is observed with increasing
$h/p$
, which is the reason for the trend in
$K$
observed in figure 4c. The integral of the pressure profiles in figure 5b exhibits a decreasing trend in response to the increasing velocity. It should also be noted that the velocity profiles become more uniform as
$h/p$
increases, and that an apparent self-similarity is present in both the velocity and pressure profiles for
$h/p \leq \sim 2$
. This fact has been recognised by earlier authors (Montes, Reference Montes1998; Castro-Orgaz & Hager, Reference Castro-Orgaz and Hager2017), who have shown that an analytical solution from potential flow theory is available for the coupled velocity and pressure profiles in the critical section for
$h/p$
values in this range. However, further experimental investigation is required to validate these preliminary findings.
Figure 5. (a) Streamwise velocity profiles in the critical section. (b) Pressure profiles in the critical section.
$E$
is the minimum specific energy head (see figure 1).
Figure 6. Dimensionless contour plots of the vertical pressure gradient for data summarised in table 2.
To help explain the trends in
$C_{c}$
observed in figure 4d, figure 6 plots the contours of the dimensionless vertical pressure gradient. Upon examination of (
1.1), it can be seen for the convergent, steady and two-dimensional flow cases analysed here, the greatest influence on the vertical pressure gradient (
$\partial P/ \partial z$
) will come from the vertical inertia (
$u\frac {\partial w}{\partial x} + w \frac {\partial w}{\partial z}$
) terms. In the absence of vertical inertia, a hydrostatic pressure condition is present. When the sum of the vertical inertia terms is equal to the gravitational acceleration, the nappe falls as a zero-pressure-gradient (ZPG) jet. The transition from essentially 1-D hydrostatic open-channel flow to a free-falling jet can be seen in figure 6. Upstream of the weir and away from any boundaries, the velocity distribution is uniform and horizontal streamlines result in a hydrostatic pressure distribution. If the obstacle is large enough, flow separation occurs due to an adverse pressure gradient, creating a recirculation zone on the upstream corner formed by the channel bed and the weir. The boundary layer eventually reattaches along the vertical weir face and then the fluid accelerates upwards toward the crest with significant vertical inertia. This, along with approach stagnation at the weir face, results in a zone of reinforced hydrostatic pressure just below the crest (
$-(\partial P/\partial z)/\gamma \gt 1$
). As the flow springs from the crest, gravitational effects work against the upward vertical inertia along the lower nappe surface, causing it to decelerate and create an adverse pressure gradient zone (
$-(\partial P/\partial z)/\gamma \lt 0$
). On the upper nappe surface, draw-down effects caused by the weight of the falling jet reduce the hydrostatic pressure gradient (
$-(0\lt \partial P/\partial z)/\gamma \lt 1$
). These two opposing gradients on the upper and lower surfaces of the nappe work against the remaining pressure until the flow becomes a free-falling jet with atmospheric pressure throughout its whole thickness.
Observing the trend in
$\varepsilon /h$
with
$h/p$
shown in Table 2 reveals that at small values of
$h/p$
, contraction effects along the lower nappe surface are significant due to the relatively large amount of upward vertical inertia along the weir face that must be navigated around the sharp crest. This results in local losses, causing
$C_{v}\lt 1$
. As
$h/p$
increases towards the free overfall, the contraction is initially reduced as streamwise inertia increases and the streamlines in the approach flow become more horizontal. However, increased draw-down in the free surface (
$\delta /h$
) occurs in the sill regime due to the insufficiency of the sill to cause any appreciable flow separation and stagnation near the crest, which in the weir regime, allowed for the retainment of supportive pressure below the nappe. These combined effects explain the observed non-monotonic trend in
${ Fr}_{h}$
with
$h/p$
.
5.3 The weir and sill regimes
With an understanding now established of the flow dynamics occurring within the wall overflow and free overfall limits, we turn to considerations for flow measurement. It should be noted that entrance into the sill regime comes with increased practical complications. As
$h/p$
grows, the relative influence of the turbulent boundary layer on the flow field also increases (Rouse, Reference Rouse1964), and increased mixing caused by turbulence results in a greater energy dissipation rate and friction slope. Eventually, the friction slope becomes appreciable enough that the free surface can no longer be assumed to asymptotically approach a horizontal line. Rouse (Reference Rouse1936) recognised this and noted that for the free overfall case, the critical depth,
$D_{c}$
(which for the free overfall is equivalent to the gauging depth), could be consistently located at an upstream distance from the crest equal to
$4D_{c}$
. However, it should be noted that this guideline must be met through successive iteration because
$q$
is not known a priori, and thus neither is
$D_{c}$
.
Figure 7 shows how the friction slope increases with
$h/p$
. Above an
$h/p$
value of
${\sim}5$
, small deviations in the location of the gauging section measurement can lead to variations in the estimation of
${ Fr}_{h}$
greater than 2 %. This is consistent with the findings of other authors (Kindsvater & Carter, Reference Kindsvater and Carter1957; Sinclair et al., Reference Sinclair, Venayagamoorthy and Gates2022). Additionally, as
${ Fr}_{h}$
approaches unity, surface waves become more significant, further complicating a reliable measurement of
$h$
. It can also be seen that for values of
$h/p\lt {\sim} 8$
, the bed pressure in the approach to the crest increases or at least is held constant, whereas for flow cases above this threshold, the bed pressure drastically begins to decrease towards the crest. This is evidence of flow separation and stagnation being present for
$h/p\lt \sim 8$
, and absent otherwise.
Figure 7. Pressure head on the channel bottom made dimensionless by the flow depth at the gauging section (located at a distance of
$4h$
upstream of the crest), indicating changes in the flow behaviour in the approach to the weir crest with varying
$h/p$
. Dash-dotted lines represent variation in the dimensionless flow depth that indicate a
$\pm 2\,\%$
variation in the estimation of
${ Fr}_{h}$
.
The previous studies of Kandaswamy & Rouse (Reference Kandaswamy and Rouse1957) and Ramamurthy et al. (Reference Ramamurthy, Tim and Rao1987) implicitly suggest that the practical application range of weir flows with regards to
$h/p$
is unlimited, so long as a sufficiently large drop in the bed elevation is present to allow for full ventilation of the overflowing nappe. Ramamurthy et al. (Reference Ramamurthy, Tim and Rao1987) identified two regimes for these types of flows: for
$h/p \leq 10$
, the flow was considered weir flow; and for
$h/p \gt 10$
, the flow was considered sill flow. The maximum discharge capacity was found to be at the transition point of
$h/p=10$
. However, defining this threshold based solely on
$h/p$
lacks a physical justification. To investigate this transition, we plot in figure 8 the data summarised in table 1, showing both
${ Fr}_{h}$
and
$h/p$
against the channel Froude number,
${ Fr}$
.
Figure 8. Dependency of (a)
${ Fr}_{h}$
and (b)
$h/p$
with the channel Froude number,
${ Fr}$
. The black trend line indicates a logistic relationship between the two variables. Respective limits for the self-similar, ideal operating, transition and sill regimes are given.
6. Limitations and practical considerations
Figure 8 reveals that
${ Fr}_{h}$
and
$h/p$
exhibit a logistic relationship with the channel Froude number,
${ Fr} = q/ (\sqrt {g}D^{3/2})$
, until a distinct break point at approximately
${ Fr}=0.95$
. From examining figure 8b, this value is equivalent to
$h/p\simeq 8$
. These trends reveal that the weir flow regime exists for up to 95 % of the critical
${ Fr}$
, after which it enters into an unstable regime characteristic of sill flow. However, figure 7 reveals practical limitations of accurately measuring
$h$
above an
$h/p$
value of 5. We therefore suggest that the ideal operating regime for weir flows be considered as
${ Fr}\lt 0.75;\; h/p\lt 5$
, with a transition regime between high-inertia weir flows and sill flows existing in the range
$0.75\leq { Fr} \lt 0.95$
. A practical upper limit of
$h/p=5$
has been suggested by other authors (Kandaswamy & Rouse, Reference Kandaswamy and Rouse1957). However, in pursuit of even further precision, one could also suggest that weir flows used in practice should remain below the threshold of
$h/p=2$
, equivalent to
${ Fr} \simeq 0.4$
, so that the velocity and pressure profiles within the nappe critical section remain self-similar – as was observed in figure 5. This practical upper limit of
$h/p=2$
for accurate flow measurement is one that can be found in several previous reliable studies (Kindsvater & Carter, Reference Kindsvater and Carter1957; Castex, Reference Castex1969; Ackers et al., Reference Ackers, White, Perkins and Harrison1978). Furthermore, other authors have suggested an upper limit of
${ Fr}=0.5$
for flow measurement using long-throated flumes to avoid wave action on the free surface that hinders an accurate measurement of
$h$
(Bureau of Reclamation, 2001; Bos et al., Reference Bos, Replogle and Clemmens1984). Finally, it should be noted that the sill regime (
${ Fr}\gt 0.95$
) does not represent a realistic condition for flow measurement due to the breakdown in the monotonic trend between
$h/p$
and
${ Fr}_{h}$
that renders the weir-discharge equation unusable, along with the difficulty of defining the gauging depth due to a considerable friction slope and free-surface fluctuations.
7. Conclusion
Dimensional analysis has been well-leveraged in the past to develop functional and intuitive engineering solutions to classical hydraulics problems, such as for pressure losses in closed conduit flow. The weir discharge equation has thus far relied upon a derivation taken from its treatment as an orifice problem, but without much consideration for the relevant physical mechanisms at play. This has lead to ambiguity in a physical explanation of the weir discharge coefficient,
$C_{d}$
. We show how dimensional analysis using Buckingham
$\Pi$
theorem provides fundamental insights on the traditional weir-discharge equation, namely that
$C_{d}$
is most fundamentally a ratio of inertial to gravitational effects, represented by the dimensionless weir Froude number,
${ Fr}_{h}$
, equal to
$(2\sqrt {2}/3) C_{d}$
. Its behaviour is governed by the combined effects of local friction losses, upstream flow inertia and contraction in the overflowing nappe. The contribution of flow inertia to the overall discharge characteristics was found to be the most significant factor, with the non-monotonic trend in the contraction coefficient helping to delineate a transition between weir and sill flows. We investigate a physical explanation for the transition from weir to sill flow, finding that the growth rate of
${ Fr}_{h}$
with
$h/p$
is linear over a large range of the subcritical flow regime, until the channel Froude number,
${ Fr} = q/\left (\sqrt {g}D^{3/2}\right )$
, is approximately equal to 0.95. Above this threshold, the flow becomes unstable as it reaches criticality, with the stagnation pressure at the base of the weir diminishing so that increased free-surface draw-down of the overflowing nappe causes a marked increase in contraction and a reduction in discharge capacity. A practical upper limit of the ideal operating regime for weir flows is defined as
${ Fr}=0.75$
due to the increased effect of boundary friction and free-surface slope above this threshold, while
${ Fr}\leq 0.4$
is found to be a threshold for self-similarity in the pressure and velocity profiles.
The results of this study can help inform practicing hydraulic engineers working with weirs on (1) a physical basis for the discharge capacity of a structure, and (2) what effects are primary and secondary in determining
$C_{d}$
. A clear delineation between weir and sill flows helps ensure that flow measurement practises are completed in the most reliable regime. Together, these insights will help in providing engineers with practical guidance on the design, calibration and operation of overflow structures to achieve accurate and dependable discharge measurement and regulation.
Acknowledgements
The authors gratefully acknowledge the help of Céline Berni of the Hydraulics and Hydromorphology Lab at INRAE-Lyon, and Marie Rastello of LEGI-Grenoble for their help in conducting experiments related to this work. The authors also thank the anonymous reviewers for their constructive feedback and insights on the initial draft of this manuscript.
Data availability statement
Data summarised in tables 1 and 2 are currently available at https://datadryad.org/dataset/ doi:10.5061/dryad.pg4f4qs2n.
Funding Statement
This work was supported by the Colorado Agricultural Experiment Station under Grant No. COL00424.
Competing interests
The authors declare no conflict of interest.
Open access
