Generally, roofs are the best candidates for rainwater harvesting. In this context, the correct evaluation of the quantity and quality of runoff from roofs is essential to effectively design rainwater harvesting systems. This study aims to evaluate the performance of a kinematic wave based numerical model in simulating runoff on sloping roofs, by comparing the numerical results with the ones obtained from laboratory rainfall simulations on a real-scale *Lusa* ceramic tile roof. For all studied slopes, simulated discharge hydrographs had a good adjust to observed ones. Coefficient of determination and Nash–Sutcliffe efficiency values were close to 1.0. Particularly, peak discharges, times to peak, peak durations and runoff volumes were very well simulated.

## INTRODUCTION

Rainwater harvesting is a water saving alternative strategy that presents many advantages and can provide solutions to address major water resources problems, such as fresh water scarcity, urban stream degradation and flooding (e.g. Leusbrock *et al.* 2015). In recent years, these problems have become global challenges, due to climatic change, population growth and increasing urbanisation.

Generally, roofs are the first to come into contact with rainwater; thus, they are the best candidates for rainwater harvesting. In this context, the correct evaluation of roof runoff quantity and quality is essential to effectively design rainwater harvesting systems (e.g. Farreny *et al.* 2011). Despite this, many studies usually focus on the qualitative aspects (e.g. Lye 2009; Egodawatta *et al.* 2012) in detriment of the quantitative aspects. Quantifying roof runoff depends not only on the climatic conditions (e.g. rainfall) but also on the roof type. In the Mediterranean region, ceramic tile roofs are a well-known trademark (Farreny *et al.* 2011). In Portugal, one of the most used is the *Lusa* ceramic tile (Sousa *et al.* 1998).

A reliable rainfall–runoff numerical model can provide an efficient and economic tool to predict beforehand possible impacts of climatic changes and land use changes (e.g. increasing urbanisation), allowing the development of alternative water management policies and water saving techniques (e.g. Borris *et al.* 2013; Jung *et al.* 2015). Despite the large amount of studies, rainfall–runoff modelling still involves several difficulties, given the high spatial and temporal rainfall variability (e.g. Notaro *et al.* 2013).

For decades, kinematic wave modelling has been used as a reliable tool in water resources (see the reviews in Singh 1996, 2001). Despite its wide acceptance, several authors continue investigating kinematic wave accuracy and applicability in a wide range of runoff studies (e.g. de Lima & Singh 2002; Deng *et al.* 2005; Isidoro & de Lima 2013; Abrantes *et al.* 2015; Mallari *et al.* 2015).

Laboratory and field studies using rainfall simulators have been widely used to investigate rainfall–runoff processes (e.g. Egodawatta & Goonetilleke 2008; Arguelles *et al.* 2013; de Lima *et al.* 2013; Montenegro *et al.* 2013; Ries *et al.* 2014). These studies enabled a detailed exploration and systematic replication of a large range of hydrologic conditions, such as rainfall spatial and temporal characteristics, providing for a fast way to obtain precise and consistent data that can be used to calibrate and validate numerical models (e.g. Deng *et al.* 2005; Abrantes *et al.* 2015; Mallari *et al.* 2015).

This study aims to evaluate the performance of a kinematic wave based numerical model in simulating runoff on sloping roofs, by comparing the numerical results with the ones obtained from laboratory rainfall simulations on a real-scale *Lusa* ceramic tile roof.

## MATERIALS AND METHODS

### Laboratory installation

*et al.*(2013), Montenegro

*et al.*(2013), Carvalho

*et al.*(2014) and Isidoro & de Lima (2015)), a

*Lusa*ceramic tile roof and a continuous discharge measuring system.

The rainfall simulator (Figure 1) has a pressurised hydraulic system comprised of: (i) a steady downward-oriented full-cone nozzle (1/4-HH-14W FullJet from Spraying Systems Co., USA), with 3.58 mm orifice diameter, positioned 2.2 m above the geometric centre of the roof surface; (ii) a hydraulic system attached just in front of the nozzle to eliminate pressure fluctuations (more details in Isidoro & de Lima 2015); and (iii) a submerged pump (76.2 mm SQ from Grundfos Holding A/S, Denmark), installed in a constant head reservoir supplied with tap water. This system allowed a steady operating pressure of 1.4 bar at the nozzle, producing a total discharge of 7.5 L min^{−1}, with a spray angle of 120 °.

The ceramic tile roof (Figure 1) consisted on an assembly of seven rows by seven lines of *UM Lusa* ceramic tiles (from Umbelino Monteiro SA, Portugal), covering a total area of 4.2 m^{2} (1.5 × 2.8 m^{2}), supported by a metal structure with adjustable longitudinal slope. According to the manufacturer, *UM Lusa* ceramic tiles' dimensions are 0.455 × 0.254 m^{2}, weight 3.15 kg, have low permeability, high resistance to frost and are produced according to Portuguese technical standards. Manning's roughness coefficient for ceramic tiles was considered 58.8 m^{1/3} s^{−1} (Chow 1959).

Experiments were conducted only considering the five most central tile rows of the roof, positioned at the slopes of 3.5, 23.0 and 40.5%. The two most peripheral rows worked as buffer zones in order to compensate for water ejected outside the other rows due to splash.

### Rainfall characteristics

Raindrop mean diameter and mean velocity, calculated from measurements with a laser precipitation monitor (Thies Laser Distrometer from Adolf Thies GmbH & Co., Germany), were, respectively, 0.52 mm and 1.41 m s^{−1} (more details in Carvalho *et al.* (2014)).

### Tiles absorption intensity

*UM Lusa*ceramic tiles absorption intensity, for the three defined slopes. The tiles, pre-dried in an oven during 24 h at 105 °C, were placed over high precision weighting devices and were subjected to simulated rainfall during 30 min. Monitoring the combined weight of the tiles and absorbed water, it was possible to characterise the tiles absorption intensity temporal evolution, as shown in Figure 3. These tests were repeated three times for each defined slope.

### Experimental procedure

Rainfall simulations were conducted to produce 15 discharge hydrographs for the five most central rows of the roof positioned at the defined slopes of 3.5, 32.0 and 40.5%. Each simulation comprised one constant rainfall event, with the duration of 1 min, approximately equal to the time of concentration (time to peak) of the roof positioned at the lowest slope of 3.5% (the highest time of concentration of the three slopes). All simulations were repeated three times to assure their reliability.

At the beginning of the experiments, the roof was wetted for 5 min. After that, rainfall events were simulated consecutively with a 5–10 min interval between them. Thus, the tiles' initial moisture and absorption condition were approximately the same for all simulations.

Discharge hydrographs at the roof downstream end were obtained using the continuous discharge measuring system (Figure 1). This system consisted of plastic recipients placed over weighting devices with a precision of 0.01 g (one for each row), to monitor runoff during rainfall events and runoff recessions, in a total of 3 min, with a temporal resolution of 1 s. Runoff was conveyed to the plastic recipients using funnels, minimising runoff losses and the potential turbulence caused by falling water, thus obtaining smoother hydrographs.

## KINEMATIC WAVE MODELLING OF ROOF RUNOFF

*h*(

*x*,

*t*) is the runoff depth [L],

*Q*(

*x*,

*t*) is the runoff discharge [L

^{3}T

^{−1}],

*q*(

*x*,

*t*) is the rainfall excess intensity [L T

^{−1}],

*W*is the runoff width [L],

*t*is the simulation time [T] and

*x*is distance along the flow direction [L].

*S*) equals the friction slope (

*S*) and using existing open channel flow friction equations, we can express the runoff discharge at any point and time as function of the runoff depth, as follows: where

_{f}*α*is the kinematic wave resistance parameter [L

^{1/3}T

^{−1}] and

*β*is the Bakmeteff dimensionless exponent [-]. Since runoff was assumed as turbulent, Manning's formula can be used to estimate

*α*and

*β*, as and , where

*k*is the Manning's roughness coefficient [L

_{M}^{1/3}T

^{−1}] and

*S*is the bed slope [L L

^{−1}].

*q*) was estimated as the rainfall intensity (

*p*) subtracted by the tiles absorption intensity (

*a*), as follows: where

*p*is the rainfall intensity [L T

^{−1}] and

*a*is the tiles absorption intensity [L T

^{−1}].

*t*) adjusted by an initial absorption time (

*t*) that represents the tiles absorption intensity at the beginning of the simulation, as follows: where

_{a}*t*is the initial absorption time [T],

_{a}*A*is a constant of the power curve [L T

^{−2}] and

*B*is an dimensionless exponent of the power curve [-]. The power function was adjusted from the absorption tests (Figure 3). Since the tiles' initial moisture and absorption condition were approximately the same for all simulations, simulated hydrographs were obtained considering only three initial absorption times (one for each slope), representing the tiles' absorption intensity at the beginning of the simulation. Therefore, according to the rainfall characteristics (e.g. spatial pattern) and the absorption characteristics of the tiles (e.g. initial moisture condition), rainfall excess intensity varied in space and time.

*i*identifies the time (

*t*),

*j*identifies the space (

*x*),

*Δt*is the time discretisation [T] and

*Δx*is the space discretisation [L].

*L*is the total length of the surface along the flow direction [L].

*r*

^{2}), the Nash–Sutcliffe efficiency (

*NSE*) and the Nash–Sutcliffe efficiency with logarithmic values (

*lnNSE*), comparing observed to simulated discharge at the downstream boundary, as follows (e.g. Krause

*et al.*2005): where r

^{2}is the Determination coefficient [-], NSE is the Nash-Sutcliffe efficiency [-], lnNSE is the Nash-Sutcliffe efficiency with logarithmic values [-], is the observed discharge at the downstream boundary and instant

*i*[L

^{3}T

^{−1}], is the simulated discharge at the downstream boundary and instant

*i*[L

^{3}T

^{−1}], is the average observed discharge at the downstream boundary [L

^{3}T

^{−1}], is the average simulated discharge at the downstream boundary [L

^{3}T

^{−1}] and

*n*is the total runoff discharge data. The range of

*r*

^{2}lies between 0 and 1. A value of 0 means no correlation at all whereas a value of 1 means that the dispersion of the simulation is equal to that of the observation. The range of

*NSE*and

*lnNSE*lies between −∞ and 1. An efficiency of 1 means a perfect fit and an efficiency lower than 0 indicates that the mean value of the observed time series would have been a better predictor than the model.

## RESULTS AND DISCUSSION

### Experimental results

Considering the standard deviation, both runoff initiations, peak discharges and runoff volumes were very similar for all roof slopes, being, respectively, 13.3 ± 0.8 s, 7.2 ± 0.5 mL s^{−1} and 424.8 ± 20.2 mL. For all roof slopes, runoff coefficient of the roof was approximately 0.87. Furthermore, small differences on these parameters can be attributed to differences of rainfall intensity on each slope. Nevertheless, the gentle slope of 3.5% clearly produced a delay on the rising and recession limbs, when compared to the other slopes. Time to peak was 57.4 s, 23.0 s and 20.2 s, respectively, for the slopes of 3.5%, 33.0% and 40.5%.

### Numerical results

Coefficient of determination (*r*^{2}), *NSE* and *lnNSE*, comparing the observed to simulated discharge (Equations (10)–(12)) are also shown in Figure 5. The numerical model had a good performance on simulating runoff for all three slopes, since in r* ^{2}* and

*NSE*values were close to 1.0. However,

*lnNSE*values were very low. While

*NSE*strongly overestimates peak discharge values whereas lower values are neglected, with

*lnNSE*the influence of the low discharge values is increased in comparison to the peak ones. Therefore, it is clear that the numerical model better estimated peak discharge values than the lower one. This happened due to the already noticed differences in time to runoff.

On average, the best model performance was obtained for the middle slope of 23.0%, with an *r*^{2} and an *NSE* of 0.98 and an *lnNSE* of 0.46. The worst adjustment was obtained for the highest slope of 40.5%, with an *r*^{2} and an *NSE* of 0.96 and an *lnNSE* of 0.42. However, differences in *r ^{2}*,

*NSE*and

*lnNSE*values for all three slopes were negligible.

## CONCLUSIONS

A kinematic wave based numerical model was used to runoff on a *Lusa* ceramic tile roof, with good adjustment for all studied slopes. Coefficient of determination and NSE values were close to 1.0. Simulated peak discharges, times to peak and peak durations were close to the observations from the rainfall experiments. However, lower discharge values at the beginning of runoff were not simulated with a good performance.

Future work should include other types of ceramic tiles (e.g. *Marseille*, *Canudo*), other types of roofs, initial moisture conditions and absorption conditions and rainfall intensities (e.g. low rainfall intensity below absorption intensity).

## ACKNOWLEDGEMENTS

The first author acknowledges CNPq, Brazil for the financial support through the Post Doctoral Grant 206872/2014-3. The second and third authors acknowledge FCT, Portugal for the financial support through the Doctoral Grant SFRH/BD/103300/2014 and the Project PTDC/ECM/HID/4259/2014. The authors also acknowledge the contribution of Antonio Mousinho Fernandes and Gabriel Silveira for their help during the experimental work.

The laboratory experiments described in this study were conducted at the Laboratory of Hydraulics, Water Resources and Environment of the Department of Civil Engineering of the Faculty of Science and Technology of the University of Coimbra, Portugal.