Estimating time parameters of overland flow on impervious surfaces by the particle tracking method

Abstract

Travel time and time of concentration T_c are important time parameters in hydrological designs. Although T_c is the time for the runoff to travel to the outlet from the most remote part of the catchment, most researchers have used an indirect method such as hydrograph analysis to estimate T_c. A quasi two-dimensional diffusion wave model with particle tracking for overland flow was developed to determine the travel time, and validated for runoff discharges, velocities, and depths. Travel times for 85%, 95% and 100% of particles arrival at the outlet of impervious surfaces (i.e. T_t85, T_t95, and T_t100) were determined for 530 model runs. The correlations between these travel times and T_c estimated from hydrograph analysis showed a significant agreement between T_c and T_t85. All the travel times showed nonlinear relationships with the input variables (plot length, slope, roughness coefficient, and effective rainfall intensity) but showed linear relationships with each other.

Editor D. Koutsoyiannis; Associate editor S. Grimaldi

Résumé

Le temps de transfert et le temps de concentration T_c sont des paramètres de temps importants pour l’ingénierie hydrologique. Bien que T_c soit défini comme le temps mis par le ruissellement pour aller de la partie la plus éloignée du bassin versant jusqu’à l’exutoire, la plupart des chercheurs ont utilisé une méthode indirecte, telle que l’analyse de l’hydrogramme, pour estimer T_c. Un modèle d’écoulement de surface par onde diffusante en quasi-deux dimensions avec suivi de particules a été développé pour déterminer le temps de transfert, et validé pour les débits, les vitesses et les hauteurs d’écoulement. Les temps de transfert pour 85%, 95% et 100% d’arrivées de particules à l’exutoire de surfaces imperméables (c’est-à-dire T_t85, T_t95, et T_t100) ont été déterminés pour 530 exécutions du modèle. Les corrélations entre ces temps de transfert et les valeurs de T_c estimées par analyse d’hydrogramme ont montré un bon accord entre T_c et T_t85. Tous les temps de transfert ont montré des relations non-linéaires avec les variables d’entrée (longueur de la parcelle, pente, coefficient de rugosité, et intensité de la pluie efficace) mais ont montré des relations linéaires entre eux.

Key words:

• travel time
• overland flow
• time of concentration
• impervious surface
• runoff
• particle tracking
• model
• surface hydrology

Mots clefs:

• temps de transfert
• ruissellement de surface
• temps de concentration
• surface imperméable
• ruissellement
• suivi de particules
• modèle
• hydrologie de surface

INTRODUCTION

The problem of estimating travel velocities and thus travel time from digital elevation models (DEMs) is still an open topic. Maidment (1993) was one of the first to relate travel velocities to geomorphic properties, such as contributing area and slope. However, to specify travel speed for each particle along the travel path was challenging and Maidment (1993) used a DEM for overland flow routing based on constant velocity. Similarly, Olivera and Maidment (1999) developed a raster-based, spatially distributed routing technique based on a first-passage time response function. Smedt et al. (2000) and Liu et al. (2003) proposed a flow routing method along flow paths determined by the topography depending upon slopes, flow velocities, and dissipation characteristics along the flow paths. Luo and Harlin (2003) derived simplified theoretical travel times based on watershed hypsometry ignoring friction and stated that the theoretical travel time is not equal to the real measured time-to-peak discharge. Admiraal et al. (2004) used a particle tracking routine in their two-dimensional flow model to release particles (water parcels) at the inlet and computed their paths and travel times through the wetland from the velocity field. Noto and La Loggia (2007) computed the travel time of two watersheds implementing a cell-to-cell flow path through the landscape determined from a DEM. Cleveland et al. (2008) used a particle tracking model based on a DEM to derive synthetic unit hydrographs for 126 Texas watersheds. The velocity of each particle was calculated using an equation similar to Manning’s equation with resistance coefficient and mean flow depth as watershed calibration parameters. Grimaldi et al. (2010) improved the width function in an instantaneous unit hydrograph based on travel time function and tested it in 11 case studies. Froehlich (2011) developed an overland flow travel time formula based on rainfall intensity–duration relationships using the kinematic wave approximation. Niri et al. (2012) derived a semi-analytical travel time formula for time-area model applications. Travel time is often closely linked to another time parameter, time of concentration T_c commonly used in hydrological designs and analyses (Kent 1972, Pilgrim 1976, McCuen et al. 1984, Saghafian and Julien 1995). Most of the empirical formulas developed for T_c are based on estimation of T_c from observed or simulated hydrographs.

Even though T_c is commonly defined as the time it takes runoff to travel from the most distant point along a hydraulic pathway in the watershed to the outlet, many researchers have used an indirect method, such as hydrograph analysis, to estimate T_c instead of measuring or determining travel time for particles to reach the outlet. Kull and Feldman (1998) assumed the travel time of runoff from each cell in a watershed is proportional to T_c scaled by the ratio of travel length of the cell over the maximum travel length. This method is implemented in HEC-GEOHMS (Doan 2000). However, this means that the average runoff velocity travelling from any point to the outlet is assumed uniform and constant (Saghafian et al. 2002). This assumption of uniform velocity gives an incorrect estimate of both travel time and T_c. Since overland flow is both unsteady and non-uniform, the spatial variability in flow characteristics both across and down the slope makes it difficult to obtain reliable estimates of all the velocities in the flow (Abrahams et al. 1986).

An alternative approach to obtain the travel time is to measure the time required for the passage of dye along a prespecified distance. Emmett (1970), Dunne and Dietrich (1980), Roels (1984), Govers (1992), and Rouhipour et al. (1999) recorded the travel time of the front of the dye. However, they measured the travel time for the fastest fronts of flow while neglecting the slower-moving parts of the flow. Consequently, their measurements of surface velocity were overestimated (Emmett 1970, Dunne and Dietrich 1980). Moreover, the larger the fraction in a flow that is outside the measured fast flowing front, the more severe will be the overestimation of the velocity (Abrahams et al. 1986). Similarly, coarse space and time resolutions of such data have limited their usability (Grimaldi et al. 2012). Similarly, while there is an increasing availability of high resolution DEM for river basins, there is often a lack of adequate hydrological data for small and ungauged catchments (Grimaldi et al. 2010). As a result, travel time estimation for such small and poorly monitored catchments is a topic of increasing interest.

Hence, a direct approach for determining the travel time of both fast and slow moving particles to reach the outlet of impervious surfaces is presented in this paper. Even so, the travel of particles in a watershed may occur on the surface or below it or in a combination of both; however, only the travel of particles on the surface (impervious) is considered in this study. A quasi two-dimensional diffusion wave model with particle tracking (DWMPT) was developed from the diffusion hydrodynamic model (DHM) (Hromadka II and Yen 1986) and validated extensively using observed discharge, flow depth and velocity data from published studies. The travel times of all particles, including both slow moving particles at the upstream and fast moving particles at the downstream, are consolidated into representative time parameters. Based on the percentage of particles reaching the outlet, travel times for 85%, 95% and 100% of particles arrival at the outlet are defined as T_t85, T_t95 and T_t100, respectively. A total of 530 DWMPT runs were performed to determine travel times on impervious surfaces using diverse combinations of the four physically-based input variables: length (L), slope (S_o), Manning’s surface roughness coefficient (n), and effective rainfall intensity (i). The dataset of travel times was generated using a particle tracking method with simulated velocities from the validated DWMPT for conducting a parametric study. Finally, relationships between these travel times and input variables were developed for impervious overland flow planes, and relationships between these travel times were investigated. Relationships between T_c estimated from hydrograph analysis and these travel times were developed and presented.

DIFFUSION WAVE MODEL WITH PARTICLE TRACKING

Hydrological models are numerical models, which attempt to mimic the physical processes within a watershed. Both kinematic and diffusion wave models have been used to simulate these physical processes of surface water movement (Kazezyilmaz-Alhan and Medina Jr. 2007, López-Barrera et al. 2012) in hydrological–hydraulic models. This study deals with numerical simulations of overland flows to investigate direct estimation of travel times of particles on impervious surfaces. Overland flow has been simulated using one- and two-dimensional (1D or 2D) kinematic or diffusion wave models (Henderson and Wooding 1964, Woolhiser and Liggett 1967, Singh 1976, Yen and Chow 1983, Abbott et al. 1986, Chen and Wong 1993, Wong 1996, Jia et al. 2001, Ivanov et al. 2004) and dynamic wave models (Morgali and Linsley 1965, Yeh et al. 1998, Su and Fang 2004). The kinematic wave model is frequently used for the development of T_c formulas (Wong 2005). McCuen and Spiess (1995) suggested that the kinematic wave assumption should be limited to the kinematic wave number  < 100. Hromadka II and Yen (1986) developed the quasi-2D diffusion hydrodynamic model (DHM) to incorporate the pressure effects neglected by the kinematic wave approximation. The diffusion wave approximation is generally accurate for most of the overland flow conditions (Singh and Aravamuthan 1995, Moramarco and Singh 2002, Singh et al. 2005). It is understood that DHM is inaccurate for cases in which the inertial terms play prominent roles such as when the slope of the surface is small (Yeh et al. 1998). Studying overland flows on low-sloped planes has been presented elsewhere (KC et al. 2014) and is beyond the scope of this study.

In this study, a quasi-2D diffusion wave model (DWM), that includes the particle tracking (PT) capability called DWMPT, was developed by modifying the quasi-2D DHM (Hromadka II and Yen 1986) for simulating overland flow on impervious surfaces/planes. The quasi-2D DWMPT solves the momentum and continuity equations (Akan and Yen 1984). Details of the governing equations and numerical solution method used in the DHM are described in Hromadka II and Yen (1986). A number of modifications have been made in DWMPT over DHM. DWMPT has efficient numerical stability criteria using the Courant number, checks for mass and momentum conservation errors, and incorporates detailed rainfall loss models (fractional loss model for impervious surface), and a sub-module for determining time of concentration from the outflow hydrograph. These modifications and improvements are described in detail in KC et al. (2014). The 1D particle tracking model (PTM), a new sub-module of DWMPT, is presented below.

Particle tracking model

The fundamental idea of PTM is to introduce discrete particles (or parcels) in the simulation domain that move according to some kinematics. Various approaches of particle tracking are employed in tracking the flow and pollutants in different fields of computational fluid dynamics. A particle-in-cell approach for particle tracking in fluid dynamics was used by Harlow (1963) and Amsden (1966). Since then, a number of different approaches have been used for tracking particles like the random-walk method (Ahlstron et al. 1977, Prickett et al. 1981), method of characteristics (Konikow and Bredehoeft 1978) and marker-and-cell method (Harlow and Welch 1966, Taylor 1983). Similarly, there are many measurement studies of particle tracking that often involve the use of tracers such as dyes (Emmett 1970, Dunne and Dietrich 1980, Roels 1984, Abrahams et al. 1986), electrolytes (Planchon et al. 2005, Tingwu et al. 2005, Lei et al. 2010), magnetic materials (Ventura et al. 2001, Hu et al. 2011), isotopes (Gardner and Dunn III 1964, Berman et al. 2009), buoyant and fluorescent objects (Bradley et al. 2002, Dunkerley 2003, Meselhe et al. 2004, Tauro et al. 2010, 2012a). Among these, particle tracking velocimetry (PTV) and particle image velocimetry (PIV) are well suited for measurements of surface velocity of particles in physical models (Admiraal et al. 2004). The credit for initial development of PTV and PIV goes to Chang et al. (1985), Liu et al. (1991), and Landreth et al. (2004) who contributed primarily with their laboratory experimental studies of turbulent flows. In particular, PIV is a powerful technique to perform velocity measurements for a large number of flows (Adrian 1991, Raffel et al. 1998) and has been used in velocity measurements for shallow flows (Weitbrecht et al. 2002, Meselhe et al. 2004, Tauro et al. 2012b). However, development of physically-based regression equations, to calculate travel times using particle tracking in hydrology for small and ungauged catchments, has not been attempted.

In this study, the travel time of each particle in the experimental watershed is computed using the PTM sub-module that uses flow velocities simulated by DWMPT at each time step (Fig. 1). Over each computational time step, the particle travels over a distance determined by the product of the appropriate velocities and the time interval. The PTM is integrated within the framework of the DWMPT model to simulate the travel time required by the particle at the farthest cell (or designated cells) to arrive at the outlet. Thus, the travel time estimated can be used to define T_c for the watershed because it is the time required by the farthest particle to reach the outlet. The algorithm used to calculate the travel time using PTM is summarized below:

1. Flow velocities such as v_N and v_S in Fig. 1 at the inter-cell faces are calculated by the DWMPT by solving the flow equations. The PTM sub-module is activated at the desired time (at the beginning of the rainfall event for impervious surfaces). A single particle at each cell of the discretized domain is initialized.

2. Check the location of the particle; if the particle is at the outlet cell, it is removed from the domain; and PTM counts number of particles arriving at the outlet. Find the direction of steepest gradient to determine the direction of flow. The direction of flow can also be designated in the input file if the travel time at a defined path (e.g. along a stream) has to be determined.

3. Find the first estimate of the particle velocity v₁ (Fig. 1) by linear interpolation. Based on the current position of the particle, its velocity is interpolated using modelled velocities of the two opposite inter-cell faces such as v_N and v_S in Fig. 1 in the direction of flow.

4. Calculate the distance ∆d that is travelled by the particle using the velocity v₁ during the time interval ∆t.

5. Find the second estimate of particle velocity v₂ at the new location d + ∆d by linear interpolation. Using an average of the velocities at current (v₁) and new location (v₂), a new estimate of ∆d is calculated. This process is iterated until the error between the old and new estimate of ∆d is reduced to a pre-specified tolerance. The Courant condition (C_r ≤ 1) is implemented in DWMPT, and the particle does not travel more than a cell length during a computational time interval.

6. Check whether or not the particle crosses the cell boundary. If the particle travels to the next cell, a new velocity is interpolated using the new particle cell location. The distance travelled for the remaining time interval with the new velocity is calculated similar to the above steps.

7. The travel path length and travel time of the particle are updated along with the particle position.

8. The above steps are iterated for the next time step with updated velocities simulated by DWMPT until all the particles in the discretized domain reach the outlet.

Fig. 1 Two-dimensional DWMPT finite difference grids illustrating the movement of particle i from the old location d at time t to new location d + ∆d at time t + ∆t. The grid sizes Δx and Δy are the same for DWMPT.

Before PTM was used to determine the particle travel time, DWMPT was validated using observed velocity, flow depth, and discharge measurements published in the literature.

MODEL VALIDATION

Model validation using discharge hydrographs

The DWMPT was first validated using observed discharge hydrographs published by Chow (1967) and Ben-Zvi (1984), who conducted experiments in the Watershed Experimentation System of the University of Illinois. The experiments were conducted on a 12.2 × 12.2 m (40 × 40 ft) masonite and aluminum-surfaced watershed. The rainfall intensities in the experiments varied from 51 to 305 mm h^-1 (2–12 in h^-1), and the rainfall durations from 0.5 to 8 min. The longitudinal slope of the watersheds varied from 0.5% to 3% and cross slope was fixed at 1%. Areas exposed to the rainfall had nine different configurations (Ben-Zvi 1984); the configuration called Lot 3, which had rainfall in the entire watershed, was modelled using DWMPT in this study. Simulated discharge hydrographs at the outlet for four rainfall events on different slopes, selected from Lot 3 data on an aluminum surface from Ben-Zvi (1984) experiments, were compared with observed ones (Fig. 2). Observed and simulated hydrographs along with rainfall hyetographs for a longitudinal slope of 0.5%, 1.0%, 1.5%, and 2.0% are shown in Fig. 2(a)–(d), respectively. The hyetographs with durations of 8 min of the experiments presented in Fig. 2(a)–(d) are uniform rainfall intensities of 215 mm h^-1 (8.45 in h^-1), 291 mm h^-1 (11.46 in h^-1), 220 mm h^-1 (8.65 in h^-1) and 221 mm h^-1 (8.72 in h^-1), respectively. DWMPT simulated hydrographs were generated using calibrated Manning’s roughness coefficient of 0.02 for an aluminum surface. The agreement between observed and simulated peak discharges and hydrographs is excellent (Fig. 2 and Table 1).

Fig. 2 Observed effective rainfall hyetographs and observed and simulated hydrographs for an aluminum surface, 12.2 m long and 12.2 m wide, with a cross slope of 1% and longitudinal slope of: (a) 0.5%, (b) 1.0%, (c) 1.5%, and (d) 2.0%. Observed data presented above are from Lot 3 of Ben-Zvi (1984).

Table 1 Input variables, peak discharges Q_p estimated using the Rational method, observed data, and modelled using DWMPT, and corresponding model performance parameters calculated between observed and modelled discharge hydrographs for impervious overland flow planes.

L (m)	So (%)	n	I (mm h^-1)	Qp (×10^-3 m^3 s^-1)			Ns	RMSE (×10^-3 m^3 s^-1)
				Rational method	Observed	Modelled
152.4 ^a	2.0	0.013	189	0.086	0.086	0.086	0.98	0.116
152.4 ^a	0.5	0.011	50	0.023	0.023	0.023	0.98	0.031
152.4 ^a	0.5	0.035	51	0.023	0.023	0.023	0.99	0.024
152.4 ^a	0.5	0.011	22	0.010	0.010	0.010	0.95	0.024
12.2 ^b	0.5	0.020	215	0.313	0.338	0.313	0.94	0.895
12.2 ^b	1.0	0.020	291	0.425	0.439	0.425	0.97	0.831
12.2 ^b	1.5	0.020	220	0.321	0.334	0.320	0.98	0.497
12.2 ^b	2.0	0.020	221	0.323	0.331	0.323	0.97	0.733

Note: Input variables for the experimental overland flow planes are L: length (m), S_o: slope (%), n: Manning’s roughness coefficient (-), and i: effective rainfall intensity (mm h^-1). Model performance parameters are Ns: Nash-Sutcliffe coefficient and RMSE: root mean square error.

^aExperimental data from Yu and McNown (1963).

^bExperimental data from Ben-Zvi (1984).

Yu and McNown (1963) measured runoff and depth hydrographs for different combinations of slope, roughness, and rainfall intensity at the airfield drainage system in Santa Monica, California. Table 1 lists four rainfall events from Yu and McNown (1963), and four events from Ben-Zvi (1984), along with input variables (L in m, S_o in percent, n, and i in mm h^-1), peak discharge, Q_p (in m³ s^-1) estimated using the Rational method, observed and simulated Q_p values, and model performance parameters. The goodness of fit between modelled and observed quantities was measured using the Nash-Sutcliffe coefficient, Ns (Legates and McCabe 1999) and root mean square error, RMSE. For a hydrograph simulation, a good agreement between simulated and measured discharges is achieved when Ns exceeds 0.7 (Bennis and Crobeddu 2007). The average Ns was 0.97 (range: 0.94–0.99), and average RMSE was 0.04 × 10^-3 m³ s^-1 (range: 0.02–0.09 × 10^-3 m³ s^-1) for eight simulated hydrographs (Table 1). The Q_p estimates from the Rational method agreed well with measured and simulated Q_p values (Table 1). The average error between the Rational method Q_p and observed Q_p was 2.5% (range: 0.2–7.4%), and average error between simulated Q_p and observed Q_p was also 2.5% (range: 0.1–7.5%). These statistics indicate close agreement between measured and simulated hydrographs.

Model validation using depth hydrographs

The DWMPT was also validated using observed depth hydrographs from Yu and McNown (1963). They measured observed depth to an accuracy of a few ten-thousandths of a foot. Figure 3(a) and (b) show simulated and observed depth hydrographs from concrete surfaces of 152.4 m (500 ft) by 0.3 m (1 ft) with a slope of 2% and 0.5%, respectively, and Fig. 3(c) for a turf surface of the same dimension with a slope of 0.5%. Depth hydrographs are reported at 142.3 m (467 ft) from the upstream boundary in Fig. 3(a) and (c) but at 101.5 m (333 ft) for Fig. 3(b). The hyetographs for the experiments in Fig. 3(a), (b) and (c) are effective rainfall intensities of 189 mm h^-1 (7.44 in h^-1) with duration of 8 min, 49 mm h^-1 (1.93 in h^-1) with duration of 16 min, and 49 mm h^-1 (1.93 in h^-1) with duration of 28 min, respectively. Depth hydrographs simulated using DWMPT matched well with the observed ones (Fig. 3). The average Ns was 0.96, and average RMSE was 0.63 mm for three simulated depth hydrographs and indicate close agreement with measured depth hydrographs.

Fig. 3 Observed effective rainfall hyetographs and observed and simulated runoff depth hydrographs for a plot, 152.4 m long and 0.3 m wide with a slope of: (a) 2% at 142.3 m from upstream, (b) 0.5% at 101.5 m from upstream, and (c) 0.5% at 142.3 m from upstream. Observed data presented in (a) and (b) are for concrete surfaces, and (c) for a turf surface from Yu and McNown (1963).

Model validation using velocity observations

The DWMPT was also validated using calculated and observed velocity data from Yu and McNown (1963), Izzard and Augustine (1943), and Mügler et al. (2011). Yu and McNown (1963) used an ogee-shaped weir to measure runoff at the outlet. Since the critical flow condition was created for the measurement of the observed discharge hydrograph, the calculated velocity can be obtained from the observed discharge hydrograph using the critical flow equation. Four example comparisons using calculated velocity data from Yu and McNown (1963) are shown in Fig. 4(a)–(d). Calculated and simulated runoff velocity hydrographs at the outlet of a 152.4 m (500 ft) long and 0.3 m (1 ft) wide concrete surface with slopes of 2%, 0.5% and 0.5% are shown in Fig. 4(a), (b) and (d), respectively, and for a turf surface with a slope of 0.5% in Fig. 4(c). The hyetographs for the experiments presented in Fig. 4(a)–(d) were uniform effective rainfall intensities of 189 mm h^-1 (7.44 in h^-1) with a duration of 8 min, 50 mm h^-1 (1.98 in h^-1) with duration 16 min, 51 mm h^-1 (2.0 in h^-1) with duration 28 min, and 22 mm h^-1 (0.85 in h^-1) with duration 22 min, respectively. Even though the turf surface is a pervious surface, it was simulated by DWMPT as an impervious surface. The rainfall in the experiment was applied for a longer time to make the rainfall loss uniform. The rainfall input to the DWMPT was the effective rainfall intensity. DWMPT simulated velocity hydrographs were generated using calibrated Manning’s roughness coefficients of 0.011–0.013 for concrete and 0.035 for turf surfaces (Table 2). The calculated velocities compared well with the simulated velocities at the outlet having the critical flow condition (Fig. 4).

Fig. 4 Observed effective rainfall hyetographs and observed calculated and simulated runoff velocity hydrographs at the outlet for a plot, 152.4 m long and 0.3 m wide, with a slope of: (a) 2%, (b) 0.5%, (c) 0.5%, and (d) 0.5%. Observed and calculated data presented in (a), (b) and (d) are for concrete surfaces, and (c) for turf surface from Yu and McNown (1963) based on the critical flow condition at the outlet.

Table 2 Input variables, model performance parameters for simulated velocity hydrographs, time of concentration T_c, travel times (T_t85, T_t95 and T_t100) for impervious overland flow planes.

L (m)	So (%)	n	I (mm h^-1)	Ns	RMSE (cm s^-1)	Tc (min)			Travel time (min)
						Experiment	Model	KC et al. (2014)	Tt85	Tt95	Tt100
152.4 ^a	2.0	0.013	189	0.98	1.6	4.6	4.1	4.1	4.7	5.4	6.2
152.4 ^a	0.5	0.011	50	0.97	1.1	11.7	10.8	10.2	12.0	13.9	16.2
152.4 ^a	0.5	0.035	51	0.97	1.2	22.6	21.3	21.5	24.4	28.1	33.1
152.4 ^a	0.5	0.011	22	0.96	1.1	16.9	14.9	14.2	16.8	19.5	22.7
3.7 ^b	2	0.013	49	0.66	2.5	3.2	2.0	0.9	0.9	1.0	1.3
21.9 ^b	0.1	0.013	46	0.96	1.6	8.0	8.1	7.3	7.0	8.2	9.6
21.9 ^b	0.1	0.013	94	0.97	1.6	6.3	6.5	5.6	5.4	6.4	7.6
21.9 ^b	0.1	0.013	98	0.96	2.2	6.7	6.5	5.5	5.2	6.3	7.5

Note: Input variables for the experimental overland flow planes are L: length in m, S_o: slope (%), n: Manning’s roughness coefficient (-), and i: effective rainfall intensity in (mm h^-1). Travel times for 85%, 95% and 100% of particles arriving at the outlet are T_t85, T_t95, and T_t100, respectively. Model performance parameters are Ns: Nash-Sutcliffe coefficient and RMSE: root mean square error.

^aExperimental data from Yu and McNown (1963).

^bExperimental data from Izzard and Augustine (1943).

Izzard and Augustine (1943) collected rainfall and runoff data from paved and turf surfaces for the Public Roads Administration in 1942. Runoff was measured by a 0.4 ft wide HS flume with head measured by a point gauge. The calculated velocity hydrograph was obtained from the observed discharge hydrograph using the critical flow equation. Table 2 lists four rainfall events from Yu and McNown (1963) and four events from Izzard and Augustine (1943), along with input variables L, S_o, n and i; and model performance parameters (Ns and RMSE) developed between simulated and calculated velocity hydrographs. The average Ns was 0.89 (range: 0.66–0.98), and average RMSE was 2.0 cm s^-1 (range: 1.1–2.5 cm s^-1) for eight simulated velocity hydrographs (Table 2). These statistics indicate a strong agreement between measured and simulated velocity hydrographs.

Mügler et al. (2011) collected velocity measurements to compare different roughness models in a 10 m (32.8 ft) long by 4 m (13.1 ft) wide sandy soil surface with a 1% longitudinal slope. The rainfall simulator provided a constant average intensity of 69 mm h^-1 (Esteves et al. 2000). Rainfall was applied for more than 2 h to maintain steady runoff before the flow velocity measurements were performed. The steady state discharge was 0.5 × 10^-3 m³ s^-1 (0.018 ft³ s^-1). The flow velocities were measured at 60 monitoring points on the surface using salt velocity gauge technology, an automated, miniaturized device based on the inverse modelling of the propagation of a salt plume (Tatard et al. 2008). The observed flow velocities ranged from 0.006 to 0.27 m s^-1 and provided a valuable dataset for comparison with DWMPT simulated velocities. Figure 5 shows observed vs simulated velocities at 60 locations within the plot. The Ns for velocity data pairs is 0.76, and RMSE is 0.03 m s^-1. All the above statistics, derived from velocity observations from Yu and McNown (1963), Izzard and Augustine (1943), and Mügler et al. (2011), and simulations by DWMPT, indicate that the flow model can predict runoff velocities with an acceptable accuracy. It can, therefore, be used for particle tracking to determine travel times of the particles.

Fig. 5 Observed and simulated velocities for a plot of sandy soil surface, 10 m long and 4 m wide, with a slope of 1.0% and constant average rainfall intensity of 70 mm h^-1. Observed data are from Mügler et al. (2011).

COMPUTATION OF TRAVEL TIMES AND TIME OF CONCENTRATION

Flow velocity of overland flow

Previous studies (Abrahams et al. 1986, Maidment 1993, Kull and Feldman 1998, Grimaldi et al. 2012) have acknowledged the difficulty of measuring the runoff speed of overland flow in the whole study domain with a fine spatial resolution. Hence, some researchers have assumed constant velocity to compute travel time of the overland flow (Maidment 1993, Kull and Feldman 1998). However, the flow characteristics (such as flow depth in Fig. 2 and flow velocity in Figs 4, 5 and 6) of the overland flows are spatially varied (non-uniform) and can be unsteady before reaching equilibrium condition and after the cessation of the rainfall event. Figure 6 shows simulated velocity distributions from the upstream boundary to the downstream outlet at different simulation times for a concrete surface of 152.4 m (500 ft) by 0.3 m (1 ft) with a slope of 2% under a constant effective rainfall intensity of 189 mm h^-1 (7.44 in h^-1) with duration of 8 min. The length of the plot was discretized into 500 cells, and the x-axis represents the distance of each cell (grid) centre from the outlet (152.2 m for the upstream boundary cell and 0.2 m for the outlet cell). Before surface runoff reaches equilibrium condition (time t = 1, 2 and 3 min), flow velocity increases from the upstream cell (near zero velocity) towards downstream cells up to a certain distance and then becomes uniform. Therefore, before reaching equilibrium, there is a downstream portion with the front of fast moving particles and an upstream portion with slow moving particles trailing behind. At about 4 min, the uniform flow portion disappears, which indicates that the fast moving particles (front) from upstream have arrived at the outlet. From 4 to 8 min, the velocity distribution remains more or less the same and is spatially varied. When the rainfall stops at 8 min (Fig. 6), runoff velocities at all cells begin to decrease with time. The velocity distributions after the cessation of the rainfall event (t > 8 min) are spatially varied over the whole plot length and are different from the distributions prior to equilibrium (t < 4 min). In comparison to the velocity distributions for a steady rainfall event, velocity distributions of overland flow under unsteady rainfall intensities would be highly unsteady and spatially varied and are not graphically illustrated here.

Fig. 6 Simulated runoff velocity at different times (t = 1, 2, 3, 4, 8, 9 and 10 min) at each cell from upstream to downstream for a plot, 152.4 m long and 0.3 m wide, with a slope of 2%. The length of the plot is discretized into 500 cells.

Travel times

Based on spatially varied velocity distributions (for example, Fig. 6), pollutant particles released at different locations on the plot would have different travel times reaching the outlet. The particles near the upstream boundary of the plot travel slowly until they travel to a certain distance away where they gradually start picking up velocity. Therefore, the travel time of a particle released or washed-off by the runoff near the upstream boundary would be significantly longer than the travel time of other particles released at the other middle locations of the plot. Hence, to correctly estimate the overall travel time of runoff/pollutant particles on a plot under a rainfall event, both the travel times of the slow-moving particles at the upstream boundary area and the fast-moving particles at the moving front should be considered. In this study, one particle is assigned for each computational cell and initiated when the rainfall begins, for example, the simulation domain of the plot for Fig. 6 has 500 cells with 500 virtual particles used for PTM, i.e. 500 particles were tracked by DWMPT to determine the percentage of particles exiting the outlet in order to study various characteristic travel times for overland flow.

Including these slower particles near the upstream boundary makes the overall travel time of the plot longer. A similar argument also applies to discharge when it approaches equilibrium under a constant rainfall. It takes much longer time for discharge to reach equilibrium when the rate of change of discharge near equilibrium is very small (KC et al. 2014). Hence, to avoid this sensitivity of travel time to computational equilibrium (McCuen 2009), T_c is usually computed at the time less than 100% contribution of the watershed to peak discharge (Su and Fang 2004, Wong 2005, KC et al. 2014). Following similar arguments, characteristic travel times for 85%, 95% and 100% of particles arrival at the outlet were computed and were defined as T_t85, T_t95 and T_t100, respectively, based on the time when a specified percentage of particles have arrived at the outlet of the plot. It should be noted that the authors intended to choose some percentages less than 100%; however, the exact selection of T_t85 and T_t95 was arbitrary.

Time of concentration

Another commonly used hydrological time parameter (McCuen 2009), closely linked to travel time, is T_c. Despite its frequent use, T_c lacks a universally accepted definition and method to quantify it. Several definitions of T_c can be found in the literature along with several related estimation methods for each definition. Grimaldi et al. (2012), using their case studies, have demonstrated that available approaches for the estimation of T_c may yield numerical predictions that differ from each other by up to 500%.

McCuen (2009) divided multiple definitions of T_c into two broad categories; travel time for runoff particles reaching the watershed outlet and lag time relating to the distribution of rainfall and runoff. McCuen then divided the latter definition into six more computational methods for estimating T_c. However, the estimation of T_c by the latter definition requires knowledge of the distribution of the rainfall hyetograph and runoff hydrograph, which are not available for ungauged watersheds. Hence, T_c defined as the travel time of runoff has often been used in many previous studies and also in this study. Even though researchers as early as Mulvany (1851) and Kuichling (1889) to McCuen (2009) defined T_c similarly, as the time it takes for runoff to travel from the most distant point along a hydraulic pathway in the watershed to the outlet, their computational implementation of estimating T_c differed. An indirect method such as hydrograph analysis is often used to estimate T_c and then develop T_c regression formulas with respect to certain independent variables.

Several researchers (Kuichling 1889, Hicks 1942, Muzik 1974) estimated T_c as the time until discharge from a given drainage area reaches equilibrium. Izzard and Hicks (1946) defined T_c from the beginning of rainfall until runoff reaches 97% of the rainfall input rate. Wong (2005) considered T_c as the time from the beginning of effective rainfall to the time when flow reaches 95% of the equilibrium discharge. Su and Fang (2004) and KC et al. (2014) determined T_c as the time from the beginning of effective rainfall to the time when the flow reaches 98% of the peak discharge. The same method is also used in this study to compute T_c.

Figure 7 shows an example of how three characteristic travel times (T_t85, T_t95 and T_t100) and T_c were computed from either the percent of particles arrival at the outlet or from the outflow hydrograph. Both observed and simulated hydrographs are shown in Fig. 7 for an asphalt pavement, 21.9 m long and 1.83 m wide with a slope of 0.1% and a rainfall of 98.3 mm h^-1 for 11 min (Izzard and Augustine 1943). The occurrence times of 98% of peak discharges (Q_p) and 100% Q_p in the discharge hydrograph are shown in Fig. 7 with three characteristic travel times (T_t85, T_t95, and T_t100). Table 2 lists three T_c estimates (in min) using observed and simulated hydrographs and an empirical equation from KC et al. (2014) for non-low slopes, and travel times T_t85, T_t95, and T_t100 (in min) derived from DWMPT for four rainfall events from Yu and McNown (1963) and four events from Izzard and Augustine (1943). Three T_c estimates agree well with each other. The estimates of travel times T_t85, T_t95, and T_t100 are different from T_c because travel times and T_c are computed using two separate methods. Figure 7 shows that T_t100 occurs when the discharge is near the equilibrium discharge and travel times T_t95 and T_t85 occur when the discharge approaches the equilibrium discharge. The rate of change of discharge over time near equilibrium is very small as shown in Fig. 7, however, it also depends on the input parameters. The discharge approaches the equilibrium faster when the plot length is shorter, the slope is steeper, the surface roughness is smaller, and the rainfall intensity is larger. Figure 7 as an example and Table 2 for eight experiments show that T_c may correlate with T_t85 and T_t95, which will be further explored below.

Fig. 7 Observed and simulated hydrographs, percentage of total number of particles arrival at the outlet, and time of concentration (T_c), and travel times for 85% (T_t85), 95% (T_t95) and 100% (T_t100) particles arrival for an asphalt pavement, 21.9 m long and 1.83 m wide, with a slope of 0.1% and an effective rainfall of 98.3 mm h^-1. The location occurrence times of 98% of peak discharge (Q_p) and 100% of Q_p in the discharge hydrograph are also shown. Observed hydrographs are from Izzard and Augustine (1943).

In this study, the travel times are computed based on the percentage of particles arrival at the outlet using PTM coupled with DWM velocity simulation, which is to directly compute the travel time of particles through the simulation domain. However, T_c is computed from the discharge hydrograph as the time from the beginning of effective rainfall to the time when the flow reaches 98% of the peak discharge. Therefore, T_c directly depends on the flow depth of a single cell at the outlet. However, the travel times T_t85, T_t95 and T_t100 are consolidated time parameters synthesizing the arrival times of most or all particles at the outlet cell and are more representative time parameters to describe the movement of runoff/pollutant particles inside the simulation domain. To develop regression equations for estimating travel times (T_t85, T_t95 and T_t100), based on physically-based input parameters, provides a simple and easy method to determine the overall travel time characteristics of the overland flows. The parametric study to find the relationships between the travel times and the input variables (L, S_o, n, and i) is provided in the next section.

PARAMETRIC STUDY OF THE TRAVEL TIME OF OVERLAND FLOW

Several authors (Morgali and Linsley 1965, Woolhiser and Liggett 1967, Su and Fang 2004, KC et al. 2014) have derived estimation formulas using the four physically-based input variables: L, S_o, n, and i. Hence, the same input variables are chosen for parametric study computing travel times of overland flows. More than 530 separate values of each time parameter (T_c, T_t85, T_t95, and T_t100) were determined from numerical experiments using DWMPT by varying input variables, L, S_o, n, and i to extend the dataset available for regression analysis. The input variable L was varied from 5 to 305 m (16–1000 ft), S_o from 0.1% to 10%, n from 0.01 to 0.80, and i from 2.5 to 254 mm h^-1 (0.1–10.0 in h^-1). The numerical experiments were simulated holding the three variables constant and varying the fourth one by 10–20%.

A generalized power relation, equation (1) was chosen for developing the regression equations:

(1)

where L is in m; S_o is in m m^-1 (-); i is in mm h^-1; and C₁, k₁, k₂, k₃, and k₄ are regression parameters. Equation (1) was log-transformed and a nonlinear regression was used to estimate parameter values. The resulting equation for T_t85 is:

(2)

where T_t85 is in min, and other variables are as previously defined. Regression results are presented in Table 3. Statistical results indicate that the input variables L, S_o, n, and i have a high level of significance with the p value < 0.0001 (Table 3) and are critical variables in the determination of T_t85. The regression parameters (C₁, k₁, k₂, k₃, and k₄) have small standard errors and small ranges of variation at the 95% confidence interval (Table 3). Predicted values from equation (2) compare well with the estimates from DWMPT numerical experiments (Fig. 8).

Fig. 8 Travel time for 85% of particles arrival at the outlet (T_t85) predicted using regression equation (2) vs T_t85 derived from numerical experiments using DWMPT for impervious overland flow planes.

Table 3 Parameter estimates of the regression equation (2) for the independent variables of travel time for 85% of particles arrival at the outlet (T_t85) for impervious overland flow surfaces.

Parameter	Parameter estimate	95% confidence limits		Standard error	t value	p value
In(C2)	2.225	2.209	2.240	0.008	284.42	<0.0001
k1 for L	0.599	0.597	0.602	0.001	473.66	<0.0001
k2 for So	−0.303	−0.304	−0.301	0.001	−331.74	<0.0001
k3 for n	0.609	0.607	0.611	0.001	699.21	<0.0001
k4 for i	−0.399	−0.401	−0.397	0.001	−381.35	<0.0001

Similar regression analysis resulted following equations for T_t95:

(3)

and for T_t100:

(4)

The estimation equations for T_t85, T_t95, and T_t100 along with the statistical parameters: coefficient of determination R², RMSE, and p values are summarized in Table 4. The R² and RMSE were developed against respective time parameter (T_t85, T_t95, and T_t100) data generated from 530 DWMPT model runs for the parametric study. The estimation equations have an excellent average R² of 0.99 and low average RMSE of 1.06 min. The p values reported herein were developed between respective time parameters and all four input variables (L, S_o, n, i) and are <0.0001, showing the high level of statistical significance for each estimation equation (Table 4). The regression equations for three travel time parameters have the same or almost the same exponents for each of the input variables L, S_o, n, and i (Table 4). Travel times are directly proportional (nonlinearly) to L and n, and inversely proportional (nonlinearly) to S_o and i. Saghafian et al. (2002) also showed that for spatially uniform, yet temporally variable excess rainfall intensity, the travel time is inversely proportional to the rainfall intensity raised to a power of 0.4. Designers can use the estimating equations (Table 4) to determine the overall travel times of 85%, 95%, and 100% of particles arriving at the outlet of an impervious plot by simply using physically based input parameters.

Table 4  Statistical error parameters for regression equations of three time parameters (T_t85, T_t95, and T_t100) developed for impervious overland flow surfaces.

Function, equation	Formula	R^2	RMSE ^a (min)	p value ^b
Tt85 = f (L, So, n, i), equation (2)		0.999	1.1	<0.0001
Tt95 = f (L, So, n, i), equation (3)		0.998	1.0	<0.0001
Tt100 = f (L, So, n, i), equation (4)		0.996	1.1	<0.0001

^aStatistical parameter R² and RMSE were developed against respective time parameters (T_t85, T_t95, and T_t100) generated from 530 DWMPT model runs for the parametric study.

^bThe p values reported herein were developed between respective time parameters and all four input variables (L, S_o, n, i) in each regression equation.

The time parameters T_c, T_t95, and T_t100 show a linear variation with T_t85 (Fig. 9). This is because the regression equation for T_c developed by KC et al. (2014) has similar exponents to four input variables as the equations in Table 5. Hence, the time parameters T_c, T_t95, and T_t100 can also be expressed in linear relationships with T_t85 (Table 5). The equation for T_c in terms of T_t85 is:

(5)

Fig. 9 Time of concentration (T_c), travel times for 95% (T_t95) and 100% (T_t100) of particle arrival at the outlet vs travel time for 85% of particle arrival at the outlet (T_t85) developed from numerical experiments using DWMPT for impervious overland flow planes. Three linear regression lines between three time parameters (T_c, T_t95, and T_t100) and T_t85 are also shown.

Table 5 Statistical error parameters for linear regression equations between three time parameters (T_c, T_t95, and T_t100) and T_t85 developed for impervious overland flow surfaces.

Equation	Formula	R^2	RMSE ^a (min)
Equation (5)		0.949	7.83
Equation (6)		0.998	1.55
Equation (7)		0.996	7.97

^aStatistical parameter R² and RMSE were developed against respective time parameters (T_c, T_t95, and T_t100) derived from DWMPT and predicted using the corresponding equation for 530 data points for the parametric study.

Similarly the equations for T_t95 and T_t100 in terms of T_t85 are:

(6)

and

(7)

Equation (5) shows that T_t85 estimated using the particle-tracking method could be used to compute T_c with acceptable accuracy. Equations (6) and (7) show T_t95 and T_t100 are about 17% and 43% larger than T_t85, respectively. The equations for T_c, T_t95, and T_t100 in terms of T_t85 along with R² and RMSE are given in Table 5. The R² and RMSE were developed against respective time parameters (T_c, T_t95, and T_t100) from DWMPT and predicted using the corresponding equation (Table 5) for 530 data points. The linear equations have an excellent average R² of 0.98 and average RMSE of 5.78 min (Table 5). Estimations of three travel times were useful to check the inter-relationships between these travel times. From Tables 4 and 5, it can be seen that T_t95 can be approximately expressed as an average of T_t85 and T_t100. Even though all the travel times (T_t85, T_t95, and T_t100) show the nonlinear relationships with the input variables L, S_o, n, and i, the characteristic time parameters showed linear relationships between each other (T_t85, T_t95, and T_t100) and T_c).

SUMMARY AND CONCLUSIONS

A quasi two-dimensional diffusion wave model with particle tracking, DWMPT, was developed from the diffusion hydrodynamic model DHM (Hromadka II and Yen 1986) for calculating the travel time of all the particles in the discretized domain of the impervious overland surface. A comprehensive validation of the model was performed using discharge, velocity and depth data from four different sources (Izzard and Augustine 1943, Yu and McNown 1963, Ben-Zvi 1984, Mügler et al. 2011). The travel times of all particles in the domain consisting of both slow moving particles at the upstream and fast moving particles at the downstream are considered in DWMPT, and consolidated representative time parameters T_t85, T_t95 and T_t100 are extracted. Based on the percentage of particles that have arrived at the outlet, travel times for 85%, 95% and 100% particles arrival at the outlet are defined as T_t85, T_t95 and T_t100, respectively. A total of 530 DWMPT runs were performed to determine travel times on impervious surfaces using diverse combinations of the four physically-based input variables L, S_o, n and i, and were used to develop estimation equations for T_t85, T_t95 and T_t100. All the travel times (T_t85, T_t95 and T_t100) showed nonlinear relationships with the input variables (L, S_o, n, and i) in equations (2), (3), and (4); and are directly proportional to L and n (both raised approximately to a power 0.6) and inversely proportional to S_o (raised approximately to a power 0.3) and i (raised approximately to a power 0.4). Similarly, T_c was estimated from hydrograph analysis as the time when outflow discharge reaches 98% of Q_p, and T_t85 showed a close correlation with T_c as shown in equation (5) and can be used to estimate T_c. Even though all the travel times showed nonlinear relationships with the input variables, they showed linear relationships with each other in equations (6) and (7), and T_t95 can be approximately expressed as an average of T_t85 and T_t100.

Acknowledgements

The authors wish to express their thanks to the Texas Department of Transportation (TxDOT) and its members for their guidance and support for the study; and to the Civil Engineering Department at Auburn University for teaching assistantship support. The authors thank Dr A. Ben-Zvi and Dr C. Mugler for providing data; their data were used in Figs 2 and 5. The authors also thank Dr David B. Thompson and Dr Theodore Cleveland for their inputs to the study. The anonymous reviewers are thanked for their useful comments.

Additional information

Funding

This study was partially supported by TxDOT Research Project 0–6382.