- Email: [email protected]
Adding sediment transport to the integrated hydrology model (InHM): Development and testing
Advances in Water Resources 29 (2006) 930–943 www.elsevier.com/locate/advwatres
Adding sediment transport to the integrated hydrology model (InHM): Development and testing Christopher S. Heppner a, Qihua Ran a, Joel E. VanderKwaak b, Keith Loague a
a,*
Department of Geological and Environmental Sciences, Stanford University, Stanford, CA 94305-2115, USA b 3DGeo Development Inc., Santa Clara, CA 95054, USA Received 31 May 2005; received in revised form 6 August 2005; accepted 13 August 2005 Available online 13 October 2005
Abstract The addition of a sediment transport algorithm to the comprehensive hydrologic-response model known as the Integrated Hydrology Model (InHM) is discussed. The first test of the sediment transport version of InHM is reported, using field data from a series of erosion experiments conducted by Gabet and Dunne [E.J. Gabet, T. Dunne, Sediment detachment by rain power, Water Resour Res 39 (2003) 1002]. The performance of the sediment transport component of InHM, in both calibration and validation phases, is judged to be successful, based upon quantitative statistical criteria. The ability to simulate sub-plot-scale interactions between surface water hydrology and rain-induced sediment transport with InHM is demonstrated. Sensitivity analysis reveals that the rainfall intensity exponent has a substantial impact on simulated sediment discharge. Future work, related to both testing InHM and much needed field experiments for model parameterization, is discussed. 2005 Elsevier Ltd. All rights reserved. Keywords: Hydrologic-response simulation; Sediment transport simulation; Rainsplash erosion; Hydraulic erosion; Plot-scale experiments; Model performance; Calibration; Validation; Sensitivity analysis
1. Introduction Understanding the relationships between near-surface hydrologic response and sediment transport is of interest to those in the fields of ecology, engineering, geomorphology, hydrology, and land-use planning. Dating back to the pioneering efforts of Robert Horton [22] there has been considerable work focused on identifying the linkages between hydrology and erosion (see the reviews of Merritt et al. [37], Gerits et al. [18], Kirkby [29], and Lane et al. [31]). During the International Hydrologic Decade, there was a heroic effort by the Agricultural Research Service (ARS) within the US Department of Agriculture (USDA) to make detailed *
Corresponding author. Tel.: +1 650 723 3090; fax: +1 650 725 0979. E-mail address: [email protected] (K. Loague). 0309-1708/$ - see front matter 2005 Elsevier Ltd. All rights reserved. doi:10.1016/j.advwatres.2005.08.003
measurements of both runoff and sediment yields for several unit-source watersheds [47]. More recently, many field-based studies have investigated hydrologically-driven erosion processes related to, for example, cumulative watershed effects [44,55], fire [23,25,45], and roads [36,53,54]. It is our opinion that it is currently not necessary (or even possible) to track, through observation or simulation, every drop of water or particle of sediment (even at relatively small scales) for problems with a landscape evolution component. We do, however, feel that it is very important to include the best characterization of hydrologic response, in terms of both process and dimension, when attempting to address transport-limited sediment transport [34]. Table 1 lists a few selected models that have been developed to consider both hydrologic response and sediment transport. Perusal of Table 1 shows that, in terms of process representation
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943
931
Table 1 Dimensionality of selected models that consider both hydrologic response and sediment transport, in comparison to the dimensionality of InHM Reference
Acronym
Year
Hydrologic response Subsurface
Surface
Channel
Overland
[30] [2] [27] [39] [21] [52] [7] [51]a [38] [20] This study
CREAMS ANSWERS IALLUVIAL STARS SEDICOUP KINEROS2 FLUVIAL-12 SHESED EUROSEM CH3D-SED InHM
1980 1980 1982 1987 1990 1990 1993 1996 1998 1999 2005
1D na na na na na na 3D na na 3D
1D 2D 1D 2D 1D 1D 2D 2D 1D 3D 2D
na 2D 1D 2D 1D 1D 2D 1D 1D 3D 2D
1D 2D na na na 1D na 2D 1D na 2D
a
Sediment transport
The hydrologic-response model driving SHESED is the European Hydrology System SHE [1].
and dimensionality, not all models are created equally. Table 2 compares the characteristics of two well-known hydrologic-response models that have sediment transport algorithms with the Integrated Hydrology Model (InHM), which is the focus of the effort reported herein. It is well known that hydrologically-driven erosion models must faithfully simulate surface responses (e.g., Horton overland flow and channel flow). However, in many settings (e.g., wetlands and steep convex hillslopes) the contributions of subsurface flow (e.g., Dunne overland flow and groundwater discharge) to the surface can be significant. Hydrologically-driven erosion models that neglect subsurface contributions will without question, for many situations, misrepresent near-surface hydrologic response and, subsequently, sediment transport. Inspection of Tables 1 and 2 clearly shows that InHM provides a firm hydrologic-response foundation upon which to build a sediment transport model. It is well beyond the scope of the effort reported here to explore fully all of the nuances encountered in the realm of physics-based sediment transport modeling. Rather, the primary objective of this study was to develop a sediment transport algorithm for InHM. The InHM sediment transport algorithm facilitates the simulation of both erosion and deposition. Built into the development is consideration for multiple species sediment transport, with detachment from rainsplash and/ or hydraulic erosion driven by spatially variable surface
water depths and velocities. The sediment transport version of InHM is tested using the field data from the rainsplash experiments (i.e., raindrop detachment with transport by flow [28]) of Gabet and Dunne [16].
2. Integrated hydrology model The comprehensive InHM was designed, in the spirit of the Freeze and Harlan [14] blueprint, to simulate quantitatively, in a fully-coupled approach, 3D variably-saturated flow and solute transport in porous media and macropores and 2D flow and solute transport over the surface and in open channels. The important and innovative characteristics of InHM, including no a priori assumption of a specific hydrologic-response mechanism (i.e., Horton or Dunne overland flow, discharge from subsurface stormflow or groundwater), are discussed elsewhere [33–35,48,49]. InHM has been successfully employed for catchment scale, event-based rainfall-runoff simulation [33,34,49] and for solute transport simulations [26,48]. It should be noted that Panday and Huyakorn [40] have reported a comprehensive physics-based hydrologic-response model that is comparable to InHM. The underlying InHM equations employed in this study are presented below. For the effort reported here, the macropore and solute transport components of
Table 2 Characteristics of the hydrology components of two well-known sediment transport models, in comparison to the characteristics of InHM Reference
[52] [51] This study
Acronym
KINEROS2 SHESED InHM
Year
1990 1996 2005
Subsurface
1D, U, A 1D, U, N/R; 2D, S, N/G 3D, U/S, N/R
Surface Overland
Channel
1D, N/K 2D, N/DW 2D, N/DW
1D, N/K 1D, N/DW 2D, N/DW
Coupling
Process
SQ SQ FO
H H, D, SS, GW H, D, SS, GW
U (unsaturated); S (saturated); U/S (unsaturated/saturated); N (numerical solution); A (analytical solution); R (Richards equation); G (groundwater flow equation); DW (diffusion wave); K (kinematic wave); SQ (sequential); FO (first-order); H (Horton overland flow); D (Dunne overland flow); SS (subsurface stormflow); GW (groundwater).
932
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943
InHM are not employed. For hydrologic response, only the basic equations of subsurface and surface water movement are provided; the reader requiring more detail is referred to VanderKwaak [48]. For sediment transport, each component of the algorithm is described. A simplified flow chart for the sediment transport version of InHM is given in Fig. 1.
tween the subsurface and surface continua (T1), / is porosity (–), Sw is water saturation (–), t is time (T), fa is the area fraction associated with each continuum (–), and fv is the volume fraction associated with each continuum. The Darcy flux is given by q g ~ q ¼ k rw w ~ krðw þ zÞ ð2Þ lw
2.1. Hydrologic response
where krw is the relative permeability (–), qw is the density of water (ML3), g is the gravitational acceleration (LT2), lw is the dynamic viscosity of water (ML1T1), ~ k is the intrinsic permeability vector (L2), z is the elevation head (L), and w is the pressure head (L). The transient flow of water on the land surface is estimated by the diffusion wave approximation of the depth-integrated shallow water equations. Such 2D
Subsurface flow is calculated by r f a~ q qb qe ¼ f v
o/S w ot
ð1Þ
where ~ q is the Darcy flux (LT1), qb is a specified rate source/sink (T1), qe is the rate of water exchange be-
Fig. 1. Simplified flow chart for InHM with sediment transport. Note, the superscript plus signs (+) refer to model components included in InHM without sediment transport.
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943
surface flow is conceptualized as a second continuum that interacts with the underlying variably saturated porous medium through a thin soil layer of thickness as (L). Assuming a negligible influence of inertial forces and a shallow depth of water, ws (L), the conservation of water on the land surface is described by ~ r wmobile qs as qb as qe ¼ s
oðS ws hs þ wstore Þ s ot
ð3Þ
where ~ qs is the surface water velocity (LT1), qb is the source/sink rate (i.e., rainfall/evaporation) (T1), qe is the surface–subsurface water exchange rate (T1), and hs is the average height of non-discretized surface microtopography (L); the superscripts mobile and store refer to water exceeding and held in depression storage, respectively [48]. Surface water velocities are calculated utilizing a two-dimensional form of the ManningÕs equation given by ~ qs ¼
ðwmobile Þ s ~ nU1=2
2=3
rðws þ zÞ
ð4Þ
where ~ n is the ManningÕs surface roughness tensor (TL1/3) and U is the friction (or energy) slope (–). The linkages between (1) and (3) are through first-order, physically-based flux relationships driven by pressure head gradients. The flow equations are discretized in space using the control volume finite-element method, which combines the geometric flexibility of finite elements with the local conservation characteristics of control volumes. Each coupled system of non-linear equations in an InHM simulation is solved implicitly using Newton iteration [48]. Efficient and robust iterative sparse matrix methods are used to solve the large sparse Jacobian systems [48]. InHM provides the spatio-temporally variable surface water depths and velocities required for simulating sediment transport, as described in the next section.
component and the terms in brackets are the sources and sinks of sediment concentration. The net erosion rate is calculated by esed ¼ es þ eh
ð6Þ 3 1
where es is the rainsplash erosion rate (L T ) and eh is the hydraulic erosion rate (L3T1). The rainsplash erosion rate in (6) is estimated as follows ( cf F ðws ÞðcosðhÞ rÞb A3D q > 0 es ¼ ð7Þ 0 q<0 where cf is the rainsplash coefficient ((TL1)b1), b is the rainfall intensity exponent (–), F(ws) is a factor representing the reduction in splash erosion with increasing surface water depths (–), h is the angle of the element from horizontal (–), r is the rainfall rate (LT1), A3D is the three-dimensional area associated with the node (L2), and q is the sum of rainfall rate and infiltration rate (LT1). The F(ws) function in (7) is estimated by [52] F ðws Þ ¼ expðcw ws Þ
eh ¼ ased ðC sedmax C sed Þ
ð9Þ
where ased is the hydraulic erosion transfer coefficient (L3T1) and C sedmax is the concentration at equilibrium transport capacity (L3L3). The C sedmax term in (9) is calculated by [10] C sedmax ¼ 0:05
where Csed is volumetric sediment concentration (L3L3), esed is the volumetric rate of soil erosion and/or deposition (L3T1), qbsj represents the rate of water added/removed via the jth boundary condition (L3T1), C sedj is the sediment concentration of the water added/removed via the jth boundary condition (L3L3), BC is the total number of boundary conditions, and Vw is the volume of water at the node (L3). On the right hand side of (5) the first term is the advective transport
ð8Þ
where cw is a parameter representing the damping effectiveness of surface water (L1). Prior to the accumulation of surface water and for deep flows F(ws) equals, respectively, one and zero. When multiple sediment species are considered, the total rainsplash erosion rate is apportioned to each species by scaling by the source fraction (ri) of that species. The hydraulic erosion rate in (6) is conceptualized as a kinetic transfer process [46] and is estimated as follows [52]
2.2. Sediment transport Depth-integrated multiple species sediment transport, restricted to the surface continuum, for each particle species is calculated by [3] ! BC X oC sed 1 b ¼ r ~ qs C sed þ esed þ qsj C sedj ð5Þ Vw ot j¼1
933
~ qs q3 g2 Dsed ws ðcsed 1Þ2
ð10Þ
where q* is the local shear velocity (LT1), Dsed is particle diameter (L), and csed is the particle specific gravity (i.e., the ratio of particle density to fluid density) (–). The shear velocity is calculated by pffiffiffiffiffiffiffiffiffiffiffiffiffi q ¼ gws rz ð11Þ where $z is the local surface slope (–), given by rz ¼ abs½tanðhÞ
ð12Þ
The ased term in (9) is defined for erosion (i.e., C sedmax > C sed Þ and deposition (i.e., C sedmax < C sed Þ as C sed > C sedmax 2vsed n ased ¼ A ð13Þ u~ qs ws v C sed < C sedmax
934
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943
where A is the area associated with the node (L2), vsed is the particle settling velocity (calculated using StokeÕs Law for the terminal velocity of spheres in calm fluid) (LT1), n is a coefficient [46] related to turbulence in the surface water due to raindrop impact (–), u is an erodibility coefficient [13] related to surface properties and texture (L1), and v is the particle erodibility factor (ranging from zero to one) (–). During erosion ased is a function of water depth, water velocity, surface properties (e.g., roughness or texture) and, in the case of multiple species, the relative source fraction and the armoring effect of any immobile species. During deposition ased is a function of particle settling velocity and the raindropinduced turbulence. Alternative transport capacity formulae [12,42,46] could easily be implemented with InHM. For each sediment species, denoted by the subscript i, the v term in (13) is calculated by vi ¼
ri 1k ð1 rimmobile Þ i
ð14Þ
where ri is the species source fraction (–), 1i is the cohesion coefficient (ranging from zero to one for soils with considerable or no cohesion, respectively) (–), rimmobile is the fraction of immobile species (–), and k is the armoring factor (–). First, erodibility (in a binary sense) is calculated with a critical shear stress test. A species is erodible when the critical shear stress (scr) for incipient motion is exceeded, and unerodible when this is not the case. The relationship for scr (ML1T2) used here is [32] scr ¼ ds ðcs cw ÞDsed
ð15Þ
where ds is a constant [19] (equal to 0.047 for most flow conditions) (–), cs is the specific weight of the sediment
(ML2T2), and cw is the specific weight of water (ML2T2). The shear stress is given by ð16Þ
s ¼ cw ws S
where S is the energy slope (–). The rimmobile term in (14) is calculated as rimmobile ¼
k mx X
ð17Þ
ri
i¼k sm
where ksm is the index of the smallest immobile particle species and kmx is the index of the largest particle species. The source fractions for all species must sum to one. The k term in (14), which represents the effect of immobile species on the erodibility of mobile species, is calculated by k ¼ max½0; ð1 arimmobile Þ
ð18Þ
where a is an empirical shadow factor (–), with a value greater than one. The non-linear equations describing transport of multiple sediment species on the land surface are solved implicitly using Newton iteration. The use of appropriate numerical methods ensures mass conservation and accurate monotone solutions [48].
3. The Gabet and Dunne experiments The careful rainsplash plot experiments reported by Gabet and Dunne [16,17] were selected to begin testing the sediment transport version of InHM. The Gabet and Dunne [16] fieldwork was carried out at the University of California Sedgwick Reserve, near Santa Barbara. The information required to excite the InHM sediment transport algorithm is given in Table 3. The
Table 3 Parameters and variables required for simulating sediment transport with InHM Symbol
Definition (units)
Source, this study
Dsed qsed r $z r ~ qs ws S vsed scr cw b n u f a cf
Species average grain diameter (L) Particle density (ML3) Species source fraction (–) Surface slope (–) Rainfall rate (LT1) Water velocity (LT1) Water depth (L) Energy slope (–) Species settling velocity (LT1) Species critical shear stress (ML1T2) Rainsplash depth dampening factor (L1) Rainfall intensity exponent (–) Rain-induced turbulence coefficient (–) Surface erodibility coefficient (–) Sediment cohesion coefficient (–) Shadow factor (–) Rainsplash coefficient ((TL1)b1)
0.00001 m [16] 2650 kg m3 1.0 Table 4 [16] Table 4 [16] Flow simulation, InHM Flow simulation, InHM Flow simulation, InHM StokeÕs Law (15) 600 m1 [52] 1.6 [16] 0.25 [46] NAa NAa NAa Calibration, Table 4
a
Not required for this study.
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943
precious data provided by Gabet and Dunne [16] (i.e., steady-state runoff and sediment concentrations) does not include spatio-temporally variable information that is typically employed to check (compare against) InHM simulated near-surface hydrologic response [33,35,49]. The missing data include (i) surface water depths, (ii) water table depths, (iii) total head values, and (iv) soilwater content values. However, the site information reported by Gabet and Dunne [16] is sufficient to parameterize InHM for the simulations reported herein. For example, for the different plots the slopes were between 4 and 17, the groundcover (vegetation) ranged between 18% and 94%, the soils were silty clay loams, and hydrologic response was only by the Horton overland flow mechanism. Details of the experimental set up for the sprinkler experiments (e.g., simulator operation, nozzle size) are provided by Gabet and Dunne [16].
4. Methods A split sample (calibration/validation) approach, grounded to the Gabet and Dunne [16] experimental data, was used in this study to evaluate the sediment transport version of InHM. It is important to emphasize that the simulations reported here consider erosion for a single species by rainsplash only, despite the ability of InHM to address more complex problems. The chronology of the model testing effort employed in this study can be summarized by four phases: (i) boundaryvalue problem (BVP) development, (ii) hydrologicresponse simulations (calibration), (iii) sediment transport simulations (calibration), and (iv) sediment transport simulations (validation). The four phases are described below. 4.1. Phase I, Boundary-value problem Each of the Gabet and Dunne [16] plots had the same surface dimensions (6.0 m by 2.4 m). Fig. 2 shows the finite-element mesh used for each of the BVPs in this study. The BVPs for the different experimental plots had different slopes, as summarized in Table 4. The dimensions of the mesh are 6.0 m long, 2.4 m wide, and 3.0 m deep. The vertical nodal spacing (Dz) in the mesh varies from 0.05 to 0.25 m; the horizontal nodal spacing (Dx, Dy) is 0.2 m. The total number of nodes/ elements in the mesh are, respectively, 9009 and 16,128. With reference to the lettered coordinate locations identified in Fig. 2, the boundary conditions for each face of each BVP are as follows: impermeable (A–E–G–C, B–F–H–D, C–G–H–D), local sink (A–B– D–C, E–F–H–G), and flux (A–E–F–B). Relative to the near-surface soil-hydraulic properties, the entire subsurface domain was considered to be homogeneous and isotropic; the porosity and compressibility were set at 0.46
935
Z
E
Y
X
F
A
H B
G
C D 1m Fig. 2. Boundary-value problem used to represent the Gabet and Dunne [16] experiments for the InHM simulations of hydrologic response and sediment transport.
(i.e., average value for silty clay loam soils [6,43]) and 1 · 107 ms2 kg1 (i.e., mid range value for clay [15]), respectively. The soil-water retention and permeability functions, based upon the van Genuchten [50] approach and characteristic information on silty clay loam soils from Carsel and Parrish [6] and Rawls et al. [43], are given in Fig. 3. For each InHM simulation the mobile water depth (i.e., depression storage) was set to 0.0005 m and the height of microtopography was set at 0.01 m. The slope and groundcover percentage for each of the experiments included in this study are given in Table 4. For each experiment, the ManningÕs roughness coefficient was calculated using [16] n ¼ 0:053e2:7Cv
ð19Þ
where Cv is the groundcover percentage (expressed as a decimal, 70% = 0.70). The entire subsurface domain was unsaturated throughout each of the simulations reported herein (i.e., no water table within 3 m of the surface), with surface runoff resulting from Horton overland flow. The initial conditions, taken to be the same for each experiment, were set by upstream and downstream local head boundary conditions so that the initial water table was 5.0 m below the land surface. The validity of the assumed initial conditions (i.e., water table position) was
936
Table 4 Summary of field experiment characteristics, parameterization, and results from hydrologic response/sediment transport simulations with InHM Rainfall ratea,b (105 m/s)
Slopea ()
Groundcovera (%)
ManningÕs roughness coefficient (–)
Calibrated saturated hydraulic conductivity (105 m/s)
Observed steady-state water dischargea (104 m3/s)
Simulated steady-state water discharge (104 m3/s)
Best-fit rainsplash coefficient (–)
Estimated rainsplash coefficient (–)
Observed steady-state sediment dischargea (104 kg/s)
Simulated steady-state sediment discharge (104 kg/s)
2/1 3/1 5/2 6/2 7/2 8/2 12/3 14/4 15/4 16/5 17/5 18/5 19/6 20/6 21/6 22/7 23/7 24/7 25/8 26/8 27/8 28/8 29/8 30/8 31/8 32/9 33/9 35/9 36/10 37/10 38/10 39/10 40/10
3.47c 1.81 1.67 3.36 1.72 3.25 3.31 1.39 1.67 2.06 1.53 1.53 1.25 3.06 3.47 2.92 2.69 2.03 1.83 2.28 3.06 1.89 2.81 2.03 3.22 1.86 3.08 3.14 1.86 3.69 3.33 3.39 3.39
10 10 17 17 17 17 14 4 4 17 17 17 13 13 13 9 9 9 13 13 13 13 13 13 13 13 13 13 5 5 5 5 5
94 88 72 72 72 72 68 98 98 82 74 66 70 70 56 74 64 68 64 62 62 39 35 32 32 34 34 34 86 76 52 45 46
0.671 0.570 0.370 0.370 0.370 0.370 0.332 0.747 0.747 0.485 0.391 0.315 0.351 0.351 0.240 0.391 0.298 0.332 0.298 0.283 0.283 0.152 0.136 0.126 0.126 0.133 0.133 0.133 0.540 0.413 0.216 0.179 0.184
1.10d 1.10 0.66 0.66 0.66 0.66 1.30 1.13 1.13 0.10 0.10 0.10 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 0.70 0.70 0.70 0.80 0.80 0.80 0.80 0.80
3.36e 1.10 1.44 3.34 1.44 3.34 2.42 0.36 0.65 2.50 1.99 1.99 0.36 2.40 2.86 2.50 2.50 1.66 1.25 1.49 2.86 1.44 2.86 1.44 3.36 1.61 3.34 3.34 1.66 4.01 3.34 4.01 3.34
3.36 0.96 1.37 3.78 1.45 3.62 2.79 0.34 0.73 2.73 1.98 1.98 0.33 2.88 3.49 2.73 2.40 1.44 1.14 1.78 2.89 1.24 2.53 1.44 3.16 1.63 3.41 3.48 1.49 4.14 3.65 3.74 3.74
1.09 0.16h 0.88 2.27 0.80 3.14h 1.55 0.09h 0.12h 0.89 0.92 0.82 0.81 2.84h 1.73h 2.29 1.66 0.46h 1.43h 1.05 4.93 1.63h 5.87 1.90 4.89h 1.54 5.40 7.86h 0.67 3.49h 5.84 6.10h 5.35
0.61
1.34f 0.11 1.18 7.37 1.14 9.84 4.05 0.02 0.04 1.80 1.57 1.57 0.34 7.01 6.66 4.09 3.56 0.51 1.86 2.03 14.51 3.51 23.62 4.40 25.17 3.73 28.62 42.37 0.53 5.73 16.68 21.04 18.25
0.75g 0.11i 1.43g 3.44g 1.51g 9.85i 6.75g 0.02i 0.04i 3.08g 3.65g 5.28g 1.04g 7.00i 6.63i 3.91g 6.09g 0.60i 1.84i 5.91g 8.87g 3.05i 20.33g 12.31g 25.31i 12.60g 27.22g 42.29i 0.84g 5.69i 11.16g 21.04i 14.52g
a b c d e f g h i
Gabet and Dunne [16]. Irrigation rates from Gabet and Dunne [16] are taken as rainfall rates. 3.47 · 105. 1.10 · 105. 3.36 · 104. 1.34 · 104. Validation simulation. Value used in calibration. Calibration simulation.
2.28 2.28 2.28 2.58
1.52 2.13 2.74 2.43
2.13 2.89
3.04 3.04 5.09 5.32 5.17 5.17 1.22 3.80 4.26
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943
Experiment number/ plot numbera
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943 0.5 Soil-water content Relative permeability
0.8
0.3
0.6
0.2
0.4
0.1
0.2
0
Relative permeability (-)
Soil-water content (-)
4.3. Phase III, Sediment transport simulations— calibration
1
0.4
0 0
-1
-2
-3
-4
937
-5
Pressure head (m) Fig. 3. Characteristic curves used to represent the Gabet and Dunne [16] experiments for the InHM simulations of hydrologic response and sediment transport.
confirmed through steady-state hillslope-scale long-profile InHM simulations. 4.2. Phase II, Hydrologic-response simulations The focus of this effort was not to test the hydrologicresponse capability of InHM, which has already been reported for more complete data sets [33,35,49]. Hydrologic-response simulations with InHM were conducted for each of the 44 experiments reported by Gabet and Dunne [16]. For each of the 11 plots the saturated hydraulic conductivity was adjusted, until the best match was achieved between the observed and simulated steady-state discharges. Each of the InHM simulations successfully captured the Horton overland flow response reported by Gabet and Dunne [16]. Eleven of the Gabet and Dunne [16] experiments were removed from further consideration in this study for two reasons: • It was not possible to identify experiments 4, 41, 42, 43, and 44 (from Table 1 of Gabet and Dunne [16]) in Fig. 5 of Gabet and Dunne [16]. Therefore, these five experiments were not included in this study. • It was not possible to obtain satisfactory matches between observed and simulated discharges for experiments 1, 9, 10, 11, 13, and 34 (from Gabet and Dunne [16]) during the hydrology calibration. Greater runoff was observed for these experiments (from the same plots) for lower irrigation intensities, which is most likely due to variations (undocumented) in initial conditions. Therefore, these six experiments were not included in this study to avoid further calibration. Removing 11 of the 44 Gabet and Dunne [16] experiments, and subsequently one of the 11 plots, yields 33 experiments, conducted on ten plots, for further consideration in this study.
The sediment discharge for each experiment was calculated by multiplying the InHM simulated concentration and water discharge. Four steps were employed to calibrate the sediment transport version of InHM using the Gabet and Dunne [16] experimental data. The first step was to rerun the calibrated InHM (Phase II) for each of the 33 experiments listed in Table 4. For each InHM simulation, the rainsplash coefficient in (7) was adjusted until a best-fit was achieved between the observed and simulated sediment discharges. The second step was to develop a relationship (multiple regression) between the best-fit rainsplash coefficients and the groundcovers and slopes reported for the ten plots. The third step was to select 13 of the 33 experiments, representing the entire range of slopes and groundcovers, as calibration experiments for the first part of a split sample test (see Table 4). The fourth and final step was to develop a relationship (linear regression) for the 13 experiments selected in step 3 to estimate rainsplash coefficients from whichever characteristic (groundcover or slope) in step 2 provides the greatest correlation to the best-fit rainsplash coefficients. 4.4. Phase IV, Sediment transport simulations— validation The InHM sediment transport validation simulations, for the 20 experiments in the split sample (see Table 4), build upon the hydrologic response and sediment transport calibration simulations in Phases II and III, respectively. The rainsplash coefficients used in the validation simulations were estimated based upon the relationship established in step 4 of Phase III. 5. Model performance evaluation The two measures of model performance [24,41] used to test InHM are modeling efficiency (EF) and the relative root mean square error (RRMSE). The mathematical expressions for EF and RRMSE are " #, n n n X X X 2 2 2 EF ¼ ðOi OÞ ðP i Oi Þ ðOi OÞ i¼1
i¼1
" RRMSE ¼
n X 2 ðP i Oi Þ =n i¼1
i¼1
#0:5
ð20Þ 100 O
ð21Þ
where Pi are the predicted values, Oi are the observed values, n is the number of samples, and O is the mean of the observed data. When the observed and simulated values are identical the EF and RRMSE statistics are 1.0 and 0.0, respectively.
938
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943
6. Results Simulated steady-state water discharge (10-4 m3/s)
5
6.1. Hydrologic-response simulations The rainfall rates, ManningÕs roughness coefficients, saturated hydraulic conductivity values, and observed/ simulated water discharges for each of the 33 experiments included in this study are provided in Table 4. Results from the InHM hydrologic-response calibration simulations are summarized in Table 5 and Fig. 4. The model performance statistics given in Table 5 and the observed versus simulated discharges illustrated in Fig. 4 show that InHM was able, as should be expected in a calibration mode, to capture the steady-state hydrologic response (i.e., Horton overland flow) for each of the 33 experiments (i.e., EF of 0.94).
Plot 1 Plot 2 Plot 3 Plot 4 Plot 5 Plot 6 Plot 7 Plot 8 Plot 9 Plot 10
4
3
2
1
0 0
6.2. Sediment transport simulations—calibration
1
2
3
4 -4
5 3
Observed steady-state water discharge (10 m /s)
The best-fit rainsplash coefficients and observed/simulated steady-state sediment discharges for each of the 33 experiments and the subset of 13 calibration experiments are given in Table 4. Results from the InHM sediment transport calibration simulations are summarized in Figs. 5 and 6 and Table 5. Fig. 5a shows the relationship between the best-fit and multiple regression estimates of the rainsplash coefficient for the 33 experiments. The correlation between the best-fit rainsplash coefficients and the slopes and groundcovers for
Fig. 4. InHM simulated hydrologic response results for the Gabet and Dunne [16] experiments. The solid line is 1 to 1; the dashed lines are plus/minus 10% of the 1 to 1 line.
the 33 experiments are shown in Fig. 5b and c, respectively. By comparing the R2 values in Fig. 5 it is obvious, for this study, that groundcover provides a greater correlation to the best-fit rainsplash coefficient than slope. Fig. 6 shows the relationship between the best-fit rainsplash coefficient and groundcover for the subset of 13
Table 5 Performance results for InHM simulations of hydrologic response and sediment transport for the Gabet and Dunne [16] experiments Process
Number of experiments/ simulationsa
Calibration EF (–)
RRMSE (%)
Hydrologic response Sediment transport Sediment transport Sediment transport
33 33 13 20
0.94 1.00 1.00
10.6 1.5 1.5
Validation EF (–)
RRMSE (%)
0.78
54.8
EF is modeling efficiency, (20); RRMSE is the relative root mean square error(21). a See Table 4.
(a) 10
(b) 10
2
R = 0.51
Cfbest (-)
(c) 10
2
R = 0.02
R = 0.45
8
8
8
6
6
6
4
4
4
2
2
2
0
0 0
2
4
2
6
Cfreg (-)
8
10
0 0
5
10
15
Slope (degrees)
20
0
25
50
75
100
Groundcover (%)
Fig. 5. Best-fit rainsplash coefficient relationships from 33 InHM simulations for the Gabet and Dunne [16] experiments. (a) Best-fit to multiple regression (slope and groundcover). (b) Best-fit to slope. (c) Best-fit to groundcover. The solid line in (a) is 1 to 1; the dashed lines in (b) and (c) are best-fit regression lines.
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943
calibration experiments. The model performance statistics given in Table 5 show that InHM performed well, in a calibration mode, for the entire set of experiments and for the subset of 13 experiments (i.e., EF of 1.00). 6.3. Sediment transport simulations—validation The regression estimated rainsplash coefficients and observed/simulated steady-state sediment discharges for the 20 validation experiments are given in Table 4. Results from the InHM sediment transport validation simulations are summarized in Table 5 and Fig. 7. The model performance statistics given in Table 5 show that
Best-fit rainsplash coefficient (-)
9 8 2
7
R = 0.54
6 5 4 3 2 1 0 20
30
40
50
60
70
80
90
100
Groundcover (%) Fig. 6. Best-fit regression relationship (solid line) for the rainsplash coefficient relative to groundcover based upon InHM simulations for 13 of the Gabet and Dunne [16] experiments (i.e., the subset of calibration experiments).
Simulated steady-state sediment discharge (10-4 kg/s)
30
25
20
15
10
5
0 0
5
10
15
20
25
30
Observed steady-state sediment discharge (10-4 kg/s) Fig. 7. Comparison of observed and InHM simulated sediment discharge rates for the Gabet and Dunne [16] experiments. The solid line is 1 to 1.
939
InHM performed reasonably well in a validation mode for the subset of 20 experiments (i.e., EF of 0.78), especially considering the relatively simple method of estimating the rainsplash coefficient. Fig. 7 shows that sediment discharges are overestimated and underestimated, respectively, for the smaller and larger observed values.
7. Discussion 7.1. InHM as a diagnostic or concept-development tool Fig. 8 illustrates the type of surface water (depth, discharge, runoff rate) and sediment (concentration, discharge, erosion rate) information that can be generated with InHM. The steady-state snapshot results shown in Fig. 8 are from a longitudinal transect of 33 equallyspaced monitoring (simulated) points for experiment 21 (see Table 4). Relative to the capabilities of the model, experiment 21 is relatively simple application of InHM in that only sheet flow is considered (i.e., no flow convergence related to topography; essentially a 1D BVP) and that the only hydrologic response and sediment detachment processes are, respectively, Horton overland flow and rainsplash. The two principal variables extracted from the experiment 21 results are surface water depth and sediment concentration. The surface water velocity and discharge values shown in Fig. 8 were calculated from the depth information using ManningÕs equation. Inspection of Fig. 8 shows that discharge increases linearly downslope, which is expected as the rainfall rate for this experiment (i.e., all experiments in this study) exceeds the saturated hydraulic conductivity, Ks (LT1), of the soil (i.e., Horton by definition). The difference in discharge per unit area from one point to the next point downslope equals the difference between the rainfall rate and the saturated hydraulic conductivity. The relationship between saturated hydraulic conductivity and permeability (used in InHM) is Ks = kqwg/lw. In Fig. 8 the longitudinal sediment concentration profile decreases in magnitude with distance downslope. There are two reasons for this decrease: (i) rainsplash erosion rate is a negative function of the water depth, which increases downslope; (ii) transported sediment becomes diluted in the greater volumes of water as total discharge increases downslope. Hydraulic erosion is minimal in this study [16]. If hydraulic erosion were not minimal, the increased downslope depths and velocities would contribute to greater concentrations. The sediment discharge rates shown in Fig. 8 are the product of water discharge and sediment concentration. The erosion rates shown in Fig. 8 were calculated as the difference in sediment discharge from one monitoring point to the next downstream monitoring point divided by
0.003
2
0.0004
7
Water depth Discharge Runoff rate Sediment concentration Sediment discharge Erosion rate
6
5
4 0.002 3
2 0.001
0.0003
0.0002
0.0001
3
Water depth (m), Sediment discharge (kg/s)
0.004
1
0
0 0
1
2
3
4
5
Discharge (m /s), Runoff rate (m/s), Erosion rate (kg/m s)
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943
Sediment concentration (kg/m3)
940
0
6
Distance from downstream boundary (m) Fig. 8. InHM simulation results (i.e., water depth, discharge, runoff rate, sediment concentration, sediment discharge, and erosion rate) along a longitudinal transect of plot number 6 for experiment number 21 (see Table 4).
the area. Similarly, the runoff rate was calculated as the difference in water discharge per unit width from one monitoring point to the next downstream monitoring point. Inspection of Fig. 8 shows, for experiment 21, that the erosion rate decreases downslope even though the total sediment discharge increases. The slight spatial oscillations in the runoff rate and erosion rate curves can be traced to variations in the elemental geometry and nodal connectivities for the nodes along the centerline of the plot, resulting in small variations in calculated discharge. The diagnostic and concept-development capabilities of InHM are most useful for hydrologic response/sediment transport BVPs that are more challenging than the one discussed here. Related to hydrologic response, rather than just focusing on a Horton mechanism response, the contributions from both the Horton and Dunne mechanisms, as-well-as subsurface stormflow and groundwater, can all be considered simultaneously with InHM. Related to sediment transport, the contributions from rainsplash (as reported here) and hydraulic
erosion can be considered individually or simultaneously with InHM. 7.2. Sensitivity analysis The sensitivity of the InHM simulated hydrologic response/sediment transport results reported herein, for reasonable uncertainties in the Gabet and Dunne [16] experimental data, was investigated for experiment 21 by individually increasing/decreasing (i.e., relatively small plus/minus 10% changes) the values of six parameters and conducting 12 new simulations. The six parameters considered in the sensitivity analysis are (i) the rainfall rate, (ii) the saturated hydraulic conductivity, (iii) the ManningÕs roughness coefficient, (iv) the rainsplash coefficient, (v) the rainfall intensity exponent, and (vi) the rainsplash depth dampening factor. The base case (experiment 21) parameter values are given in Tables 3 and 4. Table 6 presents the results from the 12 new simulations as percent differences relative to the base-case simulation for both water and sediment
Table 6 Sensitivity analysis results for InHM simulations of hydrologic response and sediment transport for experiment 21 Parameter
Percent differences (relative to the base case) for increased/decreased parameter values Water discharge
Rainfall rate, r Saturated hydraulic conductivity, Ks ManningÕs roughness coefficient, n Rainsplash coefficient, cf Rainfall intensity exponent, b Rainsplash depth dampening factor, cw
Sediment discharge
10% increase
10% decrease
10% increase
10% decrease
14 4 0 na na na
14 4 0 na na na
12 3 7 10 81 14
12 3 8 10 417 16
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943
discharges. Perusal of the results in Table 6 leads to the following comments: • Increasing/decreasing the rainfall rate resulted in relatively small and symmetrical changes in the simulated water and sediment discharges (i.e., the 10% increase and decrease in r lead to 14% greater and smaller, respectively, water discharge and 12% greater and smaller, respectively, sediment discharges). The disproportionate rise in water discharge is related to the effect of depression storage, which makes the actual depth of runoff less than the total depth. As depression storage is full in both the base case and increased rainfall rate scenarios, the additional rain in the increased rainfall rate scenario adds only to the mobile water. The sensitivity of sediment discharge to the rainfall rate is complicated by negative feedbacks (e.g., higher rainfall rates would seem to produce higher rainsplash erosion rates, but increased water depths reduce the rainsplash effect). • Increasing/decreasing the saturated hydraulic conductivity resulted in relatively small and symmetrical changes in the simulated water and sediment discharges (i.e., the 10% increase and decrease in Ks lead to 4% smaller and greater, respectively, water discharges and 3% smaller and greater, respectively, sediment discharges). In general, greater infiltration rates and smaller discharges (for the Horton case) result from increased Ks values. There is a negative feedback between water depth and the rainsplash erosion rate. • Increasing/decreasing the ManningÕs roughness coefficient had no effect on the simulated water discharge but resulted in relatively small and symmetrical changes in the simulated sediment discharges (i.e., the 10% increase and decrease in n lead to 7% smaller and 8% greater, respectively, sediment discharges). Differences in the depth–velocity relationship coupled with the depth dependence of the rainsplash erosion rate produce the changes in the sediment discharges. • Increasing/decreasing the rainsplash coefficient resulted in directly proportional changes in the simulated sediment discharge (i.e., the 10% increase and decrease in cf lead to 10% greater and smaller sediment discharges, respectively). • Increasing/decreasing the rainfall intensity exponent resulted in large and non-symmetrical changes in the simulated sediment discharge (i.e., the 10% increase and decrease in b lead to 81% smaller and 417% greater, respectively, sediment discharges). Small shifts in b produce large changes in the rainsplash erosion rate [see (7)], emphasizing the need for continued efforts to thoroughly quantify the relationship between rainsplash detachment rate and rainfall intensity [16].
941
• Increasing/decreasing the rainsplash depth dampening factor resulted in relatively small and symmetrical changes in the simulated sediment discharge (i.e., the 10% increase and decrease in cw lead to 14% smaller and 16% greater, respectively, sediment discharges). Again, the negative feedback between water depth and the rainsplash erosion rate is apparent.
7.3. Future work The use of the Gabet and Dunne [16] experiments to test the rainsplash erosion component of the sediment transport version of InHM is but a first step. Future efforts will test the hydraulic erosion component of InHM. It will also be important to test the sediment transport version of InHM for combined hillslope and channel systems. The first and second authors, respectively, are currently using InHM for long-term landscape evolution problems at the catchment (i.e., the 0.1 km2 R-5 catchment near Chickasha, OK) and regional (i.e., the 117 km2 Hawaiian island of KahoÕolawe) scales. The measure and model approach pioneered by Robert Horton (see [4,5]) is the best protocol in 21st century hydrology. In general, at the catchment scale, computer resources do not limit the use of InHM for multidimensional non-linear problems. There are, however, significant information shortfalls that must be met before process-based hydrologic response/landscape evolution simulation with InHM can be considered routine, even in the research arena. The information shortfall in hydrology/geomorphology is not news to those practicing the measure and model approach, especially those who have focused on the assessment of process-based models. There have been several calls for detailed long-duration field monitoring programs at selected experimental catchments (see [8,9,11]). The information needed to rigorously test and implement InHM for hydrologically-driven landscape evolution problems includes the data needed to excite both the hydrology [33,35,49] and sediment transport (see Table 3 and Fig. 1) components of the model, as-well-as observations to compare against the simulated results. It is also important to understand the interactions between the different types of data (e.g., groundcover and cohesion). The information needs to be acquired from field experiments that consider the spatio-temporal variability of near-surface parameters (which may not be time invariant [34]) and hydrologic/sediment transport responses.
8. Summary If the objective of a modeling study is to estimate surface erosion based upon simulated runoff, the full range
942
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943
of near-surface hydrologic-response possibilities needs to be considered. InHM, as developed by VanderKwaak [48], is a comprehensive physics-based hydrologic-response model, without an a priori assumption for a specific hydrologic-response mechanism (e.g., Horton overland flow), providing a firm foundation to simulate sediment transport. We have developed a spatially distributed sediment transport algorithm for InHM and done preliminary testing of the rainsplash component. The sediment transport version of InHM has been designed (like the original version of InHM [48]) to meld measurements with modeling. We have demonstrated the spatial capabilities of the new model (see Fig. 8) and stress the need for additional spatio-temporal measurements to move forward with process-based hydrologically-driven erosion research.
[12] [13]
[14] [15] [16] [17]
[18]
[19]
Acknowledgements [20]
The effort reported here was supported by National Science Foundation grant EAR-0438749. Our conversations with Bob Street on sediment transport modeling have been very helpful. The thoughtful comments from Steve Burges and four anonymous reviewers on an earlier version of this paper are greatly appreciated.
[21]
[22]
[23]
References [24] [1] Abbott MB, Bathurst JC, Cunge JA, OÕConnell PE, Rasmussen JL. An introduction to the European Hydrology System SHE, 2: Structure of a physically-based, distributed modeling system. J Hydrol 1986;87:61–77. [2] Beasley DB, Huggins LF, Monke EJ. ANSWERS: a model for watershed planning. Trans Am Soc Agric Eng 1980;23:938–44. [3] Bennett JP. Concepts of mathematical modeling of sediment yield. Water Resour Res 1974;10:485–92. [4] Beven K, Robert E. HortonÕs perceptual model of infiltration processes. Hydrol Process 2004;18:3447–60. [5] Beven K, Robert E. Horton and the abrupt rises of ground water. Hydrol Process 2004;18:3687–96. [6] Carsel RF, Parrish RS. Developing joint probability distributions of soil water retention characteristics. Water Resour Res 1988;24:755–69. [7] Chang HH. Fluvial-12, Mathematical model for erodible channels (userÕs manual). San Diego, CA: San Diego State University; 1993. [8] Dunne T. Relation of field studies and modeling in the prediction of storm runoff. J Hydrol 1983;65:25–48. [9] Dunne T. Wolman lecture: hydrologic science . . . in landscapes . . . on a planet . . . in the future. Hydrologic sciences: taking stock and looking ahead. Washington, DC: National Academy Press; 1998. p. 10–43. [10] Engelund F, Hansen E. A monograph on sediment transport in alluvial streams. Copenhagen, Denmark: Technical University of Denmark; 1967. [11] Entekhabi D, Asrar GR, Betts AK, Beven KJ, Bras RL, Duffy CJ, et al. An agenda for land surface hydrology research and a call
[25]
[26]
[27]
[28] [29]
[30]
[31]
[32]
for the second international hydrological decade. Bull Am Meteorol Soc 1999;80:2043–58. Ferro V. Evaluating overland flow sediment transport capacity. Hydrol Process 1998;12:1895–910. Foster GR. Modeling the erosion process. In: Haan, Johnson, Brakensiek, editors. Hydrologic modeling of small watersheds. St. Joseph, MI: American Society of Agricultural Engineers; 1982 [chapter 8]. Freeze RA, Harlan RL. Blueprint for a physically-based digitally simulated, hydrologic response model. J Hydrol 1969;9:237–58. Freeze RA, Cherry JA. Groundwater. Englewood Cliffs, NJ: Prentice-Hall; 1979. Gabet EJ, Dunne T. Sediment detachment by rain power. Water Resour Res 2003;39:1002. doi:10.1029/2001WR000656. Gabet EJ, Dunne T. Correction to sediment detachment by rain power. Water Resour Res 2004;40. doi:10.1029/2004WR003422 [W08901]. Gerits JJP, de Lima JLMP, van den Broek TMW. Overland flow and erosion. In: Anderson, Burt, editors. Process studies in hillslope hydrology. Chichester, West Sussex, England: John Wiley and Sons Ltd.; 1990 [chapter 6]. Gessler J. The beginning of bedload movement of mixtures investigated as natural armouring in channels. Reprint, Laboratory of Hydraulics and Water Resources, CIT, 1965. Gessler D, Hall B, Spasojevic M, Holly F, Pourtaheri H, Raphelt N. Application of 3D mobile bed, hydrodynamic model. J Hydraul Eng 1999;125:737–49. Holly F, Rahuel J-L. New numerical/physical framework for mobile bed modeling, Part 1: Numerical and physical principles. J Hydraul Res 1990;28:401–16. Horton RE. Erosional development of streams and their drainage basins: hydrophysical approach to quantitative morphology. Bull Geol Soc Am 1945:275–370. Inbar M, Tamir M, Wittenberg L. Runoff and erosion processes after a forest fire in Mount Carmel, a Mediterranean area. Geomorphology 1998;24:17–33. James LD, Burges SJ. Selection, calibration and testing of hydrologic models. In: Haan, Johnson, Brakensiek, editors. Hydrologic modeling of small watersheds. St. Joseph, MI: American Society of Agricultural Engineers; 1982 [chapter 14]. Johansen MP, Hakonson TE, Breshears DD. Post-fire and erosion from rainfall simulation: contrasting forests with shrublands and grasslands. Hydrol Process 2001;15:2953–65. Jones JP, Sudicky EA, Brookfield AE, Park Y-J. Quantifying groundwater contributions to streamflow generation: the reliability of tracer-based hydrograph separation techniques. Water Resour Res, in review. Karim MF, Kennedy JF. IALLUVIAL: a computer-based flowand sediment-routing for alluvial streams and its application to the Missouri River. Iowa Institute of Hydraulic Research, The University of Iowa, Iowa City, IA, Report No. 250, 1982. Kinnell PIA. Raindrop-impact-induced erosion processes and prediction: a review. Hydrol Process 2005;19:2815–44. Kirkby MJ. Implications for sediment transport. In: Kirkby, editor. Hillslope hydrology. Chichester, West Sussex, England: John Wiley and Sons Ltd.; 1978 [chapter 9]. Knisel WG. CREAMS: a field scale model for chemicals, runoff and erosion from agricultural management systems. US Department of Agriculture, Conservation Research Report No. 26, 1980. Lane LJ, Shirley ED, Singh VP. Modeling erosion on hillslopes. In: Anderson, editor. Modelling geomorphological systems. Chichester, West Sussex, England: John Wiley and Sons Ltd; 1988 [chapter 10]. Li RM. Water and sediment routing from watersheds. In: Shen, editor. Modeling of rivers. New York, NY: Wiley; 1979 [chapter 9].
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943 [33] Loague K, VanderKwaak JE. Simulating hydrologic response for the R-5 catchment: comparison of two models and the impact of the roads. Hydrol Process 2002;16:1015–32. [34] Loague K, VanderKwaak JE. Physics-based hydrologic response simulation: platinum bridge, 1958 Edsel, or useful tool. Hydrol Process 2005;18:2949–56. [35] Loague K, Heppner CS, Abrams RH, VanderKwaak JE, Carr AE, Ebel BA. Further testing of the integrated hydrology model (InHM): event-based simulations for a small rangeland catchment located near Chickasha, Oklahoma. Hydrol Process 2005;19: 1373–98. [36] Luce CH, Black TA. Sediment production from forest roads in western Oregon. Water Resour Res 1999;35:2561–70. [37] Merritt WS, Letcher RA, Jakeman AJ. A review of erosion and sediment transport models. Environ Model Software 2003;18: 761–99. [38] Morgan RPC, Quinton JN, Smith RE, Govers G, Poesen JWA, Auerswald K, et al. The European soil erosion model (EUROSEM): a dynamic approach for predicted sediment transport from fields and small catchments. Earth Surf Process Landforms 1998;23:527–44. [39] Orvis C, Randle R. STARS: sediment transport and river simulation model. Denver, CO: Technical Service Center, US Bureau of Reclamation; 1987. [40] Panday S, Huyakorn PS. A fully coupled physically-based spatially-distributed model for evaluating surface/subsurface flow. Adv Water Resour 2004;27:361–82. [41] Pebesma EJ, Switzer P, Loague K. Error analysis for the evaluation of model performance: Rainfall-runoff time series. Hydrol Process 2005;19:1529–48. [42] Prosser IP, Rustomji P. Sediment transport capacity relations for overland flow. Prog Phys Geography 2000;24:179–93. [43] Rawls WJ, Gish TJ, Brakensiek DL. Estimating soil water retention from soil physical properties and characteristics. Adv Soil Sci 1991;16:213–34. [44] Reid LM. Research and cumulative watershed effects. PSWGTW-141, Pacific Southwest Research Station, Forest Service, US Department of Agriculture, 1993. 118p.
943
[45] Roering JJ, Gerber M. Fire and the evolution of steep, soil-mantled landscapes. Geology 2005;33:349–52. doi:10.1130/G21260.1. [46] Tayfur G. Applicability of sediment transport capacity models for nonsteady state erosion from steep slopes. J Hydrol Eng 2002; 7:252–9. [47] US Department of Agricultural, Agricultural Research Service (USDA-ARS). Hydrology, erosion, and water quality studies in the southern great plains research watershed, Southwestern Oklahoma, 1961–1978. Agricultural Reviews and Manuals, ARM-S-29, 1983. [48] VanderKwaak JE. Numerical simulation of flow and chemical transport in integrated surface–subsurface hydrologic systems. PhD dissertation, Department of Earth Sciences, University of Waterloo, Ontario, Canada, 1999. [49] VanderKwaak JE, Loague K. Hydrologic-response simulations for the R-5 catchment with a comprehensive physics-based model. Water Resour Res 2001;37:999–1013. [50] van Genuchten MT. A closed-form equation for predicting the hydraulic conductivity of unsaturated soils. Soil Sci Soc Am J 1980;44:892–8. [51] Wicks JM, Bathurst JC. SHESED: a physically based, distributed erosion and sediment yield component for the SHE hydrological modeling system. J Hydrol 1996;175:213–38. [52] Woolhiser DA, Smith RE, Goodrich DC. KINEROS, a kinematic runoff and erosion model:Documentation and user manual. Agricultural Research Service, US Department of Agriculture, Report No. ARS-77, 1990. [53] Ziegler AD, Giambelluca TW, Southerland RA. Improved method for modeling sediment transport on unpaved roads using KINEROS2 and dynamic erodibility. Hydrol Process 2002;16: 3079–89. [54] Ziegler AD, Southerland RA, Giambelluca TW. Acceleration of Horton overland flow and erosion by footpaths in an upland agricultural watershed in northern Thailand. Geomorphology 2001;41:249–62. [55] Ziemer RR, Lewis J, Rice RM, Lisle TE. Modeling the cumulative watershed effects of forest management strategies. J Environ Quality 1991;20:36–42.
Adding sediment transport to the integrated hydrology model (InHM): Development and testing Christopher S. Heppner a, Qihua Ran a, Joel E. VanderKwaak b, Keith Loague a
a,*
Department of Geological and Environmental Sciences, Stanford University, Stanford, CA 94305-2115, USA b 3DGeo Development Inc., Santa Clara, CA 95054, USA Received 31 May 2005; received in revised form 6 August 2005; accepted 13 August 2005 Available online 13 October 2005
Abstract The addition of a sediment transport algorithm to the comprehensive hydrologic-response model known as the Integrated Hydrology Model (InHM) is discussed. The first test of the sediment transport version of InHM is reported, using field data from a series of erosion experiments conducted by Gabet and Dunne [E.J. Gabet, T. Dunne, Sediment detachment by rain power, Water Resour Res 39 (2003) 1002]. The performance of the sediment transport component of InHM, in both calibration and validation phases, is judged to be successful, based upon quantitative statistical criteria. The ability to simulate sub-plot-scale interactions between surface water hydrology and rain-induced sediment transport with InHM is demonstrated. Sensitivity analysis reveals that the rainfall intensity exponent has a substantial impact on simulated sediment discharge. Future work, related to both testing InHM and much needed field experiments for model parameterization, is discussed. 2005 Elsevier Ltd. All rights reserved. Keywords: Hydrologic-response simulation; Sediment transport simulation; Rainsplash erosion; Hydraulic erosion; Plot-scale experiments; Model performance; Calibration; Validation; Sensitivity analysis
1. Introduction Understanding the relationships between near-surface hydrologic response and sediment transport is of interest to those in the fields of ecology, engineering, geomorphology, hydrology, and land-use planning. Dating back to the pioneering efforts of Robert Horton [22] there has been considerable work focused on identifying the linkages between hydrology and erosion (see the reviews of Merritt et al. [37], Gerits et al. [18], Kirkby [29], and Lane et al. [31]). During the International Hydrologic Decade, there was a heroic effort by the Agricultural Research Service (ARS) within the US Department of Agriculture (USDA) to make detailed *
Corresponding author. Tel.: +1 650 723 3090; fax: +1 650 725 0979. E-mail address: [email protected] (K. Loague). 0309-1708/$ - see front matter 2005 Elsevier Ltd. All rights reserved. doi:10.1016/j.advwatres.2005.08.003
measurements of both runoff and sediment yields for several unit-source watersheds [47]. More recently, many field-based studies have investigated hydrologically-driven erosion processes related to, for example, cumulative watershed effects [44,55], fire [23,25,45], and roads [36,53,54]. It is our opinion that it is currently not necessary (or even possible) to track, through observation or simulation, every drop of water or particle of sediment (even at relatively small scales) for problems with a landscape evolution component. We do, however, feel that it is very important to include the best characterization of hydrologic response, in terms of both process and dimension, when attempting to address transport-limited sediment transport [34]. Table 1 lists a few selected models that have been developed to consider both hydrologic response and sediment transport. Perusal of Table 1 shows that, in terms of process representation
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943
931
Table 1 Dimensionality of selected models that consider both hydrologic response and sediment transport, in comparison to the dimensionality of InHM Reference
Acronym
Year
Hydrologic response Subsurface
Surface
Channel
Overland
[30] [2] [27] [39] [21] [52] [7] [51]a [38] [20] This study
CREAMS ANSWERS IALLUVIAL STARS SEDICOUP KINEROS2 FLUVIAL-12 SHESED EUROSEM CH3D-SED InHM
1980 1980 1982 1987 1990 1990 1993 1996 1998 1999 2005
1D na na na na na na 3D na na 3D
1D 2D 1D 2D 1D 1D 2D 2D 1D 3D 2D
na 2D 1D 2D 1D 1D 2D 1D 1D 3D 2D
1D 2D na na na 1D na 2D 1D na 2D
a
Sediment transport
The hydrologic-response model driving SHESED is the European Hydrology System SHE [1].
and dimensionality, not all models are created equally. Table 2 compares the characteristics of two well-known hydrologic-response models that have sediment transport algorithms with the Integrated Hydrology Model (InHM), which is the focus of the effort reported herein. It is well known that hydrologically-driven erosion models must faithfully simulate surface responses (e.g., Horton overland flow and channel flow). However, in many settings (e.g., wetlands and steep convex hillslopes) the contributions of subsurface flow (e.g., Dunne overland flow and groundwater discharge) to the surface can be significant. Hydrologically-driven erosion models that neglect subsurface contributions will without question, for many situations, misrepresent near-surface hydrologic response and, subsequently, sediment transport. Inspection of Tables 1 and 2 clearly shows that InHM provides a firm hydrologic-response foundation upon which to build a sediment transport model. It is well beyond the scope of the effort reported here to explore fully all of the nuances encountered in the realm of physics-based sediment transport modeling. Rather, the primary objective of this study was to develop a sediment transport algorithm for InHM. The InHM sediment transport algorithm facilitates the simulation of both erosion and deposition. Built into the development is consideration for multiple species sediment transport, with detachment from rainsplash and/ or hydraulic erosion driven by spatially variable surface
water depths and velocities. The sediment transport version of InHM is tested using the field data from the rainsplash experiments (i.e., raindrop detachment with transport by flow [28]) of Gabet and Dunne [16].
2. Integrated hydrology model The comprehensive InHM was designed, in the spirit of the Freeze and Harlan [14] blueprint, to simulate quantitatively, in a fully-coupled approach, 3D variably-saturated flow and solute transport in porous media and macropores and 2D flow and solute transport over the surface and in open channels. The important and innovative characteristics of InHM, including no a priori assumption of a specific hydrologic-response mechanism (i.e., Horton or Dunne overland flow, discharge from subsurface stormflow or groundwater), are discussed elsewhere [33–35,48,49]. InHM has been successfully employed for catchment scale, event-based rainfall-runoff simulation [33,34,49] and for solute transport simulations [26,48]. It should be noted that Panday and Huyakorn [40] have reported a comprehensive physics-based hydrologic-response model that is comparable to InHM. The underlying InHM equations employed in this study are presented below. For the effort reported here, the macropore and solute transport components of
Table 2 Characteristics of the hydrology components of two well-known sediment transport models, in comparison to the characteristics of InHM Reference
[52] [51] This study
Acronym
KINEROS2 SHESED InHM
Year
1990 1996 2005
Subsurface
1D, U, A 1D, U, N/R; 2D, S, N/G 3D, U/S, N/R
Surface Overland
Channel
1D, N/K 2D, N/DW 2D, N/DW
1D, N/K 1D, N/DW 2D, N/DW
Coupling
Process
SQ SQ FO
H H, D, SS, GW H, D, SS, GW
U (unsaturated); S (saturated); U/S (unsaturated/saturated); N (numerical solution); A (analytical solution); R (Richards equation); G (groundwater flow equation); DW (diffusion wave); K (kinematic wave); SQ (sequential); FO (first-order); H (Horton overland flow); D (Dunne overland flow); SS (subsurface stormflow); GW (groundwater).
932
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943
InHM are not employed. For hydrologic response, only the basic equations of subsurface and surface water movement are provided; the reader requiring more detail is referred to VanderKwaak [48]. For sediment transport, each component of the algorithm is described. A simplified flow chart for the sediment transport version of InHM is given in Fig. 1.
tween the subsurface and surface continua (T1), / is porosity (–), Sw is water saturation (–), t is time (T), fa is the area fraction associated with each continuum (–), and fv is the volume fraction associated with each continuum. The Darcy flux is given by q g ~ q ¼ k rw w ~ krðw þ zÞ ð2Þ lw
2.1. Hydrologic response
where krw is the relative permeability (–), qw is the density of water (ML3), g is the gravitational acceleration (LT2), lw is the dynamic viscosity of water (ML1T1), ~ k is the intrinsic permeability vector (L2), z is the elevation head (L), and w is the pressure head (L). The transient flow of water on the land surface is estimated by the diffusion wave approximation of the depth-integrated shallow water equations. Such 2D
Subsurface flow is calculated by r f a~ q qb qe ¼ f v
o/S w ot
ð1Þ
where ~ q is the Darcy flux (LT1), qb is a specified rate source/sink (T1), qe is the rate of water exchange be-
Fig. 1. Simplified flow chart for InHM with sediment transport. Note, the superscript plus signs (+) refer to model components included in InHM without sediment transport.
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943
surface flow is conceptualized as a second continuum that interacts with the underlying variably saturated porous medium through a thin soil layer of thickness as (L). Assuming a negligible influence of inertial forces and a shallow depth of water, ws (L), the conservation of water on the land surface is described by ~ r wmobile qs as qb as qe ¼ s
oðS ws hs þ wstore Þ s ot
ð3Þ
where ~ qs is the surface water velocity (LT1), qb is the source/sink rate (i.e., rainfall/evaporation) (T1), qe is the surface–subsurface water exchange rate (T1), and hs is the average height of non-discretized surface microtopography (L); the superscripts mobile and store refer to water exceeding and held in depression storage, respectively [48]. Surface water velocities are calculated utilizing a two-dimensional form of the ManningÕs equation given by ~ qs ¼
ðwmobile Þ s ~ nU1=2
2=3
rðws þ zÞ
ð4Þ
where ~ n is the ManningÕs surface roughness tensor (TL1/3) and U is the friction (or energy) slope (–). The linkages between (1) and (3) are through first-order, physically-based flux relationships driven by pressure head gradients. The flow equations are discretized in space using the control volume finite-element method, which combines the geometric flexibility of finite elements with the local conservation characteristics of control volumes. Each coupled system of non-linear equations in an InHM simulation is solved implicitly using Newton iteration [48]. Efficient and robust iterative sparse matrix methods are used to solve the large sparse Jacobian systems [48]. InHM provides the spatio-temporally variable surface water depths and velocities required for simulating sediment transport, as described in the next section.
component and the terms in brackets are the sources and sinks of sediment concentration. The net erosion rate is calculated by esed ¼ es þ eh
ð6Þ 3 1
where es is the rainsplash erosion rate (L T ) and eh is the hydraulic erosion rate (L3T1). The rainsplash erosion rate in (6) is estimated as follows ( cf F ðws ÞðcosðhÞ rÞb A3D q > 0 es ¼ ð7Þ 0 q<0 where cf is the rainsplash coefficient ((TL1)b1), b is the rainfall intensity exponent (–), F(ws) is a factor representing the reduction in splash erosion with increasing surface water depths (–), h is the angle of the element from horizontal (–), r is the rainfall rate (LT1), A3D is the three-dimensional area associated with the node (L2), and q is the sum of rainfall rate and infiltration rate (LT1). The F(ws) function in (7) is estimated by [52] F ðws Þ ¼ expðcw ws Þ
eh ¼ ased ðC sedmax C sed Þ
ð9Þ
where ased is the hydraulic erosion transfer coefficient (L3T1) and C sedmax is the concentration at equilibrium transport capacity (L3L3). The C sedmax term in (9) is calculated by [10] C sedmax ¼ 0:05
where Csed is volumetric sediment concentration (L3L3), esed is the volumetric rate of soil erosion and/or deposition (L3T1), qbsj represents the rate of water added/removed via the jth boundary condition (L3T1), C sedj is the sediment concentration of the water added/removed via the jth boundary condition (L3L3), BC is the total number of boundary conditions, and Vw is the volume of water at the node (L3). On the right hand side of (5) the first term is the advective transport
ð8Þ
where cw is a parameter representing the damping effectiveness of surface water (L1). Prior to the accumulation of surface water and for deep flows F(ws) equals, respectively, one and zero. When multiple sediment species are considered, the total rainsplash erosion rate is apportioned to each species by scaling by the source fraction (ri) of that species. The hydraulic erosion rate in (6) is conceptualized as a kinetic transfer process [46] and is estimated as follows [52]
2.2. Sediment transport Depth-integrated multiple species sediment transport, restricted to the surface continuum, for each particle species is calculated by [3] ! BC X oC sed 1 b ¼ r ~ qs C sed þ esed þ qsj C sedj ð5Þ Vw ot j¼1
933
~ qs q3 g2 Dsed ws ðcsed 1Þ2
ð10Þ
where q* is the local shear velocity (LT1), Dsed is particle diameter (L), and csed is the particle specific gravity (i.e., the ratio of particle density to fluid density) (–). The shear velocity is calculated by pffiffiffiffiffiffiffiffiffiffiffiffiffi q ¼ gws rz ð11Þ where $z is the local surface slope (–), given by rz ¼ abs½tanðhÞ
ð12Þ
The ased term in (9) is defined for erosion (i.e., C sedmax > C sed Þ and deposition (i.e., C sedmax < C sed Þ as C sed > C sedmax 2vsed n ased ¼ A ð13Þ u~ qs ws v C sed < C sedmax
934
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943
where A is the area associated with the node (L2), vsed is the particle settling velocity (calculated using StokeÕs Law for the terminal velocity of spheres in calm fluid) (LT1), n is a coefficient [46] related to turbulence in the surface water due to raindrop impact (–), u is an erodibility coefficient [13] related to surface properties and texture (L1), and v is the particle erodibility factor (ranging from zero to one) (–). During erosion ased is a function of water depth, water velocity, surface properties (e.g., roughness or texture) and, in the case of multiple species, the relative source fraction and the armoring effect of any immobile species. During deposition ased is a function of particle settling velocity and the raindropinduced turbulence. Alternative transport capacity formulae [12,42,46] could easily be implemented with InHM. For each sediment species, denoted by the subscript i, the v term in (13) is calculated by vi ¼
ri 1k ð1 rimmobile Þ i
ð14Þ
where ri is the species source fraction (–), 1i is the cohesion coefficient (ranging from zero to one for soils with considerable or no cohesion, respectively) (–), rimmobile is the fraction of immobile species (–), and k is the armoring factor (–). First, erodibility (in a binary sense) is calculated with a critical shear stress test. A species is erodible when the critical shear stress (scr) for incipient motion is exceeded, and unerodible when this is not the case. The relationship for scr (ML1T2) used here is [32] scr ¼ ds ðcs cw ÞDsed
ð15Þ
where ds is a constant [19] (equal to 0.047 for most flow conditions) (–), cs is the specific weight of the sediment
(ML2T2), and cw is the specific weight of water (ML2T2). The shear stress is given by ð16Þ
s ¼ cw ws S
where S is the energy slope (–). The rimmobile term in (14) is calculated as rimmobile ¼
k mx X
ð17Þ
ri
i¼k sm
where ksm is the index of the smallest immobile particle species and kmx is the index of the largest particle species. The source fractions for all species must sum to one. The k term in (14), which represents the effect of immobile species on the erodibility of mobile species, is calculated by k ¼ max½0; ð1 arimmobile Þ
ð18Þ
where a is an empirical shadow factor (–), with a value greater than one. The non-linear equations describing transport of multiple sediment species on the land surface are solved implicitly using Newton iteration. The use of appropriate numerical methods ensures mass conservation and accurate monotone solutions [48].
3. The Gabet and Dunne experiments The careful rainsplash plot experiments reported by Gabet and Dunne [16,17] were selected to begin testing the sediment transport version of InHM. The Gabet and Dunne [16] fieldwork was carried out at the University of California Sedgwick Reserve, near Santa Barbara. The information required to excite the InHM sediment transport algorithm is given in Table 3. The
Table 3 Parameters and variables required for simulating sediment transport with InHM Symbol
Definition (units)
Source, this study
Dsed qsed r $z r ~ qs ws S vsed scr cw b n u f a cf
Species average grain diameter (L) Particle density (ML3) Species source fraction (–) Surface slope (–) Rainfall rate (LT1) Water velocity (LT1) Water depth (L) Energy slope (–) Species settling velocity (LT1) Species critical shear stress (ML1T2) Rainsplash depth dampening factor (L1) Rainfall intensity exponent (–) Rain-induced turbulence coefficient (–) Surface erodibility coefficient (–) Sediment cohesion coefficient (–) Shadow factor (–) Rainsplash coefficient ((TL1)b1)
0.00001 m [16] 2650 kg m3 1.0 Table 4 [16] Table 4 [16] Flow simulation, InHM Flow simulation, InHM Flow simulation, InHM StokeÕs Law (15) 600 m1 [52] 1.6 [16] 0.25 [46] NAa NAa NAa Calibration, Table 4
a
Not required for this study.
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943
precious data provided by Gabet and Dunne [16] (i.e., steady-state runoff and sediment concentrations) does not include spatio-temporally variable information that is typically employed to check (compare against) InHM simulated near-surface hydrologic response [33,35,49]. The missing data include (i) surface water depths, (ii) water table depths, (iii) total head values, and (iv) soilwater content values. However, the site information reported by Gabet and Dunne [16] is sufficient to parameterize InHM for the simulations reported herein. For example, for the different plots the slopes were between 4 and 17, the groundcover (vegetation) ranged between 18% and 94%, the soils were silty clay loams, and hydrologic response was only by the Horton overland flow mechanism. Details of the experimental set up for the sprinkler experiments (e.g., simulator operation, nozzle size) are provided by Gabet and Dunne [16].
4. Methods A split sample (calibration/validation) approach, grounded to the Gabet and Dunne [16] experimental data, was used in this study to evaluate the sediment transport version of InHM. It is important to emphasize that the simulations reported here consider erosion for a single species by rainsplash only, despite the ability of InHM to address more complex problems. The chronology of the model testing effort employed in this study can be summarized by four phases: (i) boundaryvalue problem (BVP) development, (ii) hydrologicresponse simulations (calibration), (iii) sediment transport simulations (calibration), and (iv) sediment transport simulations (validation). The four phases are described below. 4.1. Phase I, Boundary-value problem Each of the Gabet and Dunne [16] plots had the same surface dimensions (6.0 m by 2.4 m). Fig. 2 shows the finite-element mesh used for each of the BVPs in this study. The BVPs for the different experimental plots had different slopes, as summarized in Table 4. The dimensions of the mesh are 6.0 m long, 2.4 m wide, and 3.0 m deep. The vertical nodal spacing (Dz) in the mesh varies from 0.05 to 0.25 m; the horizontal nodal spacing (Dx, Dy) is 0.2 m. The total number of nodes/ elements in the mesh are, respectively, 9009 and 16,128. With reference to the lettered coordinate locations identified in Fig. 2, the boundary conditions for each face of each BVP are as follows: impermeable (A–E–G–C, B–F–H–D, C–G–H–D), local sink (A–B– D–C, E–F–H–G), and flux (A–E–F–B). Relative to the near-surface soil-hydraulic properties, the entire subsurface domain was considered to be homogeneous and isotropic; the porosity and compressibility were set at 0.46
935
Z
E
Y
X
F
A
H B
G
C D 1m Fig. 2. Boundary-value problem used to represent the Gabet and Dunne [16] experiments for the InHM simulations of hydrologic response and sediment transport.
(i.e., average value for silty clay loam soils [6,43]) and 1 · 107 ms2 kg1 (i.e., mid range value for clay [15]), respectively. The soil-water retention and permeability functions, based upon the van Genuchten [50] approach and characteristic information on silty clay loam soils from Carsel and Parrish [6] and Rawls et al. [43], are given in Fig. 3. For each InHM simulation the mobile water depth (i.e., depression storage) was set to 0.0005 m and the height of microtopography was set at 0.01 m. The slope and groundcover percentage for each of the experiments included in this study are given in Table 4. For each experiment, the ManningÕs roughness coefficient was calculated using [16] n ¼ 0:053e2:7Cv
ð19Þ
where Cv is the groundcover percentage (expressed as a decimal, 70% = 0.70). The entire subsurface domain was unsaturated throughout each of the simulations reported herein (i.e., no water table within 3 m of the surface), with surface runoff resulting from Horton overland flow. The initial conditions, taken to be the same for each experiment, were set by upstream and downstream local head boundary conditions so that the initial water table was 5.0 m below the land surface. The validity of the assumed initial conditions (i.e., water table position) was
936
Table 4 Summary of field experiment characteristics, parameterization, and results from hydrologic response/sediment transport simulations with InHM Rainfall ratea,b (105 m/s)
Slopea ()
Groundcovera (%)
ManningÕs roughness coefficient (–)
Calibrated saturated hydraulic conductivity (105 m/s)
Observed steady-state water dischargea (104 m3/s)
Simulated steady-state water discharge (104 m3/s)
Best-fit rainsplash coefficient (–)
Estimated rainsplash coefficient (–)
Observed steady-state sediment dischargea (104 kg/s)
Simulated steady-state sediment discharge (104 kg/s)
2/1 3/1 5/2 6/2 7/2 8/2 12/3 14/4 15/4 16/5 17/5 18/5 19/6 20/6 21/6 22/7 23/7 24/7 25/8 26/8 27/8 28/8 29/8 30/8 31/8 32/9 33/9 35/9 36/10 37/10 38/10 39/10 40/10
3.47c 1.81 1.67 3.36 1.72 3.25 3.31 1.39 1.67 2.06 1.53 1.53 1.25 3.06 3.47 2.92 2.69 2.03 1.83 2.28 3.06 1.89 2.81 2.03 3.22 1.86 3.08 3.14 1.86 3.69 3.33 3.39 3.39
10 10 17 17 17 17 14 4 4 17 17 17 13 13 13 9 9 9 13 13 13 13 13 13 13 13 13 13 5 5 5 5 5
94 88 72 72 72 72 68 98 98 82 74 66 70 70 56 74 64 68 64 62 62 39 35 32 32 34 34 34 86 76 52 45 46
0.671 0.570 0.370 0.370 0.370 0.370 0.332 0.747 0.747 0.485 0.391 0.315 0.351 0.351 0.240 0.391 0.298 0.332 0.298 0.283 0.283 0.152 0.136 0.126 0.126 0.133 0.133 0.133 0.540 0.413 0.216 0.179 0.184
1.10d 1.10 0.66 0.66 0.66 0.66 1.30 1.13 1.13 0.10 0.10 0.10 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 1.00 0.70 0.70 0.70 0.80 0.80 0.80 0.80 0.80
3.36e 1.10 1.44 3.34 1.44 3.34 2.42 0.36 0.65 2.50 1.99 1.99 0.36 2.40 2.86 2.50 2.50 1.66 1.25 1.49 2.86 1.44 2.86 1.44 3.36 1.61 3.34 3.34 1.66 4.01 3.34 4.01 3.34
3.36 0.96 1.37 3.78 1.45 3.62 2.79 0.34 0.73 2.73 1.98 1.98 0.33 2.88 3.49 2.73 2.40 1.44 1.14 1.78 2.89 1.24 2.53 1.44 3.16 1.63 3.41 3.48 1.49 4.14 3.65 3.74 3.74
1.09 0.16h 0.88 2.27 0.80 3.14h 1.55 0.09h 0.12h 0.89 0.92 0.82 0.81 2.84h 1.73h 2.29 1.66 0.46h 1.43h 1.05 4.93 1.63h 5.87 1.90 4.89h 1.54 5.40 7.86h 0.67 3.49h 5.84 6.10h 5.35
0.61
1.34f 0.11 1.18 7.37 1.14 9.84 4.05 0.02 0.04 1.80 1.57 1.57 0.34 7.01 6.66 4.09 3.56 0.51 1.86 2.03 14.51 3.51 23.62 4.40 25.17 3.73 28.62 42.37 0.53 5.73 16.68 21.04 18.25
0.75g 0.11i 1.43g 3.44g 1.51g 9.85i 6.75g 0.02i 0.04i 3.08g 3.65g 5.28g 1.04g 7.00i 6.63i 3.91g 6.09g 0.60i 1.84i 5.91g 8.87g 3.05i 20.33g 12.31g 25.31i 12.60g 27.22g 42.29i 0.84g 5.69i 11.16g 21.04i 14.52g
a b c d e f g h i
Gabet and Dunne [16]. Irrigation rates from Gabet and Dunne [16] are taken as rainfall rates. 3.47 · 105. 1.10 · 105. 3.36 · 104. 1.34 · 104. Validation simulation. Value used in calibration. Calibration simulation.
2.28 2.28 2.28 2.58
1.52 2.13 2.74 2.43
2.13 2.89
3.04 3.04 5.09 5.32 5.17 5.17 1.22 3.80 4.26
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943
Experiment number/ plot numbera
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943 0.5 Soil-water content Relative permeability
0.8
0.3
0.6
0.2
0.4
0.1
0.2
0
Relative permeability (-)
Soil-water content (-)
4.3. Phase III, Sediment transport simulations— calibration
1
0.4
0 0
-1
-2
-3
-4
937
-5
Pressure head (m) Fig. 3. Characteristic curves used to represent the Gabet and Dunne [16] experiments for the InHM simulations of hydrologic response and sediment transport.
confirmed through steady-state hillslope-scale long-profile InHM simulations. 4.2. Phase II, Hydrologic-response simulations The focus of this effort was not to test the hydrologicresponse capability of InHM, which has already been reported for more complete data sets [33,35,49]. Hydrologic-response simulations with InHM were conducted for each of the 44 experiments reported by Gabet and Dunne [16]. For each of the 11 plots the saturated hydraulic conductivity was adjusted, until the best match was achieved between the observed and simulated steady-state discharges. Each of the InHM simulations successfully captured the Horton overland flow response reported by Gabet and Dunne [16]. Eleven of the Gabet and Dunne [16] experiments were removed from further consideration in this study for two reasons: • It was not possible to identify experiments 4, 41, 42, 43, and 44 (from Table 1 of Gabet and Dunne [16]) in Fig. 5 of Gabet and Dunne [16]. Therefore, these five experiments were not included in this study. • It was not possible to obtain satisfactory matches between observed and simulated discharges for experiments 1, 9, 10, 11, 13, and 34 (from Gabet and Dunne [16]) during the hydrology calibration. Greater runoff was observed for these experiments (from the same plots) for lower irrigation intensities, which is most likely due to variations (undocumented) in initial conditions. Therefore, these six experiments were not included in this study to avoid further calibration. Removing 11 of the 44 Gabet and Dunne [16] experiments, and subsequently one of the 11 plots, yields 33 experiments, conducted on ten plots, for further consideration in this study.
The sediment discharge for each experiment was calculated by multiplying the InHM simulated concentration and water discharge. Four steps were employed to calibrate the sediment transport version of InHM using the Gabet and Dunne [16] experimental data. The first step was to rerun the calibrated InHM (Phase II) for each of the 33 experiments listed in Table 4. For each InHM simulation, the rainsplash coefficient in (7) was adjusted until a best-fit was achieved between the observed and simulated sediment discharges. The second step was to develop a relationship (multiple regression) between the best-fit rainsplash coefficients and the groundcovers and slopes reported for the ten plots. The third step was to select 13 of the 33 experiments, representing the entire range of slopes and groundcovers, as calibration experiments for the first part of a split sample test (see Table 4). The fourth and final step was to develop a relationship (linear regression) for the 13 experiments selected in step 3 to estimate rainsplash coefficients from whichever characteristic (groundcover or slope) in step 2 provides the greatest correlation to the best-fit rainsplash coefficients. 4.4. Phase IV, Sediment transport simulations— validation The InHM sediment transport validation simulations, for the 20 experiments in the split sample (see Table 4), build upon the hydrologic response and sediment transport calibration simulations in Phases II and III, respectively. The rainsplash coefficients used in the validation simulations were estimated based upon the relationship established in step 4 of Phase III. 5. Model performance evaluation The two measures of model performance [24,41] used to test InHM are modeling efficiency (EF) and the relative root mean square error (RRMSE). The mathematical expressions for EF and RRMSE are " #, n n n X X X 2 2 2 EF ¼ ðOi OÞ ðP i Oi Þ ðOi OÞ i¼1
i¼1
" RRMSE ¼
n X 2 ðP i Oi Þ =n i¼1
i¼1
#0:5
ð20Þ 100 O
ð21Þ
where Pi are the predicted values, Oi are the observed values, n is the number of samples, and O is the mean of the observed data. When the observed and simulated values are identical the EF and RRMSE statistics are 1.0 and 0.0, respectively.
938
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943
6. Results Simulated steady-state water discharge (10-4 m3/s)
5
6.1. Hydrologic-response simulations The rainfall rates, ManningÕs roughness coefficients, saturated hydraulic conductivity values, and observed/ simulated water discharges for each of the 33 experiments included in this study are provided in Table 4. Results from the InHM hydrologic-response calibration simulations are summarized in Table 5 and Fig. 4. The model performance statistics given in Table 5 and the observed versus simulated discharges illustrated in Fig. 4 show that InHM was able, as should be expected in a calibration mode, to capture the steady-state hydrologic response (i.e., Horton overland flow) for each of the 33 experiments (i.e., EF of 0.94).
Plot 1 Plot 2 Plot 3 Plot 4 Plot 5 Plot 6 Plot 7 Plot 8 Plot 9 Plot 10
4
3
2
1
0 0
6.2. Sediment transport simulations—calibration
1
2
3
4 -4
5 3
Observed steady-state water discharge (10 m /s)
The best-fit rainsplash coefficients and observed/simulated steady-state sediment discharges for each of the 33 experiments and the subset of 13 calibration experiments are given in Table 4. Results from the InHM sediment transport calibration simulations are summarized in Figs. 5 and 6 and Table 5. Fig. 5a shows the relationship between the best-fit and multiple regression estimates of the rainsplash coefficient for the 33 experiments. The correlation between the best-fit rainsplash coefficients and the slopes and groundcovers for
Fig. 4. InHM simulated hydrologic response results for the Gabet and Dunne [16] experiments. The solid line is 1 to 1; the dashed lines are plus/minus 10% of the 1 to 1 line.
the 33 experiments are shown in Fig. 5b and c, respectively. By comparing the R2 values in Fig. 5 it is obvious, for this study, that groundcover provides a greater correlation to the best-fit rainsplash coefficient than slope. Fig. 6 shows the relationship between the best-fit rainsplash coefficient and groundcover for the subset of 13
Table 5 Performance results for InHM simulations of hydrologic response and sediment transport for the Gabet and Dunne [16] experiments Process
Number of experiments/ simulationsa
Calibration EF (–)
RRMSE (%)
Hydrologic response Sediment transport Sediment transport Sediment transport
33 33 13 20
0.94 1.00 1.00
10.6 1.5 1.5
Validation EF (–)
RRMSE (%)
0.78
54.8
EF is modeling efficiency, (20); RRMSE is the relative root mean square error(21). a See Table 4.
(a) 10
(b) 10
2
R = 0.51
Cfbest (-)
(c) 10
2
R = 0.02
R = 0.45
8
8
8
6
6
6
4
4
4
2
2
2
0
0 0
2
4
2
6
Cfreg (-)
8
10
0 0
5
10
15
Slope (degrees)
20
0
25
50
75
100
Groundcover (%)
Fig. 5. Best-fit rainsplash coefficient relationships from 33 InHM simulations for the Gabet and Dunne [16] experiments. (a) Best-fit to multiple regression (slope and groundcover). (b) Best-fit to slope. (c) Best-fit to groundcover. The solid line in (a) is 1 to 1; the dashed lines in (b) and (c) are best-fit regression lines.
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943
calibration experiments. The model performance statistics given in Table 5 show that InHM performed well, in a calibration mode, for the entire set of experiments and for the subset of 13 experiments (i.e., EF of 1.00). 6.3. Sediment transport simulations—validation The regression estimated rainsplash coefficients and observed/simulated steady-state sediment discharges for the 20 validation experiments are given in Table 4. Results from the InHM sediment transport validation simulations are summarized in Table 5 and Fig. 7. The model performance statistics given in Table 5 show that
Best-fit rainsplash coefficient (-)
9 8 2
7
R = 0.54
6 5 4 3 2 1 0 20
30
40
50
60
70
80
90
100
Groundcover (%) Fig. 6. Best-fit regression relationship (solid line) for the rainsplash coefficient relative to groundcover based upon InHM simulations for 13 of the Gabet and Dunne [16] experiments (i.e., the subset of calibration experiments).
Simulated steady-state sediment discharge (10-4 kg/s)
30
25
20
15
10
5
0 0
5
10
15
20
25
30
Observed steady-state sediment discharge (10-4 kg/s) Fig. 7. Comparison of observed and InHM simulated sediment discharge rates for the Gabet and Dunne [16] experiments. The solid line is 1 to 1.
939
InHM performed reasonably well in a validation mode for the subset of 20 experiments (i.e., EF of 0.78), especially considering the relatively simple method of estimating the rainsplash coefficient. Fig. 7 shows that sediment discharges are overestimated and underestimated, respectively, for the smaller and larger observed values.
7. Discussion 7.1. InHM as a diagnostic or concept-development tool Fig. 8 illustrates the type of surface water (depth, discharge, runoff rate) and sediment (concentration, discharge, erosion rate) information that can be generated with InHM. The steady-state snapshot results shown in Fig. 8 are from a longitudinal transect of 33 equallyspaced monitoring (simulated) points for experiment 21 (see Table 4). Relative to the capabilities of the model, experiment 21 is relatively simple application of InHM in that only sheet flow is considered (i.e., no flow convergence related to topography; essentially a 1D BVP) and that the only hydrologic response and sediment detachment processes are, respectively, Horton overland flow and rainsplash. The two principal variables extracted from the experiment 21 results are surface water depth and sediment concentration. The surface water velocity and discharge values shown in Fig. 8 were calculated from the depth information using ManningÕs equation. Inspection of Fig. 8 shows that discharge increases linearly downslope, which is expected as the rainfall rate for this experiment (i.e., all experiments in this study) exceeds the saturated hydraulic conductivity, Ks (LT1), of the soil (i.e., Horton by definition). The difference in discharge per unit area from one point to the next point downslope equals the difference between the rainfall rate and the saturated hydraulic conductivity. The relationship between saturated hydraulic conductivity and permeability (used in InHM) is Ks = kqwg/lw. In Fig. 8 the longitudinal sediment concentration profile decreases in magnitude with distance downslope. There are two reasons for this decrease: (i) rainsplash erosion rate is a negative function of the water depth, which increases downslope; (ii) transported sediment becomes diluted in the greater volumes of water as total discharge increases downslope. Hydraulic erosion is minimal in this study [16]. If hydraulic erosion were not minimal, the increased downslope depths and velocities would contribute to greater concentrations. The sediment discharge rates shown in Fig. 8 are the product of water discharge and sediment concentration. The erosion rates shown in Fig. 8 were calculated as the difference in sediment discharge from one monitoring point to the next downstream monitoring point divided by
0.003
2
0.0004
7
Water depth Discharge Runoff rate Sediment concentration Sediment discharge Erosion rate
6
5
4 0.002 3
2 0.001
0.0003
0.0002
0.0001
3
Water depth (m), Sediment discharge (kg/s)
0.004
1
0
0 0
1
2
3
4
5
Discharge (m /s), Runoff rate (m/s), Erosion rate (kg/m s)
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943
Sediment concentration (kg/m3)
940
0
6
Distance from downstream boundary (m) Fig. 8. InHM simulation results (i.e., water depth, discharge, runoff rate, sediment concentration, sediment discharge, and erosion rate) along a longitudinal transect of plot number 6 for experiment number 21 (see Table 4).
the area. Similarly, the runoff rate was calculated as the difference in water discharge per unit width from one monitoring point to the next downstream monitoring point. Inspection of Fig. 8 shows, for experiment 21, that the erosion rate decreases downslope even though the total sediment discharge increases. The slight spatial oscillations in the runoff rate and erosion rate curves can be traced to variations in the elemental geometry and nodal connectivities for the nodes along the centerline of the plot, resulting in small variations in calculated discharge. The diagnostic and concept-development capabilities of InHM are most useful for hydrologic response/sediment transport BVPs that are more challenging than the one discussed here. Related to hydrologic response, rather than just focusing on a Horton mechanism response, the contributions from both the Horton and Dunne mechanisms, as-well-as subsurface stormflow and groundwater, can all be considered simultaneously with InHM. Related to sediment transport, the contributions from rainsplash (as reported here) and hydraulic
erosion can be considered individually or simultaneously with InHM. 7.2. Sensitivity analysis The sensitivity of the InHM simulated hydrologic response/sediment transport results reported herein, for reasonable uncertainties in the Gabet and Dunne [16] experimental data, was investigated for experiment 21 by individually increasing/decreasing (i.e., relatively small plus/minus 10% changes) the values of six parameters and conducting 12 new simulations. The six parameters considered in the sensitivity analysis are (i) the rainfall rate, (ii) the saturated hydraulic conductivity, (iii) the ManningÕs roughness coefficient, (iv) the rainsplash coefficient, (v) the rainfall intensity exponent, and (vi) the rainsplash depth dampening factor. The base case (experiment 21) parameter values are given in Tables 3 and 4. Table 6 presents the results from the 12 new simulations as percent differences relative to the base-case simulation for both water and sediment
Table 6 Sensitivity analysis results for InHM simulations of hydrologic response and sediment transport for experiment 21 Parameter
Percent differences (relative to the base case) for increased/decreased parameter values Water discharge
Rainfall rate, r Saturated hydraulic conductivity, Ks ManningÕs roughness coefficient, n Rainsplash coefficient, cf Rainfall intensity exponent, b Rainsplash depth dampening factor, cw
Sediment discharge
10% increase
10% decrease
10% increase
10% decrease
14 4 0 na na na
14 4 0 na na na
12 3 7 10 81 14
12 3 8 10 417 16
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943
discharges. Perusal of the results in Table 6 leads to the following comments: • Increasing/decreasing the rainfall rate resulted in relatively small and symmetrical changes in the simulated water and sediment discharges (i.e., the 10% increase and decrease in r lead to 14% greater and smaller, respectively, water discharge and 12% greater and smaller, respectively, sediment discharges). The disproportionate rise in water discharge is related to the effect of depression storage, which makes the actual depth of runoff less than the total depth. As depression storage is full in both the base case and increased rainfall rate scenarios, the additional rain in the increased rainfall rate scenario adds only to the mobile water. The sensitivity of sediment discharge to the rainfall rate is complicated by negative feedbacks (e.g., higher rainfall rates would seem to produce higher rainsplash erosion rates, but increased water depths reduce the rainsplash effect). • Increasing/decreasing the saturated hydraulic conductivity resulted in relatively small and symmetrical changes in the simulated water and sediment discharges (i.e., the 10% increase and decrease in Ks lead to 4% smaller and greater, respectively, water discharges and 3% smaller and greater, respectively, sediment discharges). In general, greater infiltration rates and smaller discharges (for the Horton case) result from increased Ks values. There is a negative feedback between water depth and the rainsplash erosion rate. • Increasing/decreasing the ManningÕs roughness coefficient had no effect on the simulated water discharge but resulted in relatively small and symmetrical changes in the simulated sediment discharges (i.e., the 10% increase and decrease in n lead to 7% smaller and 8% greater, respectively, sediment discharges). Differences in the depth–velocity relationship coupled with the depth dependence of the rainsplash erosion rate produce the changes in the sediment discharges. • Increasing/decreasing the rainsplash coefficient resulted in directly proportional changes in the simulated sediment discharge (i.e., the 10% increase and decrease in cf lead to 10% greater and smaller sediment discharges, respectively). • Increasing/decreasing the rainfall intensity exponent resulted in large and non-symmetrical changes in the simulated sediment discharge (i.e., the 10% increase and decrease in b lead to 81% smaller and 417% greater, respectively, sediment discharges). Small shifts in b produce large changes in the rainsplash erosion rate [see (7)], emphasizing the need for continued efforts to thoroughly quantify the relationship between rainsplash detachment rate and rainfall intensity [16].
941
• Increasing/decreasing the rainsplash depth dampening factor resulted in relatively small and symmetrical changes in the simulated sediment discharge (i.e., the 10% increase and decrease in cw lead to 14% smaller and 16% greater, respectively, sediment discharges). Again, the negative feedback between water depth and the rainsplash erosion rate is apparent.
7.3. Future work The use of the Gabet and Dunne [16] experiments to test the rainsplash erosion component of the sediment transport version of InHM is but a first step. Future efforts will test the hydraulic erosion component of InHM. It will also be important to test the sediment transport version of InHM for combined hillslope and channel systems. The first and second authors, respectively, are currently using InHM for long-term landscape evolution problems at the catchment (i.e., the 0.1 km2 R-5 catchment near Chickasha, OK) and regional (i.e., the 117 km2 Hawaiian island of KahoÕolawe) scales. The measure and model approach pioneered by Robert Horton (see [4,5]) is the best protocol in 21st century hydrology. In general, at the catchment scale, computer resources do not limit the use of InHM for multidimensional non-linear problems. There are, however, significant information shortfalls that must be met before process-based hydrologic response/landscape evolution simulation with InHM can be considered routine, even in the research arena. The information shortfall in hydrology/geomorphology is not news to those practicing the measure and model approach, especially those who have focused on the assessment of process-based models. There have been several calls for detailed long-duration field monitoring programs at selected experimental catchments (see [8,9,11]). The information needed to rigorously test and implement InHM for hydrologically-driven landscape evolution problems includes the data needed to excite both the hydrology [33,35,49] and sediment transport (see Table 3 and Fig. 1) components of the model, as-well-as observations to compare against the simulated results. It is also important to understand the interactions between the different types of data (e.g., groundcover and cohesion). The information needs to be acquired from field experiments that consider the spatio-temporal variability of near-surface parameters (which may not be time invariant [34]) and hydrologic/sediment transport responses.
8. Summary If the objective of a modeling study is to estimate surface erosion based upon simulated runoff, the full range
942
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943
of near-surface hydrologic-response possibilities needs to be considered. InHM, as developed by VanderKwaak [48], is a comprehensive physics-based hydrologic-response model, without an a priori assumption for a specific hydrologic-response mechanism (e.g., Horton overland flow), providing a firm foundation to simulate sediment transport. We have developed a spatially distributed sediment transport algorithm for InHM and done preliminary testing of the rainsplash component. The sediment transport version of InHM has been designed (like the original version of InHM [48]) to meld measurements with modeling. We have demonstrated the spatial capabilities of the new model (see Fig. 8) and stress the need for additional spatio-temporal measurements to move forward with process-based hydrologically-driven erosion research.
[12] [13]
[14] [15] [16] [17]
[18]
[19]
Acknowledgements [20]
The effort reported here was supported by National Science Foundation grant EAR-0438749. Our conversations with Bob Street on sediment transport modeling have been very helpful. The thoughtful comments from Steve Burges and four anonymous reviewers on an earlier version of this paper are greatly appreciated.
[21]
[22]
[23]
References [24] [1] Abbott MB, Bathurst JC, Cunge JA, OÕConnell PE, Rasmussen JL. An introduction to the European Hydrology System SHE, 2: Structure of a physically-based, distributed modeling system. J Hydrol 1986;87:61–77. [2] Beasley DB, Huggins LF, Monke EJ. ANSWERS: a model for watershed planning. Trans Am Soc Agric Eng 1980;23:938–44. [3] Bennett JP. Concepts of mathematical modeling of sediment yield. Water Resour Res 1974;10:485–92. [4] Beven K, Robert E. HortonÕs perceptual model of infiltration processes. Hydrol Process 2004;18:3447–60. [5] Beven K, Robert E. Horton and the abrupt rises of ground water. Hydrol Process 2004;18:3687–96. [6] Carsel RF, Parrish RS. Developing joint probability distributions of soil water retention characteristics. Water Resour Res 1988;24:755–69. [7] Chang HH. Fluvial-12, Mathematical model for erodible channels (userÕs manual). San Diego, CA: San Diego State University; 1993. [8] Dunne T. Relation of field studies and modeling in the prediction of storm runoff. J Hydrol 1983;65:25–48. [9] Dunne T. Wolman lecture: hydrologic science . . . in landscapes . . . on a planet . . . in the future. Hydrologic sciences: taking stock and looking ahead. Washington, DC: National Academy Press; 1998. p. 10–43. [10] Engelund F, Hansen E. A monograph on sediment transport in alluvial streams. Copenhagen, Denmark: Technical University of Denmark; 1967. [11] Entekhabi D, Asrar GR, Betts AK, Beven KJ, Bras RL, Duffy CJ, et al. An agenda for land surface hydrology research and a call
[25]
[26]
[27]
[28] [29]
[30]
[31]
[32]
for the second international hydrological decade. Bull Am Meteorol Soc 1999;80:2043–58. Ferro V. Evaluating overland flow sediment transport capacity. Hydrol Process 1998;12:1895–910. Foster GR. Modeling the erosion process. In: Haan, Johnson, Brakensiek, editors. Hydrologic modeling of small watersheds. St. Joseph, MI: American Society of Agricultural Engineers; 1982 [chapter 8]. Freeze RA, Harlan RL. Blueprint for a physically-based digitally simulated, hydrologic response model. J Hydrol 1969;9:237–58. Freeze RA, Cherry JA. Groundwater. Englewood Cliffs, NJ: Prentice-Hall; 1979. Gabet EJ, Dunne T. Sediment detachment by rain power. Water Resour Res 2003;39:1002. doi:10.1029/2001WR000656. Gabet EJ, Dunne T. Correction to sediment detachment by rain power. Water Resour Res 2004;40. doi:10.1029/2004WR003422 [W08901]. Gerits JJP, de Lima JLMP, van den Broek TMW. Overland flow and erosion. In: Anderson, Burt, editors. Process studies in hillslope hydrology. Chichester, West Sussex, England: John Wiley and Sons Ltd.; 1990 [chapter 6]. Gessler J. The beginning of bedload movement of mixtures investigated as natural armouring in channels. Reprint, Laboratory of Hydraulics and Water Resources, CIT, 1965. Gessler D, Hall B, Spasojevic M, Holly F, Pourtaheri H, Raphelt N. Application of 3D mobile bed, hydrodynamic model. J Hydraul Eng 1999;125:737–49. Holly F, Rahuel J-L. New numerical/physical framework for mobile bed modeling, Part 1: Numerical and physical principles. J Hydraul Res 1990;28:401–16. Horton RE. Erosional development of streams and their drainage basins: hydrophysical approach to quantitative morphology. Bull Geol Soc Am 1945:275–370. Inbar M, Tamir M, Wittenberg L. Runoff and erosion processes after a forest fire in Mount Carmel, a Mediterranean area. Geomorphology 1998;24:17–33. James LD, Burges SJ. Selection, calibration and testing of hydrologic models. In: Haan, Johnson, Brakensiek, editors. Hydrologic modeling of small watersheds. St. Joseph, MI: American Society of Agricultural Engineers; 1982 [chapter 14]. Johansen MP, Hakonson TE, Breshears DD. Post-fire and erosion from rainfall simulation: contrasting forests with shrublands and grasslands. Hydrol Process 2001;15:2953–65. Jones JP, Sudicky EA, Brookfield AE, Park Y-J. Quantifying groundwater contributions to streamflow generation: the reliability of tracer-based hydrograph separation techniques. Water Resour Res, in review. Karim MF, Kennedy JF. IALLUVIAL: a computer-based flowand sediment-routing for alluvial streams and its application to the Missouri River. Iowa Institute of Hydraulic Research, The University of Iowa, Iowa City, IA, Report No. 250, 1982. Kinnell PIA. Raindrop-impact-induced erosion processes and prediction: a review. Hydrol Process 2005;19:2815–44. Kirkby MJ. Implications for sediment transport. In: Kirkby, editor. Hillslope hydrology. Chichester, West Sussex, England: John Wiley and Sons Ltd.; 1978 [chapter 9]. Knisel WG. CREAMS: a field scale model for chemicals, runoff and erosion from agricultural management systems. US Department of Agriculture, Conservation Research Report No. 26, 1980. Lane LJ, Shirley ED, Singh VP. Modeling erosion on hillslopes. In: Anderson, editor. Modelling geomorphological systems. Chichester, West Sussex, England: John Wiley and Sons Ltd; 1988 [chapter 10]. Li RM. Water and sediment routing from watersheds. In: Shen, editor. Modeling of rivers. New York, NY: Wiley; 1979 [chapter 9].
C.S. Heppner et al. / Advances in Water Resources 29 (2006) 930–943 [33] Loague K, VanderKwaak JE. Simulating hydrologic response for the R-5 catchment: comparison of two models and the impact of the roads. Hydrol Process 2002;16:1015–32. [34] Loague K, VanderKwaak JE. Physics-based hydrologic response simulation: platinum bridge, 1958 Edsel, or useful tool. Hydrol Process 2005;18:2949–56. [35] Loague K, Heppner CS, Abrams RH, VanderKwaak JE, Carr AE, Ebel BA. Further testing of the integrated hydrology model (InHM): event-based simulations for a small rangeland catchment located near Chickasha, Oklahoma. Hydrol Process 2005;19: 1373–98. [36] Luce CH, Black TA. Sediment production from forest roads in western Oregon. Water Resour Res 1999;35:2561–70. [37] Merritt WS, Letcher RA, Jakeman AJ. A review of erosion and sediment transport models. Environ Model Software 2003;18: 761–99. [38] Morgan RPC, Quinton JN, Smith RE, Govers G, Poesen JWA, Auerswald K, et al. The European soil erosion model (EUROSEM): a dynamic approach for predicted sediment transport from fields and small catchments. Earth Surf Process Landforms 1998;23:527–44. [39] Orvis C, Randle R. STARS: sediment transport and river simulation model. Denver, CO: Technical Service Center, US Bureau of Reclamation; 1987. [40] Panday S, Huyakorn PS. A fully coupled physically-based spatially-distributed model for evaluating surface/subsurface flow. Adv Water Resour 2004;27:361–82. [41] Pebesma EJ, Switzer P, Loague K. Error analysis for the evaluation of model performance: Rainfall-runoff time series. Hydrol Process 2005;19:1529–48. [42] Prosser IP, Rustomji P. Sediment transport capacity relations for overland flow. Prog Phys Geography 2000;24:179–93. [43] Rawls WJ, Gish TJ, Brakensiek DL. Estimating soil water retention from soil physical properties and characteristics. Adv Soil Sci 1991;16:213–34. [44] Reid LM. Research and cumulative watershed effects. PSWGTW-141, Pacific Southwest Research Station, Forest Service, US Department of Agriculture, 1993. 118p.
943
[45] Roering JJ, Gerber M. Fire and the evolution of steep, soil-mantled landscapes. Geology 2005;33:349–52. doi:10.1130/G21260.1. [46] Tayfur G. Applicability of sediment transport capacity models for nonsteady state erosion from steep slopes. J Hydrol Eng 2002; 7:252–9. [47] US Department of Agricultural, Agricultural Research Service (USDA-ARS). Hydrology, erosion, and water quality studies in the southern great plains research watershed, Southwestern Oklahoma, 1961–1978. Agricultural Reviews and Manuals, ARM-S-29, 1983. [48] VanderKwaak JE. Numerical simulation of flow and chemical transport in integrated surface–subsurface hydrologic systems. PhD dissertation, Department of Earth Sciences, University of Waterloo, Ontario, Canada, 1999. [49] VanderKwaak JE, Loague K. Hydrologic-response simulations for the R-5 catchment with a comprehensive physics-based model. Water Resour Res 2001;37:999–1013. [50] van Genuchten MT. A closed-form equation for predicting the hydraulic conductivity of unsaturated soils. Soil Sci Soc Am J 1980;44:892–8. [51] Wicks JM, Bathurst JC. SHESED: a physically based, distributed erosion and sediment yield component for the SHE hydrological modeling system. J Hydrol 1996;175:213–38. [52] Woolhiser DA, Smith RE, Goodrich DC. KINEROS, a kinematic runoff and erosion model:Documentation and user manual. Agricultural Research Service, US Department of Agriculture, Report No. ARS-77, 1990. [53] Ziegler AD, Giambelluca TW, Southerland RA. Improved method for modeling sediment transport on unpaved roads using KINEROS2 and dynamic erodibility. Hydrol Process 2002;16: 3079–89. [54] Ziegler AD, Southerland RA, Giambelluca TW. Acceleration of Horton overland flow and erosion by footpaths in an upland agricultural watershed in northern Thailand. Geomorphology 2001;41:249–62. [55] Ziemer RR, Lewis J, Rice RM, Lisle TE. Modeling the cumulative watershed effects of forest management strategies. J Environ Quality 1991;20:36–42.