The stochastic runoff-runon process: Extending its analysis to a finite hillslope

Journal of Hydrology xxx (2016) xxx–xxx

Contents lists available at ScienceDirect

Journal of Hydrology journal homepage: www.elsevier.com/locate/jhydrol

The stochastic runoff-runon process: Extending its analysis to a finite hillslope O.D. Jones a,⇑, P.N.J. Lane b, G.J. Sheridan b a b

School of Mathematics and Statistics, University of Melbourne, VIC 3010, Australia School of Ecosystem and Forest Sciences, University of Melbourne, VIC 3010, Australia.

a r t i c l e

i n f o

Article history: Received 31 July 2014 Received in revised form 20 June 2016 Accepted 25 June 2016 Available online xxxx This manuscript was handled by K. Georgakakos, Editor-in-Chief, with the assistance of Venkat Lakshmi, Associate Editor Keywords: Infiltration excess runoff Overland flow Stochastic runoff-runon process Single server queue Finite hillslope

a b s t r a c t The stochastic runoff-runon process models the volume of infiltration excess runoff from a hillslope via the overland flow path. Spatial variability is represented in the model by the spatial distribution of rainfall and infiltration, and their ‘‘correlation scale”, that is, the scale at which the spatial correlation of rainfall and infiltration become negligible. Notably, the process can produce runoff even when the mean rainfall rate is less than the mean infiltration rate, and it displays a gradual increase in net runoff as the rainfall rate increases. In this paper we present a number of contributions to the analysis of the stochastic runoff-runon process. Firstly we illustrate the suitability of the process by fitting it to experimental data. Next we extend previous asymptotic analyses to include the cases where the mean rainfall rate equals or exceeds the mean infiltration rate, and then use Monte Carlo simulation to explore the range of parameters for which the asymptotic limit gives a good approximation on finite hillslopes. Finally we use this to obtain an equation for the mean net runoff, consistent with our asymptotic results but providing an excellent approximation for finite hillslopes. Our function uses a single parameter to capture spatial variability, and varying this parameter gives us a family of curves which interpolate between known upper and lower bounds for the mean net runoff. Ó 2016 Elsevier B.V. All rights reserved.

1. Introduction The volume of catchment discharge that reaches a stream via the overland flow path is critical for water quality prediction, because it is via this pathway that most particulate pollutants are generated and transported to the stream channel, via surface erosion processes. Two of the key properties determining this volume are the rainfall rate and the infiltration rate. In natural systems both these rates are variable in both space and in time. Suppose that our hillslope is divided into cells. If the rainfall rate exceeds the infiltration rate in a given cell, then the excess will flow overland to the next cell downhill. Thus the water flowing into a cell is given by the sum of the rainfall and runon from the cell above. Any excess, after infiltration is taken into account, becomes runoff. The resulting system is highly non-linear, because runoff is truncated below at zero. Nahar (2003) showed that for soils with moderate to high mean saturated conductivity relative to rainfall rate, the runoff-runon process plays an important part in determining the total overland discharge for a hillslope. These ⇑ Corresponding author. E-mail address: [email protected] (O.D. Jones).

conditions are typical in temperate forests, where saturated conductivity values are usually high, and are common in many other landscapes for the majority of rainfall events (Dunkerley, 2008). Because of the complexity of the problem, models that incorporate both spatial and temporal variability have, to date, been analysed using numerical simulation methods. Our interest is in analytic solutions. The most common simplification made in this context is to neglect spatial variability and model rainfall and infiltration as a function of time only. This can be attributed to the early development of analytical expressions for the temporal change in infiltration rate at a point (Green and Ampt, 1911). For catchment scale predictions these point-scale results have generally been scaled up by optimizing the infiltration parameters using catchment or hillslope runoff time-series data. As a result of this scaling process, the parameters lose their physical meaning (e.g. see discussion by Grayson et al., 1992). A recent alternative is the stochastic runoff-runon process introduced by Jones et al. (2009) and developed in Jones et al. (2013) and Harel and Mouche (2013, 2014). The stochastic runoff-runon process allows for spatial variability but assumes temporal stationarity. It does however admit analytic asymptotic solutions, with parameters that retain their physical meaning. In this paper we pay particular attention

http://dx.doi.org/10.1016/j.jhydrol.2016.06.056 0022-1694/Ó 2016 Elsevier B.V. All rights reserved.

Please cite this article in press as: Jones, O.D., et al. The stochastic runoff-runon process: Extending its analysis to a finite hillslope. J. Hydrol. (2016), http:// dx.doi.org/10.1016/j.jhydrol.2016.06.056

2

O.D. Jones et al. / Journal of Hydrology xxx (2016) xxx–xxx

to how the model behaves on finite hillslopes, when the previously obtained asymptotic solutions are not available. Hawkins and Cundy (1987) were the first to propose an analytic solution to the net runoff generation problem incorporating variability in the spatial dimension. Hawkins and Cundy showed that for an area with constant rainfall and spatially variable infiltration there exist maximum and minimum curves relating the net runoff rate to the precipitation rate; see Fig. 1. The curves are derived by arranging the point infiltration values from largest to smallest, or vice versa. The true (but generally unknown) function relating precipitation rate to net runoff rate must lie between these enveloping curves. Assuming that the distribution of infiltration rates has an exponential density with mean mi and that the rainfall is constant with rate mp (both in mm/h), the minimum net runoff rate is mp
mi , for mp P mi , and the maximum is

mp
mi ð1
e
mp =mi Þ:

ð1Þ

Note that the function depends on precipitation rate and not time, as the system is assumed to be in temporal equilibrium. We also see that runoff is generated even when the precipitation rate is lower than the average infiltration rate, and it increases gradually as the precipitation rate increases. These are characteristic consequences of including spatial variability. The Hawkins and Cundy model has not received widespread attention, despite the fact that Yu and others have reported considerable success using the maximum net runoff curve as the basis of a runoff model at the plot scale (Yu et al., 1997, 1998; Yu, 1999; Fentie et al., 2002; Kandel et al., 2005). This approach was found to perform better than the time-variant, space-invariant Green and Ampt (1911) model for the prediction of infiltration excess runoff at the plot scale (Yu, 1999). One of the main contributions of the present paper is the derivation of a family of curves for the mean net runoff, which smoothly interpolate between the upper and lower bounds of Hawkins and Cundy; see Section 5.1, Eq. (19). As a consequence of modelling runoff from each cell, the stochastic runoff-runon process can be used to analyse hillslope connectivity, whereby adjacent cells are connected if there is runoff from one to the other. In this context the process has been called the stochastic runoff connectivity (SRC) process (Sheridan et al., 2009a,b). Harel and Mouche (2014) also consider connectivity through the lens of the stochastic runoff-runon process, extending the model to include lateral diffusion of runoff. We do not consider connectivity explicitly in this paper, though note that in Section 3.2 we do draw conclusions about the proportion of the hillslope that contributes to net runoff.

The structure of the paper is as follows. In Section 2 we give a definition of the (time-stationary one-dimensional) stochastic runoff-runon process, and then fit it to some experimental data using maximum likelihood. The fit is quite good, lending credence to the model. In Section 3 we present an analysis of the asymptotic properties of the model, as the length of the hillslope tends to infinity. Three regimes emerge, depending on whether the mean rainfall rate is less than, equal to, or greater than the mean infiltration rate. We refer to these regimes as subcritical, critical and supercritical respectively. Then in Section 4 we use Monte Carlo simulation to quantify the scales at which asymptotic results can be used to approximate runoff behaviour on finite hillslopes. In Section 5 we develop a function for the mean runoff rate on a finite hillslope, as a function of the mean rainfall rate, which spans all three regimes, subcritical, critical and supercritical. We compare our function for the mean runoff rate to that of Hawkins and Cundy, and show that as you increase the spatial variability of the infiltration rate, our function smoothly transitions from their lower bound to their upper bound. A discussion and summary of our results are given in Sections 6 and 7. 2. The stochastic runoff-runon process The stochastic runoff-runon process is a stochastic timeinvariant model for the flow of infiltration-excess runoff down a planar hillslope. We model the hillslope as a series of parallel and independent strips perpendicular to the bottom edge. Each strip is broken up into a line of blocks or cells of equal size, and we suppose that within each block the rate of rainfall and infiltration are fixed. The flow of runoff from one cell to the next downslope can be considered a stochastic process spatially indexed by the position of the cell down the strip. The runoff-runon model is constructed in two steps. Firstly we consider the runoff generation down a single strip of land from the top to the bottom of the hillslope, perpendicular to the contours, with a random arrangement of rainfall and infiltration capacity down its length. Analysis of this component of the model draws on queuing theory. Next, we consider the properties of the aggregated output from these strips. Analysis of this component of the model depends on the central limit theorem. We can use the central limit theorem because we assume that runoff is confined within strips, so that the net runoff from individual strips is independent. This follows from our assumption that the hillslope is ‘‘planar”, that is, the contours are parallel.

30 10

20

mean infiltration rate 50 mm/h maximum net runoff rate minimum net runoff rate

0

mean net runoff rate mm/h

40

2.1. Single strip runoff-runon model

0

20

40

60

80

rainfall rate mm/h Fig. 1. Bounds on the relationship between precipitation rate (mm/h) and net runoff rate (mm/h), shown for the case of a mean infiltration rate of 50 mm/h. Modified from Hawkins and Cundy (1987).

We consider a single strip of land, width lx , divided into blocks of length ly . Number the blocks 1; 2; . . . ; n, starting at the top of the slope. Let X k be the rate at which water runs from block k to block k þ 1, that is, the flow from block k to k þ 1, in m3 h
1. Let pk be the precipitation (rainfall) rate and ik the infiltration rate for block k, (both are fluxes, measured in mm h
1), assumed to be constant over time. If pk > ik then on average runoff builds up down the length of the block, conversely if pk < ik then on average the runoff declines. Let the depth of water at the end of block k be dk and let its speed be v k , then the volume of water leaving the block per unit time is lx dk v k ¼ X k . If v k were constant, then we would have dk / ly . Let P k ¼ lx ly pk =1000 be the flow of rain falling on block k, and let Ik ¼ lx ly ik =1000 be the maximum flow of water absorbed by block k (both in m3 h
1). Here we assume Pk represents incident rainfall if there is no canopy or over-story, or through-fall if there is an over-storey. If we assume that there is no significant runoff onto

Please cite this article in press as: Jones, O.D., et al. The stochastic runoff-runon process: Extending its analysis to a finite hillslope. J. Hydrol. (2016), http:// dx.doi.org/10.1016/j.jhydrol.2016.06.056

3

O.D. Jones et al. / Journal of Hydrology xxx (2016) xxx–xxx

our given strip from neighboring strips, then we have the runoffrunon equation:


Xk ¼

maxf0; P1
I1 g;

k¼1

maxf0; X k
1 þ Pk
Ik g; k > 1

:

ð2Þ

Eq. (2) can also be found in stochastic queuing theory, where it gives the waiting time in a single server first-in first-out (FIFO) queue. If we let P k be the (random) service time for customer k and let Ik be the (random) inter-arrival time between customers k and k þ 1, then X k is the (random) waiting time for customer k þ 1, that is, the time between arriving and service commencing. From the queuing literature we have results on the asymptotic distribution of X k as k ! 1. Most of the standard results in this area can be found, for example, in Asmussen (2003) and Prabhu (1998) or Kleinrock (1975a,b). To obtain analytic results we need to make a number of assumptions about the random variables fIk gnk¼1 and fP k gnk¼1 . Firstly we suppose that they are time invariant, that is pk and ik , and thus P k and Ik , are independent of t. Next we suppose that they are mutually independent and identically distributed (iid) sequences. For the infiltration this implies that there is no spatial correlation in the infiltration rate at the scale of blocks used in the model. Small-scale spatial correlation in the water infiltration rate ik has been observed (Loague and Gander, 1990; Buttle and House, 1997), thus the validity of this assumption depends on the scales lx and ly : we need lx and ly large enough that the correlation between Ik and Ikþ1 is negligible. Similarly we are assuming that there is no spatial correlation in the rainfall values at the scale of blocks used in the model. As for the infiltration, small-scale spatial correlation in the rainfall rate pk has been observed (Loustau et al., 1992; Whelan and Anderson, 1996), so to justify this assumption we again need lx and ly to be relatively large, compared to the rainfall correlation scale. Independence between the P k and Ik is clearly perfectly reasonable. A more subtle assumption of the model is that the infiltration rate ik is independent of surface water depth dk . This assumption means that the effect of the infiltration rate in cell k appears only via Ik . For illustrative simulations of the stochastic runoff-runon process the reader is referred to Jones et al. (2013) and Harel and Mouche (2013). In particular, Harel and Mouche consider the effect of different infiltration distributions on the asymptotic distribution of runoff. 2.2. Field experiment We demonstrate the suitability of the stochastic runoff-runon process by fitting it to experimental data from a rainfall simulator experiment conducted in September 2005 in the East Kiewa Research catchment in NE Victoria, Australia. The Eucalyptus forested catchment is described by Lane et al. (2006) and was partially burnt by wildfire in early 2005. The rainfall simulator delivers raindrops from three oscillating fan-shaped sprays at 1 m centres sweeping intermittently across the plot surface and is based on an original design by Bubenzer and Meyer (1965), which was modified by Loch (1989). Water was supplied to the rainfall simulator from a 1000 L tank filled at nearby forest streams. Rainfall was applied at nominal rates ranging from 100, 150 and 250 mm h
1 to plots with a slope of 25° 1.5 m wide with minimal cross slope and plot lengths of 0.125, 0.25, 0.5, 1.0, and 2.0 m. Plots were bordered using steel sheets hammered into the ground, and overland flow was collected at the downslope end of the plot in stainless steel troughs, sealed to the soil surface with petroleum jelly. Timed runoff samples were collected in 500 mL plastic bottles for measurement of the runoff rate, with several samples collected for each combination of rainfall rate and plot length to establish

the steady-state runoff rate for a given set of conditions. The rainfall and runoff data are given in Table 1. For this data we have 15 independent observations, each corresponding to a single strip. The data are given as fluxes, that is mm h
1, and we will work with them as such. To get flows in m3 h
1 they just need to be scaled by lx ly =1000. For each strip the rainfall fpk gnk¼1 are constant, and we used a gamma ða; bÞ distrin bution for the infiltrations fik gk¼1 , where the shape parameter a and the rate parameter b are fixed over all 15 strips. Suppose that pk ¼ r for a particular strip, then the distribution of xn ¼ 1000X n =ðlx ly Þ, the runoff flux from the nth block, is a mixed distribution with an atom at 0 and continuous on ð0; nrÞ. Let g be n the gamma density function for the fik gk¼1 , then from (2) we obtain the distribution function F n of xn as follows

F 1 ðx; rÞ ¼

8 > <0 R > :

1 r
x

1

x<0 gðyÞdy 0 6 x < r xPr

8 0 > > > < R 1 F ðy
r þ x; rÞgðyÞdy n
1 F n ðx; rÞ ¼ Rr
x 1 > F ðy
r þ x; rÞgðyÞdy n
1 > > : 0 1

x<0 06x
ð3Þ

x P nr

We fit a and b by maximum likelihood. Let nðiÞ be the number of blocks for the ith observation, rðiÞ the rainfall flux, and xðiÞ the observed value of xnðiÞ . Then the log likelihood is

lða; b; nð1Þ; rð1Þ; xð1Þ; . . . ; nð15Þ; rð15Þ; xð15ÞÞ X X ¼ log F nðiÞ ð0; rðiÞÞ þ log F 0nðiÞ ðxðiÞ; rðiÞÞ: i:xðiÞ¼0

ð4Þ

i:xðiÞ>0

Note that we get separate terms for the discrete and continuous parts, and that in the second part we use the derivative of F n ðx; rÞ w.r.t. x. In our case the first term is empty. The F n were calculated numerically, using Simpson’s method. Note that gðyÞ behaves like ya
1 at 0, so if a P 1 then g is bounded and the integrals are all well behaved. However for a < 1 g is unbounded, in which case we used R F ðy
r þ x; rÞgðyÞdy  F n
1 ðx
r; rÞGðÞ, the approximation 0 n
1 where G is the integral of g (that is, the gamma distribution function). To calculate F 0n ðx; rÞ for x > 0 one can differentiate the RHS of (3) to obtain a similar equation for F 0n ðx; rÞ in terms of F 0n
1 ð; rÞ and F n
1 ð0; rÞ. However if a < 1 then F 0n inherits singularities from g and as n increases this quickly becomes unmanageable. Our solu-

Table 1 A summary of runoff-rate data from rainfall simulation experiments on 1.5 m wide plots conducted at a range of rainfall intensities and for different slope lengths. Strip length (m) 0.125 0.125 0.125 0.25 0.25 0.25 0.5 0.5 0.5 1 1 1 2 2 2

Num blocks n

Rainfall r (mm/h)

Final runoff xn (mm/h)

1 1 1 2 2 2 4 4 4 8 8 8 16 16 16

95.31 159.19 269.86 95.31 159.19 237.89 92.82 203.68 235.69 97.32 194.76 267.18 96.85 187.98 282.52

61.64 142.56 240.00 88.26 206.06 279.06 109.72 276.12 437.20 136.24 339.28 520.40 150.08 321.12 550.56

Please cite this article in press as: Jones, O.D., et al. The stochastic runoff-runon process: Extending its analysis to a finite hillslope. J. Hydrol. (2016), http:// dx.doi.org/10.1016/j.jhydrol.2016.06.056

In this section we first consider a single strip and look at the distribution of X n as n ! 1. We then consider the runoff from a hillslope made up of m strips, as m ! 1. Some of these results appeared previously in Jones et al. (2009, 2013), but are included here as part of a larger picture that includes supercritical and critical as well as the previously considered subcritical regimes. 3.1. Asymptotic properties of the stochastic runoff-runon model for a single strip Let mP ¼ EP k , mI ¼ EIk and q ¼ mP =mI (q is known as the ‘‘traffic intensity” in queuing theory). The asymptotic behaviour of X n as n ! 1 changes as q is greater, equal to, or less than 1. In what follows we write P and I for a generic Pk and Ik respectively. The quantities in Eqs. (5)–(9) below are all flows, measured in m3 h
1 and with mean and standard deviation proportional to lx ly . 3.1.1. Case 1: rainfall dominates infiltration q > 1 (supercritical) In this case the mean and variance of the runoff both grow linearly. If r2 ¼ VarP þ VarI < 1 then (Asmussen, 2003 Chapter X Proposition 1.2, Prabhu, 1998 Chapter 1 Theorem 10)

where Z  Nð0; 1Þ. Note that mP
mI / lx ly .

15

-

-

10

runoff

-

-

-

-

----0.0

0.5

1.0

--

-- -

1.5

2.0

-

-

--

--

2.5

rainfall Fig. 2. Observed values of runoff flux xn plotted against rainfall flux r for a strip of length n. Observed values are given by the circles, 90% prediction intervals are given by the vertical bars, and medians are given by the horizontal ticks. Note that the axes have been scaled by a factor of 100.

3.1.2. Case 2: balanced infiltration and rainfall q ¼ 1 (critical) pffiffiffi In this case the mean runoff grows like n, but the variance grows like n. If r2 ¼ VarP þ VarI < 1 then (Asmussen, 2003 Chapter X Proposition 1.2)

Xn pffiffiffi ! jrZj in distribution n where Z  Nð0; 1Þ. Note that jrZj has mean ance r2 ð1
2=pÞ.

3. Asymptotic distribution of runoff

Xn ! mP
mI with probability 1 n X n
nðmP
mI Þ pffiffiffi ! rZ in distribution n

n 1 2 4 8 16

5

tion in this case is just to use a numerical derivative of F n . Maximising (4) numerically we obtained a ¼ 0:97 and b ¼ 0:0056, corresponding to a point infiltration capacity with mean 173 mm h
1 and standard deviation 176 mm h
1. To illustrate the fit we calculated 90% prediction intervals for each observation, which are plotted with the observations in Fig. 2. A good fit is indicated, with 14 of the 15 observations falling within these prediction intervals. However, we do note that from the fitted model we would expect some zero observations. This discrepancy can in fact be seen as a discretisation effect, due to treating the infiltration as constant over each block. Taking smaller blocks would reduce this effect, but would introduce the problem of spatial dependence between blocks. We remark that since a ¼ 0:97, the fitted infiltration distribution is very close to the exponential. This observation is consistent with independently collected infiltration data from similar soils. Nyman (2007) reports measurements of the spatial distribution of saturated conductivity K sat (mm h
1) in freshly burnt soils, with water repellence and ash deposits, that show an exponential distribution. Elsewhere K sat has been measured for many different soils and is widely reported as log-normal (Neilsen et al., 1973; Sharma et al., 1980; Brakensiek et al., 1981; Luxmoore, 1981; Watson and Luxmoore, 1986; Wilson and Luxmoore, 1988; Loague and Gander, 1990; Price, 1994; Loague and Kyriakidis, 1997). In this experiment the rainfall does not vary spatially, and this is presumably appropriate for open areas at scales relevant to infiltration excess overland flow. However, the partial interception of rainfall by vegetation (such as a forest canopy) can alter the incident rainfall distribution, and this throughfall has been reported as normally distributed. The coefficient of variation of throughfall appears to decrease exponentially with the total rainfall for the storm, from around 40 to 50% for storms of < 5 mm, rapidly asymptoting to about 8% for storms of > 5 mm (Llorens et al., 1997; Carlyle-Moses et al., 2004).

20

O.D. Jones et al. / Journal of Hydrology xxx (2016) xxx–xxx

0

4

ð6Þ pffiffiffiffiffiffiffiffiffi

r 2=p / lx ly and vari-

3.1.3. Case 3: infiltration dominates rainfall q < 1 (subcritical) If q < 1 then as n ! 1, X n converges in distribution to some limit X. That is, the cumulative distribution function (cdf) F n of X n converges to the cdf F of some random variable X. We write F n ðxÞ ¼ PðX n 6 xÞ ! FðxÞ ¼ PðX 6 xÞ. Formally, imagine a strip of infinite length, and choose a random location to inspect the level of runoff from one block to the next. The (randomly) sampled value will have the distribution (or law) F. The form of F is not known in general, however in the special case where P  expðlÞ and I  expðkÞ, where l ¼ 1=mP , k ¼ 1=mI and q ¼ k=l, we have the following exact result

FðxÞ ¼ 1
qe
ðl
kÞx :

ð7Þ

That is, PðX ¼ 0Þ ¼ 1
q and XjX > 0  expðl
kÞ. If I  expðkÞ but the distribution of P is allowed to be general, then the famous Pollaczek-Khintchine formula gives the characteristic function of X. In theory we can invert the characteristic function to obtain the distribution F of X, however in practice it is not usually possible to do this analytically and we have to resort to numerical inversion techniques. None-the-less we can obtain the mean and variance of X directly from the characteristic function (Kleinrock, 1975b Vol. II, S2.3):

EX ¼

  E P2k

2ðmI
mP Þ

ð5Þ VarX ¼ ðEXÞ2 þ

  E P3k

3ðmI
mP Þ

Please cite this article in press as: Jones, O.D., et al. The stochastic runoff-runon process: Extending its analysis to a finite hillslope. J. Hydrol. (2016), http:// dx.doi.org/10.1016/j.jhydrol.2016.06.056

5

O.D. Jones et al. / Journal of Hydrology xxx (2016) xxx–xxx

Fig. 3. Monte Carlo estimates and approximations of the limiting mean net runoff. In both left and right plots the rainfall intensity is (truncated) normal with mean q and coefficent of variation (CV) equal to 0.1. In the left plot the infiltration is gamma distributed with mean 1 and CV equal to 0.75, 1.0, and 1.25. In the right plot the infiltration is log-normal distributed, again with mean 1 and CV 0.75, 1.0, and 1.25. In each plot the Monte Carlo estimates of the mean are given with error bars corresponding to 2 standard deviations. The dashed lines are the approximation of Marchal (1976) (Eq. (9)), and the solid green lines are the approximation of Kramer and Lagenbach-Belz (1976) (Eq. (8)). Note that the results have been scaled by 2ð1
qÞ. Also note that when the CV for the infiltration is 1, the two approximations coincide. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

If we allow the infiltration to have a general (non-negative) distribution, then we don’t have any exact results for the mean or variance of X. There are however bounds available (Daley et al., 1992; Rao and Feldman, 2001) and also some approximations (Kramer and Lagenbach-Belz, 1976; Marchal, 1976; Bhat, 1993; Whitt, 1993). To assess the usefulness of these bounds/ approximations, we compared them to Monte Carlo estimates obtained using both a gamma and a log-normal distribution for the infiltration, with parameters suggested by the data given in Section 2.2, and a truncated normal distribution for the rainfall. We found the best fits were given by the mean approximation of Kramer and Lagenbach-Belz (1976) (see Bhat, 1993), and the variance approximation of Bhat (1993) Eq. (6). Both approximations are exact when the infiltration is exponential. Let cP and cI be the coefficient of variation of P and I respectively, then we have

8 2
2ð1
qÞð1
c2 Þ >   > < 3q c2 þc2 I 2 2 2 mP cI þ cP g ð I PÞ EX  where log g ¼ 2 > 2ðmI
mP Þ > :
ð1
2qÞðcI2
1Þ cI þ4cP

 VarX 

c2I m2P þ r2P 2ðmI
mP Þ

þ

A summary of our assessments of the formulas given in Eqs. (8) and (9) is given in Fig. 3 (for the mean) and Fig. 4 (for the standard deviation).

3.2. Asymptotic properties of the stochastic runoff-runon model at the hillslope scale Eqs. (5)–(9) describe the asymptotic properties of a single strip. We consider now the aggregated asymptotic properties of a tilted plane, consisting of many parallel linear strips. A plane can be represented as a collection of adjacent strips placed perpendicular to a downslope boundary. The proposed plane model depends on the central limit theorem and assumes that adjacent strips are independent and that there are no lateral inflows or outflows from a strip. Let X ðiÞ be the runoff flow from the ith strip. By the central limit

if

c2I

<1

if

c2I

P1

theorem, the sum Z ¼ X ð1Þ þ X ð2Þ þ    þ X ðmÞ of m independent and

2

  EðP
mP Þ3 þ m3P 3r4I
mI EðI
mI Þ3 =m4I þ 3c2I mP r2P

identically distributed (iid) random variables X ðiÞ (of any distribution) with mean l and finite variance r2 , is approximately normally distributed as m becomes large, with mean ml and variance mr2 . In our case we suppose that we have a hillslope made up of m independent strips of width lx and length Ly :¼ nly . Let Lx ¼ mlx be the width of the hillslope. The distribution of the total runoff Z at the lower boundary of the will be

þ

3ðmI
mP Þ ð8Þ

Z¼

m 2 X   l 2 r X ðiÞ  N mlX ; mr2X ¼ N Lx ly X ;Lx lx ly 2 X2 lx ly lx ly i¼1

!
1

in m3 h ;

ð10Þ

The variance approximation above is in fact a minor modification of that of Bhat: the factor 3r4I
mI EðI
mI Þ3 has been replaced by its positive part. This was found to improve the approximation when the skewness of I is large. A simpler approximation to the mean, which none-the-less does quite well, was given by Marchal (1976) (see Kleinrock, 1975b Vol II, S2.3). It also is exact when the infiltration is exponential.

EX 

1 þ c2P



r2I þ r2P

ð1=qÞ þ c2P 2ðmI
mP Þ 2

ð9Þ

where lX and r2X are the mean and variance of the runoff flow, given in Eqs. (5)–(9). As before we take cases. In all cases we will assume that r2 ¼ VarP þ VarI < 1, and we stress that these approximations are valid as the number of blocks per strip n ! 1 and the number of strips m ! 1. The capacity to represent the net runoff from a tilted plane simply from the mean and variance of the runoff from a single strip is fortunate, because the entire distribution of the runoff rate for a strip is not available in general. Fig. 5 illustrates lx , ly , Lx , Ly , n and m in the context of a hillslope divided into strips and cells.

Please cite this article in press as: Jones, O.D., et al. The stochastic runoff-runon process: Extending its analysis to a finite hillslope. J. Hydrol. (2016), http:// dx.doi.org/10.1016/j.jhydrol.2016.06.056

6

O.D. Jones et al. / Journal of Hydrology xxx (2016) xxx–xxx

Fig. 4. Monte Carlo estimates and approximations of the standard deviation of the limiting runoff. In both left and right plots the rainfall intensity is (truncated) normal with mean q and coefficent of variation (CV) equal to 0.1. In the left plot the infiltration is gamma distributed with mean 1 and CV equal to 0.75, 1.0, and 1.25. In the right plot the infiltration is log-normal distributed, again with mean 1 and CV 0.75, 1.0, and 1.25. In each plot the Monte Carlo estimates are given with error bars corresponding to 2 standard deviations. The dashed lines are the approximation of Whitt (1993) and the solid green lines are the approximation of Bhat (1993) (Eq. (8)). Note that the results have been scaled by 2ð1
qÞ. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

lx

pffiffiffiffiffi Lx Ly . Thus there is a strip along the bottom of the hillslope that is contributing to the generation of net runoff, but the height of this strip depends on the length of the slope.

X1 X2 Ly

3.2.3. Case 3: infiltration dominates rainfall q < 1 (subcritical) Write lX and r2X for the (approximate) mean and variance of X given by Eq. (8). From (10) we have

ly

Z¼ (1)

(m)

Xn

X

! m X l lx 2 r2 X ðiÞ  Lx N ly X ; ly 2 X2 : lx ly Lx lx ly i¼1

X

As before

Lx Fig. 5. A hillslope of width Lx and length Ly divided into m strips, each divided into n cells of width lx and length ly . X ðiÞ is the runoff from strip i and X j is the runoff from cell j in a given strip.

3.2.1. Case 1: rainfall dominates infiltration q > 1 (supercritical) We have m X mP
mI r2 Z¼ X ðiÞ  N Lx Ly ; Lx lx Ly ly 2 2 lx ly lx ly i¼1 ! mP
mI lx ly r2 ; ¼ Lx Ly N lx ly Lx Ly l2x l2y

Here ðmP
mI Þ=ðlx ly Þ and



!

ð11Þ

ties, and we see that the mean runoff is proportional to the area of the hillslope, Lx Ly . That is, practically all of the hillslope is contributing to the generation of runoff. 3.2.2. Case 2: balanced infiltration and rainfall q ¼ 1 (critical) We have

Here

pffiffiffiffiffiffiffiffiffi





lX =ðlx ly Þ and r2X = l2x l2y are dimensionless quantities. We

see that in this case the mean runoff is proportional to the width Lx of the hillslope, but the length Ly plays no part. That is, runoff into the stream is coming from a strip whose height depends on the rainfall and infiltration distributions, but not the length of the slope. In all three cases above, we can think of the length ly as a system parameter that quantifies the spatial correlation scale of rainfall intensity and infiltration capacity, in the direction of flow. It should be just large enough that the P k and Ik appear to be uncorrelated. It is particularly important when q ¼ 1, when the mean runoff scales pffiffiffiffi like ly , and when q < 1, in which case the runoff scales linearly with ly . 4. Finite hillslopes



r2 = l2x l2y are both dimensionless quanti-

! qffiffiffiffiffiffiffiffi pffiffiffiffiffiffiffiffiffi m X r 2=p r2 ð1
2=pÞ Z¼ X ðiÞ  N Lx Ly ly ; Lx lx Ly ly 2 2 lx ly lx ly i¼1 ! p ffiffiffiffiffiffiffiffi ffi qffiffiffiffiffi qffiffiffiffi r 2=p l r2 ð1
2=pÞ x ¼ Lx Ly N ; ly ly 2 2 Lx lx ly lx ly



ð13Þ

The runoff distributions (11)–(13) are all limiting results as the slope length and slope width go to infinity. It is natural to ask at what scale they make good approximations. The answer depends on q ¼ mP =mI , m ¼ Lx =lx and n ¼ Ly =ly . 4.1. Monte Carlo analysis for a single strip

ð12Þ



r 2=p=ðlx ly Þ and r2 ð1
2=pÞ= l2x l2y are dimensionless quan-

tities. In this (rather special) case the mean runoff is proportional to

We consider first a single strip. We tested the sensitivity of the analytic model to n, the number of cells in the strip, by using numerical simulations of net runoff X n in the case of exponentially distributed rainfall and infiltration capacity. For a variety of q and n values, 100,000 realisations of the net runoff were calculated using Eq. (2), and the mean and standard deviation of the net runoff estimated. The estimates were then compared to Eqs. (5)–(7). As in Section 3 we take cases.

Please cite this article in press as: Jones, O.D., et al. The stochastic runoff-runon process: Extending its analysis to a finite hillslope. J. Hydrol. (2016), http:// dx.doi.org/10.1016/j.jhydrol.2016.06.056

7

O.D. Jones et al. / Journal of Hydrology xxx (2016) xxx–xxx

4.1.1. Case 1: rainfall dominates infiltration q > 1 (supercritical) Assuming exponential rainfall and infiltration, we have from Eq. (5) that the expected net runoff grows like nðmP
mI Þ and the qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi  ffi standard deviation like n m2P þ m2I . Fig. 6 gives estimates of the mean and standard deviation of X n for n between 10 and 100 and q between 1 and 3. We see that we approach the limit faster when q is larger. For example, for q P 2 Eq. (5) gives a very good approximation for n as small as 20. However we observe that as q approaches 1 the approximation deteriorates. 4.1.2. Case 2: balanced infiltration and rainfall q ¼ 1 (critical) From Eq. (6) we have that, for exponential rainfall and infiltraqffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi   tion, the mean net runoff grows like 2 m2P þ m2I n=p and the qffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi  standard deviation like m2P þ m2I nð1
2=pÞ. In Fig. 7 we give estimates of the mean and standard deviation of X n for n from 1 to 1000. We see that the standard deviation converges very quickly to its asymptotic value, but the mean takes a very long time to converge, as n ! 1. This is a reflection of the fact that there is a phasechange in the behaviour of the system as we move from q > 1 to q < 1. 4.1.3. Case 3: infiltration dominates rainfall q < 1 (subcritical) As we have exponential rainfall and infiltration, the net runoff has the limiting distribution given by Eq. (7). From this we have PðX ¼ 0Þ ¼ 1
mP =mI , EX ¼ m2P =ðmI
mP Þ and VarX ¼ m3P ð2mI
mP Þ=ðmI
mP Þ2 . Our simulation estimates of the zero probability, mean and standard deviation of X n are given in Fig. 8. As in the case q > 1, we see that the level of agreement between the numerical simulation of a finite strip of n pixels, and the asymptotic distribution given by Eq. (7), is sensitive to both the value of q and the strip length n. For example, for q 6 0:6 the asymptotic result gives a very good approximation for values of n as small as 20, but as q increases to 1 the quality of the approximation deteriorates rapidly. Taking all three cases together, we see that the asymptotic approximations given by Eqs. (5)–(7) are only reliable for values of q away from 1. We deal with this problem in Section 5, where we give an expression for the mean net runoff that provides an excellent approximation for all q. 4.2. Monte Carlo analysis for a hillslope

5. Mean runoff on a finite hillslope In this section we derive an approximation for the mean net runoff from a single strip of finite length. To get the mean net runoff for a hillslope of m strips we just multiply by m. We will consider the special case Ik  expðkÞ, so that mI ¼ 1=k and c2I ¼ 1. Fix n and define

  f mP ; n; mI ; c2P :¼ EX n As mP ! 0 the number of pixels connected to the lower boundary decreases, and in the limit we do not notice that the strip has finite length. Thus we can treat the strip as if it has infinite length, so we have from Eq. (8)

f





    m2P 1 þ c2P m2P 1 þ c2P :   2ðmI
mP Þ 2mI

  f mP ; n; mI ; c2P  nðmP
mI Þ: With these two equations we have asymptotic representations for f at 0 and 1. We would like to find a function to interpolate these. Clearly f is increasing in mP and n, and bounded above by     0 f mP ; 1; mI ; c2P ¼ m2P 1 þ c2P =ð2ðmI
mP ÞÞ. We will assume that f

1.5

2.0

2.5

mean rainfall

2.0

n = 100 50 20 10 mean infiltration = 1

1.0

mean infiltration = 1

2.5

3.0

0.5

1.0

100

1.5

stdev net runoff/sqrt(n)

1.5

50

1.0

mP ; n; mI ; c2P

As mP ! 1 we have from Eq. (5)

n = 10 20

0.0

mean net runoff/n

2.0

Eqs. (11)–(13) assume that there are sufficient strips m across the plane that the distribution of the sum of the runoff flows from

many strips will be approximately normal. In other words, m must be large enough for the central limit theorem to apply. When q P 1 the net runoff from each strip is already approximately normal, so the net hillslope runoff will be approximately normal even when m is small. However when q < 1 the net runoff from each strip can be quite far from normal (c.f. Eq. (7)), which raises the question how many strips are required for the asymptotic distribution to give a good approximation. We again used simulation to assess this question. Assuming exponentially distributed rainfall and infiltration capacity, and sufficiently long strip lengths, the net runoff distribution for a single strip is given by F in Eq. (7). Simulating m independent realizations of X from this distribution and then summing gives us a realization of Z, as per Eq. (10). In Fig. 9 we estimate the distribution of Z for various values of q and m. We see that if q is close to 1, then normality is approached rapidly and approximately 40 strips are sufficient. However as q decreases, more strips are required for the sum to appear normal. The simulations suggest that approximately m ¼ 300 are required when q ¼ 0:1.

3.0

1.0

1.5

2.0

2.5

3.0

mean rainfall

Fig. 6. Estimates of the mean (left) and standard deviation (right) of the net runoff for various values of the mean rainfall mP and strip length n. The mean infiltration is set to pffiffiffi 1, so that q ¼ mP . The mean is scaled by n and the standard deviation by n, and the asymptotic limit as n ! 1 is given by the dashed line in each case. The ordering of the lines corresponds to the ordering of n in each case.

Please cite this article in press as: Jones, O.D., et al. The stochastic runoff-runon process: Extending its analysis to a finite hillslope. J. Hydrol. (2016), http:// dx.doi.org/10.1016/j.jhydrol.2016.06.056

O.D. Jones et al. / Journal of Hydrology xxx (2016) xxx–xxx

1

2

3

4

5

6

1.010 1.005

stdev net runoff (scaled) 0

1.000

0.9 0.8 0.7 0.6 0.5 0.4

mean net runoff (scaled)

1.0

8

7

0

1

2

log n

3

4

5

6

7

log n

n=5 n = 10 n = 20 n = 50

6 4

n = 20 n = 10 n=5 n=2 n=1

0

0.0

mean infiltration = 1

mean net runoff

n=2

n = 50

mean infiltration = 1

2

0.6 0.4

n=1

0.2

P(net runoff = 0)

0.8

8

1.0

Fig. 7. ffiEstimates of the mean (left) andpstandard deviation (right) of the net runoff for mP ¼ mI ¼ 1 (q ¼ 1) and various values of the strip length n. The mean is scaled by pffiffiffiffiffiffiffiffiffiffiffi ffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi ffi 4n=p and the standard deviation by 2nð1
2=pÞ, and the asymptotic limit as n ! 1 is given by the dashed line in each case.

0.0

0.2

0.4

0.6

0.8

1.0

0.0

0.2

0.4

0.6

0.8

1.0

mean rainfall

6

n = 50

4

mean infiltration = 1

n = 20 n = 10

2

stdev net runoff

8

mean rainfall

0

n=5 n=2 n=1

0.0

0.2

0.4

0.6

0.8

1.0

mean rainfall Fig. 8. Estimates of the probability of no net runoff, the mean net runoff, and the standard deviation of the net runoff, for various values of the mean rainfall mP and strip length n. The mean infiltration is set to 1, so that q ¼ mP . In each case the asymptotic value, for large n, is given by the dashed line.

is an increasing function. Let g be our model for f , then we suppose that g 0 has the following sigmoid shape

h ðh
dÞ expðbðx
cÞÞ
ad
d¼ 1 þ a expð
bðx
cÞÞ a þ expðbðx
cÞÞ h gðxÞ ¼ logða þ expðbðx
cÞÞÞ
dx b

g 0 ðxÞ ¼

0

know about f . Using the equations f ð0Þ ¼ 0, f ð0Þ ¼ 0 and   00 f ð0Þ ¼ 1 þ c2P =mI we get

h logða þ expð
bcÞÞ ¼ 0 ) a þ e
bc ¼ 1 b h h
d g 0 ð0Þ ¼
d ¼ 0 ) ebc ¼ ) hð1
aÞ ¼ d 1 þ a expðbcÞ ad   hab expðbcÞ ¼ 1 þ c2P =mI g 00 ð0Þ ¼ 2 ð1 þ a expðbcÞÞ

gð0Þ ¼

ð14Þ Clearly there are other curves that can give us a sigmoid shape for g 0 . We chose this form for its tractability; it is the same sigmoid used for logistic regression. By differentiating it can easily be verified that g is as given. To fit g we match its asymptotic behaviour at 0 and 1 to what we

Using f ðxÞ  nðx
mI Þ as x ! 1, we get

nðx
mI Þ ¼ hðx
cÞ
dx ) h
d ¼ n and nmI ¼ ch:

Please cite this article in press as: Jones, O.D., et al. The stochastic runoff-runon process: Extending its analysis to a finite hillslope. J. Hydrol. (2016), http:// dx.doi.org/10.1016/j.jhydrol.2016.06.056

O.D. Jones et al. / Journal of Hydrology xxx (2016) xxx–xxx

9

Fig. 9. Relative frequency histograms showing the effect of q (the ratio of the expected value of rainfall intensity to the expected value of infiltration capacity), and the number of strips m, on the shape of the distribution of the net hillslope runoff rate from a tilted plane.

Solving these we get

mean-squared error between the line g and the estimated f . The fit is remarkably good, for small and large values of n. To facilitate fitting the curve to observed net runoff data, we reparameterise it. Recall the second Einstein function E2 ðxÞ ¼ x=ðe x
1Þ, then from (15) we have

d ¼ nð1
aÞ=a h ¼ n=a c ¼ amI b ¼
logð1
aÞ=ðamI Þ

a ¼ E
1 2

where a is the solution to (for a 2 ð0; 1Þ)




1
a 1 þ c2P logð1
aÞ ¼ : a n

ð15Þ

Note that the restriction a 2 ð0; 1Þ implies 1 þ < n. Putting the constants back into (14) we get, for a ¼
logð1
aÞ, c2P

nmI xnð1
aÞ logða þ e
ðx
amI Þ logð1
aÞ=ðamI Þ Þ
a
logð1
aÞ  

   x 1 þ c2P nmI ¼ logð1
e
a þ exp x 1 þ c2P ea =ðnmI Þ
a Þ


gðxÞ ¼

a

a

ð16Þ

A simulation experiment was run to compare the interpolating curve g to f ¼ EX n . We used exponential infiltration with a mean of 10, and gamma distributed rainfall with a mean ranging from 0 pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi to 20, and coefficient of variation 0:0625. The results are given in Fig. 10, where each pane takes a different value of n. Here the dots are Monte Carlo estimates of f ðmP Þ, the solid line is gðmP Þ where the constant a was chosen by solving Eq. (15), and the dotted line is also gðmP Þ but using a value of a chosen to minimise the

 1 þ c2P : n

ð17Þ

Since cP and n ¼ Ly =ly are dimensionless, so is a. We put   b ¼ ly 1 þ c2P and use it to parameterize the spatial correlation scale and the variability of the rainfall. We have a ¼ E
1 2 ðb=Ly Þ, so our model for the mean net runoff flow from a slope of width Lx and length Ly becomes

RðmP ; mI ; bÞ ¼ mgðmP Þ

mI 1 mP b=Ly logð1
e
a þ expfea ðb=Ly ÞðmP =mI Þ
agÞ
¼ Lx Ly lx ly a lx ly a ð18Þ We call this the stochastic runoff mean curve. We see that the spatial variability is expressed entirely through the (dimensionless) argument b=Ly . Thus the relative effect of spatial variability (on the mean at least), decreases as the slope length increases. Note that the requirement 1 þ c2P < n translates to 0 < b=Ly < 1.

Please cite this article in press as: Jones, O.D., et al. The stochastic runoff-runon process: Extending its analysis to a finite hillslope. J. Hydrol. (2016), http:// dx.doi.org/10.1016/j.jhydrol.2016.06.056

10

O.D. Jones et al. / Journal of Hydrology xxx (2016) xxx–xxx

m_I = 10 n = 20 c2_P = 0.0625

150 50 0 0

5

10

15

20

0

5

10

15

20

m_P

m_I = 10 n = 40 c2_P = 0.0625

m_I = 10 n = 80 c2_P = 0.0625

0

0

200

400

200

E runoff

600

300

800

400

m_P

100

E runoff

100

E runoff

60 40 0

20

E runoff

80

100

200

m_I = 10 n = 10 c2_P = 0.0625

0

5

10

15

20

0

5

10

m_P

15

20

m_P

  Fig. 10. f mP ; n; mI ; c2P by simulation (black circles) with two approximations: gðmP Þ using a as given in Eq. (15) (solid blue line), and gðmP Þ using a chosen to minimise a least squares fit of the curve to the simulated values (red dashed line). (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

1

a

  c log 1
e
a þ expfea cðmp =mi Þ
ag
mp

a

ð19Þ

where c ¼ b=Ly 2 ð0; 1Þ measures spatial variability and a ¼ E
1 2 ðcÞ. c is inversely proportional to the number of cells in each strip, and increases as the rainfall and infiltration become more variable. The limit as c ! 0 is the minimum net runoff rate, namely mp
mi for mp P mi , and the limit as c ! 1 is the maximum net runoff rate mp
mi ð1
e
mp =mi Þ. Hawkins and Cundy (1987) obtained these bounds by arranging the point infiltration values in either increasing or decreasing order, as you move down the slope.
1

The results for mi ¼ 10 mm h are given in Fig. 11, for various values of c ¼ b=Ly . We see that as c increases from 0 to 1, the curve given by Eq. (19) interpolates between the lower and upper bounds of Hawkins and Cundy. This helps explain why, as noted in the introduction, the upper bound of Hawkins and Cundy has been seen to give a good approximation in certain situations. 6. Discussion In Section 2.2 the applicability of the stochastic runoff-runon model was demonstrated by fitting it to data from a field experiment. However we emphasize that for this experiment the soil was saturated before measurements were taken, so that the soil infiltration in each cell could be regarded as constant over time.

5

mi

0

We compare the stochastic runoff mean curve R given in Eq. (18) with the upper and lower bounds given by Hawkins and Cundy (Eq. (1)). Eq. (1) is given in terms of fluxes rather than flows, thus to compare them we put mi ¼ mI =ðlx ly Þ and mp ¼ mP =ðlx ly Þ, and we divide R by the area of the hillslope Lx Ly . The mean runoff per square meter is then approximately

net runoff rate 10 15

20

5.1. Comparison with Hawkins and Cundy

0

5

10

15

20

25

30

rainfall rate Fig. 11. Comparison of the net runoff rates given by Hawkins and Cundy (1987) (Eq. (1)) and the family of stochastic runoff mean curves (Eq. (18)). The mean infiltration rate is 10 mm h
1 and we take b=Ly = 0.0001, 0.01, 0.1, 0.3, 0.7 (respectively the red, pale green, dark green, blue and purple lines). The dashed lines give the upper and lower bounds of Hawkins and Cundy. (For interpretation of the references to colour in this figure legend, the reader is referred to the web version of this article.)

In general runoff depends on both spatial and temporal factors. For example the infiltration of the soil depends on how wet it is, and if the duration of the rainfall event is short compared to the travel time of flow down the hillslope, then flow paths may not be able to connect with the downslope boundary (Reaney et al., 2007; Wainwright and Parsons, 2002; Corradini et al., 1998). Julien and Moglen (1990) found that when the rainfall event duration was less than the time of concentration, spatial variability had little effect on the net runoff. Conversely, when the time of concentration was less than the rainfall event duration, it was very important. Inherent in the stochastic runoff-runon model is the

Please cite this article in press as: Jones, O.D., et al. The stochastic runoff-runon process: Extending its analysis to a finite hillslope. J. Hydrol. (2016), http:// dx.doi.org/10.1016/j.jhydrol.2016.06.056

O.D. Jones et al. / Journal of Hydrology xxx (2016) xxx–xxx

assumption that any infiltration excess flow generated in a cell has time to runoff to the next cell downslope. The asymptotic analysis in Section 3 highlights the qualitatively different behaviour of net runoff when q ¼ mP =mI < 1 and q > 1 (and when q ¼ 1, though this case is of more theoretical than practical interest). When q > 1 the whole hillslope is potentially contributing to net runoff, but when q < 1 only a strip along the stream edge contributes. Bren and Turner (1979) measured runoff from bounded forested plots 2.5, 13 and 19 m long under natural rainfall and concluded that ‘‘the volume of overland flow generated is not determined by the length of slope above the point of measurement”. Similarly, Ronan (1986) measured natural runoff from bounded forested plots 20 m long and concluded that the ‘‘volume harvested will depend on the length of the collection gutter, not on the area above the gutter”. In both these studies it is likely that q < 1 because of the typically high infiltration found in temperate forests with well-structured macro-porous forest soils. These results are consistent with other studies reporting inverse relationships between unit area volumes of overland flow, and plot lengths (Gomi et al., 2008; Bryan and Poesen, 1989; Joel et al., 2002; Parsons et al., 2006). Similarly the stochastic runoff-runon model agrees with the observation of Smith and Goodrich (2000), that the net infiltration rate of the plot (rainfall minus runoff) increases with increasing rainfall intensity, particularly when q is small. This is because increasing the rainfall intensity results in more of the available infiltration capacity within the plot being utilized. This has been observed for plot data by Yu (1999) and Yu et al. (1998). The stochastic runoff mean curve (Eqs. (18) and (19)) gives a practical tool for quantifying average net runoff. It interpolates between the upper and lower bounds of Hawkins and Cundy using   only a single additional parameter, b=Ly , where b ¼ ly 1 þ c2P . Small c ¼ b=Ly puts us close to the lower bound and large c puts us close to the upper bound. Thus we are in a position to improve upon previous applications of the Hawkins and Cundy bound to the problem of runoff/infiltration prediction at the hillslope scale (Yu et al., 1997, 1998; Yu, 1999; Fentie et al., 2002; Kandel et al., 2005). The parameter b summarizes the spatial variability of the hillslope, and we see that net runoff increases as both cP , the coefficient of variation for rainfall, and ly , the spatial correlation scale, increase. (Recall that we set cI ¼ 1 when deriving the stochastic runoff mean curve.) Our definition of ly , as that scale at which spatial correlation in rainfall and infiltration becomes negligible, is arguably imprecise, but none-the-less has a physical interpretation. The spatial correlation of rainfall intensity in real landscapes is not often measured, and the limited literature is contradictory (Loustau et al., 1992; Whelan and Anderson, 1996). When measuring spatial correlation, we necessarily discretize what is naturally a continuous system. The scale of this discretization often results from the scale of the instrument used for measurement, which in the case of rainfall and soil properties is typically of the order of 0.01 m2. Spatial correlation of saturated conductivity at this scale is reported as highly variable between soil types (Binley et al., 1989; Loague and Gander, 1990; Luxmoore, 1981; Helmers and Eisenhauer, 2006; Buttle and House, 1997). From our analysis we conclude that determination of typical spatial correlation scales for throughfall and infiltration capacity in different landscapes will contribute to the understanding of the variability of runoff responses reported in the literature.

7. Summary Jones et al. (2013) showed that the runoff-runon phenomenon across a spatially variable area could be represented conceptually and mathematically using queuing theory. However the assump-

11

tions and restrictions required to achieve the simple analytic solutions presented by Jones et al. (2013) constrained the practical application of the new approach. In particular, the approach was limited to the cases where (i) mean rainfall intensity was less than the mean infiltration capacity and (ii) the slope length was considered infinite (there was no upper boundary condition). The main contribution of the present paper is to address these practical constraints, whilst maintaining the underlying conceptual approach based on queuing theory. The result is the derivation of a family of curves for the mean net runoff, which, with the addition of a single additional parameter, smoothly interpolate between the upper and lower bounds of the rainfall-runoff functions for a spatially variable area proposed Hawkins and Cundy (1987). Acknowledgments This project was funded in part by Melbourne Water and the Department of Environment, Land, Water and Planning, Victoria. References Asmussen, 2003. Applied Probability and Queues, second ed. Springer. Bhat, B.N., 1993. Approximation for the variance of the waiting time in a GI/G/1 queue. Microelectron. Reliab. 33, 1997–2002. Binley, A., Elgy, J., Beven, K., 1989. A physically based model of heterogeneous hillslopes 1: runoff production. Water Resour. Res. 25 (6), 1219–1226. Brakensiek, D.L., Englemen, R.L., Rawls, W.J., 1981. Variation within texture classes of soil water parameters. Trans. Am. Soc. Agric. Eng. 24 (1), 335–339. Bren, L.J., Turner, A.K., 1979. Overland flow on a steep, forested infiltrating slope. Aust. J. Soil Res. 30, 43–52. Bryan, R.B., Poesen, J., 1989. Laboratory experiments on the influence of slope length on runoff, percolation and rill development. Earth Surf. Proc. Land. 14, 211–231. Bubenzer, G.D., Meyer, L.D., 1965. Simulation of rainfall and soils for laboratory research. Trans. Am. Soc. Agric. Eng. 8, 73–75. Buttle, J.M., House, D.A., 1997. Spatial variability of saturated hydraulic conductivity in shallow macroporous soils in a forested basin. J. Hydrol. 203, 127–142. Carlyle-Moses, D.E., Flores Laureano, J.S., Price, A.G., 2004. Throughfall and throughfall spatial variability in Madrean oak forest communities of northeastern Mexico. J. Hydrol. 297, 124–135. Corradini, C., Morbidelli, R., Melone, F., 1998. On the interaction between infiltration and Hortonian runoff. J. Hydrol. 204, 52–67. Daley, D.J., Ya. Kreinin, A., Trengove, C.D., 1992. Inequalities concerning the waitingtime in single-server queues: a survey. In: Bhat, U.N., Basawa, I.V. (Eds.), Queuing and Related Models. Clarendon Press. Dunkerley, D., 2008. Rain event properties in nature and in rainfall simulation experiments: a comparative review with recommendations for increasingly systematic study and reporting. Hydrol. Process. 22, 4415–4435. Fentie, B., Yu, B., Silburn, M.D., Ciesiolka, C.A.A., 2002. Evaluation of eight different methods to predict hillslope runoff rates for a grazing catchment in Australia. J. Hydrol. 261, 102–114. Gomi, T., Sidle, R.C., Miyata, S., Kosugi, K., 2008. Dynamic runoff connectivity of overland flow on steep forested hillslopes: scale effects and runoff transfer. Water Resour. Res. 44. Grayson, R.B., Moore, I.D., McMahon, T.A., 1992. Physically based hydrologic modelling 2. Is the concept realistic? Water Resour. Res. 26 (10), 2659–2666. Green, W.H., Ampt, G.A., 1911. Studies on soil physics: I. Flow of air and water through soils. J. Agric. Sci. 4, 1–24. Harel, M.-A., Mouche, E., 2013. 1-D steady state runoff production in light of queuing theory: heterogeneity, connectivity, and scale. Water Resour. Res. 49, 7973–7991. Harel, M.-A., Mouche, E., 2014. Is the connectivity function a good indicator of soil infiltrability distribution and runoff flow dimension? Earth Surf. Proc. Land. 39, 1514–1525. Hawkins, R.H., Cundy, T.W., 1987. Steady-state analysis of infiltration and overland flow for spatially varied hillslopes. Water Resour. Bull. 23 (2), 251–256. Helmers, M.J., Eisenhauer, D.E., 2006. Overland flow modeling in a vegetative filter considering non-planar topography and spatial variability of soil hydraulic properties and vegetation density. J. Hydrol. 328, 267–282. Joel, A., Nessing, I., Sequel, O., Casanova, M., 2002. Measurement of surface water runoff from plots of two different sizes. Hydrol. Process. 16, 1467–1478. Jones, O.D., Sheridan, G.J., Lane, P.N., 2013. Using queuing theory to describe steadystate runoff-runon phenomena and connectivity under spatially variable conditions. Water Resour. Res. 49 (11), 7487–7497. Jones, O.D., Sheridan, G.J., Lane, P.N., 2009. A stochastic runoff model incorporating spatial variability. In: Anderssen, R., Braddock, R., Newham, L. (Eds.), 18th World IMACS Congress and MODSIM09 International Congress on Modelling and Simulation. The Modelling and Simulation Society of Australia and New Zealand, pp. 1865–1871. Julien, P.Y., Moglen, G.E., 1990. Similarity and length scale for spatially varied overland flow. Water Resour. Res. 26 (8), 1819–1832.

Please cite this article in press as: Jones, O.D., et al. The stochastic runoff-runon process: Extending its analysis to a finite hillslope. J. Hydrol. (2016), http:// dx.doi.org/10.1016/j.jhydrol.2016.06.056

12

O.D. Jones et al. / Journal of Hydrology xxx (2016) xxx–xxx

Kandel, D.D., Western, A.W., Grayson, R.B., 2005. Scaling from process timescales to daily time steps: a distribution function approach. Water Resour. Res. 41. Kleinrock, L., 1975a. Queueing Systems: Theory, vol. I. Wiley. Kleinrock, L., 1975b. Queueing Systems: Computer Applications, vol. II. Wiley Interscience, New York. Kramer, W., Lagenbach-Belz, M., 1976. Approximate formula for the delay in the queueing system GI/G/1. In: Proceedings ITC 8, Melbourne. Lane, P.N.J., Sheridan, G.J., Noske, P.J., 2006. Changes in sediment loads and discharge from small mountain catchments following wildfire in south eastern Australia. J. Hydrol. 331, 495–510. Llorens, P., Poch, R., Latron, J., Gallart, F., 1997. Rainfall interception by a Pinus sylvestris forest patch overgrown in a Mediterranean mountainous abandoned area. I: monitoring design and results down to the event scale. J. Hydrol. 199, 331–345. Loague, K., Gander, G.A., 1990. R-5 revisited, 1, spatial variability of infiltration on a small rangeland watershed. Water Resour. Res. 26 (5), 957–971. Loague, K., Kyriakidis, P.C., 1997. Spatial and temporal variability in the R-5 infiltration data set: Déjà vu and rainfall-runoff simulations. Water Resour. Res. 33 (12), 2883–2895. Loch, R.J., 1989. Aggregate Breakdown Under Rain: its Measurement and Interpretation PhD Thesis. University of New England, Armidale, NSW. Loustau, D., Berbigier, P., Granier, A., El Hadj Moussa, F., 1992. Interception loss, throughfall and stemflow in a maritime pine stand. I: variability of throughfall and stemflow beneath the pine canopy. J. Hydrol. 138, 449–467. Luxmoore, R.J., 1981. Micro-, meso-, and macroporosity of soil. Soil Sci. Soc. Am. J. 45 (3), 671–672. Marchal, W.G., 1976. An approximate formula for waiting time in single server queues. IIE Trans. 8 (4), 473–474. Nahar, N., 2003. Influence of Run-on on Field Scale Surface and Subsurface Water and Contaminant Movement Over Spatially Variable Hillslopes PhD Thesis. Nielsen, D.R., Biggar, J.W., Erh, K.T., 1973. Spatial variability of field measured soilwater properties. Hilgardia 42, 215–259. Nyman, P., 2007. The Effect of Hydrophobicity on Hydraulic Conductivity in a Macroporous Forest Soil: Implications For Infiltration and Runoff Generation Hons Thesis. The University of Melbourne. Parsons, A.J., Brazier, R.E., Wainwright, J., Powell, D.M., 2006. Scale relationships in hillslope runoff and erosion. Earth Surf. Proc. Land. 31, 1381–1393. Prabhu, N.U., 1998. Stochastic Storage Processes: Queues, Insurance Risk, Dams and Data Communication, second ed. Springer.

Price, A.G., 1994. Measurement and variability of physical properties and soil water distribution in a forest podzol. J. Hydrol. 161, 347–364. Rao, B.V., Feldman, R.M., 2001. Approximations and bounds for the variance of steady-state waiting times in a GI/G/1 queue. Operat. Res. Lett. 28, 51–62. Reaney, S.M., Bracken, L.J., Kirkby, M.J., 2007. Use of the Connectivity of Runoff Model (CRUM) to investigate the influence of storm characteristics on runoff generation and connectivity in semi-arid areas. Hydrol. Process. 21, 894–906. Ronan, N.M., 1986. Water Supply Catchment Hydrology Research: The Hydrological Effects of Fuel Reduction Burning and Wildfire at Wallaby Creek. Melbourne and Metropolitan Board of Works. Sheridan, G.J., Jones, O.D., Lane, P.N., 2009a. Stochastic rainfall-runoff equations for quantifying runoff and pollutant connectivity between hillslopes and streams. EGU General Assembly Conference Abstracts, vol. 11, p. 3845. Sheridan, G.J., Jones, O.D., Smith, H., Cawson, J., Lane, P.N., 2009b. Stochastic runoff connectivity (SRC) equations: integration with erosion models for water quality prediction. EGU General Assembly Conference Abstracts, vol. 11, p. 13981. Sharma, M.L., Gander, G.A., Hunt, C.G., 1980. Spatial variability of infiltration in a watershed. J. Hydrol. 45, 101–122. Smith, R.E., Goodrich, D.C., 2000. Model for rainfall excess patterns on randomly heterogeneous areas. J. Hydrol. Eng. 5, 355–362. Wainwright, J., Parsons, A.J., 2002. The effect of temporal variations in rainfall on scale dependency in runoff coefficients. Water Resour. Res. 38 (12). Watson, K.W., Luxmoore, R.J., 1986. Estimating macroporosity in a forest watershed by use of a tension infiltrometer. Soil Sci. Soc. Am. J. 50, 578–582. Whelan, M.J., Anderson, J.M., 1996. Modelling spatial patterns of throughfall and interception loss in a Norway spruce (Picea abies) plantation at the plot scale. J. Hydrol. 186, 335–354. Whitt, W., 1993. Approximations for the GI/G/m queue. Prod. Operat. Manage. 2, 114–161. Wilson, G.V., Luxmoore, R.J., 1988. Infiltration, macroporosity, and mesoporosity distributions on two forested watersheds. Soil Sci. Soc. Am. J. 52, 329–335. Yu, B., 1999. A comparison of the Green-Ampt and a spatially variable infiltration model for natural storm events. Trans. Am. Soc. Agric. Eng. 42 (1), 89–97. Yu, B., Rose, C.W., Coughlan, K.J., Fentie, B., 1997. Plot-scale rainfall runoff characteristics and modelling at six sites in Australia and South East Asia. Trans. Am. Soc. Agric. Eng. 40 (5), 1295–1303. Yu, B., Cakurs, U., Rose, C.W., 1998. An assessment of methods for estimating s at the plot scale. Trans. Am. Soc. Agric. Eng. 41 (3), 653–661.

Please cite this article in press as: Jones, O.D., et al. The stochastic runoff-runon process: Extending its analysis to a finite hillslope. J. Hydrol. (2016), http:// dx.doi.org/10.1016/j.jhydrol.2016.06.056