Journal o f Hydrology, 145 (1993) 357-370
357
Elsevier Science Publishers B.V., A m s t e r d a m
[31
An application of a physically based semidistributed model to the Balquhidder catchments A.J. Robson a, P.G. Whiteheada and R.C. Johnson b ~Institute of Hydrology, Walling)CordOXIO 8BB, UK blnstitute of Hydrology, Balquhidder, Lochearnhead FKI9 8PQ, UK (Received 23 August 1991; revision accepted 9 November 1992)
ABSTRACT Robson, A.J., Whitehead, P.G. and Johnson, R.C., 1993. An application of a physically based semi-distributed model to the Balquhidder catchments. J. Hydrol., 145: 357-370. A physically based semi-distributed model, TOPMODEL, is applied to the two catchments at Balquhidder. The model uses a topographic index which highlights hydrologically significant areas within the catchments. The model is used to simulate runoff and to compare the behaviour of the two catchments. The results indicate that a large proportion of stream water is generated from saturated contributing areas (source areas); the Monachyle has higher contributions of water of this type. The results may also indicate that the hydrology of the Monachyle catchment has changed following agricultural improvement by increased drainage; a damped subsurface response is observed.
INTRODUCTION
•
T O P M O D E L (Beven and Kirkby, 1979; Beven and Wood, 1983; Beven et al., 1984) is a physically based semi-distributed model developed for use in predicting and subsequently understanding rainfall-runoff mechanisms. It provides a compromise between the complexity of fully distributed models (e.g. I H D M (Calver, 1988)) and the relative crudity of lumped models (Birkenes (Christophersen et al., 1982), Pulse (Bergstrom et al., 1985) and Hyrom (Blackie and Eeles, 1985)). Heterogeneity in catchment topography is incorporated into T O P M O D E L by means of a topographical index, and the movement of water through the catchment is founded on a simple representation of physical processes. In this paper both the semi-distributed and the physical aspects of TOPMODEL are exploited. T O P M O D E L is applied to the Monachyle and Kirkton catchments at Balquhidder, described by Johnson and Whitehead C o r r e s p o n d e n c e to: A.J. R o b s o n , Institute o f Hydrology, Wallingford OX10 8BB, U K .
0022-1694/93[$06.00
© 1993 - - Elsevier Science Publishers B.V. All rights reserved
358
A.J. ROBSON E"I- AL.
(1993), with the aims of (1) comparing the hydrological response of the two catchments, (2) relating these responses to the flow generation mechanisms operating in the catchments and (3) consideration of the effects that land-use changes may have brought about. OUTLINE THEORY OF TOPMODEL
A brief summary of the concepts and the basic equations upon which T O P M O D E L is founded are given here. Full details of the derivation of the equations and the rationale behind the model can be found elsewhere (e.g. Beven and Kirkby, 1979; Beven, 1986). T O P M O D E L represents catchment topography by means of a topographic index, In (a/tan/3), where a is the area draining through a grid square per unit length of contour and tan/3 is the average outflow gradient from the square. The index is calculated from a digital terrain map (DTM) across a grid covering the catchment. The grid must be sufficiently fine to resolve important characteristic slope formations. A high index value usually indicates a wet part of the catchment; this can arise either from a large contributing drainage area (valley bottoms or convergent hollows) or from very flat slopes (bogs on hill tops). Areas with low index values are usually drier, resulting from either steep slopes or a small contributing drainage area. Grid squares with the same index value are assumed to behave in a hydrologically similar manner. As a result of this assumption, the catchment's topography may be summarised by the distribution of the index values. This version of T O P M O D E L identifies two sources of stream water: water draining from subsurface saturated zones and water displaced from saturated and near-saturated parts of the catchment. The second of these is generated by rainfall landing on saturated contributing areas and causing rapid movement to the stream via macropore flow, overland flow or old-water displacement. It is distinguished from the subsurface saturated zone discharge because of the much faster speed of response. The saturated contributing area will both grow and decline during the course of a storm event. For any point on the hillslope the downhill flow from the saturated zone is assumed to decrease exponentially with depth to the water table and to be proportional to the local gradient: qsi
=
To
exp C-:~) tan/3i
where qs~ is the local lateral saturated flow per unit length of contour (m 2 h - t ), z~ is the local depth to the water table (m),/3 i is the local slope angle, To is the lateral transmissivity when the soil is saturated to the surface (m 2 h -1 ) and f describes the exponential decrease in transmissivity with depth (m ~).
MODELLING OF BALQUHIDDER CATCHMENTS
359
Contributions are summed over the catchment to give a total flow per unit area from the saturated zone: qs
To exp (-~I-~)
=
where 2 is the areal average of the index In (ai/tan[3i): A
2
1/A l-ln (ai/tan//i) di
=
0
5 is the mean catchment depth to the water table and A is the total catchment area. The local depth to the water table and the mean catchment depth to the water table are linked: zi
=
5 +
[2 -
ln(ai/tan
[3~)]/f
This is based upon the assumption of there being quasi-steady flow within the catchment. Knowledge of local depth to the water table allows determination of areas where the water table is at the surface (i.e. the saturated contributing area, SCA), and this is used to calculate qscA, the average flow per unit area of the catchment generated by rainfall on the saturated contributing areas. Hence at a given time t the total flow per unit area is qs, + qSCA
q~ =
At each time step the value of z is updated ready for use in the next interval: 5,+1
=
2, + (qs, -
qvt)/AO
where AO is the storage capacity of the soil as a proportion of total soil volume and qv, is the total vertical flow through the unsaturated zone down to the saturated zone; qv t is calculated by summing the local values of vertical drainage q v , . These are derived by assuming that the hydraulic conductivity has an exponential profile with depth (with the same exponential decay parameter f as for lateral flow) and that near the water table there is a unit hydraulic gradient. Local vertical flow is then given by qv,,
=
K o exp (-z'y)
where K0 is the vertical conductivity at the surface (mh-I). This allows for increasing mean residence time in the unsaturated zone with increasing depth to the water table. D I G I T A L T E R R A I N MAP A N D T O P O G R A P H I C I N D EX D I S T R I B U T I O N
Ln(a/tan/~) maps provide information which can be used to characterise
360
A.J. ROBSON ET AL.
catchment hydrological and hydrochemical behaviour (Wolock et al., 1990; Wolock et al., 1989). The maps can be used to help identify source areas; these are potentially important in the control of the chemical characteristics of the stream (e.g. the organic compounds in the Monachyle) and in sediment transport. The maps are also valuable in assessing the likely impacts of land-use changes (afforestation or drainage) and of herbicide or liming applications (Waters et al., 1991). In addition, the maps can be used in combination with soil data, if available; Beven (1986) has shown that spatial variation in the transmissivity, To, can be taken into account in a combined topographic-soil index, ln(a/Totanfl). A D T M was produced for assessment of the ln(a/tanfl) index. This map was generated at the 100m x 100m scale by interpolation from the contours of a 1:10000 Ordnance Survey map. Index values were calculated using a multidirectional routing approach as described by Quinn et al. (1991). The area, a, used in this index represents the area of catchment which drains through a given grid square. These areas are accumulated down the catchment until they reach the outlet. The algorithm used to calculate this drainage area allows a weighted proportion of the area to be passed on downhill. The weighting is based on the slope and the direction of the downhill squares (immediately adjacent squares receive a higher weighting than diagonally adjacent squares). The slope, tan/3, is an averaged value of the outflow slopes; it too is weighted for direction. Maps of the index values show that the catchments are different in terms of the spatial distribution of ln(a/tan/3) (Figs. 1 and 2). The index varies smoothly from east to west across the Kirkton catchment, increasing almost monotonically with height. High index values are concentrated along a continuous band near the stream; this band widens towards the catchment outlet. Two additional source areas are also visible in the higher parts of the catchment; one to the northeast and the other to the west near two small lochans. The spatial distribution of the ln(a/tanfl) index for the Monachyle catchment is much more irregular - - a direct result of the uneven topography and the presence of glacial deposits. A large proportion of high index values occur in the fiat upper Monachyle sub-catchment. Further areas of high In(a/tan~3) occur alongside the lower stream banks of the Monachyle: this band is narrower and more irregular than that of the Kirkton. In general, the index maps correspond well to the impressions of catchment wetness gained through field experience. For example, both catchments contain dry rocky areas on parts of their upper slopes; these agree well with the low In(a/tan~3) values shown on the maps. Similarly, many of the wetter parts of the upper slopes (near the small lochans, for example) show higher index values. The predominance of high index values in the upper
MODELLING OF BALQUHIDDER CATCHMENTS
361
LN (A/TAN e &BOrg
:
r--I
re.o-: ~.o-
p.
u-
g,l)7,0-
mR
5.0-
ms m.ow
6.0 Lo
Fig. 1. T w o - d i m e n s i o n a l
ln(a/tanfl)
index m a p for the K i r k t o n c a t c h m e n t .
Monachyle catchment ties in well with the many pockets of wet peat seen in this area. The cumulative distributions of the ln(a/tan/~) indices are compared in Fig. 3. The distributions are similar, although the Monachyle catchment has a slightly bigger proportion of the higher In(a/tan//) values. This difference means there is likely to be more flow generated in the Monachyle from the saturated contributing areas (source areas). CALIBRATION AND RESULTS T O P M O D E L was calibrated for both the Monachyle and Kirkton catchments using daily rainfall and flow records for the year June 1984-June 1985. For a daily time step the channel flow time is negligible. Three parameters were unknown and were allowed to vary in each of the optimisations. These were f, the exponential decay rate in transmissivity, K0, the
362
A.J. ROBSON ET AL.
le at Cr
LN (A/TAN E~
ABOVE t~c 19.0 11,0 I0.09.0 8.07.0~,05,0BELOW
,8) 14.0 t4~o 1~.0 IP-,.O
ILO 10.0 9.0 8.0 7.0 6.Q 5.0
Fig. 2. Two-dimensional ln(a/tanfl) index map for the Monachyle catchment.
surface conductivity, and To, the transmissivity at saturation to the surface. Optimisation was carried out using the Nelder-Mead Simplex procedure (Nash, 1979) to minimise the least-squares error between observed and modelled flow; optimised parameter values are shown in Table 1. Figures 4(a) and 4(b) and 5(a) and 5(b) show the rainfall, and the simulated and observed flow for the two catchments. In Figs. 4(c) and 5(c) the modelled flow is subdivided into subsurface flow and saturated overland flow. Both catchments are very flashy, and a relatively high proportion of the total response is made up of water from saturated contributing areas. The Kirkton catchment gives the highest proportion of subsurface flow; this catchment also gives the betteF fit (Kirkton: R 2 =- 0.83; Monachyle: R 2 = 0.78). The higher saturated overland flow for the Monachyle corresponds to the extensive peat areas in the upper parts of the catchment and the higher observed rainfall. A validation run, using the estimated parameters for 1984-1985, was
363
MODELLING OF BALQUH1DDER CATCHMENTS
1.0
K,r onS
0.8 .Q
0.6 Q. >
0.4 E "-i o
0.2
0
i
0
i
;~
i
4
i
i
6
i
i
8
i
i
i
10
i
12
i
i
14
i
i
16
i
18
In(a/tan B) value
Fig. 3. Cumulative distribution for the ln(a/tan/~) index of the Monachyle and Kirkton catchments.
undertaken for both catchments between June 1987 and June 1988 (Figs. 6 and 7). The model gives a satisfactory fit, although during winter 1987-1988 both catchments show a period where modelled and observed results agree poorly because of the heavy snow falls. A calibration run was carried out for TABLE 1 Estimated TOPMODEL parameters 1984-1985
1987-1988
Calibration
Validation
0.97 8.47 5.5 67 0.83
0.97 8.47 5.5 66 0.80
1.08 6.43 4.96 71 0.80
To K0 f
0.89 80.6 12.6
0.89 80.6 12.6
0.42 387.0 4.86
% Subsurface flow R2
50 0.78
47 0.64
51 0.68
Calibration
Kirkton To Ko f % Subsurface flow R2
Monachyle
364
A.J. R O B S O N
E T AL.
60-
(a)
40
200
i
m
m
"l
I
0
,LII]J i
m
i
~
50
.[
I
I
i
i
I
i
,
,
'
I
150
100
'
'
,
i
I
200
i
i
,
i
I
250
'
;
I
i
I
300
i
350
"~ 50' (b)
~ 25' 0
I
i
J
i
i'
50
.~ 50°
I
i
i
i"i
I
100
i
i
,
i
I
150
i
i
i
i
200
I
i
i
,
I
I
250
'
i
r
i
I
300
i
350
I1
(c)
[25 o J
0
I
I
I
I
50
I
I
i
I
l
100
I
I
J
I
I
I
t
I
I
l
150 uy~a-s 200
i
I
I
i
l
250
i
I
i
i
l
300
i
,
i
'
l
J
350
Fig. 4. Kirkton catchment 1984-1985 (calibration). (a) Rainfall; (b) observed (full line) and modelled (broken line) flow; (c) components of simulated flow--saturated contributing area flow (black) and subsurface flow (grey).
the 1987-1988 period to help assess whether any change in response has occurred between the two periods (Table 1; Fig. 8). Little improvement in fit between the validation and calibration runs is seen for the Kirkton catchment. The changes in the parameters are relatively unimportant and there is no evidence to suggest that any change in flow generation mechanisms or flow routing has occurred for this catchment. Some change in K0 occurs between the validation and calibration runs for this period; however, this parameter was found to be insensitive in the model, giving similar results over a wide range of values. The values of K0 are, in fact, rather higher than field measurements at other upland sites, and so may not be physically realistic; this is perhaps a consequence of the poor sensitivity of this parameter. For the Monachyle catchment, calibration for the 1987-1988 period results in an improved fit relative to the validation run, and there are noticable differences in the optimised parameter values. To illustrate these differences
MODELLING OF BALQUHIDDER CATCHMENTS
365
8060-
(a)
= 40-
2o0
I
I
1
I
[
I
I
50
0
I
I
I
I
I
100
I
I
I
I
I
150
I
I
I
I
i
I
300
25O
200
350
v
(b) 25
0
~
"~ 50 °- -
(c)
I
I
I
I
'•
i
I
/"1
I
I
I
I
I
50
100
150
200
250
50
100
150 xa y n a - s 200
251
I
I
I
I
I
I
I
I
I
'
350
300
25"
8 "-
0
I
0
I
I
I
I
300
I
I
I
I
I
I
350
Fig. 5. Monachyle catchment 1984-1985 (calibration). (a) Rainfall; (b) observed (full line) and modelled (broken line) flow; (c) components of simulated flow--saturated contributing area flow (black) and subsurface flow (grey).
the transmissivity profiles for both Monachyle and Kirkton are plotted in Fig. 9. DISCUSSION OF RESULTS
Physical interpretation of changes in fitted parameters are valid only if the inherent assumptions contained within the equations hold for the catchment in question. For example, it has been assumed that the transmissivity decays exponentially with depth and that quasi-steady conditions apply in the unsaturated zone. Neither of these assumptions has been explicitly tested for the Balquhidder catchments. Studies at other sites suggest that the assumption of exponential transmissivity is reasonable (Beven, 1982; Chappell et al., 1990). The second assumption may not hold where soil moisture content varies rapidly with depth; this may be a problem when water tables are shallow. In
366
A.J. ROBSON ET AL.
8060-
(a)
v 40200
'
I
0
i
,
i
50
i00
150
200
250
800
350
50
100
150
200
250
300
350
"~ 5O (b)
25 0 ,~,
,
50 °
(c)
~25. g 0
~u 0
I
0
i
I
I
|
50
I
I
t
[
[
100
I
i
l
I
~
I
I
I
I
t
150 ~yr'a--s 200
I
I
I
[
|
250
t
I
I
i
|
300
r
I
i
[
|
F
350
Fig. 6. Kirkton catchment 1987-1988 (validation). (a) Rainfall; (b) observed (full line) and modelled (broken line) flow; (c) components of simulated flow--saturated contributing area flow (black) and subsurface flow (grey).
addition, the interpretation relies on an assumption that the parameters are uniquely identifiable and without this the physical links can only be tentative. Nevertheless, given the limitations of these assumptions, the transmissivity profile may be interpreted as an indication of how quickly water moves laterally at different depths. For this application it is not expected that this profile is representative of any particular soil profile; this would require much more detailed data and division of the catchment into homogeneous subcatchments. Instead, the transmissivity profile represents the net effective transmissivity of the heterogeneous soils and soil structures (e.g. macropores, micropores, flushes and drains) contained within the catchment. For 1984-1985 the transmissivity profiles (Fig. 9) suggest that there is a net fast movement near the surface in both catchments - - but with water moving comparatively slowly at depth in the Monachyle. There is little difference for the Kirkton catchment profiles between 1984-1985 and 1987-1988, but for the
MODELLING
OF
BALQUHIDDER
367
CATCHMENTS
8060-
(a)
,~ 4 0 200
I
I
I
I
I
I
1
I
I
I
I
i
i
i
I
i
i
i
i
I
0
50
100
150
200
50_0
50
loo
150
200
i
I
i
I
I
I
I
I
I
250
I
1
I
i
I
300
I
I
350
"~ 50
(b)
0
I
I
I
I
I
I
I
I
I
|
I
I
I
I
300
250
I
I
350
~e
(c)
[25
0
~ o
I
I
I
I
I
50
I
I
I
I
I
100
I
I
I
I
I
I
I
I
I
[
150 oyna-s 200
I
I
I
I
I
250
I
I
I
I
I
300
I
I
I
I
I
I
350
Fig. 7. Monachyle catchment 1987-1988 (validation). (a) Rainfall; (b) observed (full line) and modelled (broken line) flow; (c) components of simulated flow--saturated contributing area flow (black) and subsurface flow (grey).
Monachyle in 1987-1988 there is a substantial reduction in the subsurface flow just beneath the surface and a corresponding relative increase in the subsurface flow at depth. The lack of change in the Kirkton may be an indication that the clear-felling, which involved minimal disturbance to the catchment, has not had a noticeable effect on catchment drainage. It should be noted that no change in the evapotranspiration, and thus water balance, has been observed, for the catchment as a whole, since the clear-felling. The changes in the Monachyle need not represent a physical change to the soils but could instead result from structural modifications to the drainage system affecting both flow movement and water storage. These may be a result of the land-use changes in the Monachyle catchment, in which 6% of the catchment was ploughed. Although 6% is a fairly small area to produce an observable effect, it is likely that these areas have an enhanced hydrological effect because they form a broken band across the lower slopes near the stream. The
368
A.J.
ROBSON
ET
AL.
60-
(a) "~4 0 2o-
e -],L'"' .... l
i
0
I
' i
i
i
i
|
50
i
i
i
i
I
100
,
i
i
i
i
|
150
i
i
. J.l.
i
i
I
200
i
i,,
i
i
i
I
250
,i..iv i
i
J..,
i
i
|
300
i
350
"~ 5O s
(b)
25 I
I
I
I
|
I
I
I
50
60 o
(c)
I
I
I
I
I
I
100
|
I
I
I
l
J
'
i
,
I
200
150
LL/d i
'
I
I
250
f
I
J
i
I
I
350
300
~25 ta
a, 0
f
I
0
I
I
|
50
I
I
I
I
|
100
I
I
I
f
|
I
I
I
I
|
150 JJYr'a-s 200
I
1
I
I
|
250
I
I
I
I
|
300
1
I
I
I
i
I
350
Fig, 8, Monachyle catchment 1987-1988 (calibration). (a) Rainfall; (b) observed (full line) and modelled (broken line) flow; (c) components of simulated flow--saturated contributing area flow (black) and subsurface flow (grey).
ploughed areas have medium to high index ln(a/tanfl) values (around 8-10), and the total contributing area draining through these ploughed parts is significantly higher than 6% of the catchment area. Despite the altered transmissivity profile between the two periods, the total subsurface flow from the Monachyle catchment has apparently not changed as a proportion of the total flow (Table 1). Instead, the pattern of subsurface flow emergence has changed; the modelled subsurface response is more damped in 1987-1988 than in 1984-1985. CONCLUSIONS
The results suggest the following: (1) The Monachyle has the larger proportion of area where generation of
MODELLING
OF
BALQUHIDDER
369
CATCHMENTS
1.2 .
1.0- "\
11984] 85 I 1987/88 .....
~'\
¢-
"~ 0.8 ~~, ,,~
~e
Kirkton
I
0.6 E
0.4
I.-
0.2
0
0
|
012
|
!
|
!
0.4 0.6 Depth (m)
|
i
0.8
! i
1.0
Fig. 9. Transmissivity profiles for the Kirkton and Monachyle catchments.
flow through saturation is likely; this proportion is particularly large in the upper reaches of the catchment. (2) The flow response of both catchments to rainfall includes a large contribution from saturated areas - - more for the Monachyle catchment than for the Kirkton. (3) No effects of land-use change are detected for the Kirkton catchment. For the Monachyle catchment the subsurface response seen in 1987-1988 is damped compared with that of 1984-1985 - - this may follow from the ploughing of the catchment before afforestation. ACKNOWLEDGEMENTS
The authors are grateful to Dr. Keith Beven for helpful comments and to Daniel Butterfield for his assistance with the graphics. REFERENCES Bergstrom, S., Carlsson, B., Sandberg, G. and Maxe, L., 1985. Integrated modelling of runoff, alkalinity and pH on a daily basis. Nordic Hydrol., 16: 89-104. Beven, K., 1982. On subsurface storm: predictions with simple kinematic theory for saturated and unsaturated flows. Water Resour. Res., 18(6): 1627-1633. Beven, K., 1986. Runoff production and flood frequency in catchments of order n: an alternative approach. In: V.K. Gupta, I. Rodriguez-Iturbe and E.F. Wood (Editors), Scale Problems in Hydrology. D. Reidel, Dordrecht, pp. 107-131. Beven, K.J. and Kirkby, M.J., 1979. A physically based variable contributing area model of basin hydrology. Hydrol. Sci. Bull., 24(1): 43-69.
370
A J. ROBSON ET AL.
Beven, K.J. and Wood, E.J., 1983. Catchment geomorphology and the dynamics of runoff contributing areas. J. Hydrol., 65: 139-158. Beven, K.J., Kirkby, M.J., Schoffield, N. and Tagg, A., 1984. Testing a physically based flood forecasting model TOPMODEL for three U.K. catchments. J. Hydrol. 69:119-143. Blackie, J.R. and Eeles, C.W.O,, 1985. Lumped catchment models. Cpt 11. In: M.G. Anderson and T.P. Butt (Editors), Hydrological Forecasting. Wiley, pp. 311-343. Calver, A., 1988. Calibration, sensitivity and validation of a physically based rainfall-runoff model. J Hydrol., 103: 103-105. Chapell, N.A., Ternan, J.L., Williams, A.G. and Reynolds, B., 1990. Preliminary analysis of water and solute movement beneath a coniferous hillslope in Mid-Wales, UK. J. Hydrol., 116: 201-217. Christophersen, N., Seip, H.M. and Wright, R.F., 1982. A model for streamwater chemistry at Birkenes, Norway. Water Resour. Res., 18: 977-996. Johnson, R.C. and Whitehead, P.G., 1993. An introduction to the research in the Balquhidder experimental catchments. J. Hydrol., 145: 231-238. Nash, J.C., 1979. Compact Numerical Methods for Computers. Hilger, Bristol, pp. 141-152. Quinn, P., Beven, K., Chevallier, P. and Planchon, O., 1991. The prediction of hillslope flow paths for distributed hydrological models using digital terrain models. Hydrol. Processes, 5(1): 59-80. Waters, D., Jenkins, A., Staples, T. and Donald, A., 1991. The Importance of Hydrological Source Areas in Terrestrial Liming. Journal of Institute of Water and Environmental Management, pp. 336-341. Wolock, D.M., Hornberger, G.M. and Musgrove, T.J., 1990. Topographic effects on flow path and surface water chemistry of Llyn Brianne catchments in Wales. J. Hydrol., 115: 243-259. Wolock, D.M., Hornberger, G.M., Beven, K.J. and Campbell, W.G., 1989. The relationship of catchment topography and soil hydraulic characteristics to lake alkalinity in the North Eastern U.S. Water Resour. Res., 25(5): 829-837.