|
|
Bridging the Gap between Dual-Porosity Model and Discrete Fracture
Network Model for the Analysis of Groundwater Flow in a Fractured Rock Mass
LEE, JUNG-WOO
Graduate Student, Department of Civil Engineering, Yonsei
University, Seoul, Korea,120-749
Tel: 82-2-361-2805, Fax: 82-2-364-5300, E-mail:
ljw007@netian.com
ABSTRACT
The
dual-porosity model and the discrete fracture network model are being used for
the analysis of groundwater flow in a fractured rock mass. The dual-porosity
model assumes that the medium consists of two continua at the macroscopic
level, that is, fracture and rock matrix block. Saturated water flow in the
matrix as well as in the fracture pore system is described using two separate
flow equations which are coupled by means of water transfer term. The
formulation leads to two coupled systems of nonlinear partial differential
equations. The discrete fracture network model assumes that fractures form the
paths for seeping water and that the conductivity in rock matrix is neglected.
Node equations are developed assuming that cubic law governs the water flow
through fractures. The relationship between the dual-porosity approach and the
discrete fracture network approach was sought in order to bridge the gap
between two models.
Keywords : discrete
fracture network model, dual-porosity model, cubic law
INTRODUCTION
Basically
three types of approaches are used for the analysis of groundwater flow in a
fractured rock mass. The first is the continuum approach that treats rock
formation as an equivalent porous medium. This approach is based on
well-developed theories. The second is the dual-porosity (DP) approach which
assumes that the medium consists of two continua, one associated with the
fracture networks and the other with a less permeable pore system of rock
matrix block. The third is the discrete fracture network (DFN) approach which
hypothesizes that fractures form the paths for seeping water and the rock
matrix is impervious.
The
continuum approach is good for both saturated and unsaturated flows in porous
media which are not too heterogeneous. However, if there is a strong
non-homogeniety in a porous medium due to fractures, fissures, cracks, faults,
and worm hole, the flow and solute transport can not be predicted realistically
by this method. Instead, the DP model or the DFN model should be used for the
analysis of such preferential flows and contaminant transport by them. The DP
model seems to be preferred by hydraulic engineers while the DFN model is by
geotechnical engineers and geologists. Both models have been applied to similar
problems, however, no previous efforts have been made to show distinct features
of each approach. The purposes of this study were to investigate the DP and the
DFN approaches for fractured rock mass and to bridge the gap between two
approaches.
DUAL-POROSITY MODEL
Governing
Equations
The DP concept assumes that the porous
medium can be separated into two distinct pore systems, which are considered as
homogeneous media with separate flow and solute transport properties (Gerke and
van Genuchten, 1993a). The two porous media exchange water and solutes
reciprocally in response to pressure head and concentration gradients.
Therefore macroscopically, the porous medium is characterized by two pore water
velocities, two pressure heads, and two solutes concentrations.
The flow governed by Darcy's law is
considered in both the fracture and matrix pore system, while the transfer of
water between the two pore systems is described macroscopically using a space
and time dependent exchange term. Two-dimensional water flow in two pore
systems is described by
(1a)
(1b)
where (1a) and (1b) describe the flow of
water in the fracture (subscript f)
and matrix (subscript m) pore
systems, respectively. In these equations, 
is
the total head (L), 
is
the specific storage (L-1), 
is
the hydraulic conductivity (LT-1), 
is
the fracture volume fraction, 
is
time, and 
is
the exchange term (T-1) describing the transfer of water between the
two pore systems. The exchange term is assumed to be proportional to the
difference in pressure head between the two porous media as follows (Gerke and
van Genuchten, 1993b).
, 
(2)
Where 
is
a first-order transfer coefficient for water (L-1T-1), 
is
a factor depending on the geometry of the matrix blocks, 
means the distance (L) from the center of a fictitious matrix to
the fracture boundary, and 
is an empirical coefficient. van Genuchten
and Dalton (1986) suggested that the geometry factor 
equals 3 for rectangular slab, 15 for spheres, and Gerke and van
Genuchten (1993a) proposed that the coefficient 
is
found to be 0.4. Herein, the governing equations of the DP model are solved
simultaneously by using the finite difference method.
Analytical
Solution of Dual-Porosity Model
Önder
(1998) derived an analytical solution of the head variation resulting from a
sudden rise of stream stage in a confined finite fractured aquifer bounded by a
stream on one side and by an impervious boundary on the other. In this paper, Önder's analytical solution is modified by
revising a part of transformations. If we neglect the conductivity of the
matrix, and consider the storage of the matrix and the quasi-steady water
transfer between fractures and pores of blocks, the governing equations for
one-dimensional confined fracture flow are written by
(3a)
(3b)
where subscripts 1, 2 means fracture and
matrix, respectively, 
is
the total head, 
is
the coefficient of transmissivity, 
is
the storage coefficient, and 
relates to the geometry of the fractured rock. For the analytical
solution of equation (3), integral transforms are used: the space derivative, 
may be replaced by using a finite Fourier sine transform and time
derivative, 
by
Laplace transform, which leads to these equation (4a) and (4b),
(4a)
(4b)
where 
, 
, 
, 
, 
, 
, 
=
length of the aquifer in the 
-direction, 
=
initial total head, 
=
initial head increase, 
=
Fourier transform variable, 
, 
, 
, 
, and 
.
To demonstrate the
validity of the DP model a comparison with analytical solution is presented.
For this, a confined finite fractured aquifer (800m
40m) bounded by a stream on one side and
by an impervious boundary on the other is selected. We assume that the sudden
rise of stream stage adjacent to an aquifer cause the one- dimensional
transient flow in an aquifer, and also assume that a typical fractured aquifer
has the following hydraulic properties and virtual field conditions. The
fracture hydraulic conductivity 
is
10 m/day (
= 20 m2/day), the matrix hydraulic
conductivity 
is
10-4 m/day (
= 3.8
10-3 m2/day), and
the hydraulic properties for 
is
assumed to be the same as that for matrix. The storage coefficients of fracture
and matrix are 
= 1.4
10-5 and 
= 1.4
10-4, respectively. We further assume that the medium consists of
rectangular blocks (
= 3) with an average matrix block size of 
= 5 m, 
=
0.4, and 
= 5%. Fig. 1 shows the spatially transient
changes after 
= 10 hours. It is seen that the numerical
solution agrees well with the analytical solution and the DP model describes
the preferential flow which is the non-equilibrium flow caused by
heterogeneities due to fractures.
DISCRETE
FRACTURE NETWORK MODEL
Flow
Law in a Fracture
In
the DFN model, unknown heads at the intersections of fracture network are
computed based on the laws for flow through individual fractures by setting a
condition of continuity of flow for each intersection. The energy equation or
relationship between flow rate and head loss for a single fracture can be
written as
(5)
where 
is
the flow rate, 
is
the head loss, 
is
a constant which depends on the flow regime (laminar flow 
= 1, turbulent flow 
< 1), and 
is
the conductivity coefficient. The conductivity coefficient is known to be a
function of aperture, roughness of fracture, length of segment, fluid
viscosity, and gravity (Barbosa, 1990).
Most
studies of fluid flow through open fractures are based on the analogy of flow
between parallel plates with smooth walls. For this case, the flow in the
laminar regime often referred to as the 'cubic law'. This law can be derived
from the Navier-Stokes equations under the assumptions of incompressible fluid,
steady-state flow, parabolic velocity distribution, and smooth fracture walls
(Bear et al., 1993). The cubic law is given by
(6)
where 
is
the kinematic viscosity, 
is the aperture, 
is the fracture length, and 
is
the gravity. It is seen that the hydraulic
conductivity depends heavily on the aperture.
Numerical
Solution
In
the numerical solution of the network problem, intersections are numbered and
each conducting fracture is identified by the numbers of the two ending nodes.
In order to account for the flow direction within fractures, equation (6) for
the fracture i, j is rewritten as
(7)
where 
is
the flow rate from node i to node j, 
is
the head at node i, 
is
the head at node j, and 
is
the conductivity coefficient for fracture i,
j. The continuity equation at node i
can be written as
(8)
where 
is
the flow rate from node k to node i, 
is
the number of fractures connected to node i.
The energy equations and continuity equations can be combined into such node
equations (Kirchoff rule) as
(i
= 1, 2, ... , 
) (9)
where 
denotes the number of nodes with unknown head (
).
Time rate of fluid mass
storage change term is added to RHS of equation (8) to expand node equations
for transient flow. The law of conservation of mass for transient flow requires
that the net rate at which fluid enters a control volume is equal to the time
rate of change of fluid mass storage within the control volume. That is,
(10)
where 
is
the specific storage, 
is
the time step, and 
is
the control volume. Substituting the energy equation into equation (10), node
equations for transient flow is written by
(i =
1, 2, ... , 
) (11)
where
is an approximated head at 
-th time level. The solution of these
equations may be obtained by succesive iterations by using Newton-Raphson
technique.
SIMULATION EXAMPLE
Both DP and DFN models are applied to a
dam-seepage problem. Transient flow resulting from the rise of upstream stage
is analyzed by assuming that the region is surrounded by the impermeable
boundaries as shown in Fig. 2. For the DP model, values of parameters such as 
=
1%, 
= 0.82 m/sec, 
= 4.1
10-4 m/sec, 
=
4.1
10-6 m/sec, and 
=
0.002 are used. In estimating water transfer term, values of 
= 5 m, 
=
3, and 
= 0.4 are used. For the application of the
DFN model, a regular fracture net of equi-distant square grids are used. The
aperture is assumed to be 1 mm and the spacing between fractures is 10 m. In
order to make the situation same, the conductivity of the fracture in the DP
model is calculated from aperture in the DFN model. Fig. 3 and Fig. 4 show
computed equi-potential lines after 1 hour by using the DP and DFN models,
respectively. It is observed that the heads computed by the DP model are
apparently higher than the heads by the DFN model especially at the upstream
site.
RELATIONSHIP BETWEEN TWO MODELS
Each
model uses different parameters. That is, the DP model requires such parameters
as hydraulic conductivity, water exchange term, and fracture volume fraction,
and the DFN model needs aperture and fracture length. Simple sensitivity
analysis have revealed that important parameters are fracture volume fraction
in the DP model and aperture in DFN model, both of which are in the length
dimension (Lee, 1997). In the light of this, it can be deduced that the discrepancy
in the length scale in two models may cause different results.
In
order to establish a relationship between two length scales of both models,
aperture is varied until two head distributions from the DP model and the DFN
model are the same. In the previous problem, when aperture about 2.1 mm is used
in the DFN model, the computed head distributions are the same as ones from the
DP model. The results are given in Fig. 5. A value of fracture volume fraction
is now converted to aperture by 
. That is, 1 % of fracture volume fraction
corresponds to aperture of 0.05 m. Therefore it is found that the length scale
of the DP model is 10-102 times greater than the length scale of the
DFN model.
CONCLUSIONS
The
dual-porosity model and the discrete fracture network model are used for the
analysis of preferential flow in the fractured rock mass. Objectives of both
models are the same. However, each model involves different parameters, and is
being used by different group of people. In this paper, attempts have been made
to investigate distinct features of two models and to build a bridge between
two models. It is found that important parameters are fracture volume fraction
in the dual-porosity model and aperture in the discrete fracture network model.
By solving a dam seepage problem, it is shown that the length scale of the dual
porosity model is 10-102 times greater than the length scale of the
discrete fracture network model.
REFERENCES
Barbosa, R. (1990). "Discrete Element Models for
Granular Materials and Rock Masses." Ph.D. thesis, Department of Civil
Engineering, The University of Illinois at Urbana-Champaign.
Bear, J., Tsang, C.F., and Marsily, G. de. (1993). Flow and contaminant transport in fractured
rock, Academic press, Inc.
Gerke, H.H. and van Genuchten, M.T. (1993a). "A
dual-porosity model for simulating the preferential movement of water and
solutes in structured porous media." Water Resour. Res., 29(2), 305-319.
Gerke, H.H. and van Genuchten, M.T. (1993b).
"Evaluation of a first-order water transfer term for variably saturated
dual-porosity flow models." Water Resour. Res., 29(4), 1225-1238.
Lee, J-.W. (1997). "A comparative study on the
dual-porosity model and the discrete fracture network model for the analysis of
groundwater flow in a fractured rock mass." M.E. thesis, Department of Civil
Engineering, University of Yonsei, Seoul, Korea.
Önder, H. (1998). "One-dimensional transient flow
in a finite fractured aquifer system." Hydrol. Sci. J., 43(2), 243-265.
van Genuchten, M.T. and Dalton, F.N.
(1986). "Models for simulating salt movement in aggregated field
soils." Geoderma, 38, 165-183.
|
|
|
|
Fig.1 Comparison
of numerical solution and analytical solution |
Fig. 2 Dam-seepage flow region |
|
|
|
Fig. 3
Equipotential line after 1 hour (DP model) |
|
|
|
Fig. 4
Equipotential line after 1 hour (DFN model) |
Fig. 5 Comparative head distribution
