MODIFICATION of STREAM FLOW ROUTING for BANK STORAGE
Mohamed M. Hantush,1 Associate Member, ASCE, Morihiro Harada,2 and Miguel A. Marino,3
Honorary Member, ASCE
Abstract
Bank storage is a process in which volumes of water are temporally retained by alluvial
stream banks during flood events, and gradually released as baseflows. This process has
implications on ground-water resource management, reducing flood peaks, and sustaining
riparian vegetation. In this paper, analytical solutions are developed that describe ground water-
surface water interactions and the impact of bank storage on the attenuation of stream channel
discharges. In effect, the stream flow routing Muskingum method is modified for bank storage.
The analysis is based on one-dimensional lateral groundwater flow in semi-infinite homogeneous
unconfined aquifers, which are in hydraulic contact with streams through semipervious bed
sediments. Storage in the stream reach, according to the Muskingum method, is assumed to be
proportional to the reach inflow and outflow rates. The stream-reach impulse response and unit
step response functions are modified for bank storage by applying the method of Laplace
transformation to the coupled stream-aquifer system. Results indicate that increasing
conductivity of the bank reduces the stream impulse response function at earlier times but with
greater and more persisting values later. In contrast, greater aquifer conductivity decreases the
unit step response at earlier time but with a diminishing effect at later times. The stream flows
are routed for a hypothetical asymmetric flood hydrograph, and the results show increasingly
attenuated and delayed peaks and extended tailing, with greater values of the hydraulic
conductivity. The simulated stream losses to bank storage and subsequent baseflows are
significant in typical alluvial sediments. Ground-water discharges are highly dependent on the
retardation coefficient.
Introduction
The problem of stream-aquifer interactions is important to watershed management efforts
aiming at mitigating hazardous flood events and optimizing surface water and ground-water
resources, and has significant ecological implications. For example, urbanization increases the
fraction of impervious lands in watersheds, which reduces ground-water recharge and increases
surface runoff and, thus, the potential for greater stream flows during flood events. Aquifers may
provide a temporary relief for increased stream flows through bank storage and, in effect, may
reduce and delay flood peaks. The temporally stored volumes of water in stream banks provide
1 Hydrologist, Subsurface Protection and Remediation Division, National Risk Management Research Laboratory,
ORD, U.S. EPA, 919 Kerr Research Dr., Ada, OK 74820. Corresponding Author. Email:
hantush.mohamed@epa.gov
2 Associate Professor, Department of Civil Engineering, Meijo University, Tenpaku, Nagoya 468-8502, Japan.
Email: harada@meijo-u.ac.jp
3 Professor, Department of Land, Air and Water Resources, and the Department of Civil & Environmental
Engineering, University of California, Davis, CA 95616. Email: mamarino@ucdavis.edu
1
-------
moisture needed to sustain riparian vegetation, and when released gradually from storage sustain
aquatic organisms during baseflow periods.
The lateral stream losses to bank sediments during a flood stage and the subsequent
gradual releases of the volume of waters in storage to sustain baseflows between flood events is a
process called bank storage (Todd, 1955). This problem has been the topic of numerous papers in
the literature; it has been analyzed numerically by accounting for streams as time-varying
boundary conditions (Yeh, 1970; Marino, 1975), and by modeling the stream-aquifer system as
two interacting dynamic systems (Finder and Sauer, 1971; Zitta and Wiggert, 1971). Although
more restrictive, practical analytical solutions were also developed for solving the linearized
Boussinesq groundwater flow equation subject to fluctuating stream stages (Cooper and
Rorabaugh, 1963; Marino, 1973; Moench et ai, 1974; Govindaraju and KoeHiker, 1994; and
Zlotnik and Huang, 1999). These analytical solutions, however, do not consider the simultaneous
interactions, which are inherent to a dynamically coupled stream-aquifer system; rather, their use
is conditioned on a priori knowledge of the stream-stage fluctuations. Thus, they cannot be used
for routing stream flows directly, unless iterative procedures are implemented, and their utility to
watershed planning and management (e.g., evaluating the impact of watershed landscape
changes) is thereby questionable. Hunt (1990) derived an analytical solution, which couples a
linearized approximation of the kinematic wave equation in open channels to aquifer flow, and
presented a dynamic model for simulating stream-aquifer interactions. The solution, however,
required iterating between the groundwater solution and the flood-routing problem. Harada el a!.
(2000) introduced the concept of impulse response function to relate stream outflows to bank
storage, by coupling a linear reservoir model to a semi-infinite aquifer in perfect hydraulic
connection with the aquifer. The presented analytical solution, however, is limited to narrow
streams and aquifer conductivities that are much greater than typically encountered in alluvial
sediments. This paper presents an analytical solution to a dynamic stream-aquifer system, by
modifying the stream routing Muskingum method (Chow et ai, 1988) for bank storage in
alluvial aquifers of semi-infinite extent and separated from the streams by semi-pervious bed
deposits.
Stream Flow
The Muskingum method is a widely used method for hydrologic routing in streams
channels. In this method, the volume storage, S(t), in a stream channel is expressed as a
combination of a wedge storage, r\ % [/(/) - 0(0], and a prizm storage, r\ 0(f) (Chow et ai,
1988),
s(t)=n[^i(t) + (i-i>)0(t)} (i)
in which S(t) = channel storage [L3]; I(t) = inflow rate relative to the initial discahrge [L3/T]; 0(t)
= ouflow rate relative to the initial discharge [L3/T]; T] = storage time constant for the reach [T];
and 4 = a weighting factor that varies from 0 to 0.5. In natural streams, E, varies from 0 to 0.3 and
averages about 0.2. The storage time constant, rj, is approximately the kinematic wave travel
time through the reach, and can be estimated as the time interval between the inflow and outflow
peaks. If ^ = 0, S(t) = r| 0(t). This storage-discharge relationship is commonly used in level pool
routing.
2
-------
The storage in a stream reach of length L penetrating an aquifer may be described by the
continuity equation (Fig. 1):
*^± = I(t)-0(t)-2Q(t) (2)
dt
in which Q(t) = half the lateral flow rate to or out of the aquifer [L3/T] integrated over the reach
length. In the following section, we present the boundary-value problem which describes the
ground water-surface water interactions, Q(t). It is assumed that S(t) is related to the average
stream stage fluctuations, H(t),
S(t) = LH(t) (3)
Initial water-table
position
0(t)
Fig. 1 Schematic diagram of the stream-aquifer system.
in which H(t) = stream-stage fluctuations, relative to the initial equilibrium stage (Fig. 1),
averaged over the stream-reach length, L. This assumption is intended to simplify the analysis.
Also, it is assumed that the stream stage is initially at equilibrium with the water table in the
aquifer and, thus, flow in the stream is initially uniform, i.e., 1(0) = <9(0) = 0.
Ground-water Flow at the Interface
Ground-water flow in a semi-infinite homogeneous unconfined aquifer may be described
by the linearized Boussinesq partial differential equation (Fig. 1):
dh(x,t) = Dd2h(x,t)
dt dx2
3
-------
and subject to the initial and boundary conditions:
h(x,0) = 0 (5)
Tdh(0,t) _pK,[H(t)-h(0,t)\
dx b
h(oo,t) = 0 (7)
in which h(x,t) = water-table fluctuations relative to the initial equilibrium position [LI; D = T/S
= aquifer diffusivity [L2/T]; T = K h0 = average transmissivity of the aquifer [L /T]; K =
hydraulic conductivity [L/T]; h0 = average water-table elevation above the base of the aquifer
(Fig. 1) [L]; S = specific yield of the unconfined aquifer; P = half the wetted perimeter of the
stream [L]; K' = the hydraulic conductivity of the stream bed [L/T]; and b = thickness of the
stream bed [L], In Eq. (4) the water-table fluctuations are assumed to be small relative to the
average saturated thickness of the aquifer. Also, the third-type boundary condition (6) assumes
that compressibility of the aquifer below the stream can be ignored. Thus, ground-water flow
below the stream is at quasi-steady state, which may be valid only for aquifers deeply incised by
the streams.
Stream losses (or baseflow) across a channel length L is given by
Q(t) = PLK'^H^~hi0,t^ (8)
b
In the following section, we use the above equations to derive the impulse response and unit step
input response functions for streams modified for stream-aquifer interactions.
Impulse Response Functions
It can bee seen that the substitution of the Muskingum relationship (1) for S(t) in (2)
results in a linear system, which relates 0(t) to I(t) and Q(t). Thus, 0(t) is uniquely characterized
by its impulse response function and can be described by the convolution integral:
<9(0 = \u(t-x)I(x)dx (9)
o
in which u(t) = impulse response function [T'1], which describes the temporal variations of the
outflow from the stream reach due to an instantaneous input of unit amount at t = 0 at the
upstream inflow boundary. In the Laplace transform domain, Eq. (9) is described by
0(p) = u(p)7(p) (10)
where the Laplace transform of a function j{t) is defined by f(p)=\f(t)e~p'dt. The
o
application of the Laplace transformation to Eqs. (1-8) and comparison with (10) should yield
4
-------
7 = __1 1 + R^p/D
11 p ~ (1 -^)2 r| (R/jD)pV2+p + y(R/jD)4^+ [1/(^(1-^]'^ U)
where 7 = 1/(r|(l - %)) + (2 P/W)(K'/b), and R = Tb/PK is the retardation coefficient [L]
modified for partial penetration. In the specific case of P ~ h0, we have R - K bj K', which is
equivalent to a stream completely penetrating the aquifer. The inverse Laplace transform of (11)
is given by
1 1 °T \ + r47Td ez' , 5 c, ,
u(0~ i J 7=— 7= 7= r i^Z 8(/) (12)
(1 - Q r| 2ni a_,-» (./? /-Jd)z + z + y(R /-Jd)Jz + [l /(r|(l - ^)] 1-^
in which 8(/) is the Dirac delta function. Noting that the denominator in the integrand (we refer
to the integrand as F(z)) is a multiple-valued function of the complex variable z, then the integral
is evaluated on the complex plane by introducing a branch cut along the negative real axis (Fig.
2) with the argument of the principal branch of Vz defined from -n to n (z = re'e and
-k
(l-£) T1 7t W o A(^) 1-^
where
A(y) = R2y2(y-y2)2+D[\,(T](\-^)-y2^ (15)
Harada et al. (2000) obtained a closed-form solution for u(t) for the specific case of level pool
routing (^ = 0), perfect hydraulic connection with aquifer (R = 0), and when SKhQ/W2 > 1 /r|.
This latter condition, however, is limited to narrow channel widths and aquifer conductivities
much greater than those encountered in natural alluvial sediments. Figure 3 shows the effect of
the hydraulic conductivity on the impulse response function (the Dirac-delta contribution is not
shown). Notice that as K increases, u(t) shows increasingly smaller values at earlier time t but
increasingly greater values at later time t and extended tailing. This implies that bank storage
regulates stream discharges by attenuating their peaks due to initial lateral losses to the aquifer
and continuously supplying the stream outflows by gradual releases from bank storage through
5
-------
a +1 °°
z = x + i y
F, : z = u2eni, *fz = ue
T2 : z = u2e n>, -Jz •
¦ ue
a-i*°
Fig. 2 Contour integration for the Laplace inverse transformation of the multiple-valued
function F(z).
baseflows. The Laplace transform of Q(() can be shown to be
T
Q(P)
JpTd
(l-s)^ (R/jD)Pm+p + y(R/Jd)Jp+ [1/^(1-^
Hp)
(16)
From this equation one can deduce that the impulse response function that describes ground
water-surface water intearctions is given by the Laplace inverse transform of the product of the
first two terms on the right-hand side,
« (/)
t* 1 a-H<~
- I
¦4z!D
271 / Ji~{RlJD)zvl +z + y(/?/Vd)Vz + [1/(tiO-^)]
dz (17)
Similarly, the integration of the complex-valued function is carried along the contour lines
shown in Fig. 2, and can be shown to be:
dy (18)
2T4D 1 e~y/{ 2
u (t) = — — y
and Q(r) is given by the convolution integral:
6
-------
3
r] =0.4 hr
hQ~ 20 m
P = 20m
W=2Q m
2
K = 40 m/hr ;
o
o.o
0.5
2.0
t (hrs)
Fig, 3 Impulse response function for K = 0 (no-bank storage), 2, and 40 m/hr.
Q(t) = \u(t- x)I(x) dx (19)
o
Figures 4(a)-4(b) display the stream outflow and ground-water flow hydrographs in response to a
hypothetical asymmetric flood-inflow hydrograph and assumed hydraulic parameters. The inflow
hydrograph is assumed to be of the type proposed by Cooper and Rorabaugh (1963), 1(f) = N /"
e~8' [l-cos(G)/)], when 0 < t < tc, and 1(f) = 0, t > tc. The routed outflow hydrograph, (f is
the peak of the inflow flood event), is obtained by numerical integration of (9) with the impulse
response function, u(t), given by (14 and 15), for K = 0, 2, and 40 m/hr. As the aquifer hydraulic
conductivity increases, the outflow hydrograph displays greater attenuation and a delayed peak
outflow, but with extended tailing; that is, greater baseflow rates proceed after the end of the
flood event. Figure 4(b) shows ground-water flux at the interface, 2 x Q(t)/I , for K = 2 m/hr and
R = 0, 5, 20, and 50 m, and K = 40 m/hr and R = 0 m. It is clear that bank storage has greater
effect on ground water-surface water interactions as K increases. The hydraulic conductivity K =
40 m/hr may be encountered in predominantly gravelly-coarse sand sediments. Also, note that
increased retardation coefficient results in delayed but decreasing baseflow rates. Q{t)H' is
calculated by numerical integration of (19) with u (t) given by (18).
Unit Step Response Function
If the stream inflow rate increases from 0 to 1 at time t = 0, and continues indefinitely at
that rate, then the outflow rate in response to this unit step increase is given by the unit step
response function g(t) (Chow et al,, 1988). It is given by the integration of the impulse response
function,
7
-------
0(t) (K=0 m/hr)
0(t) (K = 2.0 m/hr)
0(t) (K = 40 m/hr)
R = 0 m
S = 0.2
P = 20m
W=20m
£=0.15
rj =0.4 /ir
/i = 20 m
t (hrs)
O)
g(t) —> 1, and shows less variation with K. This
implies that the impact of bank storage diminishes after large time, assuming that the flood
event persists.
The unit step response function g(t) is useful in routing stream flows when the inflow
hydrograph is of a general form and measured at discrete points rather than continuously in time,
as shown in Fig. 6. In this case, the outflow hydrograph (0„ = 0{t„)) at a discrete point in time,
t„, can be obtained by breaking the convolution integral in (9) into summation of integrals over
8
-------
eo
K = 0 m/hr
K= 2 m/hr
K = 40 m/hr
t (hrs)
Fig. 5 Unit step response function g(/) as a function of time (A^ = 0, 2, 40 m/hr).
time increments A tm = tm - tm./, between successive discrete measurements of /(/), Im, 0
-------
analytic methodology is presented, which modifies the Muskingum hydrologic routing method
for bank storage in streams cutting through alluvial aquifers. Integral expressions are obtained
for the impulse response functions, which allow for routing continuous-time inflow hydrographs
in streams, and estimating lateral losses to the surrounding bank sediments, and subsequent
baseflows. An analytical expression is presented for the unit step response function, which is
useful for routing general and discrete inflow hydrographs. Simulation results shQwed that with
K typical to alluvial sediments, bank storage can significantly reduce flood peaks and sustain
baseflows. The results also indicated that ground-water fluxes decrease significantly with
increasing retardation coefficient, thus minimizing the effect on regulating stream outflows.
Notice'. This paper has been reviewed in accordance with the U.S. Environmental Protection
Agency's peer and administrative review policies and approved for presentation and publication.
This research reported in this paper is part of Project 4940-H of the Agricultural Experiment
Station of the University of California, Davis.
References
Chow, V. T., D. R. Maidment, and L. W. Mays. 1988. Applied Hydrology. McGraw-Hill, Inc.
Cooper, H. H., and M. I. Rorabaugh. 1963. Ground-water Movements and bank storage due to
flood stages in surface streams. U.S. Geol. Surv. Water Supply Pap. 1536-J, 343-366.
Govindaraju, R. S. and J. K. Koelliker. 1994. Applicability of linearized Boussinesq equation
for modeling bank storage under uncertain aquifer parameters. J. Hydro!. 157, 349-366.
Harada, M., M. M. Hantush, and M. A. Marino. 2000. Hydraulic analysis on stream-aquifer
interaction by storage function models. To appear in the Proceedings of 1AHR International
Symposium 2000 on Groundwater.
Hall, F. R. and A. F. Moench. 1972. Application of the convolution equation to stream aquifer
relations. Water Resour. Res., 8(2): 487-493.
Hunt, B. 1990. An approximation for the bank storage effect. Water Resour. Res., 26(11), 2769-
2775.
Marino, M. A. 1973. Water table fluctuation in semipervious stream-unconfined aquifer systems.
J. Hydro!., 19, 43-52.
Marino, M. A. 1975. Digital simulation model of aquifer response to stream stage fluctuation. J.
Hydro!., 25, 51-58.
Moench, A. F., V. B. Sauer, and M. E. Jennings. 1974. Modification of routed streamflow by
channel loss and base flow. Water Resour. Res., 10(5), 963-968.
Pinder, G. F. and S. P. Sauer. 1971. Numerical simulation of flood wave modification due to
bank storage effects. Water Resour. Res., 7(1), 63-70.
Todd, D. K. 1955. Ground-water in relation to a flooding stream. Proc. Amer. Soc. Civil Engrs.,
81, pp. 1-20, separate 628.
Zitta, V. L. and J. M. Wiggert. 1971. Flood routing in channels with bank storage. 7(5), 1341-
1345.
Ych, W.W-G. 1970. Nonsteady flow to a surface reservoir. ASCE J. Hydraul., 96(3),609-618.
Zlotnik, V. A. and H. Huang. 1999. Effect of shallow penetration and streambed sediments on
aquifer response to stream stage fluctuations (analytical model). Ground Water, 37(4), 599-
10
-------
NBMRL-ADA-00221
TECHNICAL REPORT DATA
1
1. REPORT NO.
EPA/600/A-01/052
2.
4. TITLE AND SUBTITLE
Modification of Stream Flow Routing for Bank Storage
5. REPORT DATE
6, PERFORMING ORGANIZATION CODE
7. AUTHOR(S>
'Kohamed M. Hantush
'Morihiro Harada
'Miguel A. Marifto
8. PERFORMING ORGANIZATION REPORT NO.
9. PERFORMING ORGANIZATION NAME AND ADDRESS
lUSEPA, ORD, NRMRL, SPRD, P.O. BOX 1198, ADA, OK 74820
!Dept of Civil Eng., Meijo University, Tenpaku, Nagoya 468-8502, Japan
10. PROGRAM ELEMENT NO.
'Sept of Land, Air, & Water Resources and Dept of Civil & Environ. Eng
University of California, Davis, CA 9S616
11. CONTRACT/GRANT NO.
In-House
12. SPONSORING AGENCY NAME AND ADDRESS
U.S. EPA
NATIONAL RISK MANAGEMENT RESEARCH LABORATORY
SUBSURFACE PROTECTION AND REMEDIATION DIVISION
P.O. BOX 1198? ADA, OK 74820
13. TYPE OF REPORT AND PERIOD COVERED
Symposium Paper
14. SPONSORING AGENCY CODE
EPA/600/15
15. SUPPLEMENTARY NOTES
PROJECT OFFICER: Mohamed M. Hantush 580-436-8531
16. ABSTRACT
Bank storage is a process in which volumes of water are temporally retained by alluvial stream banks during
flood events, and gradually released as baseflows. This process haB implications on ground-water resource
management, reducing flood peaks, and sustaining riparian vegetation. In this paper, analytical solutions are
developed that describe ground water-surface water interactions and the impact of bank storage on the
attenuation of stream channel discharges. In effect, the stream flow routing Muskingum method is modified for
bank storage. The analysis is based on one-dimensional lateral groundwater flow in semi-infinite homogenous
unconfined aquifers, which are in hydraulic contact with streams through oemipervious bed sediments. Storage
in the stream reach, according to the Muskingum method, is assumed to be proportional to the reach inflow and
outflow rates. The stream-reach impulse response and unit step response functions are modified for bank
storage by applying the method of Laplace transformation to the coupled stream-aquifer system. Results
indicate that increasing conductivity of the bank reduces the stream impulse response function at earlier times
but with greater and more persisting values later. In contrast, greater aquifer conductivity decreases the
unit step response at earlier time but with a diminishing effect at later times. The stream flows are routed
for a hypothetical asymmetric flood hydrograph, and the results show increasingly attenuated and delayed peaks
and extended tailing, with greater values of the hydraulic conductivity. The simulated stream losses to bank
storage and subsequent baseflows are significant in typical alluvial sediments. Ground-water discharges are
highly dependent on the retardation coefficient.
17 .
KEY WORDS AND DOCUMENT ANALYSIS
A. DESCRIPTORS
B. IDENTIFIERS/OPEN ENDED TERMS
C. COSATI FIELD, GROUP
18. DISTRIBUTION STATEMENT
RELEASE TO PUBLIC
19. SECURITY CLASS(THIS REPORT)
UNCLASSIFIED
21. NO. OF PAGES
10
20. SECURITY CLASS(THIS PAGE)
UNCLASSIFIED
22. PRICE
EPA FORM 2220-1 (REV.4-77) PREVIOUS EDITION IS OBSOLETE
------- |