Caviedes 2012
Caviedes 2012
Caviedes 2012
Journal of Hydrology
journal homepage: www.elsevier.com/locate/jhydrol
a r t i c l e i n f o s u m m a r y
Article history: Hydrological simulation of rain-runoff processes is often performed with lumped models which rely on
Received 6 September 2011 calibration to generate storm hydrographs and study catchment response to rain. In this paper, a distrib-
Received in revised form 6 February 2012 uted, physically-based numerical model is used for runoff simulation in a mountain catchment. This
Accepted 1 April 2012
approach offers two advantages. The first is that by using shallow-water equations for runoff flow, there
Available online 16 April 2012
This manuscript was handled by
is less freedom to calibrate routing parameters (as compared to, for example, synthetic hydrograph meth-
Konstantine P. Georgakakos, Editor-in-Chief, ods). The second, is that spatial distributions of water depth and velocity can be obtained. Furthermore,
with the assistance of Ioannis K. Tsanis, interactions among the various hydrological processes can be modeled in a physically-based approach
Associate Editor which may depend on transient and spatially distributed factors. On the other hand, the undertaken
numerical approach relies on accurate terrain representation and mesh selection, which also affects
Keywords: significantly the computational cost of the simulations. Hence, we investigate the response of a gauged
Finite volumes catchment with this distributed approach. The methodology consists of analyzing the effects that the
Shallow-water equations mesh has on the simulations by using a range of meshes. Next, friction is applied to the model and the
Mesh generation response to variations and interaction with the mesh is studied. Finally, a first approach with the well-
SCS-CN known SCS Curve Number method is studied to evaluate its behavior when coupled with a shallow-water
model for runoff flow. The results show that mesh selection is of great importance, since it may affect the
results in a magnitude as large as physical factors, such as friction. Furthermore, results proved to be less
sensitive to roughness spatial distribution than to mesh properties. Finally, the results indicate that SCS-
CN may not be suitable for simulating hydrological processes together with a shallow-water model.
Ó 2012 Elsevier B.V. All rights reserved.
1. Introduction Experimental data for the catchment exist and are used to validate
numerical results. The objective is to understand the interactions
Runoff processes at the catchment scale have been modeled of various model parameters with the numerical model and com-
extensively with empirical and lumped methods as to generate prehend the effect that such interactions have on runoff simula-
the hydrograph response of storm events for various applications. tion. In the matter of hydrologic and hydraulic parameters
It is nevertheless possible to make use of two-dimensional shal- subject to the choice of modelers, it is necessary to characterize
low-water models to perform distributed and physically-based the response of the model to each parameter in order to provide
runoff analysis at the catchment scale. This is not done in many in- parameter sensitivity information and a grasp of the effects of
stances because of possibly large computational time as well as the uncertainty of parameter selection. Manning’s roughness coeffi-
resolution of available spatial data and difficulties when it comes cient and Soil Conservation Service Curve Number (SCS-CN)
to parameter calibration. Shallow water models allow to obtain parameters are those of particular interest in this work.
the same type of results as lumped models, namely runoff volume The solution of the shallow-water equations is done by means
and outflow hydrographs, but also spatial distribution of water as of a first order, explicit, upwind finite volume numerical scheme
well as runoff flow vectors. Hence, results are physically-based and (Murillo and Garcia-Navarro, 2010; Murillo et al., 2007). It is
therefore arguably better, but more complete, in the sense that widely known that finite volume schemes are heavily dependant
more information is obtained. on mesh properties such as cell shape (triangular, rectangular,
In this paper, a 2D shallow-water model is used to simulate etc.) and size as well as the general mesh structure (Leveque,
the response of a Pyrenean catchment during a storm event. 2002). It is quite clear that as the number of cells increases, CPU
time increases. This becomes an important issue since the very
⇑ Corresponding author. Address: Ed. Torres Quevedo, María de Luna 3, CP 50015, complex topography of mountain catchments requires a large
Zaragoza, Spain. Tel.: +34 976761881; fax: +34 976761882. number of cells. Furthermore, the explicit nature of the scheme im-
E-mail address: [email protected] (D. Caviedes-Voullième). poses a time step constraint so as to satisfy the CFL condition
0022-1694/$ - see front matter Ó 2012 Elsevier B.V. All rights reserved.
http://dx.doi.org/10.1016/j.jhydrol.2012.04.006
40 D. Caviedes-Voullième et al. / Journal of Hydrology 448–449 (2012) 39–59
(Murillo et al., 2009). The time step is smaller as finer meshes are as geological description and vegetation maps which are the basis
used. From another perspective, the finest mesh resolution for a gi- for hydrological characterization of the catchment.
ven problem is given by the available topography data, and is sel- From the experimental perspective, the response of the Arnás
dom chosen by the modeler. However, cell shape and mesh catchment to storms has been described as fast, and the shape of
structure are modeler choices. Hence, it appears important to ex- the hydrograph has been found to reflect the precipitation struc-
ploit such control in order to obtain a mesh that will allow to ade- ture (Seeger et al., 2004). From the analysis of experimental data
quately represent the complex topography and the hydraulic it has been concluded that the response to large storms shows sea-
phenomena of interest at the same time that a minimal number sonal variation, with highest runoff generation in springtime, while
of cells is obtained, thus generating the largest possible cells and in autumn and summer runoff generation is strongly limited (Lana-
only refining locally, where the topography demands it. In this pa- Renault et al., 2007). Some of the studied storms with moderate
per, the effects of mesh resolution, combined with square and tri- precipitation (total precipitation larger than 10 mm) did not gener-
angular cells are examined. This is further studied by using ate runoff and, in general, high storm-runoff variability is observed.
uniformly refined triangular meshes versus terrain-adapted lo- Furthermore, Lana-Renault et al. (2007) conclude that the rainfall
cally-refined triangular meshes. magnitude and the type of rain influence the response, while max-
The methodology for this work consists of separating the anal- imum rainfall intensity was found to correlate poorly. They also
ysis for each of the mentioned issues (i.e., mesh effects, SCS-CN, conclude that it is likely that runoff is generated mainly by infiltra-
Manning’s n) with a spectrum of cases, as well as combining sce- tion excess during the dry season (Hortonian flow), while satura-
narios in which interactions between the parameters can be seen. tion excess appears to be the primary runoff generation process
With such methodology it is possible to assess the magnitude of during the wet season. Garcia-Ruiz et al. (2005) characterize spring
the variations and to which parameter or combination of parame- and autumn as wet season, winter as an intermediate period and
ters they can be attributed. summer as the dry season and also argue that a large spatial vari-
ability of the runoff-contributing areas exists in the catchment.
2. Catchment description Both studies (Garcia-Ruiz et al., 2005; Lana-Renault et al., 2007)
conclude that the water table depth has limited influence on catch-
The Arnás catchment is located in the northern Spain Pyrenees, ment response to storms.
in the Borau valley, and has a surface of 2.84 km2, ranging in alti- From the numerical simulation point-of-view, López-Barrera et
tude from around 900 to 1340 m.a.s.l. The catchment has been al. (2011) simulated the Arnás catchment using a distributed ap-
gauged since 1995 which allows a vast recollection of hydrological proach, based on a 2D diffusive wave model for surface flow, and
data. The primary channel in the catchment is the Arnás ravine both Horton and Green-Ampt models for infiltration, finding diffi-
which flows into the Lubierre river. Flow measurements have been culties in reproducing outflow hydrographs. López-Barrera (2011)
performed at the outlet on the Arnás ravine (see Fig. 1) with 5 min also used the SCS-CN method to estimate infiltration losses in
frequency. Meteorological data have been collected by a meteoro- the Arnás catchment, not with a distributed modeling approach
logical station at the outlet location, and rainfall has been regis- (as the aforementioned 2D diffusive wave model), but with a
tered by a tipping bucket rain gauge with a frequency of 5 min. lumped approach, namely the SCS Unit Hydrograph method. Ser-
Geologically, the catchment lies over Eocene flysch formations rano-Pacheco (2009) simulated the catchment with a distributed
and has suffered land use and coverage changes in recent decades, approach, but using the complete 2D shallow-water equations
generating a mixed vegetation cover which ranges from forest for surface flow and Horton’s model for infiltration. Both studies
patches, dense and open shrubs, grassland cover and bare land conclude that infiltration plays a major role in the catchment re-
(Garcia-Ruiz et al., 2005). Previous studies (Garcia-Ruiz et al., sponse and very large precipitation losses occur because of infiltra-
2005; Lana-Renault, 2007; Lana-Renault et al., 2007; López-Barre- tion. Therefore, infiltration modeling plays a major role in the
ra, 2011; Serrano-Pacheco, 2009) have provided information such numerical response.
France
42°38'24" N
0°35'03" W
Spain
3. Surface flow model classified and numerically dealt with as belonging to the family
of hyperbolic systems. The mathematical properties of (1) include
In this paper, runoff is simulated using the 2D shallow-water the existence of a Jacobian matrix, Jn , of the flux normal to a
equation as shown in Eq. (1). In this section the 2D formulation direction given by the unit vector n, E n ¼ Fnx þ Gny , (Murillo
of the shallow-water equations is presented as depth averaged and Garcia-Navarro, 2010) defined as
equations expressing volume conservation and momentum
@E n @F @G
conservation Jn ¼ ¼ nx þ ny ð7Þ
@U @U @U
@U @FðUÞ @GðUÞ
þ þ ¼SþH ð1Þ The presented equations are solved by an explicit, first order, up-
@t @x @y
wind, finite volume scheme described in Murillo and Garcia-Navar-
where ro (2010).
U ¼ ðh; qx ; qy ÞT ð2Þ
4. Simulation cases
are the conserved variables with h representing depth (m) and
qx ¼ hu and qy ¼ hv the unit discharges (m2/s), with u and v m/s A wide range of cases was proposed with the goal to observe the
the depth averaged components of the velocity vector u along the different results of simulations with varying degrees of detail and
x and y coordinates respectively. The fluxes are given by complexity. The only investigated event in this paper (referred to
T !T as Event 1) occurred on 13/10/2005.
q2 1 2 qx qy qx qy q2y 1 2 Fig. 2 the simplest setup was to assume a no-friction, no-infil-
F ¼ qx ; x þ gh ; ; G¼ qy ; ; þ gh ð3Þ
h 2 h h h 2 tration case (F0-I0) – an ideal case of course – but which provides
interesting information. By contrasting results for such case with
where g (m/s2) is the acceleration of the gravity. The source terms of different meshes, the effect of mesh selection on the solution of
the system are split in two kind of terms. The term S is defined as the shallow-water equations for this particular problem can be
found, mass conservation can be easily verified and a lower limit
sb;x sb;y T for time-to-peak can be obtained, as well as a maximum peak out-
S ¼ 0; ; ð4Þ
q q flow. All meshes reported in Table 1 were used for this case.
with sb;x , sb;y (N/m2) the bed shear stress in the x and y direction The following series of simulations are those of friction effects.
respectively, with q (kg/m3) the density of the fluid. The friction A uniform average Manning roughness coefficient was used for the
losses are written in terms of the Manning’s roughness coefficient entire catchment obtained by averaging the roughness map re-
n (s m1/3): ported by Serrano-Pacheco (2009). This allows to observe the ef-
pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 4=3 pffiffiffiffiffiffiffiffiffiffiffiffiffiffiffiffi 4=3 fects of friction in the same mesh, and the interrelated effects
Sfx ¼ n2 u u2 þ v 2 =h ; Sfy ¼ n2 v u2 þ v 2 =h ð5Þ between mesh shape and resolution and friction. The next step
was to simulate using the roughness map instead of the uniform
The term H is defined as
average value, to assess if spatial variability (both real and artificial
T due to mesh properties) of the roughness coefficients was a signif-
@z @z
H ¼ 0; gh ; gh ð6Þ icant factor.
@x @y
The next set of simulations deals with infiltration, by means of
and expresses the variation of the pressure force along the bottom the SCS-CN method. Experimental data provide values of outlet
in the x and y direction respectively, formulated in terms of the bed runoff volume and effective precipitation. With these data, lumped
slopes of the bottom level z (m). estimates with SCS-CN can be made to calibrate CN and the initial
System (1) is time dependent, nonlinear, and contains source abstraction ratio a. With the lumped calibration as a starting point,
terms. Under the hypothesis of dominant advection it can be distributed SCS-CN simulations can be performed to obtain a new
0 0.200
Accumulated Precipitation (mm) in 5 minutes
0.175
0.5
0.150
1
0.125
Flow (m3/s)
1.5
0.100
2
0.075
2.5
0.050
3
0.025
3.5 0.000
0 50 100 150 200 250 300 350 400 450 500 550
Time (min)
7 8
6 7
Percentage of cells (%)
2
2
1 1
0 0
10 12 14 16 18 20 22 24 26 28 30 32 34 36 38 40 42 0 100 200 300 400 500 600 700 800 900 1000 1100 1200 1300 1400 1500 1600 1700
10 16
9
14
8
Percentage of cells (%)
5 8
4
6
3
4
2
2
1
0 0
0 100 200 300 400 500 600 700 800 900 1000 1100 1200 1300 0 100 200 300 400 500 600 700 800 900 1000 1100 1200 1300
2)
Cell area (m Cell area (m2)
dependency on mesh selection (Leveque, 2002), in terms of both its Consider a hydrologic analysis point of view where runoff vol-
structure and resolution. Furthermore, it is expected that the afore- ume and hydrograph properties are chief interests (ignoring for
mentioned properties of the mesh should have an effect on the now the spatial distribution). Even while considering perfect repre-
solution because of its implications in the quality of the represen- sentation of infiltration and roughness parameters, mesh resolu-
tation of hydrological parameters. Hence, in order to examine the tion impacts the resulting outflow hydrograph and runoff
significance of mesh selection, several meshes were generated for volume. The latter may be affected when using coarse meshes be-
the catchment. cause high frequency variations in the topography may not be
High mesh resolution is, in terms of accuracy, always favorable. sampled. Depending on the averaging method to produce a repre-
However, in most cases mesh resolution is limited either by avail- sentative elevation of that cell from finer spatial data, local maxima
able spatial data or by the fact that it implies an important cost of or minima may be lost or generated (López-Barrera, 2011; Merw-
computational resources. Cell size and, in consequence, mesh res- ade, 2009), thus changing the available volume for surface ponding
olution, is not bounded by physical ranges (Hardy et al., 1999) and (surface storage). Similar effects may be observed when high reso-
hence is subject to the modeler choice and experience, using lution bathymetry data is not available for a river, although high
numerical convergence as a mesh quality index, and invoking the resolution topographic data may be available for the floodplain.
idea that if a finer mesh yields very similar results as a slightly The outflow hydrograph can also be affected by poor topography
coarser mesh, then the adequate resolution has been found. When representation, since volume may be artificially lost or generated
high resolution spatial data are available, by selecting an appropri- by such mechanisms. This is also true for physically-based, distrib-
ate cell size larger than the available data, the mesh resolution can uted models, and can be further evaluated by using the inundation
be appropriately selected to obtain satisfactory results without zone as a benchmarking result (Casas et al., 2006; Cook and Merw-
excessive computational cost. But other issues play a role in this ade, 2009; Hardy et al., 1999; Horritt et al., 2006). But furthermore,
selection. Complex topography such as mountainous terrains or these errors in topographic representation affect flow computa-
small details in large terrains such as levees, ditches, roads, rail- tions since terrain slopes may be poorly represented and cell shape
roads and buildings may require high resolution meshes to cor- and mesh structure may introduce numerical viscosity of – a priori
rectly represent topographic features and small scale processes – unknown magnitude. The effects of mesh resolution and struc-
which can sometimes impact large scale flow behavior (Horritt ture can be considered part of the larger problem of appropriate
and Bates, 2001; Murillo and Garcia-Navarro, 2004; Schubert representation of spatial data, which also includes issues such as
et al., 2008). uncertainty of elevation data (Wilson and Atkinson, 2005). Thus,
44 D. Caviedes-Voullième et al. / Journal of Hydrology 448–449 (2012) 39–59
although it is not the intention of this paper, the conclusions which of the effects that uncertainty in elevation measurements or repre-
can be gathered from introducing variations in elevation represen- sentation may have on flow results. The methodology followed
tation because of mesh selection can also aid in the understanding herein generates meshes that in one way or another may represent
D. Caviedes-Voullième et al. / Journal of Hydrology 448–449 (2012) 39–59 45
1.04
Precipitated volume
and maximum cell area, thus TU-10-40 is an unstructured triangu-
Total simulated volume lar mesh with a minimum cell area of 10 m2 and maximum area of
Area factor
1.03 40 m2. The minimum area is used as an indicator of the smallest
cell, which affects the CFL condition. Locally refined meshes are
Normalized Volume
1.02 named with their total cell number, e.g, TU-LR-16666 is a locally
refined mesh with 16,666 cells.
1.01 TU meshes are more or less homogeneous, in the sense that cell
area histograms show little dispersion as shown in Fig. 4a. This is
1 not a hindrance for result quality, but the large number of triangles
that it implies when cells are small, is a problem for computational
0.99
time. Hence a local mesh refinement algorithm was implemented,
based on the gradient of terrain slope. Note that the local refine-
ment procedure described herein is not (strictly) an adaptive mesh
0.98
refinement process since refinement criteria is based only on ter-
SS
SS
SS
SS
TS
TS
TS
TS
TU
TU 0
LR 00
5
10
15
20
-2
5
10
15
20
16
18
rain data, not flow computations. Hence, the term local refinement
26
-6
-1
24
(LR) is preferred. The LR algorithm is as follows, and is illustrated in
Fig. 6. Comparison of volumes at the end of the simulation. Fig. 3. The original DEM raster data – shown in Fig. 3a – with 5 m
resolution are processed to obtain a slope vector for each pixel as
shown in Fig. 3b. Then, the gradient of the slope vector field is com-
poorly the real topography (this is, that topography which has been puted thus generating a new second rank tensor field with such
measured accurately). But the errors which are artificially gener- gradients. The magnitude of the vector field in each pixel is then
ated are similar in nature to those generated by poorly measuring computed as in Fig. 3c, and categories are defined to assign a gra-
the real terrain surface. dient index to a particular gradient range so that the user may
Roughness is generally represented by means of Manning’s choose to have as many categories as needed. Gradient indexes
coefficient. Determination of the appropriate values may be done are stored in a new raster map with resolution equal to that of
through reference data and improved by calibration. Nevertheless, the DEM raster as in Fig. 3d. A triangulation is generated using
the calibration process may be affected by other model parameters, the domain boundary with Triangle using maximum area and min-
in particular mesh resolution and structure, and the resulting cal- imum angle constrains. The next step is to map gradient indices
ibrated roughness values may be mesh-dependent, thus compen- into the triangular mesh. This is done by identifying the maximum
sating numerical errors introduced by the mesh (Hardy et al., gradient inside a triangle, an mapping onto that triangle the index
1999; Horritt and Bates, 2001). Therefore, it appears necessary to within a gradient range classification, as in Fig. 3e. When all trian-
evaluate the effects of different spatial representations of rough- gles have an index, new maximum area constraints are imposed on
ness since different permutations produce different runoff flow the triangles. The constraints are related to the gradient classes
properties (Hardy et al., 1999). For any mesh, the simplest proce- that were previously defined. This information is then passed to
dure is to assign a domain uniform average roughness value. On Triangle and a triangulation refinement process is done which uses
the other hand, it is possible to represent roughness values for the new area restrictions. Thus, a triangulation which has more tri-
smaller areas, down to cell size, thus representing each cell with angles in regions with large slope gradients is obtained. In Fig. 3f
its own particular roughness. It is therefore interesting to observe the new triangulation is shown together with the index raster,
the response of different meshes with the same roughness repre- where it can be seen how the indices govern the triangle genera-
sentation and vice versa since sensitivity to the roughness value tion. If the mesh is still too coarse the algorithm can be repeated
may change with mesh resolution (Casas et al., 2006), which is rea- until results are satisfactory. Finally, the elevation is assigned to
sonable since roughness parameters intend to represent energy the LR triangulation as the arithmetic mean elevation of those pix-
loss phenomena which can be affected by mesh resolution and els which are contained or intersected by each triangle, as shown
model scale (Hunter et al., 2007). in Fig. 3g. It is possible that more complex interpolation techniques
(area weighted, inverse distance weighted, Krigin) would yield bet-
ter discrete terrain surfaces, but this has not been explored in this
6. Mesh generation work. According to Murillo and Garcia-Navarro (2004) the initial
mesh which is to be refined also has an impact on the results.
This section is devoted to test the sensitivity of the simulation The sensitivity of the results to the first triangulation was not
to different mesh configurations. The base data for all meshes investigated in this work either. The resulting LR meshes have a
was a raster DEM with 5 m resolution. From this, both square (S) much wider cell area histogram than non refined TU meshes as
and triangular (T) meshes were generated. Triangular meshes were shown in Fig. 4. Although there is a large proportion of cells in
further subdivided into structured (TS) and unstructured meshes the smaller range, note that the maximum cell area values are sig-
(TU), and the latter were subdivided into homogeneous (TU) or lo- nificantly greater than the average cell area values.
cally-refined meshes (TU-LR). SS mesh generation was always per- Table 1 shows properties of all meshes. The area factor is a mea-
formed so that complete pixels were included inside square cells sure of the catchment area error generated by using coarser
(cell edge length as a multiple of pixel edge length). Once created, meshes, which cannot represent exactly the same catchment
the cell elevation was averaged from raster data. TS meshes were boundary as the mesh with best resolution. By coarsening the
created in the same fashion, with further division of the square mesh, the generated squares or triangles may exist in some part
cells into two right-angled triangles. TU meshes were generated outside the original boundary, hence increasing slightly the area
with Triangle (Shewchuk, 1996). The cell size was controlled by of the catchment. This is accounted for in order to assess the mag-
constraining the maximum cell area and the minimum interior an- nitude of the distortion it introduces in the computed volumes.
gle. The name codes for the meshes are as follows. For both square In order to illustrate the differences between different meshes, a
and triangular structured meshes the edge length is indicated. snapshot of them is shown in Fig. 5, showing a 3D rendering of the
Thus, TS5 is a structured triangular mesh with 5 m catheti. The terrain surface scaled 1.5 times in the vertical coordinate to facili-
names of unstructured triangular meshes indicate the minimum tate its understanding. The color scale represents the norm of slope
46 D. Caviedes-Voullième et al. / Journal of Hydrology 448–449 (2012) 39–59
22.50 22.50
100000 100000
20.00 20.00
17.50 17.50
80000 80000
Volume (m3)
Volume (m3)
15.00 15.00
Flow (m3/s)
Flow (m3/s)
Simulated outflow Simulated outflow
12.50 Outlet runoff volume 60000 12.50 Outlet runoff volume 60000
Surface storage volume Surface storage volume
10.00 Total volume 10.00 Total volume
Precipitated volume Precipitated volume
40000 40000
7.50 7.50
5.00 5.00
20000 20000
2.50 2.50
0.00 0 0.00 0
0 5000 10000 15000 20000 25000 30000 35000 0 5000 10000 15000 20000 25000 30000 35000
22.50 22.50
100000 100000
20.00 20.00
17.50 17.50
80000 80000
Volume (m3)
Volume (m3)
15.00 15.00
Flow (m3/s)
5.00 5.00
20000 20000
2.50 2.50
0.00 0 0.00 0
0 5000 10000 15000 20000 25000 30000 35000 0 5000 10000 15000 20000 25000 30000 35000
22.50 22.50
100000 100000
20.00 20.00
17.50 17.50
80000 80000
15.00 15.00
Flow (m3/s)
Flow (m3/s)
Volume (m3)
Volume (m3)
Simulated outflow Simulated outflow
12.50 Outlet runoff volume 60000 12.50 Outlet runoff volume 60000
Surface storage volume Surface storage volume
10.00 Total volume 10.00 Total volume
Precipitated volume Precipitated volume
40000 40000
7.50 7.50
5.00 5.00
20000 20000
2.50 2.50
0.00 0 0.00 0
0 5000 10000 15000 20000 25000 30000 35000 0 5000 10000 15000 20000 25000 30000 35000
Fig. 7. Hydrographs and volumes for case F0-I0-E1 with structured meshes.
vector of each cell. A good representation of the topography re- still represents fairly well the slopes in the main channel (blue1
quires that both elevation and slope are adequately modeled and zones with very little slope) as well as the very high slopes in the
this figure compares both. The figure shows a view from the south, hills. It still keeps, in particular, the less steep zone in the middle
near the outlet point. The view is in the upstream direction, i.e., of the first water divide to the right. TS15 starts to loose such reso-
northwest, as to show the valley, in particular the main channel lution, and SS20 averages the elevation so much that the slope
topography and the eastern hills with the highest slopes and high- changes significantly in the main channel and the less steep area
est slope gradients. Fig. 5a shows the SS5 mesh with elevation data of the first water divide is also lost. The next set of meshes are the
exactly as in the original raster elevation data. This is considered TU meshes shown in Fig. 5e–g. All of them appear to represent fairly
the reference and the terrain surface which best approximates
the real terrain. Fig. 5b–d shows the spatial representation by
structured meshes. Not all are shown since TS meshes are just a tri- 1
For interpretation of the references to color in Figs. 1–20, the reader is referred to
angular partition of the square cells of SS meshes. Note that TS10 the web version of this article.
D. Caviedes-Voullième et al. / Journal of Hydrology 448–449 (2012) 39–59 47
22.50 22.50
100000 100000
20.00 20.00
17.50 17.50
80000 80000
15.00 15.00
Flow (m3/s)
Flow (m3/s)
Volume (m3)
Volume (m3)
Simulated outflow Simulated outflow
12.50 Outlet runoff volume 60000 12.50 Outlet runoff volume 60000
Surface storage volume Surface storage volume
10.00 Total volume 10.00 Total volume
Precipitated volume Precipitated volume
40000 40000
7.50 7.50
5.00 5.00
20000 20000
2.50 2.50
0.00 0 0.00 0
0 5000 10000 15000 20000 25000 30000 35000 0 5000 10000 15000 20000 25000 30000 35000
Fig. 7 (continued)
well the terrain, although they introduce more variability in slope tion. Although there is no experimental result for the no-friction,
than the original SS5 representation. Nevertheless, the shapes are no-infiltration case, the result for SS5 is considered as the reference
accurate, and the extremes of the slope range appear to be well rep- solution. The response with SS5 is very quick and replicates the
resented. The last group is that of locally refined meshes. These are shape of the hyetograph, which is expected because of the friction-
shown in Fig. 5h–j. It is difficult to appreciate differences between less and impervious conditions. Note also that surface storage
them but, compared to TU and SS meshes, it can be seen that the reaches a maximum value around the same time as peak flow,
quality of the representation is as good as that of TU and sometimes when a large portion of the precipitated volume is static or moving
better that the SS or TS representations. The most important differ- on the surface. Surface storage reduces to a residual value which is
ence, however, is that LR representations require significantly less composed of flowing water and an important proportion of pond-
cells to achieve very similar topographic representation. Note that ed, static water in topographic depressions, i.e., local elevation
in all LR representations the cells near water divides, at the bottom minima. As the SS meshes become coarser the hydrograph suffers
of the valleys and in the main channel are smaller than those cells in a diffusive effect which is similar to hydrograph damping, although
the hill slopes. In the hill slopes high slopes may be present but they for SS20 some high frequency variations appear which may be
do not change significantly and hence less cells are required to rep- attributed to the coarser step-wise paths that water must travel
resent the topography correctly. to reach the outlet. In terms of time-to-peak, coarser meshes result
in larger times. In summary the hydrograph losses somewhat tem-
poral resolution, it is dampened and suffers lag. Another significant
7. Mesh effects: results and discussion result is that surface storage increases with coarser meshes. This is
due to the sum of two factors. The first is poor topographic repre-
In order to isolate the mesh effect, this section presents the re- sentation, since local minima and maxima may be poorly repre-
sults for the catchment considered frictionless and impervious. All sented with coarse meshes, which results in static, ponded water
the differences in outflow hydrographs are due to mesh resolution which cannot flow further. The second is numerical viscosity ef-
and cell properties. From the total volume and precipitated volume fects (and numerical diffusion) that are generated by coarser
data, it is clear that mass conservation is achieved at catchment le- meshes. Similar results are obtained for the TS progression. A com-
vel at all times in all cases. However, precipitated volume is parison of SS vs TS shows that SS hydrographs appears more damp-
slightly higher for the coarser SS and TS meshes, as can be seen ened than TS solutions. This is consequent with the fact that
in Fig. 6. This is expected since the increase in volume is due to rectangular meshes impose a strong directionality upon the flow,
the increase in area because of the approximation of the catchment and structured triangular meshes impose less directionality. Thus,
boundary with less points. Moreover, the increase in volume is rectangular meshes are said to be more viscous than structured tri-
identical to the increase in area due to the aforementioned effect, angular meshes.
and is in the worst case, a 3% increase in volume. Note that it only TU meshes are somewhat more complex to analyze. The three
occurs in structured meshes. TU meshes show that small dispersion of triangle sizes result in a
A summary of the results is presented in Fig. 7 for structured viscous mesh, as in the case of TU-10-40. The hydrograph is
meshes and in Fig. 8 for unstructured meshes. Both figures present dampened although high frequencies are still observed. Surface
outflow hydrographs for all meshes as well as the accumulated vol- storage is high, in the same magnitude as SS20. As the dispersion
ume from precipitation in the entire catchment. They also show grows the hydrograph is less dampened, and even higher fre-
the outlet runoff volume (total volume exiting the catchment at quencies appear. In TU-16-60 high frequency variation occur
the outlet), the surface storage volume instantly stored in the sur- around 10,000 s that are not present in such amplitude else-
face (volume which may be static or flowing, but in any case, inside where. In TU-18-100 peak flow shows high frequency variations
the domain) and the sum of both which is the total volume which which generate a significantly larger maximum instantaneous
is also compared against precipitation volume. flow. Timewise, both in TU-10-40 and TU-16-60 the second peak
From Figs. 7 and 8 the first conclusion is that of excellent mass is higher, thus changing significantly from SS5. The effect is evi-
conservation by comparing total volume and precipitated volume. dent for TU-10-40. TU-LR meshes seem to respond much better
However the distribution of the total volume is not the same for all than TU meshes. Surface storage is improved significantly and is
meshes. In this analysis the SS5 mesh is taken as the reference in the magnitude of SS5 results. These meshes have a very large
since it is the base, unprocessed topographic with highest resolu- cell size dispersion, which is consistent with the observation that
48 D. Caviedes-Voullième et al. / Journal of Hydrology 448–449 (2012) 39–59
22.50 22.50
100000 100000
20.00 20.00
17.50 17.50
80000 80000
Volume (m3)
Volume (m3)
15.00 15.00
Flow (m3/s)
Flow (m3/s)
Simulated outflow Simulated outflow
12.50 Outlet runoff volume 60000 12.50 Outlet runoff volume 60000
Surface storage volume Surface storage volume
10.00 Total volume 10.00 Total volume
Precipitated volume Precipitated volume
40000 40000
7.50 7.50
5.00 5.00
20000 20000
2.50 2.50
0.00 0 0.00 0
0 5000 10000 15000 20000 25000 30000 35000 0 5000 10000 15000 20000 25000 30000 35000
Time (s) Time (s)
22.50 22.50
100000 100000
20.00 20.00
17.50 17.50
80000 80000
Volume (m3)
Volume (m3)
15.00 15.00
Flow (m3/s)
12.50
Simulated outflow 60000 Flow (m3/s) 12.50
Simulated outflow 60000
Outlet runoff volume Outlet runoff volume
Surface storage volume Surface storage volume
10.00 Total volume 10.00 Total volume
Precipitated volume Precipitated volume
40000 40000
7.50 7.50
5.00 5.00
20000 20000
2.50 2.50
0.00 0 0.00 0
0 5000 10000 15000 20000 25000 30000 35000 0 5000 10000 15000 20000 25000 30000 35000
Time (s) Time (s)
22.50 22.50
100000 100000
20.00 20.00
17.50 17.50
80000 80000
15.00 15.00
Volume (m3)
Volume (m3)
Flow (m3/s)
Flow (m3/s)
5.00 5.00
20000 20000
2.50 2.50
0.00 0 0.00 0
0 5000 10000 15000 20000 25000 30000 35000 0 5000 10000 15000 20000 25000 30000 35000
Time (s) Time (s)
Fig. 8. Hydrographs and volumes for case F0-I0-E1 with unstructured meshes.
22.50 22.50 F0
F0 PF* = 20.38 m3/s
3
PF = 18.79 m /s 100000 TP* = 3521 s 100000
20.00 TP = 3536 s 20.00 F0
PF = 20.49 m3/s
17.50 TP = 5772 s
17.50 FM FM
PF = 15.73 m3/s 80000 PF = 15.64 m3/s
80000
TP = 4728 s TP = 4795 s
Volume (m3)
Volume (m3)
15.00 15.00
Flow (m3/s)
Flow (m3/s)
FA
FM outflow FM outflow
FA
PF = 15.67 m3/s F0 outflow PF = 15.63 m3/s
F0 outflow
12.50 TP = 4864 s FA outflow 60000 12.50 TP = 4744 s FA outflow 60000
FM Surface storage volume FM Surface storage volume
F0 surface storage volume F0 surface storage volume
10.00 FA surface storage volume
10.00 FA surface storage volume
40000 40000
7.50 7.50
5.00 5.00
20000 20000
2.50 2.50
0.00 0 0.00 0
0 5000 10000 15000 20000 25000 30000 35000 0 5000 10000 15000 20000 25000 30000 35000
22.50 22.50
F0 100000 F0 100000
3 3
20.00 PF = 17.81 m /s 20.00 PF = 18.30 m /s
TP = 3865 s TP = 3911 s
FM FM
17.50 PF = 16.28 m3/s 17.50
80000 PF = 15.89 m3/s 80000
TP = 4766 s TP = 4681 s
Volume (m3)
Volume (m3)
15.00 FA 15.00 FA
Flow (m3/s)
Flow (m3/s)
PF = 16.06 m3/s
FM outflow FM outflow
F0 outflow PF = 15.83 m3/s F0 outflow
12.50 TP = 5230 s 60000 12.50 TP = 4620 s 60000
FA outflow FA outflow
FM Surface storage volume FM Surface storage volume
10.00 F0 surface storage volume 10.00 F0 surface storage volume
FA surface storage volume FA surface storage volume
40000 40000
7.50 7.50
5.00 5.00
20000 20000
2.50 2.50
0.00 0 0.00 0
0 5000 10000 15000 20000 25000 30000 35000 0 5000 10000 15000 20000 25000 30000 35000
22.50 F0 22.50
PF = 17.51 m3/s
TP = 4500 s
100000 F0 100000
20.00 20.00 PF = 17.49 m3/s
FM TP = 4096 s
PF = 17.59 m3/s FM
17.50 TP = 4830 s 17.50
80000 PF = 15.90 m3/s 80000
FA TP = 5125 s
15.00 PF = 17.23 m3/s
15.00
Volume (m3)
Volume (m3)
Flow (m3/s)
Flow (m3/s)
FM outflow FA FM outflow
TP = 4844 s F0 outflow 3
PF = 15.87 m /s F0 outflow
12.50 FA outflow 60000 12.50 TP = 5285 s FA outflow 60000
FM Surface storage volume FM Surface storage volume
10.00 F0 surface storage volume 10.00 F0 surface storage volume
FA surface storage volume FA surface storage volume
40000 40000
7.50 7.50
5.00 5.00
20000 20000
2.50 2.50
0.00 0 0.00 0
0 5000 10000 15000 20000 25000 30000 35000 0 5000 10000 15000 20000 25000 30000 35000
is just 9% of the friction delay of SS5. This shows that LR meshes 10. Computational efficiency
outperform structured mesh coarsening and homogeneous
unstructured meshes. The wide range of simulations performed in this study are
These observations lead to the conclusion that despite mesh aimed to understand the sensitivity of results to mesh selection
coarsening may be properly done, as the mesh is coarsened numer- both by direct reasons as well as deeper interactions with other
ical friction due to the larger cell size is introduced and this effect phenomena, such as friction and its representation. However, com-
must be considered. Roughness coefficients, though properly se- putational efficiency is also of interest and must be considered
lected according to standard charts may add an excessive amount when deciding on mesh selection. Following this, a summary of
of friction to the model since the mesh may already generate sim- the most important benchmarking variables is presented in
ilar effects. It is clear then, that roughness coefficients must be cal- Fig. 12 together with CPU time for all the studied meshes. The re-
ibrated with each mesh that is to be used in order to compensate sults have been normalized by the corresponding SS5 results in
for numerical, mesh generated friction. each case for comparison.
D. Caviedes-Voullième et al. / Journal of Hydrology 448–449 (2012) 39–59 51
22.50 FM
22.50
PF = 19.83 m3/s F0
100000 PF = 17.02 m3/s 100000
20.00 TP = 5309 s 20.00 TP = 4571 s
F0
PF = 17.57 m3/s FM
17.50 TP = 4881 s 17.50 PF = 16.55 m3/s
80000 TP = 5561 s 80000
FA
Volume (m3)
Volume (m3)
15.00 15.00
Flow (m3/s)
Flow (m3/s)
PF = 17.00 m3/s FA
3
TP = 5241 s FM outflow PF = 16.44 m /s FM outflow
F0 outflow TP = 5454 s F0 outflow
12.50 60000 12.50 60000
FA outflow FA outflow
FM Surface storage volume FM Surface storage volume
10.00 F0 surface storage volume 10.00 F0 surface storage volume
FA surface storage volume 40000 FA surface storage volume 40000
7.50 7.50
5.00 5.00
20000 20000
2.50 2.50
0.00 0 0.00 0
0 5000 10000 15000 20000 25000 30000 35000 0 5000 10000 15000 20000 25000 30000 35000
Fig. 10 (continued)
Frictionless cases show only sensitivity to the mesh. It is clear the potential maximum retention, estimated by means of the Curve
from Fig. 12a that time to peak and peak surface storage grow rap- Number following
idly with coarser structured meshes. Peak outflow is mostly under-
25; 400
estimated (compared to SS5) save for TS5, TU-18-100 and all LR S¼ 254 ð9Þ
CN
meshes, with LR-31662 practically reproducing the same peak
flow, though slightly later. Time to peak surface storage is either The initial abstraction and potential storage are related by
well represented or underestimated by approximately 30%. With Ia ¼ aS ð10Þ
these results most SS and TS meshes can be discarded for inaccu-
racy as well as TU-10-40 and TU-16-60. The best choices are TS5, Inserting (10) in (8) yields
TU-18-100 and LR meshes. However, when also considering CPU
ðPt aSÞ2
time, LR-16666 and LR-22624 appear as best choices since TS5 Pe ¼ ð11Þ
P t þ Sð1 aÞ
has a higher computational cost as does LR-31662. Within LR-
16666 and LR-22624, all variables except time to peak flow are bet- where a is a coefficient sometimes called initial abstraction ratio and
ter represented by 22624 with very little additional computational frequently assumed as 0.2 from the original SCS-CN method (Ponce
effort. Thus LR-22624 is considered the best choice. and Hawkins, 1996; SCS, 1964). Nevertheless, various studies have
Fig. 12b and c confirms the observations although errors are suggested that a may be in a different range, rarely satisfying
dampened significantly because of friction. Surface storage errors a = 0.2 (Beck et al., 2009; Mishra and Singh, 1999; Mishra et al.,
are still the most significant because they depend most on mesh 2003) with particular results as low as a = 0.004 in the study of a
generation techniques. Both roughness representations confirm Greek catchment (Baltas et al., 2007) and between 0.01 and 0.154
that LR-22624 is the best choice for this terrain. in the Three Gorges area in China with an optimal value of 0.05
(Shi et al., 2009).
11. SCS Curve Number It should be noted that the SCS-CN method was not designed to
consider time (Woodward et al., 2002). It was thought as to com-
From its origins, SCS-CN was developed for agricultural water- pute the total runoff volume from a watershed basin and this
sheds (Woodward et al., 2002). Since then, it has been applied to may be a reason of why it does not adjust easily to particular tem-
a variety of watersheds and a wide range of engineering problems. poral rain patterns, or to catchments which have a significant evo-
In engineering practice it is common to have very little, if any, data lution of runoff in time. Furthermore, it has been reported that CN
for calibration of hydrologic models, in particular, infiltration data. appears to work better in agricultural watersheds, followed by
Complex, physically based solutions for infiltration can be formu- range lands and worst results are obtained in forested watersheds
lated by solving Richards’ equation. However, the required data (Woodward et al., 2002). These considerations, together with the
are complex and difficult to obtain, heterogeneity is difficult to empirical evidence that a 0.2 allow to question the validity of
deal with, and computational time is significant. Simpler models the method. The intention in this paper, is to calibrate the SCS-
such as Horton and Green-Ampt models are much less demanding CN method with two parameters, CN and a in order to simulate
in computational resources, and require significantly less field data runoff volume correctly.
to obtain results. Even so, in many applications, Horton or Green- However, the chosen methodology is not a lumped method, but
Ampt are still far too demanding of field data. Because of this, a distributed and time-advancing method. Hence, the SCS-CN
SCS-CN is still widely used in engineering practice with ungauged method is not applied once to the entire catchment. Effective pre-
catchments, and because of this widespread usage, exploring the cipitation is computed with the SCS-CN for every cell in every time
effects of such a simple and limited model is necessary. step, using the accumulated precipitation Pt since the beginning of
The SCS-CN is formulated to describe runoff per unit area of the storm to that particular time. Therefore, an accumulated effec-
catchment, or effective precipitation Pe (mm) as tive precipitation is obtained for that time. The increment in effec-
( tive precipitation is computed by subtracting the accumulated
ðPt Ia Þ2
for Pt > Ia effective precipitation of the previous time step. Further time ef-
Pe ¼ P t Ia þS ð8Þ
0 for Pt 6 Ia fects can be noticed when considering that all precipitation infil-
trates while the accumulated rainfall Pt(t) is less than the initial
where Pt (mm) is total accumulated precipitation, Ia (mm) is the ini- abstraction Ia. Ia is dependent on both CN and a values. As both val-
tial abstraction which infiltrates before runoff begins and S (mm) is ues grow total runoff volume is reduced, but also the time after the
52 D. Caviedes-Voullième et al. / Journal of Hydrology 448–449 (2012) 39–59
22.50 22.50
100000 100000
20.00 20.00 F0
F0 PF = 17.38 m3/s
PF = 16.33 m3/s TP = 5836 s
17.50 TP = 5807 s 17.50
80000 FM 80000
FM PF = 14.99 m3/s
15.00 PF = 14.33 m3/s 15.00
Volume (m3)
Volume (m3)
TP = 6107 s
Flow (m3/s)
Flow (m3/s)
TP = 6020 s FM outflow
FA F0 outflow
12.50 FA FM outflow 60000 12.50 PF = 15.09 m3/s 60000
FA outflow
PF = 14.31 m3/s F0 outflow TP = 5836 s
TP = 5969 s FA outflow FM Surface storage volume
10.00 FM Surface storage volume 10.00 F0 surface storage volume
F0 surface storage volume FA surface storage volume
FA surface storage volume 40000 40000
7.50 7.50
5.00 5.00
20000 20000
2.50 2.50
0.00 0 0.00 0
0 5000 10000 15000 20000 25000 30000 35000 0 5000 10000 15000 20000 25000 30000 35000
22.50 F0 F0
PF = 22.17 m3/s 22.50 3
PF = 23.96 m /s
100000 100000
20.00 TP = 3917 s TP = 3605 s
FM outflow 20.00
F0 outflow
17.50 FM FA outflow FM
PF = 15.36 m3/s 80000 17.50 80000
FM Surface storage volume PF = 15.63 m3/s
TP = 4742 s
Volume (m3)
Volume (m3)
15.00 F0 surface storage volume TP = 4800 s
Flow (m3/s)
Flow (m3/s)
FA FA surface storage volume 15.00 FM outflow
FA
12.50 PF = 15.39 m3/s 60000 PF = 15.59 m3/s
F0 outflow
60000
TP = 4939 s 12.50 TP = 4869 s FA outflow
FM Surface storage volume
10.00 10.00 F0 surface storage volume
FA surface storage volume
40000 40000
7.50 7.50
5.00 5.00
20000 20000
2.50 2.50
0.00 0 0.00 0
0 5000 10000 15000 20000 25000 30000 35000 0 5000 10000 15000 20000 25000 30000 35000
22.50 22.50 F0
F0
PF = 22.69 m3/s 100000 PF = 19.06 m3/s 100000
20.00 TP = 3656 s 20.00 TP = 3691 s
17.50 FM 17.50 FM
PF = 15.72 m3/s 80000 PF = 15.66 m3/s 80000
TP = 4996 s TP = 4859 s
Volume (m3)
Volume (m3)
15.00 15.00
Flow (m3/s)
Flow (m3/s)
FA FM outflow FA FM outflow
PF = 15.52 m3/s F0 outflow PF = 15.61 m3/s F0 outflow
12.50 TP = 4871 s FA outflow 60000 12.50 TP = 4769 s FA outflow
60000
FM Surface storage volume FM Surface storage volume
10.00 F0 surface storage volume 10.00 F0 surface storage volume
FA surface storage volume FA surface storage volume
40000 40000
7.50 7.50
5.00 5.00
20000 20000
2.50 2.50
0.00 0 0.00 0
0 5000 10000 15000 20000 25000 30000 35000 0 5000 10000 15000 20000 25000 30000 35000
beginning of the storm when runoff begins to appear is longer. 12. Curve Number effects: results and discussion
Hence, the outflow hydrograph can be affected in a complex way
by these parameters. In fact, it is possible that because of the tem- Traditionally, the SCS-CN method is calibrated by setting a = 0.2
poral effects of CN and a, the numerical response that could be and, with an initial guess based on land cover and use, performing
attributed to friction (both physical and numerical) is in part gen- multiple simulations so as to obtain the Curve Number which best
erated by the infiltration method. This is why it is interesting to satisfies the runoff volume. However, because of the considerations
evaluate the response to SCS-CN parameters independently from mentioned in Section 11, a was also set as a calibration parameter.
friction parameters. To address these issues, lumped estimates To establish the reference values, the experimental runoff vol-
for the catchment are performed using runoff volume as the target ume and the runoff coefficient were computed and are presented
variable, with both CN and a as calibration parameters in order to in Table 2. Total precipitation was 36 mm for the entire storm
compare calibration results with those obtained from the distrib- event, which results in a total volume of 99065.5 m3. Table 2 also
uted simulations. presents the results for eight lumped estimates using the SCS-CN
D. Caviedes-Voullième et al. / Journal of Hydrology 448–449 (2012) 39–59 53
0.5 0.01
SS
SS
SS
SS
TS
TS
TS
TS
TU
TU -40
TU -60
TU -10
TU R-
TU R-
5
10
15
20
5
10
15
20
-1
-1
-1
-L 0
-L 166
-L 226
0
M
calibrating both parameters, and it can be only thought of as a cal-
R 24
-3
16
ibration parameter with some sort of meaning which relates to
66
62
permeability.
It is clear that by calibrating a some time distribution effect
2 10
Time to peak outflow should be reflected on the effective hyetograph. In order to observe
Peak outflow
Time to peak surface storage such differences, the SCS-CN was applied in a time-distributed
Peak surface storage
CPU time fashion, though still spatially-lumped. This is, computing effective
Normalized CPU time
1.5 1
sults are shown in Fig. 13.
Note that, although all combinations of CN and a generate prac-
tically the same effective rainfall, the hyetographs can be drasti-
cally different. The two extreme cases are a = 0.2 and a = 0. Note
1 0.1 that the hyetograph shows non-zero values around 6000 s for
a = 0.2, while for a = 0 around 1000 s there exists effective precip-
itation already. All cases show similar tendency. Note that as
a ! 0, the effective hyetograph resembles the total precipitation
0.5 0.01 hyetograph, though in significantly different magnitude. Further-
more, it can be seen that hyetograph distribution is less sensitive
SS
SS
SS
SS
TS
TS
TS
TS
TU
TU -40
TU -60
TU -10
TU R-
TU R-
5
10
15
20
5
10
15
20
-1
-1
-1
-L 0
-L 166
-L 226
M
R 24
-3
62
tions might be found for a and CN that satisfy runoff volume, the
choice may affect significantly the time distribution of effective
2 10 precipitation, and hence, of runoff generation. High values of a
Time to peak outflow
Peak outflow indicate a large proportion of losses in the beginning of the storm,
Time to peak surface storage
Peak surface storage
and less incremental losses throughout the storm, while low values
CPU time of a indicate that infiltration does not occur when rainfall begins,
Normalized CPU time
Normalized Data
SS
SS
SS
TS
TS
TS
TS
TU
TU -40
TU -60
TU -10
TU R-
TU R-
10
15
20
5
10
15
20
-1
-1
-1
-L 0
-L 166
-L 226
0
62
graphs. This demands quite a lot more from the SCS-CN model as
will be shown.
Fig. 12. Comparison of benchmarking variables and CPU time for all meshes. The results from 40 different a and CN configurations did not
provide satisfactory calibration of the parameters. Neither volume,
method. The first simulation uses SCS-CN with the empirical value hydrograph peak flow or time to peak were adequately reproduced
of a = 0.2. This allows indeed for a calibration of CN which repre- by the simulations. The configurations are shown in Fig. 14. The
sents accurately the losses for event E1. However, the other seven very wide range of configurations is due to the fact that no ade-
simulations use arbitrary values of a ranging between 0.0 and 0.1 quate results were seen with the values selected based on lumped
and calibrating CN to achieve the same runoff coefficient. Note that estimates. The two best results in terms of outflow volume were
calibration is still possible and accurate values of runoff volume are obtained with CN = 52, a = 0.08 and CN = 30, a = 0.01. Note that
obtained in all cases, even with a = 0 which implies Ia = 0. Of the first case (CN = 52) is similar, though not exactly the same
course, the physical meaning of CN becomes quite diffuse when as one combination presented in the lumped estimates, though
54 D. Caviedes-Voullième et al. / Journal of Hydrology 448–449 (2012) 39–59
-5
4.5⋅10
-5 CN 65.4, α=0.2 CN 38.4, α=0.05
4.0⋅10 CN 51.5, α=0.1 CN 19.1, α=0.01
2.5⋅10-5
2.0⋅10-5
-5
1.5⋅10
-5
1.0⋅10
-6
5.0⋅10
0
0.0⋅10
-5
4.5⋅10
CN 15.2, α=0.005 CN 10.9, α=0.0005
Effective precipitation (mm)
-5
4.0⋅10 CN 11.3, α=0.001 CN 10.4, α=0.0000
-5
3.5⋅10
-5
3.0⋅10
-5
2.5⋅10
-5
2.0⋅10
-5
1.5⋅10
1.0⋅10-5
-6
5.0⋅10
0.0⋅100
0 5000 10000 15000 20000 25000 30000 35000 0 5000 10000 15000 20000 25000 30000 35000
80
smaller than in the previous case, around 2500 s, but still the sim-
70 ulated hydrograph is severely different from the experimental
results.
60
From all the simulations shown in Fig. 14, a mean value of
50 surface storage of 1657.9 m3 at the end of the simulation was
obtained, and the values ranged from zero to 2002.0 m3. For the
CN
40
LR-22624 mesh without any infiltration or friction, the surface
30 storage volume at the end of the simulation was 2026.8 m3.
Although this volume might contain moving water, an important
20
part of it is water ponded because of local topographic minima.
10 This information is significant because experimental runoff
volume, as shown in Table 2, is in the order of magnitude of the
0
0 0.01 0.02 0.03 0.04 0.05 0.06 0.07 0.08 0.09 0.1 final surface storage volume, and indeed quite close. The problem
α is that ponded volume is constituted of runoff volume from some
other area of the catchment, and SCS-CN is blind to such effects.
Fig. 14. CN–a combinations for target Vout = 1595 m3. Hence, calibration of a and CN needs to account for the fact
that the volume that is produced as runoff may not show as
outflow volume at the outlet. In other words, some distortion of
generating more effective precipitation (because of larger CN and runoff volume – as understood by SCS-CN – is introduced by the
lower a) than the second case of lumped estimates. The second topography. This can lead to the conclusion that with an effective
configuration (CN = 30) has no similar in the lumped estimates. precipitation small enough to generate a runoff volume smaller
For the CN = 52 case, both hydrograph and volume comparisons than the surface storage volume, no runoff will be observed at
are shown in Fig. 15. Note that the simulated outflow volume is the outlet. This is not necessarily true, as some part of such effec-
around 80% of the experimental outflow volume. However, the tive precipitation falls on the surface downstream of the largest
simulated peak flow is around 30% greater than the measured peak topographic local minima, as in the runoff path that such rain
flow. The rising limb of the hydrograph is severely distorted and describes, there may be no significant storage zones.
the simulated time to peak is dramatically greater. There is a lag To test this hypothesis a synthetic event consisting of an instan-
of the order of 7500 s between the simulated and experimental taneous pulse of precipitation was generated. The rainfall pulse oc-
peak flows. A great deal of the delay is due to a as can be seen in curs in the first instant of the simulation. The volume of the pulse
Fig. 13 where a value of CN = 51.5 and a = 0.1 was used, similar was set equal to the experimental runoff volume. Infiltration was
to the case currently discussed. not simulated at all, setting an impervious catchment with no fric-
Analogous results are presented in Fig. 16. In this case the sim- tion. After 32,100 s the surface storage volume was 1590 m3, and
ulated outflow volume is around 30% greater than the experimen- the runoff volume at the outlet was 2.35 m3. Results are shown
tal volume, and the peak outflow is 130% greater than the in Fig. 17. This synthetic rain case allows to see that, even when
experimental measure. The time lag between peak flows is the correct effective precipitation volume is imposed, the runoff
D. Caviedes-Voullième et al. / Journal of Hydrology 448–449 (2012) 39–59 55
0.20 1600
Simulated outflow
0.18 Experimental outflow
Simulated outlet volume 1400
Experimental outlet volume
0.16
1200
0.14
1000
Volume (m3)
Flow (m3/s) 0.12
0.10 800
0.08
600
0.06
400
0.04
200
0.02
0.00 0
0 5000 10000 15000 20000 25000 30000 35000
Time (s)
0.40 2200
2000
0.35
1800
0.30
1600
0.25 1400
Volume (m3)
Flow (m3/s)
1200
0.20
Simulated outflow 1000
Experimental outflow
Simulated outlet volume
0.15 Experimental outlet volume 800
600
0.10
400
0.05
200
0.00 0
0 5000 10000 15000 20000 25000 30000 35000
Time (s)
volume is significantly less than expected from experimental data To support this conclusion, an additional synthetic pulse event
or SCS-CN lumped estimates calibrated directly from experimental was simulated, also instantaneous, but designed to generate
data. It is clear that, in order to calibrate SCS-CN properly for its use 3200 m3. Surface volume results are shown in Fig. 18, where it
with the shallow-water simulator, the resulting effective precipita- can be seen that, once the rain pulse flows out of the catchment,
tion must be greater than the one computed from experimental the instantaneous surface storage decreases rapidly to an asymp-
data in order to account for ponded water. Clearly, this numerical totic value of around 1800 m3, 200 m3 more than expected from
issue can have a wide range of magnitudes, depending both on the previous discussion. On the other hand, although Fig. 18 shows
infiltration, rain volume and distribution as well as topography. that the simulated and experimental hydrographs are quite differ-
Hence, how much effective precipitation should be obtained from ent – which is reasonable since the hyetograph is drastically differ-
SCS-CN is difficult to know without actually simulating with the ent –, the simulated outflow volume is around 165 m3 less than the
shallow-water model at least to know how much runoff may be experimental value, which represents an error of around 10% less
lost to surface storage. Given that surface storage for this catch- volume. Both results are consistent and support the conclusions
ment is practically equal to the experimental runoff volume, the of the previous discussion that the SCS-CN effective precipitation
calibration of SCS-CN should produce an effective precipitation of volume must also account for ponded volume which cannot be ob-
around twice that volume, this is approximately 3200 m3. served in the outflow hydrograph.
56 D. Caviedes-Voullième et al. / Journal of Hydrology 448–449 (2012) 39–59
1594.0
Surface volume Because the synthetic pulse confirmed the hypothesis, addi-
1593.5
Total volume tional tests were performed to generate a rain of approximately
3200 m3, using the hyetograph for event E1 and calibrating CN
1593.0 and a for such volume.
A configuration with CN = 20.1 and a = 0 which generates a vol-
Volume (m3)
1592.5
ume in the aforementioned range was used to test the hypothesis
1592.0 together with the hyetograph. The results are shown in Fig. 19. In
terms of volume, the results are fairly satisfactory. The simulated
1591.5 volume is 3.5% smaller than the measured experimental volume.
1591.0
The simulated hydrograph still produces a larger peak flow which
is almost twice the experimental peak flow and has a delay in the
1590.5 order of 2000 s. The hydrograph is not simulated accurately.
Although a is zero, a large delay still exists which appears to be
1590.0
0 5000 10000 15000 20000 25000 30000 35000 unsolvable with the chosen methodology.
Time (s) Another result worth of mention is the one obtained with
CN = 44.4 and a = 0.05, and shown in Fig. 20. In this case the sim-
Fig. 17. Volume results for synthetic rain pulse simulation. ulated volume is around 15% less than the experimental volume.
The simulated peak flow is around 1.5 times the experimental peak
1.40 3200
Simulated outflow
Experimental outflow 3000
Simulated outlet volume
Experimental outlet volume 2800
1.20
Surface volume
2600
2400
1.00
2200
2000
Volume (m3)
Flow (m3/s)
0.80 1800
1600
0.60 1400
1200
1000
0.40
800
600
0.20
400
200
0.00 0
0 5000 10000 15000 20000 25000 30000 35000
Time (s)
0.30 1600
1400
0.25
1200
0.20
Simulated outflow 1000
Volume (m3)
Flow (m3/s)
Experimental outflow
Simulated outlet volume
Experimental outlet volume
0.15 800
600
0.10
400
0.05
200
0.00 0
0 5000 10000 15000 20000 25000 30000 35000
Time (s)
0.25 1600
Simulated outflow
Experimental outflow
Simulated outlet volume
Experimental outlet volume 1400
0.20
1200
1000
Volume (m3)
0.15
Flow (m3/s)
800
0.10
600
400
0.05
200
0.00 0
0 5000 10000 15000 20000 25000 30000 35000
Time (s)
flow and a delay in the order of 2500 s. Further simulations with chosen methodology was simply to perform a range of simulations
larger values of a resulted in larger errors both in volume and hyd- with changing parameters, without the presence of other parame-
rograph results. For CN = 68.3 and a = 0.2 the outflow volume was ters allowing to value them in their own right as well as to observe
2.57 m3 and the peak outflow was in the order of 0.0005 m3/s their interactions.
although the expected volume computed from lumped estimates The results show that mesh properties are one of the most sig-
was 3258 m3 and the simulated surface storage volume at the nificant factors. Square and right-angled triangular structured
end of the simulation was 1899 m3. meshes with different resolutions were tested, as well as homoge-
The evidence gathered from the wide range of simulations per- nous unstructured triangular meshes and locally refined unstruc-
formed leads to several conclusions concerning the use of SCS-CN tured triangular meshes. Resolution and cell size are of great
together with distributed surface flow modeling by means of shal- importance. Square structured meshes are the most numerically
low-water simulation. Firstly, SCS-CN method may generate an viscous and show flow directionality aligned with the principal
effective hyetograph which severely differs from the registered direction of the mesh. Triangular structured meshes lessen these
hyetograph. The effect observed in simulated runoff is of signifi- effects. However, high resolution structured meshes require large
cant delays in peak flow as well as hydrographs which do not computational efforts, and coarsening of structured meshes results
resemble the experimental hydrographs. Second, calibration of in a great loss of topographic representation quality, and because
SCS-CN must include not only CN but also a since the effective hye- of this, a significant loss in flow result quality. Unstructured trian-
tograph may be dramatically different with different values of a. gular meshes respond in complex ways, with mixed results, but
The third conclusion is that SCS-CN calibration cannot be per- computational time can be similar to that of structured meshes.
formed only based on the outled hydrograph and the outflow vol- Locally refined meshes by an automatic algorithm which adds cells
ume obtained from such hydrograph. Surface storage of ponded to represent correctly changes in terrain slope (slope gradient cri-
water must be considered as well. Hence, SCS-CN must be cali- teria) outperform most other meshes. By adding more and smaller
brated to generate an effective precipitation volume which is the cells were terrain complexity is high (high slope gradient) an
sum of the outlet volume and the ponded water in the catchment. accurate representation of the terrain surface is achieved, while
This volumes of surface storage may be previously obtained by keeping a low cell count, and thus, reducing the required computa-
simulating an impervious catchment and observing the volume tional effort. High precision is also achieved with locally refined
that remains ponded in the catchment after the storm. In lumped meshes.
methods, which make use of synthetic unit hydrographs, this effect Friction representation by means of Manning’s roughness is
may be masked by calibration of the unit hydrograph parameters perhaps the most commonly used representation. In this study,
but, in distributed, physically-based models, no calibration free- the effects of friction on catchment flow is studied by comparing
dom is allowed, since topography and surface flow are modeled frictionless results against results with friction. Furthermore, two
physically and explicitly. Furthermore, in this particular case there roughness representations were studied: a roughness map based
appears to be no possible calibration with SCS-CN that leads to an on soil and vegetation properties and an equivalent average rough-
adequate hydrograph or peak flow. It was only possible to calibrate ness value for the entire catchment. Both representation yielded
runoff volume. It appears that the use of shallow-water solvers for practically identical results since the distribution of roughness val-
rain-runoff transformation requires a better simulation of tempo- ues is heavily weighted near the average value. High frequency
ral distribution of infiltration and consequently effective oscillations were observed when using the roughness map, but
precipitation. their amplitude was small enough to be negligible. When compar-
ing frictionless results to results with friction the expected behav-
13. Conclusions ior is seen. Friction dampens the outflow hydrograph, delays peak
flow and smoothes oscillations. In a sense, by introducing friction,
This study was intended to understand the possible effects that mesh dependant oscillations are eliminated. Although it would ap-
different model parameter can have on simulation results. The pear that by doing so mesh selection would be less important, this
58 D. Caviedes-Voullième et al. / Journal of Hydrology 448–449 (2012) 39–59
is not so. Mesh coarsening proved to have a similar effect on the Acknowledgement
outflow hydrograph as friction. As the mesh becomes coarser, the
hydrograph is somewhat dampened and peak flow is delayed. In This work was partially supported and funded by the Spanish
fact, the delays can be larger than those generated by physical fric- Ministry of Science and Technology under research project
tion. It becomes easy to mix the effects of physical friction and CGL2008-05153-C02-02.
numerical, mesh generated friction if a frictionless analysis is not
performed. This also implies that roughness coefficients must be
calibrated with every particular mesh that is used since each mesh References
introduces its own amount of numerical friction. But it also implies
Baltas, E.A., Dervos, N.A., Mimikou, M.A., 2007. Technical note: determination of the
that a coarse mesh may be inappropriate because, on its own right
SCS initial abstraction ratio in an experimental watershed in Greece. Hydrology
without friction, it can dampen and delay a hydrograph beyond the and Earth System Sciences 11, 1825–1829.
effects of physical friction, making roughness calibration impossi- Beck, H.E., de Jeu, R.A.M., Schellekens, J., van Dijk, A.I.J.M., Bruijnzeel, L.A., 2009.
Improving Curve Number based storm runoff estimates using soil moisture
ble. On a side note, roughness discretization was performed only
proxies. IEEE Journal of Selected Topics in Applied Earth Observations and
by mapping high resolution raster data to different meshes by Remote Sensing 2, 250–259.
arithmetic averaging. It appears that, because of roughness distri- Burguete, J., Garcia-Navarro, P., Murillo, J., Garcia-Palacin, I., 2007. Analysis of the
bution in this particular case no significant error is introduced. friction term in the one-dimensional shallow-water model. Journal of Hydraulic
Engineering-ASCE 133, 1048–1063.
However, where important variations in roughness are present, Casas, A., Benito, G., Thorndycraft, V., Rico, M., 2006. The topographic data source of
more complex averaging techniques, or even considering rough- digital terrain models as a key element in the accuracy of hydraulic flood
ness variations as criteria for mesh generation (in a fashion similar modelling. Earth Surface Processes and Landforms 31, 444–456.
Cook, A., Merwade, V., 2009. Effect of topographic data, geometric configuration and
to that of slope gradient) could result in better solutions. modeling approach on flood inundation mapping. Journal of Hydrology 377,
For infiltration computations the SCS-CN method was selected 131–142.
in this study, to test its response when coupled with a physically Garcia-Ruiz, J., Arnaez, J., Begueria, S., Seeger, M., Marti-Bono, C., Regues, D., Lana-
Renault, N., White, S., 2005. Runoff generation in an intensively disturbed,
based, distributed model for surface flow. The results were poor. abandoned farmland catchment, Central Spanish Pyrenees. Catena 59, 79–92.
It was not possible to calibrate the SCS-CN parameters (CN and a) Hardy, R.J., Bates, P.D., Anderson, M.G., 1999. The importance of spatial resolution in
to achieve an accurate representation of the hydrograph. In many hydraulic models for floodplain environments. Journal of Hydrology 216, 124–
136.
cases runoff volume was not achieved, and when runoff volume Horritt, M.S., Bates, P.D., 2001. Effects of spatial resolution on a raster based model
was achieved, the hydrograph shape was poorly reproduced and of flood flow. Journal of Hydrology 253, 239–249.
peak flow was delayed far beyond the experimental results. In Horritt, M.S., Bates, P.D., Mattinson, M.J., 2006. Effects of mesh resolution and
topographic representation in 2D finite volume models of shallow water fluvial
any case, it appears necessary to calibrate a as well as CN since
flow. Journal of Hydrology 329, 306–314.
it was important effects on time distribution of infiltration. Fur- Hunter, N.M., Bates, P.D., Horritt, M.S., Wilson, M.D., 2007. Simple spatially-
thermore, when calibrating SCS-CN parameters using catchment distributed models for predicting flood inundation: a review. Geomorphology
outflow hydrographs, care must be taken to remember that the 90, 208–225.
Lana-Renault, N., 2007. Respuesta hidrológica y sedimentológica en una cuenca de
measured outflow volume is not equal to non-infiltrated volume. montaña media afectada por cambios de cubierta vegetal: la cuenca
Some of the runoff generated in a particular position of the catch- experimental de Arnás, Pirineo Central. Ph.D. thesis, Universidad de Zaragoza.
ment may become ponded elsewhere, thus never reaching the Lana-Renault, N., Latron, J., Reguees, D., 2007. Streamflow response and water-table
dynamics in a sub-Mediterranean research catchment (Central Pyrenees).
outlet. Hence, if SCS-CN is calibrating based only on outflow vol- Journal of Hydrology 347, 497–507.
ume from experimental outflow hydrographs, the simulated re- Leveque, R., 2002. Finite Volume Methods for Hyperbolic Problems. Cambridge
sults will always yield less outflow volume since some runoff University Press.
López-Barrera, D., García-Navarro, P., Brufau, P., 2011. Sources of uncertainty in the
will be lost to ponding. This has a simple solution, which is to validation of a coupled hydrological-hydraulic simulation model with sediment
measure the available ponding volume in the catchment by flood- transport. La Houille Blanche.
ing the catchment with the same storm (precipitation volume) López-Barrera, D., 2011. Desarrollo de un modelo hidrológico/hidrodinámico de
simulación de procesos de erosión y sedimentación. Ph.D. thesis, Universidad de
and allowing enough simulation time to only have static water Zaragoza.
in the catchment. In this way, ponding volume is known and Merwade, V., 2009. Effect of spatial trends on interpolation of river bathymetry.
can be added to outflow volume to calibrate SCS-CN parameters. Journal of Hydrology 371, 169–181.
Mishra, S., Singh, V., 1999. Another look at SCS-CN method. Journal of Hydrological
In this way SCS-CN will compute runoff volume, understood as
Engineering 4, 257–264.
non-infiltrated volume, which must be equal to outflow volume Mishra, S., Singh, V., Sansalone, J., Aravamuthan, V., 2003. A modified SCS-CN
and ponded volume. This behavior may become less important method: characterization and testing. Water Resources Management 17, 37–68.
as catchment scale changes. In this particular case, ponding vol- Murillo, J., Garcia-Navarro, P., 2004. A new topography adapted mesh refinement for
2D unsteady shallow water flow simulations over dry beds. In: Liong, S.-Y.,
ume is practically the same as outflow volume and errors grow Phoon, K.-K., Babovic, V. (Eds.), Proceedings of the 6th International Conference
rapidly and are very significant. In larger catchments, with less on Hydroinformatics, Singapore.
infiltration, or higher precipitated volumes and lower ponding Murillo, J., Garcia-Navarro, P., 2010. Weak solutions for partial differential equations
with source terms: application to the shallow water equations. Journal of
volumes it may not be as significant. Still, it is easy to consider. Computational Physics 229, 4327–4368.
However, the results show that even when considering this ef- Murillo, J., Garcia-Navarro, P., Burguete, J., 2009. Time step restrictions for well-
fects, resulting in correct outflow volumes, hydrographs are balanced shallow water solutions in non-zero velocity steady states.
International Journal for Numerical Methods in Fluids 60, 1351–1377.
poorly simulated. In this case, SCS-CN is not appropriate to couple Murillo, J., Garcia-Navarro, P., Burguete, J., Brufau, R., 2007. The influence of source
with a distributed model for runoff computations. In a lumped terms on stability, accuracy and conservation in two-dimensional shallow flow
approach, SCS-CN may be coupled with synthetic unit hydrograph simulation using triangular finite volumes. International Journal for Numerical
Methods in Fluids 54, 543–590.
techniques, and SCS-CN may still be valid, since it can be accu- Ponce, V.M., Hawkins, R.H., 1996. Runoff Curve Number: Has it reached maturity?
rately calibrated to simulate outflow volume, and the synthetic Journal of Hydrologic Engineering 1, 11–19.
unit hydrograph calibration parameters may have enough free- Schubert, J.E., Sanders, B.F., Smith, M.J., Wright, N.G., 2008. Unstructured mesh
generation and landcover-based resistance for hydrodynamic modeling of
dom to generate a hydrograph similar to the field measurements.
urban flooding. Advances in Water Resources 31, 1603–1621.
With a physically-based distributed model, such freedom is lost, SCS, 1964. SCS National Engineering Handbook. Soil Conservation Service, USA
and the infiltration method must accurately compute results in Department of Agriculture.
time, not only lumped volumes. Because SCS-CN appears to be Seeger, M., Errea, M., Beguería, S., Martí, C., García-Ruiz, J., Arnáez, J., 2004.
Catchment soil moisture and rainfall characteristics as determinant factors for
a poor choice, in future work, Horton and Green-Ampt methods discharge/suspended sediment hysteretic loops in a small headwater
will be tested. catchment in the Spanish Pyrenees. Journal of Hydrology 288, 299–311.
D. Caviedes-Voullième et al. / Journal of Hydrology 448–449 (2012) 39–59 59
Serrano-Pacheco, A., 2009. Simulación numérica bidimensional de procesos Wilson, M., Atkinson, P., 2005. Prediction uncertainty in floodplain elevation and its
hidrológicos e hidráulicos sobre lecho irregular deformable. Ph.D. thesis, effect on flood inundation modelling. In: Atkinson, P.M., Foody, G.M., Darby, S.E.,
Universidad de Zaragoza. Wu, F. (Eds.) Proceedings of the 7th International Conference on
Shewchuk, J.R., 1996. Triangle: engineering a 2D quality mesh generator and GeoComputation, Southampton, England, 2003, pp. 185–202.
Delaunay triangulator. In: Lin, M.C., Manocha, D. (Eds.), Applied Computational Woodward, D., Hawkings, R., Hjelmfelt, A., Van Mullen, J., Quan., Q., 2002. Curve
Geometry: Towards Geometric Engineering. Springer-Verlag, vol. 1148 of Number method: origins, applications and limitations. In: Second Federal
Lecture Notes in Computer Science, pp. 203–222. From the First ACM Interagency Hydrologic Modeling Conference, Las Vegas, Nevada.
Workshop on Applied Computational Geometry.
Shi, Z.H., Chen, L.D., Fang, N.F., Qin, D.F., Cai, C.F., 2009. Research on the SCS-CN
initial abstraction ratio using rainfall-runoff event analysis in the three gorges
area, china. Catena 77, 1–7.