A conceptual model of check dam hydraulics for gully control : efficiency , optimal spacing and relation with step-pools

There is little information in scientific literature regarding the modifications induced by check dam systems in flow regimes within restored gully reaches, despite it being a crucial issue for the design of gully restoration measures. Here, we develop a conceptual model to classify flow regimes in straight rectangular channels for initial and damfilling conditions as well as a method of estimating efficiency in order to provide design guidelines. The model integrates several previous mathematical approaches for assessing the main processes involved (hydraulic jump, impact flow, gradually varied flows). Ten main classifications of flow regimes were identified, producing similar results when compared with the IBER model. An interval for optimal energy dissipation (ODI) was observed when the steepness factor c was plotted against the design number (DN, ratio between the height and the product of slope and critical depth). The ODI was characterized by maximum energy dissipation and total influence conditions. Our findings support the hypothesis of a maximum flow resistance principle valid for a range of spacing rather than for a unique configuration. A value of c = 1 and DN ∼ 100 was found to economically meet the ODI conditions throughout the different sedimentation stages of the structure. When our model was applied using the same parameters to the range typical of step-pool systems, the predicted results fell within a similar region to that observed in field experiments. The conceptual model helps to explain the spacing frequency distribution as well as the often-cited trend to lower c for increasing slopes in step-pool systems. This reinforces the hypothesis of a close link between stable configurations of step-pool units and man-made interventions through check dams.


Introduction
A check dam is a small transverse structure designed mainly for three purposes: control water flow, conserve soil and improve land (Doolittle in Conesa-García and Lenzi, 2010).One of its most common functions is to enhance sediment deposition, reducing the bed gradient and flow velocity in order to check soil erosion within a stream, such as a gully.Although there are examples of successful projects in gully restoration using check dams (e.g.Alcali Creek Project; Heede, 1977;Weinhold, 2007), on numerous occasions faults have been reported in the performance of these structures (Heede, 1960;Iroume, 1996;Nyssen et al., 2004), including channel degradation and scouring downstream of the check dams (Porto and Gessler, 1999;Castillo et al., 2007;Conesa-García et al., 2007).
Several approaches have contributed greatly to the understanding of some of the essential processes in drop structures, such as hydraulic jump and waterfall impact (McCorqudale and Mohamed, 1994;Vischer and Hager, 1995;Chanson, 1999), allowing a precise characterization of the energy dissipation phenomena.Physically based hydraulic models have been used to evaluate flood regimes and the influence of channel geometry in ephemeral channels in arid regions (e.g.Merrit and Wohl, 2003) and may become a useful tool for contrasting the performance of conceptual models which aim to predict the free-surface water profiles in gullies controlled by hydraulic structures.Different alternatives have been proposed for determining the spacing between adjacent check dams and there is no single universally accepted criterion.The three criteria most commonly found in the literature are (a) the head-to-toe criterion, namely, the toe of the Published by Copernicus Publications on behalf of the European Geosciences Union.

C. Castillo et al.: A conceptual model of check dam hydraulics for gully control
upstream dam is at the same elevation as the top of the downstream dam (Heede, 1960); (b) the ultimate slope criterion defining an equilibrium slope for incipient sediment motion (Porto and Gessler, 1999); and (c) empirical observations of the sediment deposit gradient in restored channels (Heede, 1978;Iroume and Gayoso, 1991;Nameghi et al., 2008).It is apparent that there are great differences in these recommendations, despite the fact that spacing is a critical factor for check dam design, and this leads to an undesirable degree of technical uncertainty and potential failure of the structure.
Despite all this research, a systematic definition of the range of flow types that might occur in check dam systems has not been found by the authors in the scientific literature.In addition, there is also very limited information on how check dam design should take into account either the initial or dam-filling conditions, given that the filling period is highly variable within a specific location (Boix-Fayos et al., 2007).It is also uncertain to what extent check dams and step-pool configurations are based on the same principles and the implications that this might have on the design of gully restoration schemes.
The analogy between step-pool systems and check dam interventions has been recently recognized, to the point that check dams have been considered as the anthropogenic equivalent to step-pool sequences in steep mountain streams (Milzow, 2004).
Step-pools represent an interesting case of spontaneous, self-organized system of high stability (Chin and Phillips, 2007).Their morphological features have inspired the design criteria for artificial check dam sequences in high-gradient streams stabilization (Lenzi, 2002;Wang and Yu, 2007;Chin et al., 2009;Wohl et al. in Conesa-García and Lenzi, 2010).Furthermore, several studies (Abrahams et al., 1995;Zimmermann and Church, 2001;Lee and Ferguson, 2002) have shown that step-pool morphologies tend to maximize flow resistance, leading to minimum velocity and shear stress, which is the final cause of its stability.A number of experimental studies at the flume scale have been conducted to explore the processes of self-formation, flow regimes and flow resistance in stepped channels (Curran and Wilcock, 2005;Comiti et al., 2009).Both flume and field works have highlighted the relevance of hydraulic jumps and turbulence in the overall flow resistance (Zimmermann and Church, 2001;MacFarlane and Wohl, 2003;Comiti et al., 2009).In spite of the suggested similarities, the maximum flow resistance in step-pool series occurred with significantly shorter spacing than that recommended for gully control using check dams (Heede, 1976;Morgan, 2005).
The main aim of the work presented in this paper is to develop a conceptual model of the hydraulics of check dam systems for gully restoration in order to establish a theoretical basis for estimating their effectiveness and proposing design guidelines.For this purpose, the following specific objectives were considered: (i) to evaluate the efficiency of check dam interventions for initial and dam-filling conditions; (ii) to compare the model performance with an accepted hydrodynamic model; and (iii) to explore the similarities between check dam interventions and step-pool channels and their implications on the design of gully restoration measures.

Description of the model
The construction of check dams in a gully reach causes a flow perturbation upstream and downstream of each structure.In initial conditions (after construction, when no silting has taken place) it creates a backwater effect by increasing the water depth immediately upstream of the structure, leading to a subcritical regime (Froude number where u is the flow velocity, g the acceleration of gravity and d the flow depth).It also produces a water drop downstream of the check dam, which accelerates the flow leading to supercritical flow conditions (F > 1) at the impact zone.The spillway performs as a control section, imposing critical flow conditions (F = 1).Since, in a restored reach, subcritical flow conditions exist in the downstream sections and the regime is supercritical in the upper part, a hydraulic jump (hereafter HJ) develops in an intermediate cross section (Fig. 1a).
For dam-filling conditions, a hydraulic jump habitually occurs between the critical flow at the downstream spillway and the supercritical regime at the impact region upstream.In this case, we assumed that the top surface of the sediment wedge was a plane extending from the spillway of the downstream check dam to the toe of the upstream check dam (Fig. 1b), an assumption similar to other authors (Nyssen et al., 2004).Despite it being a geometric simplification, it presents similarities with the profile morphology reported in stepped channels (Abrahams et al., 1995;Lenzi, 2002;MacFarlane and Wohl, 2003;Curran and Wilcock, 2005).
The conceptual model of check dam hydraulics was developed to simulate the main modifications induced by check dam construction in the flow regime along restored reaches.One of the main priorities was to keep the calculations as simple as possible, while still keeping track of the relevant processes.The main simplifying assumptions were that gullies correspond to straight fixed-bed channels, uniformly sloped, presenting rectangular-shaped cross sections with clear-water flow.
The model was programmed in a standard spreadsheet implemented in a MS Excel ® file and was structured in four interrelated modules (normal flow conditions, impact flow, free-surface profiles FSP and HJ features).The lengths associated to rapidly varied flows (i.e.impact and HJ length) were explicitly considered, since their dimensions have a relevant impact on the final FSP and, therefore, on dissipation patterns and overall efficiency.The main input parameters of the model were a. unitary discharge q ranging from 0.1-1 m 2 s −1 , the typical span of values for channels, gullies and steppool reaches (Zimmerman and Church, 2001;Merrit and Wohl, 2003;Castillo, 2012); b. bed slopes ranging from 0.02 to 0.2 (m m −1 ).This is the common interval of gradient found in check dam interventions (Heede, 1978;Nameghi et al., 2008), flume experiments (Abrahams et al., 1995;Comiti et al., 2009) and step-pool systems (e.g.Zimmerman and Church, 2001;Chartrand et al., 2011); c. check dam effective height z between 0.5-1.5 m, typical of control structures in agricultural areas (Heede, 1978;Nyssen et al., 2004;Nameghi et al., 2008); d. channel roughness ranging from 0.03 to 0.06, from clean and straight to weedy winding channels (Chow et al., 1994).

Model equations
The model features a combination of mathematical approximations to the different processes involved -Manning equation of velocity in uniform flows for estimating the flow characteristics in normal conditions before check dam construction; -Rand equations for the impact flow features in straight drop structures (Chanson, 1999); -classical hydraulic jump expressions on flat or positive slopes (Vischer and Hager, 1995;Chanson, 1999) and equations of HJ on adverse slopes (McCorquodale and Mohamed, 1994); -free-surface profiles (FSP) or backwater calculations, according to Chanson (1999), derived from the continuity and energy equations.

(a) Normal flow conditions
Uniform open channel flows are characterized by a constant flow depth and mean velocity as well as by a friction slope S f in equilibrium with bed slope S, and are usually known as normal flows (Chanson, 1999).They represent the starting situation prior to check dam construction and thereby serve as a reference for comparison to evaluate the efficiency of the conservation measures.Moreover, once the intervention has been carried out, normal conditions (hereafter, NC) define the situation toward which non-uniform flow regimes tend (equilibrium between the friction losses and the gully slope).

C. Castillo et al.: A conceptual model of check dam hydraulics for gully control
In rectangular channels in which flow depth d is significantly smaller than channel width w, the hydraulic radius R can be approximated by d.Thus, the Froude number corresponding to normal conditions F n can be estimated using Manning's expression: where the subscript n indicates normal conditions for the hydraulic variables, n Manning's roughness coefficient and S the bed slope of the gully.For a given unitary discharge q, the flow depth d in a rectangular channel can be estimated using Manning's equation and again assuming R ∼ d: Solving this for d n , we obtain: Finally, using Eqs.( 1) and ( 3), an expression for F n can be calculated, which is dependent only on primary variables: Therefore, F n , d n and u n can be estimated directly from the input parameters q, S and n, allowing the direct determination of the hydraulic regime for NC.

(b) Flow regime at the impact
Rand equations, in the form provided by Chanson (1999), were applied to determine the flow characteristics at the impact for non-submerged HJ conditions: (5) where is the critical depth, z the effective height of the check dam, L i the impact length, d i the supercritical flow depth at the impact and d 2 the subcritical flow depth after the HJ.

(c) Free-surface profiles calculation
Check dam series produce gradually varied flows (GVF), in other words sub-and supercritical zones along the gully, and rapid varied flows (RVF), such as HJ and waterfall impact (Fig. 1).As a result of these transitions, flow velocity, flow depth and shear stress vary along the channel, and normal conditions (friction slope S f equal to gully slope S at all sections), are not applicable.The value for each of the hydraulic variables in a particular cross section can be calculated by iteratively solving the continuity and energy equations following the backwater-computations methodology (Chanson, 1999).This approach can be applied to the hydraulic calculations both in subcritical regimes (controlled by backwater conditions) as well as in supercritical flows (governed by upstream controls).
The FSP determination requires hydraulic calculation in both directions, from the lower dam backwards (subcritical regime imposed by the water surface elevation) and from the upper dam forward (supercritical regime after the drop).Within the FSP approach, the step method-depth calculated from distance was applied (Chanson, 1999).This method comprised the following steps: (i) definition of the control sections downstream (at the spillway, where critical conditions are reached) and upstream (where flow characteristics can be estimated at the impact zone); and (ii) application of the differential energy equation d = S − S f 1 − F 2 • s at 0.1 m intervals in both directions in order to determine the FSP (flow depth profile) and derivative variables (velocity, friction slope, shear stress) at each cross section where d is the flow depth increment and s the distance change equal to 0.1 m.

(d) Hydraulic jump
Classical HJ expressions were used to determine the main HJ characteristics in flat or positive slopes: roller length L j , the supercritical depth d 1 , the subcritical (sequent depth) d 2 , the sequent Froude number F 2 and the amount of dissipated energy H j (Chanson, 1999): In the case of negative deposition slopes in the dam-filling scenario, the HJ equations on adverse slopes were applied for estimating the roller length L j d 1 and sequent depth d 2 d 1 ratios (McCorquodale and Mohamed, 1994).For submerged HJ, or tailwater level higher than the sequent depth of the flow at the impact, the length required to achieve tailwater Hydrol.Earth Syst.Sci., 18, 1705-1721, 2014 www.hydrol-earth-syst-sci.net/18/1705/2014/ subcritical conditions was considered equal to the nonsubmerged HJ.
The HJ characteristics were estimated by graphically comparing the FSP calculated in both directions.The sequent depth condition is established at the sections where the subcritical F 2 (calculations from the spillway upwards) and the supercritical F 1 (calculations from the impact downwards) verify Eq. ( 8), being separated by a distance of L j .For this purpose, the Froude numbers at each section of the supercritical zone were transformed into their sequent values (Eq.8) considering a longitudinal offset equal to the estimated L j (Eq.11).The intersection between both curves (subcritical F 2 upwards and sequent F from F 1 with the offset) defined the location of the subcritical section of the HJ and its F 2 value.In Fig. 2, we illustrate the methodology followed for initial conditions.For the HJ on adverse slope, we applied the same methodology using the calculated HJ ratios.

Normal conditions
The areas of occurrence for sub-and supercritical regimes in normal conditions were assessed as a function of q, S and n in order to evaluate the situations in which they appear in gully networks (Eq.4).This case is relevant not only to determine the flow characteristics for the non-intervention scenario, but also to evaluate the NC influence for initial conditions after check dam installation.

Flow regime typologies after check dam construction
We have classified the possible flow regimes using four variables: (a) silting of the check dams: initial or dam-filling conditions; (b) Froude regime for NC: subcritical or supercritical; (c) type of the HJ control: normal conditions or dam influence; and (d) level of influence: partial or total influence.Here we use the term "influence" to refer to the control exerted over the HJ characteristics either by the normal conditions in the gully or the check dam downstream.The regime of influence has been classified following two criteria: (i) type of HJ control; and (ii) level of influence.
(i) Type of HJ control -Normal conditions influence (NC): for initial conditions, the HJ features are controlled by the NC (either subcritical or supercritical) when the check dams are at enough distance to avoid dam influence, allowing the establishment of a normal flow.In the dam-filling scenario, we can assess the modified NC corresponding to the slope of the sediment wedge or deposition slope.Following Abrahams et al. (1995) and Lenzi and Comiti (2003), we define the steepness factor c by analogy with step-pool terminology as where z is the effective check dam height (to the bottom of the spillway), L the check dam spacing in horizontal projection, and S the gully slope.The deposition slope S d (Heede, 1976) can be expressed as a function of the S and c: Steepness factors values below 1 produce the modified NC.For c = 1 the deposition slope is 0, the surface of the sedimentation wedge is horizontal, following the head-to-toe rule, whereas c > 1 imply negative S d .
-Dam influence: as for the initial conditions, the HJ location is controlled by the subcritical conditions imposed by the downstream check dam as a consequence of the rise in elevation of the free-surface over the spillway.As for the dam-filling situation, this influence is only exerted for c ≥ 1.In both cases, the check dam influence dominates over the NC when adjacent check dams are close enough not to allow the development of a normal flow.

(ii) Level of influence
This classification takes into consideration the HJ efficiency as a consequence of the control imposed by the tailwater level.Two cases can be considered regarding the level of influence on the HJ: -Partial hydraulic influence (PI): we termed PI that situation in which HJ occurs at a certain distance from the check dam toe and, therefore, the supercritical Froude number of the HJ downstream of the check dam (F 1 ) is smaller than the Froude number at the impact F i .Thus, there is a distance where erosive supercritical flows are established immediately after the upstream check dam before the HJ takes place and its dissipation efficiency is below the maximum.
-Total hydraulic influence (TI): TI is characterized by an HJ occurring at the toe of the check dam and a subcritical flow along the entire downstream reach between check dams.The TI threshold takes place when F i verifies the sequent depth condition for the F 2 imposed by the downstream conditions.In this situation, F i = F 1 , all the energy provided by the drop is employed to enhance the HJ performance and the efficiency is maximum.For smaller distances, the subcritical Froude number defined by the tailwater level at the toe of the upstream check dam is smaller than the sequent Froude number of F i .Thus, the HJ is submerged, producing a dissipation efficiency lower than TI threshold conditions.Energy losses in a controlled gully reach with supercritical normal conditions.At the top, the Froude number profiles are depicted to illustrate the backwater computations methodology.F i represents the Froude number at the impact; F n the Froude number for normal conditions; F 1 the Froude number at the supercritical zone; F 2 the subcritical Froude number; FSP the free-surface profile; H i the energy dissipated at the impact zone; H dis the energy dissipated due to bed friction; H j the energy dissipated at the hydraulic jump; S the gully slope; L the check dam spacing; L i the impact length; L j the roller length; L • S the difference of total head between adjacent check dams; and z the effective check dam height.

Comparing the conceptual model with IBER
We contrasted the model performance with the results provided by the hydrodynamic bi-dimensional IBER model v1.9 (GEAMA, Instituto Flumen and CIMNE, 2012).IBER is freeware for flow simulation applications featuring a hydrodynamic module based on 2D-Saint-Venant equations and a finite-volume method used for the characterization of unsteady flows and hydraulic jump formations.IBER has been recently used and validated in different hydraulic applications (González-Aguirre et al., 2012;Bladé et al., 2014).The procedure for obtaining the IBER simulation was as follows: (i) definition of the geometry of the gully channel as a triangulated irregular network of 0.1 m cells; (ii) definition of the initial and boundary conditions; (iii) creation of the model mesh for the mathematical calculations; (iv) definition of the channel roughness; and (v) calculations and extraction of results as graphics (e.g.Froude number longitudinal profiles).The 3-D input geometry reproducing the initial and dam-filling channel geometry was obtained using Matlab ® scripts (The MathWorks TM Inc., Natick, MA, USA) specifically designed by the authors for this purpose.

An energy-based approach for estimating the efficiency
In this study, a new methodology based on energy considerations is proposed to estimate the performance of the check dam construction when compared with a non-intervention scenario.The difference in total head between adjacent check dams corresponds to the L • S product.In a gully without intervention, this energy is dissipated completely through bed friction (H dis-NC ) at the wetted perimeter of the cross section at a rate given by the NC regime.Assuming uniform flow, this energy can be expressed as a function of hydraulic and geometric variables: where τ n is the shear stress for normal conditions (NC), R n the hydraulic radius for NC and γ the water specific-weight.Therefore, all the energy needs to be dissipated by exerting a drag tension over the gully bed.In this case, the friction slope is equal to the gully slope at all the cross sections.
In contrast, after the check dam construction, energy losses appear at the impact zone (H i ) and at the hydraulic jump (H j ), reducing the dissipation through bed friction (H dis ): Hydrol.Earth Syst.Sci., 18, 1705-1721, 2014 www.hydrol-earth-syst-sci.net/18/1705/2014/ where S fm is the mean friction slope, τ m is the mean shear stress and R m the mean hydraulic radius in the controlled reach.Since H dis < H dis-NC , the mean values of the hydraulic variables responsible for the erosion processes (e.g. S fm and τ m ) are reduced in the post-intervention scenario at the expense of creating local energy losses at particular locations.If this dissipation occurs under non-protected conditions, an intensification of the erosive processes might take place.
The average slope friction S fm at the corrected reach was calculated according to Eq. ( 17): The energy losses at the impact zone were estimated using Eq. ( 18): where E c is the specific energy for critical depth over the spillway and E i , d i and u i are the specific energy, flow depth and velocity at the impact zone, respectively (Eq. 6).
The energy losses at the hydraulic jump H j were estimated using Eq. ( 10) for the classical HJ and the difference of total head between d 1 and d 2 for HJ on adverse slopes.For HJ submergence conditions, the total energy dissipation H t was estimated as the difference of total head between the critical regime at the spillway (E c + z) and the subcritical conditions imposed by tailwater level at the end of the HJ roller length due to the difficulties inherent in calculating H i and H j separately for such a complex hydraulic regime.
The efficiency E s of a check dam intervention was defined as the percentage of reduction of the mean value of the friction slope S fm at the corrected reach when compared to the friction slope value for NC prior to the intervention (S fm = S):

Analysis of factors controlling check dam efficiency
An assessment of the factors' influence (q, z, S and n) on check dam efficiency was carried out for initial and damfilling conditions.This study was conducted by executing the model with increasing steepness factor c in order to characterize the efficiency for a representative sample of the regime typologies (Sect.2.2) varying one factor at a time.The mean friction slope along the reach was calculated with Eq. ( 17) and the efficiency, using Eq. ( 19).
Further, we evaluated the limits of the optimal region of c values in terms of energy dissipation in three scenarios: initial conditions with n = 0.04 as well as dam-filling conditions with n = 0.04 and n = 0.06 These cases might represent the stages of the natural evolution of a gully after a check dam intervention: dams without silting, silted dams with clean channel and finally, silted dams with vegetation occupying partially the channel.For this purpose, firstly, we determined the steepness factor value defining the TI threshold (c lo ).Next, we explored the existence of an upper limit of total influence from which a distinct drop in efficiency takes place (c up ).Finally, we studied the relationship between these c limit values and the parameter z d c • S , a variable comprising the key factors in the design of check dam series.We termed this parameter the design number (DN).

Exploring the relationships with step-pool experiments
The efficiency defined in this study presents similarities with the friction factor f = 8 used in several flume experiments, such as those by Abrahams et al. (1995) or Comiti et al. (2009), to evaluate the resistance of stepped channels.To examine more deeply this relationship, we analysed the evolution of the efficiency, f , mean velocity u m and mean Froude number F m with increasing c and varying slopes for a particular case (z = 1 m , q = 0.1 m 2 s −1 , n = 0.04).The parameter f was estimated using the gully slope S and the d m and u m values applying Eqs.(2) to (4) with the calculated S fm (Eq.17).
In addition, we carried out an exploratory analysis using the predicted optimal c limit values in Sect.2.4.2 to evaluate if these curves might help to explain some trends found in step-pool field surveys, such as the declining tendency towards smaller c for increasing slopes (Zimmerman and Church, 2001;Chartrand et al., 2011).To this aim, we selected a range of input parameters typical of step-pool systems according to the values reported in the scientific literature: S (2-20 %), q (0.1-1 m 2 s −1 ), and z (0.3-0.9 m).A total of 2000 cases were randomly generated using a uniformly distributed function for each input parameter.For each combination of input values producing a particular DN, the c limit values (c lo and c up ) were calculated using the equation for dam-filling conditions and n = 0.04, a situation considered representative of mountain streams with gravel and cobbles.We generated a random c value (c rand ) uniformly distributed within the interval c lo -c up that was plotted against the slope.Finally, we determined the frequency histogram of step spacing (L) for this random data set for comparison purposes with the results provided by Curran and Wilcock (2005) on the probability functions of L in previous step-pool field studies.

Classification of flow regime typologies
Figure 3 shows the region where sub-and supercritical flows occur in normal conditions.Supercritical flow is the predominant regime, while subcritical flow only occurs on low slopes (usually below 5 %).In addition, high hydraulic roughness or low discharges induce lower regimes and thereby promote subcritical flows.
As for the flow regimes after restoration, Fig. 4 includes a graphical depiction of the different FSP types that may develop.For initial conditions, the HJ occurrence is guaranteed as a consequence of the elevation of the water profile behind the downstream check dam.The location of the HJ is dependent on the type of HJ control.Thus, if the normal conditions are subcritical, an HJ will occur close to the impact region as soon as the flow has dissipated the excess energy and has reached the supercritical sequent depth corresponding to the subcritical depth associated with that NC (IN-SUB-NC-PI).If the NC were subcritical enough to reach the sequent depth of the impact flow, the HJ would take place at the toe of the check dam (IN-SUB-NC-TI).As for the dam-filling situation, the value of the deposition slope S d is the key factor determining the type of regime which develops, since it controls the Froude regime associated with the modified NC.If the modified NC are supercritical, an HJ will not take place since there is no downstream control to cause it (F-SUP-NHJ).Therefore, it is necessary for S d to remain sufficiently low to impose subcritical conditions (F-SUB-PI).Depending on the conditions (q, S and n factors) and the design of the intervention (z and c) for an HJ to occur may require horizontal or even negative S d (F-D-PI and F-D-TI) to be established.

Comparison with the IBER model
The results of the comparison between the conceptual model and IBER for a representative sample of FSP types are shown in Fig. 5, as well as two examples of the input geometry.These curves represent the spatial evolution of the Froude number F of the flow along the restored reach for both methods.In the same figure, the constant F corresponding to NC (for the initial situation) and also modified NC (for a damfilling scenario) are shown, in order to facilitate the understanding of FSP evolution.
Overall, the model performed well, producing comparable results with IBER.The curves are mostly coincident, especially with regard to the subcritical region, but also, importantly, in the values of F at the HJ in the supercritical region.On the other hand, there are deviations with respect to the F i at the impact zone and the location of the HJ.Higher HJ and impact lengths were predicted with our model, explaining the visible offsets between both curves at the supercritical region.The biggest differences in check dam efficiency occur when the conceptual model overestimates the Froude number at the impact (IN-SUP-CN-PI and IN-SUP-D-PI), leading to an underestimation of efficiency by the conceptual model (22.3 % instead of 36.6 % and 28.2 % compared to 44.1 %).

Assessing check dam efficiency: the optimal dissipation interval
Figure 6 shows a sample of efficiency curves as a function of c (shown in logarithmic scale).As for initial conditions (Fig. 6a-d), the efficiency curves present an irregular sigmoid shape typically featuring a linear segment at low c, a parabolic segment in the middle and an asymptotic maximum from c ∼ 1 onwards.Increasing unitary discharges q produced lower maximum efficiencies and slightly lower c thresholds (Fig. 6a).The main impact of higher z was to increase both the maximum efficiency and the threshold c (Fig. 6b).High slopes produced a lower c limit to total influence (Fig. 6c), whereas there was little difference between the roughness coefficient results (Fig. 6d).In all cases, the maximum efficiency took place around c = 1.For dam-filling conditions (Fig. 6e-h), only cases where HJ occurred were considered.Therefore, relatively high values, c > 0.8, are shown.Here, the curves presented a "plateau" shape, with rapidly increasing efficiencies for c values under 1, a plateau interval of maximum efficiencies and, from a certain c limit value onwards, decreasing efficiencies.A distinct interval of maximum efficiency was found in each curve with a lower limit c lo defined by the threshold of total influence (F i = F 1 ) and a upper limit c up by the drop of efficiency.The range of steepness factor values between c lo and c up establishes an interval of optimal spacing in terms of energy dissipation (hereafter, the optimal dissipation interval ODI).Greater ODI widths were obtained with decreasing q and S and increasing z.The roughness coefficient has an impact on the magnitude of the maximum efficiency (higher n implied lower efficiency) but no clear trend is evident regarding the interval width.In addition, total influence is achieved (c lo ) at longer spacing with lower q and S and higher z.
Overall, it was clear that in both situations high maximum efficiencies were achieved (normally over 90 %), slightly lower for dam-filling conditions.Moreover, for initial conditions, the threshold of total influence c lo was lower (0.85-1) than in the dam-filling situation (0.9-1.5).
The response in check dam efficiency results from the ratio of dissipated energy and decreasing total energy LS, as illustrated in Fig. 7 for dam-filling conditions.Along the partial influence segment, initially the growth is linear since H t remains constant (almost equal to H i since H j is negligible at this stage) but LS decreases with closer spacing.The parabolic shape (dam-influence conditions) stems from both the reduction of LS and the increase of H j losses.The total influence threshold (c lo ) defines the onset of maximum efficiency (maximum H t and low LS) which is maintained for higher c since the reduction in LS compensates for the decrease in H t (submerged jump).For initial conditions, efficiency is kept at a maximum from the total influence threshold onwards (Fig. 6a-d), since the rising of water depth due to the check dam height provokes increasingly subcritical conditions as check dams become closer (Fig. 1a).However, in the dam-filling scenario, closer spacing than the total influence threshold implied decreasingly subcritical regimes as the flow approaches to the spillway (Fig. 1b).Thus, at a certain point (c up ) the reduction of LS is not enough to counteract the decline in the energy losses at the submerged jump and efficiency tends to decrease.
Figure 8 shows the relationship between the ODI limits and DN for a representative sample of cases (n = 60).These curves illustrate more clearly the tendencies found in the analysis of factors, since DN integrates in one single variable the impact of the main input parameters.We found a potential-law correlation between c and DN for the initial conditions (see fitted equations in Fig. 8).In the dam-filling situation, the data were adjusted to a mathematical relationship similar to the flow depth -specific energy equation in hydraulics, with a pair of c values for each DN.The equation fitted for n = 0.06, plotted slightly underneath the n = 0.04 curve, in other words the roughness favoured the increase in spacing.All the situations falling inside the ODI region will verify optimal dissipation conditions (maximum efficiency and total influence).DN ∼ 100 produced the minimum c value valid during the three sedimentation stages of the structure (triple point).A minimum DN around 20 was obtained, below which total influence was not achieved.

Exploring the relationships with step-pool units
The profiles of f , efficiency, u m and F m with increasing c are shown in Fig. 9. Prior to the check dam intervention (plane bed with no steps, i.e. c = 0), f is a surrogate of the roughness coefficient n and can be calculated determining d n and v n for normal conditions.This is the minimum value of the friction factor and corresponds to the roughness resistance f rough (f grain in Comiti et al., 2009).The energy losses at the impact and HJ introduce an additional component of resistance (f step similar to f spill according to the Fig. 6.Efficiency as a function of the steepness factor c for a range of conditions.n represents the Manning roughness coefficient; q the unitary discharge; S the bed slope; and z the effective check dam height.same authors) with f = f rough + f step .Thus, the efficiency as defined in Eq. ( 19) would represent the percentage of total resistance f due to step resistance f step .The f curve presents a "peak" morphology with a marked maximum (between c = 1.5-2,Fig. 9a) unlike the efficiency plot that features a "plateau" shape, although both maxima occur at the same c value (Fig. 9b).Similarly to efficiency, u m and F m curves show fairly constant values (minimum in this case, Fig. 9c, d) at the ODI.Our results, predicting efficiencies or f step f ratios around 90-95 %, are in line with estimations from authors such as Hayward (in Zimmerman and Church, 2001) and Comiti et al. (2009).
As for the field data set, we found a good agreement between the random c-S values predicted by the model and the c-S point cloud from step-pool studies (Fig. 10).The results also show a declining tendency of c values when slope increased, derived from the shape of the ODI curve (reduced interval width for low DN corresponding to high slopes).It is worth noting that the combination of input parameters S, q and z was totally random, whereas in nature usually certain relationships between geomorphological factors are verified, such as increasing discharges or decreasing step heights for decreasing slopes (Wohl and Grodek, 1994;Chin, 1999).On the other hand, no experimental data fell inside the region of small c (between 1 and 2) and low slope values (under 4 %) predicted by the model.This fact pointed to the existence of additional processes, not examined in the present approach.
Interestingly, we found a similar distribution to that obtained by Curran and Wilcock (2005) -Fig.10, on the top right corner -in terms of overall shape (exponential Fig. 7. Profiles of the main energy parameters as a function of the steepness factor c in dam-filling conditions (q = 0.1 m 2 s −1 ; S = 20 %; n = 0.04; z = 1 m)."Impact" refers to the impact losses as the predominant dissipation process and "HJ" to the interval where the energy losses at the hydraulic jump become increasingly relevant.c lo represents the lower threshold of ODI; c up the upper threshold of ODI; H t the total energy dissipated (H i + H j ); LS the difference in total head between adjacent dams; and ODI the optimal dissipation interval.function), maximum relative frequencies (close to 0.15), maximum location (∼ 5 m) and exclusion zone (no events for L < 1 m and few under 2 m, being controlled by the c up value of the ODI curves).Consequently, the model results satisfactorily predicted not only the overall region of occurrence of step-pools, but also the relative frequencies of occurrence for step spacing.

Validity of the model
Although the hypothesis underlying the definition of the model led to major simplifications, the relative simplicity in the calculations enabled us to focus on the concepts and an overall understanding of the hydraulics.The comparison with the simulations performed by IBER was satisfactory both in the longitudinal hydraulic profiles and in overall efficiency.This suggests that it can be applied with confidence for the assessment of free-surface profiles.Despite the fact that the model proved to be successful in providing insight into the key processes involved in gully control, it does not account for many phenomena occurring in more realistic situations, such as those derived from complex gully geometries, meanders, presence of pools, stones or weeds and sedimentation dynamics, all of which may have a strong impact on flow characteristics.
Further studies should undertake more sophisticated hydraulic analysis, combined with field observations, to improve our understanding of the modifications produced by check dams in complex channel morphologies.Laboratoryflume experiments would also be extraordinarily useful to test the conclusions of our theoretical approach, since they allow the evaluation of a varied range of geometry and flow conditions.

Implications of check dam efficiency on channel stability
The study of the Froude regime for NC showed that supercritical flows are predominant in gullies.These rapid supercritical flows are associated with high shear stress and friction slope values, which are the final cause of the erosion processes that lead to the development and growth of gully networks.Their spatial extent is only controlled by the backwater effect of the downstream dam, so that, if c is low, a large part of the gully remains under supercritical conditions.If this happens, the new situation is more erosive than that prior to the intervention since the waterfall at the check dam produces an accelerated flow more intense than the former situation until it evolves to NC.In addition, HJ occurs at an intermediate, unprotected section of the reach favouring further scouring.Consequently, badly designed check dam structures might in fact be more harmful than non-intervention (Heede, 1978).
Our results support the conclusion that f and efficiency illustrate the same phenomenon (maximum energy dissipation corresponding to maximum flow resistance), but f is more  It is apparent that the head-to-toe rule is the most conservative alternative.Moreover, there are significant differences among the different authors, in absolute values as well as in the trends.Similar recommendations can be found in technical books (Coppin and Richards, 1990;Morgan, 2005).
On the other hand, we do not exactly predict horizontal deposition slopes in real situations, but recommend assuming horizontal slopes in check dam design for stability purposes in the long term.Several experiences show that nonhorizontal sediment wedges occur in restored gullies, leading on some occasions to burying of the lower part of the upstream check dam (Heede, 1960;V. O. Polyakov, personal communication, 2013).The authors' experience on check dam interventions in agricultural areas where the head-totoe criterion was applied showed that deposition slopes are frequently close to horizontal and that burying might affect not much beyond the apron, offering an additional protection for high-flow conditions to come.Apparently stable configurations under low discharge conditions may turn out to be inadequate to prevent erosion during intense storms.Cycles of degradation and aggradation have been reported along the sediment wedge throughout the relatively long lifetime of the structure (e.g.Castillo et al., 2007) and if degradation is dominant at a certain point, complete failures may occur.
Our findings highlight that c ∼ 1 is required for an effective HJ control and dissipation.Thus, it would be preferable to increase the construction costs to operate with a wider safety margin, taking into account the real risk of undercutting and reactivation of gully erosion with eventual high discharges, likely to happen during the expected lifetime of the structures (usually designed using a return period of 25 years).For this reason, even negative deposition slopes (c > 1) have been proposed in gully restoration to avoid frequent check dam failure (near 40 % check dams collapsed in 2 years in agricultural areas in northern Ethiopia, Nyssen et al., 2004).In this work we have attempted to contribute to the understanding of the basic hydraulics underlying the patterns of energy dissipation not only for academics but also for stakeholders with a technical background involved in gully control.For both audiences, not only the design criterion is important, but also the reasons why it is so significant.
Another relevant outcome is the determination of an optimal DN (between 75 and 100, approximately) to achieve the ODI conditions since the selection of the check dam height has been an issue frequently neglected in technical publications.This means that z must be related to the discharge and slope conditions.According to the principles of minimum energy expenditure along drainage networks (Rodríguez-Iturbe et al., 1992), the product S • d c may not vary substantially since relationships in the form of Q 0.5 • S = constant, u = constant and d ∼ Q 0.5 are predicted, leading to q • S = constant.Thus, in practice, for the restoration of a particular gully, an approximately constant z might be applicable along the stream profile.For instance, for z = 1 m and fixing DN = 100 (triple point), d c • S = 0.01, implying S = 10 % for q = 0.1 m 2 s −1 or S = 2 % for q ∼ 1 m 2 s −1 , which seems like reasonable magnitudes for gullies in agricultural areas with return intervals of 25 years in Southern Spain (Castillo, 2012).
Additionally, the upper threshold of the ODI region (the c up branch of the curve) defines the minimum spacing in check dam design L min and can be related with the length of the potential scour region L s .Assuming an L s equal to the classical HJ length (Eq.11) as a first approximation, we obtained a fairly uniform average ratio L min L s = 1.47 for the cases studied.Likewise, using the lower threshold c lo , an average ratio L max L s = 4.5 was estimated, where L max stands for the largest spacing to meet the ODI requirements.This is in line with the recommendation of a minimum spacing above 2 • L s as proposed by VanDine (1996).

Conclusions
The conceptual model has combined in a single framework different previous approaches related to hydraulic processes to explain the basic flow modifications that check dams produce.Its comparison with the hydrodynamic IBER model produced comparable results.Among all the possible regimes, only those situations falling inside the optimal dissipation interval (ODI) verified maximum dissipation efficiency (effectiveness requirement) and HJ control (security requirement).For check dam design, we propose the selection of c = 1 and a design number ∼ 100, taking into account the different situations likely to occur during the projected lifetime of the control structures.The c-S correlation and the spacing distribution predicted by the model fitted well the step-pool data set reported in literature.This finding highlights the parallels between check dam interventions and step-pool units as systems controlled by similar dissipation processes.Further experimental studies, either field-based or at the flume scale, are necessary to test the conclusions of our theoretical approach and to address more realistic situations regarding the complexity of gully and flow features.

Fig. 1 .
Fig. 1.Sketch of the main hydraulic features and geometric variables after check dam construction, (a) initial conditions; (b) dam-filling situation.L: spacing between adjacent check dams; LS: difference in elevation between check dams; S: Gully slope; S d : deposition slope; z: effective height of the check dam.

Fig
Fig.2.Energy losses in a controlled gully reach with supercritical normal conditions.At the top, the Froude number profiles are depicted to illustrate the backwater computations methodology.F i represents the Froude number at the impact; F n the Froude number for normal conditions; F 1 the Froude number at the supercritical zone; F 2 the subcritical Froude number; FSP the free-surface profile; H i the energy dissipated at the impact zone; H dis the energy dissipated due to bed friction; H j the energy dissipated at the hydraulic jump; S the gully slope; L the check dam spacing; L i the impact length; L j the roller length; L • S the difference of total head between adjacent check dams; and z the effective check dam height.

Fig. 3 .
Fig. 3. Froude regime in normal conditions as a function of the unitary discharge q, gully slope S, and n Manning roughness coefficient.
Figure3shows the region where sub-and supercritical flows occur in normal conditions.Supercritical flow is the predominant regime, while subcritical flow only occurs on low slopes (usually below 5 %).In addition, high hydraulic roughness or low discharges induce lower regimes and thereby promote subcritical flows.As for the flow regimes after restoration, Fig.4includes a graphical depiction of the different FSP types that may develop.For initial conditions, the HJ occurrence is guaranteed as a consequence of the elevation of the water profile behind the downstream check dam.The location of the HJ is dependent on the type of HJ control.Thus, if the normal conditions are subcritical, an HJ will occur close to the impact region as soon as the flow has dissipated the excess energy and has reached the supercritical sequent depth corresponding to the subcritical depth associated with that NC (IN-SUB-NC-PI).If the NC were subcritical enough to reach the sequent depth of the impact flow, the HJ would take place at the toe of the check dam (IN-SUB-NC-TI).Finally, downstream control by dam influence can only begin when the adjacent dams are close enough to dominate over the NC (IN-SUB-D-PI and IN-SUB-D-TI).Supercritical NC and partial influence lead to an undesirable uncertainty in the location of the HJ, since it takes place at that point where the downstream check dam produces the subcritical control (IN-SUP-NC-PI or IN-SUP-D-PI).Closer spacing is required for an effective HJ control (IN-SUP-D-TI).As for the dam-filling situation, the value of the deposition slope S d is the key factor determining the type of regime which develops, since it controls the Froude regime associated with the modified NC.If the modified NC are supercritical, an HJ will not take place since there is no downstream control to cause it (F-SUP-NHJ).Therefore, it is necessary for S d to remain sufficiently low to impose subcritical

Fig. 4 .
Fig. 4. Classification of flow regime typologies in restored gullies.d c represents the critical depth; d n the normal depth; d s the depth over the spillway; d 1 the supercritical depth at the hydraulic jump; d 2 the subcritical depth at the hydraulic jump; D the check dam influence; F the dam-filling conditions; IN the initial conditions; NC the normal conditions influence; NHJ the no hydraulic jump; PI the partial influence; SUB the subcritical normal conditions; SUP the supercritical normal conditions; TI the total influence.

Fig. 5 .
Fig. 5. Input geometry examples in IBER and comparison of Froude number F curves between the model (grey) and IBER (black) for different flow regime typologies as a function of the distance from the downstream check dam.The grey dashed line corresponds to normal conditions and the dashed-dotted line to modified normal conditions for a dam-filling scenario.c represents the steepness factor; D the check dam influence; Ef so the check dam efficiency predicted in IBER; Ef sp the check dam efficiency predicted by the conceptual model; F the dam-filling conditions; IN the initial conditions; NC the normal conditions influence; n the Manning roughness coefficient; NHJ the no hydraulic jump; PI the partial influence; q the unitary discharge; SUB the subcritical normal conditions; S the gully slope; SUP the supercritical normal conditions; TI the total influence; and z the effective check dam height.

Fig. 8 .
Fig. 8. Relationships between the steepness factor c and the design number DN defining the limits of the optimal dissipation interval (ODI).The fitted curves are depicted on the top-right corner and their equations are shown inside the main plot.c lo represents the lower threshold of ODI; c up the upper threshold of ODI; d c the critical depth; n the Manning roughness coefficient; S the bed slope; and z the effective check dam height.

Fig. 9 .
Fig.9.Relationships between the steepness factor c and the friction factor f , efficiency, mean velocity u m and mean Froude number F m for a particular case (z = 1 m , q = 0.1 m 2 s −1 , n = 0.04).

Fig. 10 .
Fig. 10.Steepness factor c -slope relationship (main plot) and spacing frequency distribution (top-right corner) for the step-pool field data set -compiled by Zimmerman and Church (2001) and Curran and Wilcock (2005) -and the random data set generated by the conceptual model.

Fig. 11 .
Fig. 11.Comparison of steepness factor c recommendations as a function of gully slope for check dam design.