SEDIMENT YIELD IN ISHITE WATERSHED USING A RUNOFF-EROSION KINEMATIC MODEL 

Celso A. G. Santos¹, Koichi Suzuki¹ and Walszon T. A. Lopes²

¹ Department of Civil and Environmental Eng., Ehime University, 3 Bunkyo-cho, Matsuyama, Ehime, 790-8577, Japan, Tel./Fax +81 89 927 9831

E-mail: celso@dpc.ehime-u.ac.jp

ksuzuki@dpc.ehime-u.ac.jp 

² Department of Civil Eng., Paraíba University, Av. Aprígio Veloso, 882, Campina Grande, Paraíba, 58109-970, Brazil, Fax: (083)310-1388

E-mail: walszon@yahoo.com 

Abstract: A kinematic distributed model named as KINEROS2 is applied for the first time in a river basin in Japan in order to simulate the runoff and sediment yield of several rainfall events during the rainy season of 1997. The physical basis of the model assures its applicability to areas with climatic differences and makes it suitable for large areas with heterogeneous reality. Thus, the model can in principle overcome many of the deficiencies of empirical and conceptual models. The results showed that the model can be used as a useful tool for such simulations in Ishite river basin as well, it can be considered a promising tool for application in other similar regions.

Keywords: erosion, surface runoff, physically-based model

1    INTRODUCTION

There are basically three categories of soil erosion models, which are empirical, conceptual or partly empirical/mixed, and physically-based model. The empirical models are based on data from field observations, mostly as standard runoff plots on uniform slopes, and are usually statistical in nature. The main limitation of this type of models is their limited applicability outside the range of conditions for which they have been developed. Their adaptation to a new environment requires a major investment of resources and time to develop the database required to drive the model. The conceptual models lie somewhere between the empirically and physically-based models, and their main limitations lie in the poor physical description of the processes which results in distortion of parameter values determined by calibration¹⁾. The physical basis of the physically-based models can, in principle, overcome many of the deficiencies of empirical and conceptual models²⁾. Therefore, a physically-based runoff-erosion model is herein applied. Since about 1970, the kinematic wave approximation has become a widely used method to simulate the movement of rainfall excess water over the land surface and through small channels. This method was applied to an arbitrary network of planes and channels in a computer model called KINGEN by Rovey et al.³⁾ This model employed a computer solution derived model for infiltration to simulate the production of runoff, and it was intended for rural or urban runoff studies using a design storm. Since then the model has been modified, with inclusion of new components, e.g., simulation of erosion and sediment transport, revision of the infiltration component and inclusion of a pond element; then, the model was called KINEROS⁴⁾. Now this model is called KINEROS2 and includes new features such as (1) the infiltration algorithm handles a two-layer soil profile and incorporates a new method, based on soil physics, to redistribute soil water during rainfall interruptions; (2) soil and sediment are characterized by a distribution of up to five particle size classes rather than a single median particle size; (3) the detention pond model accounts for seepage through the wetted area, rainfall on the pond itself, and initial storage; (4) the open channel algorithm has been extended to allow a compound cross section with an overbank level where hydraulic and infiltration parameters can differ from those in the main section; (5) base flows can be specified for open channels; (6) it automatically interpolates multi-gage rainfall intensity input to each routing element based on the spatial relationship between the element’s areal centroid and the rain gage network; and (7) time varying inflow from sources outside the model, such as measured flows or output from other models, can be injected into the model system at any point.

The following sections describe the basic equations of the KINEROS2 model and present the simulation results for the rainy season of 1997 in the Ishite river basin in Japan, which is not a large basin if compared with other basins around the world, but for the scale of Japan and considering basins without influence of urban areas, it can be considered as a large river basin.

2    THE KINEMATIC MODEL

2.1    The soil infiltration model

The model formulation allows a physically-based approximation for the redistribution of soil water, including recovery of infiltration capacity during a hiatus, and a method that determines infiltration rates following a hiatus. The soil infiltration model describes infiltration capacity f_c as a function of infiltrated depth I and it requires four basic parameters to describe the infiltration properties of a soil, which are the field effective saturated hydraulic conductivity K_s (m/s), the integral capillary drive G (m), the porosity f, and the pore size distribution index l. There is an optional parameter C_v, which describes the random variation in space of the hydraulic properties of the soil, as well as another optional parameter that allows even more explicit characterization of a soil profile, the content of large rocks, ROCK. Also, there is an event-dependent variable, which is the initial relative saturation of the upper soil layer, equal to q_i/f, in which q_i is the initial soil moisture content. The general model for infiltrability, f_c (m/s), as a function of infiltrated depth, I (m), is given by⁵⁾:

                            (1)

where B is (G + h)(q_s – q_i), combining the effects of net capillary drive, G, surface water depth, h (m), and unit storage capacity, Dq = (q_s – q_i), in which q_s is the soil moisture content at saturation. The parameter a represents the soil type; i.e., a is near 0 for a sand, in which case Eq. (1) approaches the Green-Ampt relation; and a is near 1 for a well-mixed loam, in which case Eq. (1) represents the Smith-Parlange infiltration equation. In this model, there is a soil water redistribution; i.e., it considers the hiatus between rainfall events which consist of more than one period of runoff-producing rainfall, with an intervening period during which significant drying of the soil can occur. The redistribution/reinfiltration method used in the model is described in Smith et al.⁶⁾ and Corradini et al.⁷⁾

2.2    Overland and channel flow

Overland flow can be viewed as a one-dimensional flow process in which flux is related to the unit area storage by a simple power relation:

                              (2)

where Q is discharge per unit width (m²/s), and h is the storage of water per unit area (m). Parameters a and m are related to slope, surface roughness, and flow regime, and are given by a = S^1/2/n and m = 5/3 where S is the slope, n is a Manning’s roughness coefficient for overland flow. The continuity equation for a plane is given as:

                            (3)

where t is time (s), x is the distance along the slope direction (m), and q(x,t) is the lateral inflow rate (m/s). For overland flow, Eq. (2) is substituted into Eq. (3) and the resultant equation is solved using a four-point implicit finite difference method. The continuity equation for a channel with lateral inflow is given by:

                            (4)

where A is the cross-sectional area (m²), Q is the channel discharge (m³/s), and q_c(x,t) is the net lateral inflow per unit length of channel (m²/s). The kinematic assumption is embodied in the relationship between channel discharge and cross-sectional area such that

                              (5)

where R is the hydraulic radius (m). If the Manning equation is used, a = S^1/2/n and m = 5/3. The kinematic equations for channels are solved by a four-point implicit technique similar to that for overland flow surfaces.

2.3    Upland and channel erosion

The general equation used to describe the sediment dynamics at any point along a surface flow path is a mass balance equation similar to that for kinematic water flow⁸⁾:

                    (6)

in which C_s is the sediment concentration (m³/m³), Q is the water discharge rate (m³/s), A is the cross sectional area of flow (m²), e is the rate of erosion of the soil bed (m²/s), q_s is the rate of lateral sediment inflow for channels (m²/s). For upland surfaces, e is assumed to be composed of two major components, i.e., production of eroded soil by splash of rainfall on bare soil, and hydraulic erosion (or deposition) due to the interplay between the shearing force of water on the loose soil bed and the tendency of soil particles to settle under the force of gravity. Net erosion is a sum of splash erosion rate as e_s and hydraulic erosion rate as e_h:

                                   (7)

The splash erosion rate is estimated as follows:

                                 (8)

in which r is the effective rainfall (m/s), c_f is a constant related to soil and surface properties, and  is a reduction factor representing the reduction in splash erosion caused by increasing depth of water. The parameter c_h represents the damping effectiveness of surface water, assumed to be 364.0. The hydraulic erosion rate (e_h) is estimated as being linearly dependent on the difference between the equilibrium concentration and the current sediment concentration and is given by:

                             (9)

in which C_m is the concentration at equilibrium transport capacity, C_s = C_s(x,t) is the current local sediment concentration, and c_g is a transfer rate coefficient (s^-1), which is computed as:

   if C_s £ C_m (erosion)  or        if C_s > C_m (deposition)     (10)

Where Co is the a soil cohesion coefficient, and v_s is the particle fall velocity (m/s). The model uses the transport capacity relation of Engelund and Hansen⁹⁾ modified to include the unit stream power threshold found to apply to shallow flow transport. Particle settling velocity is calculated by the following equation¹⁰⁾:

                             (11)

in which g is the gravitational acceleration (m/s²), r_s is sediment specific gravity, equal to 2.65, d is the sediment diameter (m), and C_D is the particle drag coefficient. The drag coefficient is a function of particle Reynolds number:

                           (12)

Where R_n is the particle Reynolds number, given as R_n = v_sd/n, in which n is the kinematic viscosity of water (m²/s). Settling velocity of a particle is found by solving Eqs (11), and (12) for v_s. The above series of erosion relations are applied to each of up to five particle size classes, which are used to describe a soil with a range of particle sizes. Eqs (6-12) are solved numerically at each time step used by the surface water flow equations, and for each particle size class.

The general approach to sediment transport simulation for channels is nearly the same as that for upland areas. The major difference in the equations is that splash erosion (e_s) is neglected in channel flow, and the term q_s becomes important in representing lateral inflows.

3    FIELD DATA

4    CALIBRATION AND SIMULATIONS

Considering one-layer soil profile and the fact that the soil type in Ishite river basin is sandy clay loam, several parameters were assumed as follows: G = 260 mm, f = 0.398, l = 0.32, q_s = 0.330, a = 0.85 and C_v = 3. Other parameters were set according to the basin characteristics: Manning’s n = 1.00 for the planes and 0.04 for the channel, interception = 5.0 mm, Canopy = 50%, ROCK = 30%, and K_s  = 1.5 mm/h. The initial degree of soil saturation, given by q_i/f, was optimized by trial and error to fit the measured runoff volume, and the parameter concerning the splash erosion c_f and the soil cohesion coefficient Co were optimized to fit a value coherent with the potential sediment yield in such a period. Thus, q_i = 0.080, c_f = 300,000, and Co = 0.01.

The simulated direct hydrograph and sedigraph are shown in Fig. 3, in which the sediment yield is computed just for one particle size class (0.5 mm) with total of 18,000 m³, which correspond to 27.8% of the total annual production. The computations show satisfactory agreement with the observed values. In order to quantify the sediment yield according to the particle size classes, computations were done for three particle size classes of 0.5 mm, 0.25 mm and 0.05 mm with 30%, 30% and 40%, respectively, as shown by the sedigraphs in Fig. 4.

5    CONCLUSIONS

In order to compute the sediment yield for Ishite river basin in Matsuyama city, Japan, a physically-based runoff-erosion model called KINEROS2 was used. The computations showed that the model can be used to simulate the runoff due to continuous rainfall events. The physical basis of the model allowed its applicability to such an area, which has climatic differences if compared with the basin used in the development of the model. The simulation results are satisfactory when compared with the observed hydrograph, and the sediment yield during the rainy season can be assumed to be equal to 30% of the annual sediment yield, i.e., 18,000 m³. 

Acknowledgements

The writers wish to thank Dr. Carl Unkrich of Southwest Watershed Research Center for providing the program source. The field data pro­vided by Matsuyama Construction Works Office are gratefully acknowledged.

References

[1]    Lørup. J. K. and Styczen, M.: Soil erosion modelling. In: Abbott, M.B. and Refsgaard, J.C. (eds.), Distributed Hydrological Modelling, Kluwer Academic Publishers, 93-120, 1996.

[2]    Santos, C.A.G., Watanabe, M., Suzuki, K. and Srinivasan, V.S.: Sediment yield due to heavy rainfall from a test field in Brazil and its analysis by a runoff-erosion model. Journal of Hydraulic, Coastal and Environmental Engineering, JSCE, No.586/II-42, 117-126, 1998.

[3]    Rovey, E.W., Woolhiser, D.A. and Smith, R.E.: A distributed kinematic model of upland watershed. Hydrology Paper 93, 52 pp. Colorado State University, Fort Collins, 1977.

[4]    Woolhiser, D.A., Smith, R.E. and Goodrich, D.C.: KINEROS, A kinematic runoff and erosion model: Documentation and User Manual. U.S. Department of Agriculture, Agricultural Research Service, ARS-77, 130 pp., 1990.

[5]    Parlange, J.-Y., Lisle, I., Braddock, R.D. and Smith, R.E.: The three-parameter infiltration equation. Soil Science, 133(6), 337-341, 1982.

[6]    Smith, R.E., Corradini, C. and Melone, F.: Modeling infiltration for multistorm runoff events. Water Resources Research, 29(1), 133-144, 1993.

[7]   Corradini, C., Melone, F. and Smith, R.E.: Modeling infiltration during complex rainfall sequences. Water Resources Research, 30(10), 2777-2784, 1994.

[8]    Bennett, J.P.: Concepts of mathematical modeling of sediment yield. Water Resources Research, 10(3), 485-492, 1974.

[9]    Engelund, F. and Hansen, E.: A Monograph on Sediment Transport in Alluvial Streams, Teknisk Forlag, Copenhagen, 62 pp., 1967.

[10]     Fair, G.M. and Geyer, J.C.: Water Supply and Wastewater Disposal. John Wiley and  Sons, New York, 973 pp., 1954.

Fig. 1    Map of the Ishite river basin.        Fig. 2    Observed hyetograph and hydrograph

                                             including the baseflow for 22 June to

                                             26 July 1997.

Fig. 3    Computed direct runoff and sediment yield, in which the bold line corre-sponds to the direct runoff and the thin line to the total sediment yield.

Fig. 4    Computed sediment yield, in which the old line corre-sponds to the total sediment yield and the thin lines are the sediment yield for three particle size classes (0.5 mm, 0.25 mm and 0.05 mm).to the total sediment yield.