When the capacity of the spillway of a dam is exceeded for a given flood, overtopping occurs; in such cases potentially dangerous hydrodynamic actions and scour downstream of the dam need to be foreseen. Detailed studies of jets impinging in plunge pools from overflow nappe flows are scarce. This work addresses plunge pool flows, and compares numerical results against our own experiments. The energy dissipation is larger than 75% of the impingement jet energy. Instantaneous velocities and air entrainment were obtained with the use of an Acoustic Doppler Velocimeter and optical fibre probe, respectively. Mean velocity field and turbulence kinetic energy profiles were determined. To identify the level of reliability of models, numerical simulations were carried out by using the ‘homogeneous’ model of ANSYS CFX, together with different turbulence closures. The numerical results fall fairly close to the values measured in the laboratory, and with expressions for submerged hydraulic jumps and horizontal wall jets. The observations can be well predicted for characterized profiles at a minimum distance of 0.40 m downstream from the stagnation point, horizontal velocities greater than 40% of the maximum velocity in each profile, and when the ratio of the water cushion depth to the jet thickness is lower than 20.

INTRODUCTION

There is growing consensus regarding the fact that climate change will lead to enhanced extreme flooding in certain areas of the world. This situation must be confronted with spillway capacity and special operational scenarios for large dams. If the capacity of the spillway is insufficient, the dam might be overtopped, thus generating new loading scenarios, and raising questions about potential risky hydrodynamic actions and scour downstream of the dam (Wahl et al. 2008).

Spain is the fourth country in the world according to the number of large dams – it has over 1,200. Fifty percent of these dams were built fifty years ago. In this sense, numerous dams need to be re-evaluated in their safety with respect to potential overflow, in line with what is being done in the USA (FEMA 2014).

When the rectangular jet or nappe flow occurs due to overtopping, the design considerations need to ensure that most energy is dissipated, and that there is minimal to no erosion downstream of the dam. In other words, we need to estimate the hydrodynamic actions on the bottom of the basin where the jet discharges, as a function of the characteristics of the jet (Annandale 2006).

The energy dissipation mechanisms that occur in the jet-basin structure can be grouped into the following: (a) aeration and disintegration of the jet in its fall, (b) air entrainment and diffusion of the jet into the basin, (c) impact on the basin bottom, and (d) recirculation in the basin (Figure 1).
Figure 1

Schematic of falling rectangular jets and receiving basin.

Figure 1

Schematic of falling rectangular jets and receiving basin.

Two of the variables needed to be defined in the design of the jets are the issuance conditions and the impingement conditions (Castillo 2006, 2007). The issuance conditions correspond to the flow conditions at a location where the jet leaves the spillway and starts falling freely (z=–h, where z is the vertical coordinate with origin in the crest weir, and h is the weir head). The impingement conditions correspond to the jet section before the impact with the water surface of the basin. In this location, the mean velocity, Vj, and the impingement jet thickness, Bj, must be defined. This jet thickness must include the basic thickness due to gravity Bg, and the symmetric jet lateral spreading due to turbulence and aeration effects, ξ (Castillo et al. 2015): 
formula
1
where q is the specific flow, H the fall height, and h is the energy head at the crest weir. φ = KφTu, with Tu being the turbulence intensity and Kφ an experimental parameter (1.14 for circular jets and 1.24 for the three-dimensional nappe flow case).

When the jet falls through a long-enough distance, the jet becomes fully developed (see Lb in Figure 1). Castillo et al. (2015) established different equations to calculate the jet energy dissipation in the air and in the water cushion, as a function of the Y/Bj and H/Lb ratios (where Y and H denote the depth of the water cushion at the exit and the total head, respectively, and Lb is the break-up length). During the falling, the energy dissipation is due to the air entrainment into the falling jet and the depth of water upstream of the jet. Energy dissipation in the basin by diffusion effects can only be produced if there is an effective water cushion (Y/Bj> 5.5 for the rectangular jet case (Castillo et al. 2015)).

In Figure 2, the velocity Vj and the jet thickness Bj at the impingement conditions, the core depth or minimum depth for effective water cushion and the two principal eddies that produce the dominant frequencies in the plunge pool (large scale eddies and shear layer structures) are sketched. The lowest frequencies correspond to large scale eddies that have a dimension on the order of the plunge pool depth (see Bombardelli & Gioia 2005, 2006; Gioia & Bombardelli 2005). Then, the recirculating velocity for large plunge pools is about Vr ∼ 0.035 Vj and the corresponding Strouhal number of the dominant eddies is S=fY/Vj = (VrY) (Y/Vj) ∼ 0.01 (Ervine & Falvey 1987; Falvey 1990; Ervine et al. 1997). The following dominant frequency corresponds to eddy sizes contained in half of the shear-layer width and is proportional to the entry jet velocity; then, the Strouhal number of the shear-layer eddies is equal to a constant Ss= (fsY/Vj) = K3 ∼ 0.25, and it coincides with the spread of the jet into the water cushion as shown in Figure 2 (Ervine & Falvey 1987; Ervine et al. 1997).
Figure 2

Schematic of eddy structures in effective water cushion (Y 5.5 Bj): large scale eddies size ∼Y and shear layer structures size ∼De (adapted from Ervine et al. 1997).

Figure 2

Schematic of eddy structures in effective water cushion (Y 5.5 Bj): large scale eddies size ∼Y and shear layer structures size ∼De (adapted from Ervine et al. 1997).

Within the plunge pool downstream of the impingement point, the flow resembles a flow in a submerged hydraulic jump and a wall jet. However, the situation is complicated here by the air entrainment. Several formulas have been put forward to express the horizontal velocity distribution in the vertical direction. We revisit some of these formulas later in the paper.

There are only a few well-documented references on the numerical simulations of free overflow spillways. Ho & Riddette (2010) analysed the different applications of computational fluid dynamics (CFD) code to hydraulic structures, and identified some limitations that do not allow a comprehensive analysis of the two-phase phenomena. They suggested the following principal lines of future research:

  • Validation of the air entrainment along chutes and free-falling jets.

  • Extension of models to two-fluid theories to include air entrainment and jet break-up. Validate impact pressures against experimental data.

  • Develop parametric studies with turbulence models to identify the level of reliability of the computed pressure.

Castillo et al. (2014, 2015) and Carrillo (2014) have followed the above suggestions and, specifically, they have developed laboratory research on the velocity, air, and pressure fields in the jet and the basin. They undertook experimental observations of jet velocities, and air concentrations with an optical fibre probe, and of pressures in the plunge pool with pressure transducers. Then, they applied ANSYS CFX to simulate overflow nappe impingement jets in a general way, and investigated the different turbulence closures which better represent the data. In those works, the emphasis was put on the turbulence of the falling jet, the pressure distribution near the stagnation point in the water cushion, and the horizontal distance to the stagnation point. A good agreement among numerical simulations and laboratory data was obtained.

In Castillo et al. (2014), the so-called ‘homogeneous’ theoretical model of CFX was employed. It was shown that this model is able to reproduce correctly the jet water velocity, and the averaged pressures in the plunge pool. There is always a challenge in modelling two-phase flows to discern which level of complexity is needed to represent different aspects of the flow (Bombardelli 2004; Bombardelli & Jha 2009; Jha & Bombardelli 2010). One of the objectives of this paper is to determine whether this theoretical model is sophisticated enough to represent velocities in the plunge pool.

Continuing the line of research, this work presents a systematic study which considers specific flows and water cushions in the plunge pool. New laboratory data were obtained and new three-dimensional simulations were specifically performed for this work. ANSYS CFX was again selected due to the variety of turbulence closures available in the code, the previous experience with it and, more importantly, due to the diverse two-phase flow models embedded in the package, which can allow us to expand the research further in the near and long-term future.

From an engineering point of view, knowing the parameters analysed, designers will be able to estimate the scour effects and the stability of the dam with greater certainty.

LABORATORY EXPERIMENTS

Turbulent jet experimental facility

The experimental facility was constructed at the Hydraulics Laboratory of the Universidad Politécnica de Cartagena, and was described in detail in Castillo & Carrillo (2012, 2013) and Carrillo & Castillo (2014). Here, we revisit its main features. The facility consists of a mobile mechanism which permits variation of the weir height between 1.7 and 3.5 m and flows from 10 to 150 L/s. It has an inlet channel with a length of 4.0 m and width of 0.95 m. The discharge is produced through a sharp-crested weir with a width of 0.85 m and height of 0.37 m.

The plunge pool, in which different water cushions may be simulated, is a 1.3-m high, 1.1-m wide and 3-m long methacrylate box. Turbulent kinetic energy values at the inlet channel were obtained with an acoustic Doppler velocimeter (ADV); mean velocities and air concentrations in different sections of the falling jet were acquired with optical fibre instrumentation; and instantaneous pressure values were measured with piezoresistive transducers located on the basin bottom. In addition, ADV and optical fibre probe were used in the basin to obtain velocity and air concentration profiles, respectively, downstream of the impingement point.

The flow volumetric rate was measured with a V-notch weir and was tested with ADV measurement upstream from the weir. Differences were smaller than 5%.

Figure 3 shows a picture of the experimental device in which sizable values of air concentration are apparent. Air incorporation and transport are more important for shallower water cushions.
Figure 3

Air entrainment in the plunge pool observed in the laboratory device with q = 0.082 m2/s, H = 2.19 m, and Y = 0.32 m.

Figure 3

Air entrainment in the plunge pool observed in the laboratory device with q = 0.082 m2/s, H = 2.19 m, and Y = 0.32 m.

The ADV settings

The setting characteristics of the ADV employed in these tests were determined by considering that the main objective was to characterise the turbulence. The velocity range was selected as ±4.00 m/s (the maximum available in the equipment). With this setting, the ADV was able to measure horizontal velocities up to 5.25 m/s and vertical velocities up to 1.50 m/s. Due to the fact that the Doppler equipment needs to be totally submerged into water, the first 5–6 cm of the water column could not be measured.

We used the turbulent kinetic energy measured 0.50 m upstream of the weir in the experimental facility to specify inlet boundary conditions in the numerical simulations in that location. This distance guarantees that the stream lines are horizontal upstream of the sharp weir crest (0.50 m > 5 h).

The plunge pool was surveyed in cross sections spaced 0.10 m horizontally. Four specific flows were tested (0.023, 0.037, 0.058 and 0.082 m2/s) with different water cushion depths. This covers 24 different configurations generated downstream of a rectangular free-falling jet.

In order to characterize the macro turbulence of the flow in the plunge pool, 5,000 values were recorded in each measured point by using a frequency of 10 Hz (more than 8 minutes of observation). In this way, 2006 points in the symmetrical vertical plane of the basin were obtained. As the flows are highly turbulent, the values obtained with ADV may be affected by spurious signals or ‘spikes’. Each time series must be filtered with the velocity and acceleration thresholds (Wahl 2000). Furthermore, in this particular case, the air may also affect the signal of the ADV. Frizell (2000) studied the air effect measuring concentrations varying from 0 to 3.61%. As the air concentrations increase and bubble sizes increase, correlation values drop dramatically as the acoustic signals used by the probe are absorbed and reflected by the two-phase flow mixture. Matos et al. (2002) also found that air bubbles affect the accuracy of velocity measurements taken with the ADV. However, their experimental results suggest that the ADV can provide reasonable estimates of the velocity for low air concentrations up to 8%.

Optical fibre probe

To measure the air concentration at the falling jet and at the basin, an RBI-instrumentation dual-tip probe optical fibre phase-detection instrument was used. The rise and fall of the probe signal corresponds, respectively, to the arrival and the departure of the gas phase at the tip of the sensor. The thresholding values were set to 1.0 and 2.5 V (Boes & Hager 1998). The void fraction was defined as the ratio of the total time the probe is in the gas (ΣtGi) to the experiment duration time t.

According to Stutz & Reboud (1997a, 1997b), the equipment enables measurement in water up to 20 m/s flow velocity and the relative uncertainty concerning the void fraction is estimated at about 15% of the measured value. One source of error in the estimation of air presence in the flow is due to the counting statistic of the number of air bubbles in contact with the tips of the probe (Stutz 1996). Therefore, a short duration of the sequence would contribute to an increased inaccuracy of the result. In order to evaluate the minimal measuring duration, André et al. (2005) considered the stabilization of the mean value during the measuring sequence and the quasi-stationary of the air concentration signal as statistically representative for the air concentration. Based on the sensitivity study of the probe behaviour, the authors recommend a sampling sequence of 60 s as a good compromise between accuracy and time consumption.

Boes & Hager (2003) carried out experiments with 4,000 air bubbles and sampling sequences of 30 seconds. The authors considered the accuracy of the air concentration and velocity measurements is related to the variation of the phase Nb, air to water or inversely, rather than to the sampling duration t.

Following those ideas, in this study a sampling sequence of 90 was considered. Figure 4 shows the void fraction evolution until a relative uncertainty of around 1% is reached and the bubbles number detected in the measurement.
Figure 4

Temporal convergence of the void fraction (left) and bubble number detected (right).

Figure 4

Temporal convergence of the void fraction (left) and bubble number detected (right).

Figure 5 shows the air concentration obtained by means of the optical fibre probe in different sections downstream of the jet stagnation point. The equipment measures the air-water ratio, (entrained air discharge rate to water flow rate) and, from this value, the air concentration, was calculated.
Figure 5

Air concentration in the basin for different sections downstream of the jet stagnation point. Measurements obtained by means of an optical fibre probe (q = 0.082 m2/s, H = 2.19 m, and Y = 0.32 m).

Figure 5

Air concentration in the basin for different sections downstream of the jet stagnation point. Measurements obtained by means of an optical fibre probe (q = 0.082 m2/s, H = 2.19 m, and Y = 0.32 m).

The maximum air concentration is around 12% (at a distance of 21% from the bottom) for the first sections. However, from the section 0.30 m and a distance from the bottom smaller than 70%, the air concentration is below 10%. Concentrations remain high still at the upper portion of the water depth in the basin.

Filtering the velocity records

Considering highly turbulent and aerated flow that occurs in the basin of energy dissipation, the Phase-Space Thresholding filter (Goring & Nikora 2002, modified by Castillo 2009) was used. The spikes were replaced on each record by the mean value of the twelve closest points.

This filter is based on the fact that the numerical derivative of a signal enhances its high frequency components, i.e. it enhances the spikes. The method uses the concept of a three-dimensional Poincaré map or phase-space plot in which the variable and its derivatives are plotted against each other. The points are enclosed by an ellipsoid and the points outside the ellipsoid are designated as spikes (Castillo 2009). The threshold that is usually applied arises from a theoretical result from the normal probability distribution theory which says that for n independent, identically distributed, standard, normal, random variable ξi, the expected absolute maximum is: 
formula
2a
where λU is denominated the Universal threshold (Donoho & Johnstone 1994; Goring & Nikora 2002). For a normal, random variable whose standard deviation is estimated by σ and zero mean, the expected absolute maximum is: 
formula
2b
Nonetheless, it should be noted that turbulence is not normally distributed and, therefore, the theoretical result concerning the constant applies only approximately.

The main steps of the method are as follows:

  1. Calculate surrogates for the first, , and second, , derivatives using a centered differences scheme.

  2. Calculate the standard deviations of all three variables, , , and , and thence the expected maxima using the Universal criterion

  3. Calculate the rotation angle of the principal axis of versus , using an expression for the cross correlation. Instead of using the expression by Goring & Nikora (2002): 
    formula
    3a

Castillo (2009) proposed a new relation obtained by a Gauss’ fit:

  •  
    formula
    3b
  • 4. For each pair of variables, calculate the ellipse that has maxima and minima from point 3.

    • • For versus , the major axis is and, the minor axis is .

    • • For versus , the major axis is and, the minor axis is .

    • • For versus , the major and minor axes, a and b, respectively, are the solutions of (Goring & Nikora 2002): 
      formula
      4a
       
      formula
      4b
      Castillo (2009) proposed the following system of equations instead: 
      formula
      4c
       
      formula
      4d

Figure 6 shows the original and filtered signals measured at a point located 0.40 m downstream from the jet stagnation point and 12.34% of the water depth. We can also see the three ellipsoids defined by the Universal criterion. The points outside of the ellipsoids are designated as spikes and, in this particular case, the number of points have been 690. The mean velocity of the filtered signal, 117.10 cm/s is very similar to the original signal, 124.00 cm/s; however, the standard deviation is reduced from 106.47 to 79.93 cm/s.
Figure 6

Original (a) and filtered (b) signals corresponding to the point measured to 0.40 m downstream from jet stagnation point (12.34% of the water depth); (c), (d) and (e) are the three ellipsoids defined by the Universal criterion (q = 0.082 m2/s, H = 2.19 m, and Y = 0.32 m).

Figure 6

Original (a) and filtered (b) signals corresponding to the point measured to 0.40 m downstream from jet stagnation point (12.34% of the water depth); (c), (d) and (e) are the three ellipsoids defined by the Universal criterion (q = 0.082 m2/s, H = 2.19 m, and Y = 0.32 m).

From these plots it can be concluded that the filtering attenuates the standard deviation. Further filtering affects the values of the turbulent kinetic energy.

MATHEMATICAL AND NUMERICAL MODELLING

As can be seen from Figures 3 and 5, the flow conditions in the plunge pool are such that the air concentrations are relatively elevated at the point of jet impingement and nearby areas and in the top layer of the water depth. In these areas, there is a mostly non-dilute, two-phase flow. However, as we move far from the impingement point, the flow conditions tend to become quasi-dilute. That is why we decided to solve the equations for the conservation of mass and momentum for the mixture, which may be written in compact form (ANSYS CFX Manual 2015) as: 
formula
5
where is the transported quantity, i and j are indices which range from 1 to 3, xi represents the coordinate directions (1 to 3 for x, y, z directions, respectively), and t the time. In turn, and with rk indicating the volume fraction of kth fluid, denoting the diffusion coefficient associated with the transported quantity for phase k, Np denoting the number of phases and S indicating the sources/sinks for the transported quantity (ANSYS CFX Manual 2015). In this model, phases share the same velocity field. When , S = 0, and Γ= 0, the mass conservation equation is recovered, and when , the momentum equation is recovered, with its corresponding source terms to account for the Reynolds stresses.

The theoretical model (1) comes as a result of the addition of the equations of the two phases (Drew & Passman 1999; ANSYS CFX Manual 2015). Further, . Rigorously speaking, models like this have been found to provide adequate predictions only in relatively-dilute mixtures. Jha & Bombardelli (2010) established that dilute mixtures can extend up to the range 2–5% in the context of sediment-laden flows; something similar could be assumed in bubble flows. For larger concentrations, Jha & Bombardelli (2010) found that the velocity distribution could not be well predicted relatively far from the wall with mixture models. Thus, we expect the ‘homogeneous’ model to be able to represent rather adequately those areas in which air concentrations are not that high (cf. Figure 5).

As said, the code ANSYS CFX has been used, which is based on an element-oriented, finite-volume method. It allows different types of volumes, including tetrahedral and hexahedral volumes. Solution variables are stored at the nodes (mesh vertices). More details are given in the ANSYS CFX Manual (2015).

Finally, the instantaneous values in the numerical simulations show large variations during the establishment of the quasi-steady conditions. Once the quasi-steady conditions are reached, there are small fluctuations around the mean value.

Turbulence models

In this work, some of the most usual two-equation turbulence models have been tested for the free falling jet and basin case. These models use the gradient diffusion hypothesis to relate the Reynolds stresses to the mean velocity gradients and the turbulent viscosity. The eddy viscosity hypothesis considers that eddies behave like molecules and the Boussinesq model assumes that the Reynolds stresses are proportional to the mean velocity gradients, as follows (Pope 2000): 
formula
6
with μt being the eddy viscosity or turbulent viscosity, ui is the mean square fluctuation of the velocity in i direction, is the turbulent kinetic energy and δij the Kronecker delta.

The Standard k-ɛ model is considered to be the traditional turbulence model and it is embedded in the majority of the CFD programs. The Re-Normalisation Group (RNG) k-ɛ model is based on the renormalisation group analysis of the Navier–Stokes equations (Yakhot & Orszag 1986; Yakhot & Smith 1992). The transport equations for turbulent kinetic energy and dissipation rate are similar as those for the standard model, although their respective constants are different.

The k-ω based Shear-Stress Transport (SST) model (Menter 1994) assumes that the eddy viscosity is linked to the turbulence kinetic energy, k, and the turbulent frequency, ω, as . The SST model takes into account the accuracy of the k-ω model in the near-wall region and the free stream independence of the k-ɛ model in the outer part of the boundary layer. It is considered as a hybrid model (see Rodi et al. 2012).

Free surface modelling

The free surface model assumes that each control volume contains three possible conditions:

  • rk= 0 if cell is empty (of the k-th fluid);

  • rk= 1 if cell is full (of the k-th fluid);

  • 0 < rk < 1 if cell contains the interface between the fluids.

Tracking of the interface between fluids is accomplished by the solution of the volume fraction equation (see ANSYS CFX Manual 2015).

MODEL IMPLEMENTATION IN THREE DIMENSIONS

Like in the simulation of the pressure field on the basin bottom (Castillo et al. 2014), all scenarios were obtained by a calculation time of 60 seconds, using a constant time-step of 0.05 s. The transient statistics were obtained by considering that quasi-steady-state conditions are reached after the first 20 seconds of simulation (Castillo et al. 2014). Hence, statistics were obtained from time-step 400 to time-step 1,200 (800 data in 40 s).

Boundary conditions

The model boundary conditions corresponded to the turbulence kinetic energy at the inlet obtained with ADV (located 0.50 m upstream of the weir), the upstream and downstream water levels and their hydrostatic pressure distributions. Figure 7 shows the computational domain employed.
Figure 7

Schematic of computational domain and boundary conditions, mesh blocks and water volume fraction.

Figure 7

Schematic of computational domain and boundary conditions, mesh blocks and water volume fraction.

ANSYS CFX has different treatments near the wall. ω-based turbulence models (e.g. SST) use automatic wall treatments which switch between regular wall functions (Pope 2000) and low-Reynolds wall treatment (Menter 1994). Wall functions are used when the wall adjacent vertices are in the log-law layer (y+ ≈ 20–200). The low-Reynolds wall treatment is used when the wall adjacent vertices are in the viscous sublayer (ANSYS CFX Manual 2015). Considering the wall treatment used by ANSYS CFX, the mesh sizes close to the solid boundary were smaller than in the rest of the domain. Values of y+ were smaller than 40.

For simplicity, only a half part of the model was simulated. The symmetry condition in the longitudinal plane of the plunge pool was used.

The inlet condition considers the volumetric flow rate with a normal direction to the boundary (q = 0.082 m2/s, q = 0.058 m2/s, q = 0.037 m2/s, q = 0.023 m2/s), the turbulent kinetic energy (5.1 × 10–4 m2/s2 for q = 0.082 m2/s, 3.6 × 10–4 m2/s2 for q = 0.058 m2/s, 1.9 × 10–4 m2/s2 for q = 0.037 m2/s, and 1.1 × 10–4 m2/s2 for q = 0.023 m2/s), and the water level height at the upstream deposit (2.313 m for q = 0.082 m2/s, 2.285 m for q = 0.058 m2/s, 2.263 m for q = 0.037 m2/s, 2.237 m for q = 0.023 m2/s).

The outlet condition was considered with flow normal to the boundary and hydrostatic pressure. The water level height at the outlet was modified according to the water cushion depth, Y, in the laboratory device. For all walls of the upper deposit, the weir and the dissipation basin, no slip smooth wall conditions were considered. The roughness of methacrylate was indicated in the walls. In the transverse direction, wall boundary conditions were used.

Mesh-independence tests

In Figure 8, the simulations results for the different mesh sizes (5, 7.5, 10 and 12.5 mm) in the free falling jet, obtained as a function of the vertical distance to the weir in terms of the flow velocities in the jet, are shown. Differences in velocities with the optical fibre probe measurement are smaller than 2% in all the cases (Castillo et al. 2014). In Table 1, the size of elements and the corresponding number of volumes required in the different simulations are indicated. From the analysis of Figure 8, it can be concluded that mesh-independence is reached with an element size of 10 mm. The results are in agreement with previous results obtained on pressures at the stagnation point (Castillo et al. 2014). In this way, the mesh size of 10 mm seems to be valid for the flow rates analysed.
Table 1

Number of volumes as a function of mesh size

Mesh size (mm) Number of volumes 
5.0 1,088,896 
7.5 457,810 
10.0 255,776 
12.5 134,785 
Mesh size (mm) Number of volumes 
5.0 1,088,896 
7.5 457,810 
10.0 255,776 
12.5 134,785 
Figure 8

Velocities in the falling jet as a function of the mesh size: q = 0.058 m2/s, h = 0.095 m.

Figure 8

Velocities in the falling jet as a function of the mesh size: q = 0.058 m2/s, h = 0.095 m.

Figure 9 shows a view of the free surface in the three-dimensional numerical model. We can see that the solution correctly reproduces the expected features of the jet and the plunge pool observed in the experimental facility.
Figure 9

Free surface in the falling jet for the case of a mesh of 10 mm.

Figure 9

Free surface in the falling jet for the case of a mesh of 10 mm.

Convergence criteria

To judge the convergence of iterations in the numerical solution, we monitored the residuals (Wasewar & Vijay Sarathi 2008). The solution is said to have converged in the iterations if the scaled residuals are smaller than fixed values ranging between 10–3 and 10–6. In this work, the residual values were set to 10–4 for all the variables. Each time-step required around 12 iterations to reach the convergence criterion. With this choice and for 791,354 elements (255,776 in the falling jet), the mean computational time was 7.2 × 105 seconds (≈8.3 days), using a Central Processing Unit (CPU) with sixteen processors (Intel® Xeon® E5-2699 v3 @ 2.30 GHz).

EMPIRICAL FORMULAS

In this section, we present some formulas obtained from the literature to interpret the different hydraulic parameters of plunge pools. To that end, we selected expressions for submerged hydraulic jumps and horizontal wall jets.

Velocity distribution in the plunge pool

Following Rajaratnam (1965), we can compare the velocity profiles in the forward flow if they are normalized with a velocity scale equal to the maximum velocity, Vmax, at any section, and with a length scale δl equal to the elevation y from the bottom where the local velocity V = Vmax/2, and the velocity gradient is negative (see Figure 10).
Figure 10

Velocity distribution sketch for submerged jumps (adapted from Wu & Rajaratnam 1996).

Figure 10

Velocity distribution sketch for submerged jumps (adapted from Wu & Rajaratnam 1996).

To characterize the non-dimensional velocity distribution in hydraulic jumps, some authors have proposed different expressions which we include in the Appendix. Using the results obtained by diverse researchers, Wu & Rajaratnam (1996) considered the length scale δl for the submerged jumps as a function of the distance to the beginning of the jump. They found that most of the observations for submerged jumps are contained within one standard deviation of the mean value of the wall jet, and only the data points near the end of the jump show an accelerated growth rate.

Energy dissipation in the plunge pool

In a horizontal channel, the total energy variation between the sections located upstream and downstream of the submerged hydraulic jump are, by definition: 
formula
7
where y3 and y4 are the water depths upstream and downstream of the submerged hydraulic jump generated by the jet. By using Equation (7) with the continuity equation, the energy dissipation may be obtained as (Ohtsu et al. 1990): 
formula
8
where H0 is the energy upstream of the hydraulic jump, and y0 and Fr0 the water depth and Froude number in the upstream section of the hydraulic jump. When y3=y0 and y4=y2, the free hydraulic jump expression is recovered. For the plunge pool case, H0 may be assumed as the energy upstream the weir, while y0 may be estimated with Equation (1) or the Bernoulli equation without considering the losses.

RESULTS AND DISCUSSION

Laboratory velocity profiles

Figure 11 shows the velocity profiles obtained in the laboratory, together with equations proposed by diverse authors for horizontal wall jets (Görtler 1942; Rajaratnam 1976; Ohtsu et al. 1990; Liu et al. 1998; Lin et al. 2012) (see Appendix). The overall behaviour of the observations can be predicted rather well by existing equations up to . Disagreements appear when ratio Vx/Vmax < 0.4. This seems to be related to the angle of impingement of the jet. In hydraulic jump studies, the wall jet is horizontal, while the impingement jet enters almost vertical in these tests. The higher scatter occurs when the water cushion depth is Y/Bj > 20. In this way, the self-similarity disappears when the velocity profiles are close to the stagnation point and when a very submerged condition is obtained for the hydraulic jump.
Figure 11

Results of laboratory measurements of the non-dimensional profiles of the horizontal velocity in the central vertical plane of the plunge pool (X ≥ 0.40 m). Profiles encompass all discharges and water cushion depths tested.

Figure 11

Results of laboratory measurements of the non-dimensional profiles of the horizontal velocity in the central vertical plane of the plunge pool (X ≥ 0.40 m). Profiles encompass all discharges and water cushion depths tested.

Following Wu & Rajaratnam (1996), Figure 12 shows the results of the characteristic length obtained through the plunge pool. For each horizontal velocity profile measured with the ADV, the length scale δl was obtained in each section (X distance to the stagnation point). Data have been classified as a function of the ratio water depth in the plunge pool/impingement jet thickness of each test. For ratios Y/Bj up to 20, the behaviour is similar to that found in wall jets (albeit with another slope, as expected). The values for water cushion depths Y/Bj up to 30 tend to fall within one standard deviation of the mean value. However, for larger water cushion depths, the characteristic length is higher. In this type of submerged hydraulic jump where the falling jet enters almost vertical, an equation may be obtained with the data that fall within one standard deviation of the mean value: 
formula
9
Figure 12

Characteristic length δl in submerged hydraulic jumps.

Figure 12

Characteristic length δl in submerged hydraulic jumps.

Energy dissipation in the plunge pool

Figure 13 shows the contrast between the relative energy dissipation and the Froude number at the jet impingement condition, , obtained from experiments. In addition, results coming from the use of Equation (8) have been included as a function of the ratio between the upstream water depth and the impingement jet thickness (y3/Bj). Laboratory data are well represented by Equation (8). In the laboratory device, the impingement Froude number is between 13 and 20 for the impingement jet thickness of 0.015 m plotted in Figure 14. In general, tests carried out show an energy dissipation larger than 75% of the impingement jet energy. This ratio increases when the ratio y3/Bj decreases.
Figure 13

Relative energy dissipation in the plunge pool as a function of the impingement Froude number.

Figure 13

Relative energy dissipation in the plunge pool as a function of the impingement Froude number.

Figure 14

Relative energy dissipation in the plunge pool as a function of the ratio y3/Bj for the cases Bj = 0.015 m and Fj= 13–20.

Figure 14

Relative energy dissipation in the plunge pool as a function of the ratio y3/Bj for the cases Bj = 0.015 m and Fj= 13–20.

Velocity and turbulent kinetic energy distributions

Velocities at different cross sections of the plunge pool located downstream of the stagnation point were measured with ADV. Results for the same cross sections were obtained from the CFD simulations. The velocities and turbulent kinetic energies have been made dimensionless by using the maximum horizontal values, Vmax,x and kmax, respectively (see Figure 15).
Figure 15

Horizontal velocity and turbulent kinetic energy profiles in the plunge pool downstream of the stagnation point (q = 0.082 m2/s, H = 1.993 m, Y = 0.32 m, X ≥ 0.40 m): (a) SST velocities; (b) RNG k-ε velocities; (c) standard k-ε velocities; and (d) SST turbulent kinetic energy.

Figure 15

Horizontal velocity and turbulent kinetic energy profiles in the plunge pool downstream of the stagnation point (q = 0.082 m2/s, H = 1.993 m, Y = 0.32 m, X ≥ 0.40 m): (a) SST velocities; (b) RNG k-ε velocities; (c) standard k-ε velocities; and (d) SST turbulent kinetic energy.

In general, the measured horizontal velocities, Vx, are greater than the calculated counterparts. The bellies which appear in some profiles seem to be related with the three-dimensional flow features of the flow. This level of agreement is understandable given the observed air concentrations, close to 10%. As said with this level of concentrations of the disperse phase, Jha & Bombardelli (2010) found that homogeneous-type models can only provide approximate values of velocity. In Figure 15(a)15(c) we can see that the choice of the turbulence model has a significant effect over the result of the horizontal flow profiles. The less accurate results were obtained with the RNG k-ɛ turbulence model. However, the SST and the standard k-ɛ turbulence models produce velocities fairly close to the values measured in the laboratory.

In addition to the mean velocities, the turbulent kinetic energy profiles were also compared (Figure 15(d)). Like in the mean velocity cases, the best results of turbulent kinetic energy profiles were obtained with the SST model. In general, the results from the numerical simulations show the same behaviour as the results obtained in the laboratory. However, the differences are very important close to the stagnation point, where the numerical model may not obtain accurate results due to the relatively important air entrainment into the plunge pool (Figures 3 and 5).

Figure 16 shows the non-dimensional velocity profiles obtained from the numerical simulations as well as the laboratory measurements. Results are quite similar to the profiles obtained with the above-mentioned formulae. Due to the strong recirculation with air entrainment, in the upper side of the submerged hydraulic jumps there is a bigger scatter for ratios Vx/Vmax < 0.40.
Figure 16

Velocity distribution downstream of the stagnation point (X ≥ 0.40 m) obtained through numerical simulations and laboratory experiments.

Figure 16

Velocity distribution downstream of the stagnation point (X ≥ 0.40 m) obtained through numerical simulations and laboratory experiments.

With all data, a new regression is proposed for submerged hydraulic jumps downstream of the impingement point: 
formula
10
This proposed function is the separation line between the profiles in which there is negative recirculation flow and the profiles in which the flow is moving towards downstream. For the range of flows and water cushions analyzed, the limit between both behaviours seems to be located at 0.2–0.3 m downstream of the stagnation point.
For the extreme negative and positive flow profiles, two regressions (valid for ratios Vx/Vmax < 0.40) are also proposed, respectively: 
formula
11
 
formula
12
In all cases, data collapse for ratios Vx/Vmax ≥ 0.40. Data with smaller ratios do not match with a single law. As before, the larger differences appear when big water cushions are analyzed (ratios Y/Bj > 20).

CONCLUSIONS

Observing and predicting two-phase flows in hydraulic structures is very complicated due to the rather non-dilute nature of the flow. Under non-dilute conditions, both experiments and simulations cannot be expected to lead to clean comparisons.

In this work, mean velocity and turbulent kinetic energy profiles have been analyzed in a plunge pool located downstream of a rectangular free-falling jet. The energy dissipation in a plunge pool is very high. For the test carried out, the dissipation of the impingement jet energy was between 75 and 95%. This ratio increases when the ratio water cushion/impingement jet thickness decreases.

In general, the CFD simulations provided results fairly close to the values measured in the laboratory, and to the formulas proposed by diverse authors, in spite of having used a simple two-phase flow model. ‘Homogeneous’ model seems to be able to predict rather well areas in which air concentration is not very high. However, in the highly aerated regions rather strong differences appear.

Regarding the modelling of turbulence, three closures were tested. SST turbulence model offered results closer to the laboratory measurements. The RNG k-ɛ model tends to underestimate the turbulent kinetic energy.

It was possible to propose a single mean velocity distribution law for ratios Vx/Vmax ≥0.40. For smaller values, there are necessarily diverse distribution laws.

In order to develop this work further, we plan to examine air entrainment in the stilling basin. Comparison of results with diverse CFD codes (open source and commercial ones) against data will be considered.

ACKNOWLEDGEMENTS

The first two researchers express their gratitude for the financial aid received from the Ministerio de Economía y Competitividad and the Fondo Europeo de Desarrollo Regional (FEDER), through the Natural Aeration of Dam Overtopping Free Jet Flows and its Diffusion on Dissipation Energy Basins project (BIA2011-28756-C03-02). The first author acknowledges the support of Ministerio de Educación, Cultura y Deporte of España through Estancias de Movilidad de Profesores e Investigadores Senior program (PRX 14/00367), which allowed him to develop a stay as a Visiting Scholar Researcher at the Department of Civil and Environmental Engineering of the University of California, Davis, from April to October 2015.

REFERENCES

REFERENCES
Annandale
G. W.
2006
Scour Technology. Mechanics and Engineering Practice
.
McGraw-Hill
,
New York
.
ANSYS, Inc.
2015
ANSYS CFX. Reference Guide. Release 16.0
.
ANSYS Inc., Southpoint
,
Canonsburg, PA
,
USA
.
Boes
R.
Hager
W. H.
1998
Fiber-optical experimentation in two-phase cascade flow
. In:
Proceedings of the International RCC Dams Seminar
(
Hansen
K.
, ed.).
Denver
,
CD
,
USA
.
Boes
R.
Hager
W. H.
2003
Two-phase flow characteristics of stepped spillways
.
J. Hydraul. Eng.
129
(
9
),
661
670
.
Bombardelli
F. A.
2004
Turbulence in Multi-Phase Models for Aeration Bubble Plumes
.
PhD Dissertation
,
Department of Civil and Environmental Engineering, University of Illinois at Urbana-Champaign
.
Bombardelli
F. A.
Gioia
G.
2005
Towards a theoretical model for scour phenomena
. In:
Proc. RCEM 2005, 4th IAHR Symposium on River, Coastal, and Estuarine Morphodynamics
(
Parker
G.
García
M.
, eds).
Urbana
,
IL
,
USA
,
2
, pp.
931
936
.
Carrillo
J. M.
2014
Metodología numérica y experimental para el diseño de los cuencos de disipación en el sobrevertido de presas de fábrica
.
PhD Thesis
.
Universidad Politécnica de Cartagena
,
Spain
(in Spanish)
.
Carrillo
J. M.
Castillo
L. G.
2014
Laboratory measurements and numerical simulations of overtopping nappe impingement jets
. In:
Proc. Int. Conf. 1st International Seminar on Dam Protection against Overtopping and Accidental Leakage
,
24–26 November 2014
,
Madrid, Spain
, pp.
219
230
.
Castillo
L. G.
2006
Aerated jets and pressure fluctuation in plunge pools
. In:
Proceedings of the 7th International Conference on Hydroscience and Engineering
,
Philadelphia, PA
, pp.
1
23
.
Castillo
L. G.
2007
Pressure characterisation of undeveloped and developed jets in shallow and deep pool
. In:
Proceedings of the 32nd IAHR Congress
,
1–6 July 2007
,
Venice, Italy
, pp.
645
655
.
Castillo
L. G.
2009
Measurement of velocities and characterization of some parameters inside of free and submerged hydraulic jumps
. In:
Proceedings of 33rd International Association of Hydraulic Engineering & Research Congress
,
Vancouver, Canada
.
Castillo
L. G.
Carrillo
J. M.
2012
Hydrodynamics characterization in plunge pools. Simulation with CFD methodology and validation with experimental measurements
. In:
Proc. Int. Conf. 2nd IAHR European Congress
,
27–29 June 2012
,
Munich, Germany
.
Castillo
L. G.
Carrillo
J. M.
2013
Analysis of the scale ratio in nappe flow case by means of CFD numerical simulation
. In:
Proceedings of the 2013 IAHR Congress
,
8–13 September 2013
,
Chengdu, China
.
Castillo
L. G.
Carrillo
J. M.
Sordo-Ward
A.
2014
Simulation of overflow nappe impingement jets
.
J. Hydroinform.
16
(
4
),
922
940
.
Castillo
L. G.
Carrillo
J. M.
Blázquez
A.
2015
Plunge pool mean dynamic pressures: a temporal analysis in nappe flow case
.
J. Hydraul. Res.
53
(
1
),
101
118
.
Donoho
D. L.
Johnstone
I. M.
1994
Ideal spatial adaptation by wavelet shrinkage
.
Biometrika
81
(
3
),
425
455
.
Drew
D. A.
Passman
S. L.
1999
Theory of Multicomponent Fluids. Applied Mathematical Sciences
.
Vol. 135
.
Springer
,
New York
.
Ervine
D. A.
Falvey
H. R.
1987
Behaviour of turbulent water jets in the atmosphere and plunge pools
. In:
Proceedings of the International Conference of the Institutions of Civil Engineers
,
83
, pp.
295
314
.
Ervine
D. A.
Falvey
H. T.
Withers
W. A.
1997
Pressure fluctuations on plunge pool floors
.
J. Hydraul. Res.
35
(
2
),
257
279
.
Falvey
H. T.
1990
Cavitation in Chutes and Spillways. Engineering Monograph No. 42
.
Bureau of Reclamation
,
Denver, CO
,
USA
.
FEMA
2014
Technical Manual: Overtopping Protection for Dams. Federal Emergency Management Agency. FEMA P-1014, May
.
US Department of Homeland Security
,
USA
.
Frizell
K. W.
2000
Effects of aeration on the performance of an ADV
. In:
2000 Joint Conf. on Water Resources Engineering and Water Resources Planning & Management
(
Hotchkiss
R. H.
Glade
M.
, eds).
ASCE
,
Minneapolis
,
USA
,
(CD-ROM)
.
Gioia
G.
Bombardelli
F. A.
2005
Localized turbulent flows on scouring granular beds
.
Phys. Rev. Lett.
95
,
014501
.
Goring
G.
Nikora
V.
2002
Despiking acoustic Doppler velocimeter data
.
J. Hydraul. Eng.
128
(
1
),
117
126
.
Ho
D. K. H.
Riddette
K. M.
2010
Application of computational fluid dynamics to evaluate hydraulic performance of spillways in Australia
.
Aust. J. Civil Eng.
6
(
1
),
81
104
.
Jha
S. K.
Bombardelli
F. A.
2010
Toward two-phase flow modeling of non-dilute sediment transport in open channels
.
J. Geophys. Res. Earth Surf.
115
(
F3
),
2003
2012
.
Lin
C.
Hsieh
S.-C.
Lin
I.-J.
Chang
K.-A.
Raikar
R. V.
2012
Flow property and self-similarity in steady hydraulic jumps
.
J. Exp. Fluids
53
,
1591
1616
.
Matos
J.
Frizell
K. H.
Andre
S.
Frizell
K. W.
2002
On the performance of velocity measurement techniques in air-water flows
. In:
Hydraulic Measurements and Experimental Methods Conference 2002
(
Wahl
T. L.
Pugh
C. A.
Oberg
V. A.
Vermeyen
T. B.
, eds).
ASCE
,
Estes Park, CO
,
USA
.
Ohtsu
F.
Yasuda
Y.
Awazu
S.
1990
Free and Submerged Hydraulic Jumps in Rectangular Channels. Report of the Research Institute of Science and Technology, 35
.
Nihon University
,
Tokyo
,
Japan
.
Pope
S. B.
2000
Turbulent Flows
.
Cambridge University Press
,
Cambridge
.
Rajaratnam
N.
1965
The hydraulic jump as wall jet
.
Proc. ASCE J. Hydraul. Div.
91
(
HY5
),
107
132
.
Rajaratnam
N.
1976
Turbulent Jets
.
Elsevier Scientific, Development in Water Science
,
5
,
New York
,
USA
.
Rodi
W.
Constantinescu
C.
Stoesser
T.
2012
Large-Eddy Simulation in Hydraulics, IAHR Monograph
,
CRC Press
.
Stutz
B.
1996
Analyse de la structure diphasique et instationnaire de poches de cavitation (Analysis of the two Phases and Non-Permanent Structure of Cavitation Bubbles)
.
PhD Thesis
,
Institut National Polytechnique de Grenoble
,
France
(in French)
.
Stutz
B.
Reboud
J. L.
1997a
Experiment on unsteady cavitation
.
Exp. Fluids
22
(
1997
),
191
198
.
Stutz
B.
Reboud
J. L.
1997b
Two-phase flow structure of sheet cavitation
.
Phys. Fluids
9
(
12
),
3678
3686
.
Wahl
T. L.
2000
Analyzing ADV data using WinADV
. In:
Joint Conference on Water Resources Engineering and Water Resources Planning & Management
,
July 30–August 2
,
Minneapolis, MN
.
Wahl
T. L.
Frizell
K. H.
Cohen
E. A.
2008
Computing the trajectory of free jets
.
J. Hydraul. Eng.
134
(
2
),
256
260
.
Wasewar
L.
Vijay Sarathi
J.
2008
CFD modelling and simulation of jet mixed tanks
.
Eng. Appl. Comp. Fluid Mech.
2
(
2
),
155
171
.
Wu
S.
Rajaratnam
N.
1996
Transition from hydraulic jump to open channel flow
.
J. Hydraul. Eng.
122
(
9
),
526
528
.
Yakhot
V.
Orszag
S. A.
1986
Renormalization group analysis of turbulence. I. Basic theory
.
J. Sci. Comput.
1
(
1
),
3
51
.

APPENDIX

Expressions of the velocity distribution in wall jets and hydraulic jumps which have been provided by different authors are:

Author Formula 
Görtler (1942), cited by Liu et al. (1998)   
Rajaratnam (1976)   
  
Ohtsu et al. (1990)   
Lin et al. (2012)   
where  
Author Formula 
Görtler (1942), cited by Liu et al. (1998)   
Rajaratnam (1976)   
  
Ohtsu et al. (1990)   
Lin et al. (2012)   
where