Figures have not been provided with the following write up

A LEAKY AQUIFER MODEL FOR HARD ROCK AQUIFERS

M. SEKHAR, TRANSOFT International

M.S. Mohan Kumar and K. Sridheran, Department of Civil Engineering , Indian Institute of Science, Bangalore, India.

ABSTRACT A Comparison is made between the unsteady parallel fracture model using the double porosity approach and the anisotropic aquifer – water table aquitard system. The presence of moving water table in the weathered zone (unconfined block) influenced by recharge from rainfall or other sources and leakage to the fracture zone is discussed. Effect of fracture skin on the block-fracture interface is studied for both double-porosity and aquifer-water table aquitard systems. Typical type curves are presented. It is found that the effect of both fracture skin and water table would be to have a flat draw down behavior for a large time. The applicability of the present model for parameter estimation, is illustrated on a field pumping test of 7 days duration.

INTRODUCTION

Ground water flow and storage in hard rock areas has always been a mater of great interest and importance to hydrogeologists and petroleum engineers. Current activities have been stimulated by water supply and irrigation, waste disposal and energy related concerns.

The development in well hydraulics in the 1950 's and early 1960 's provide an analytical basis for study of ground water systems, but it is generally acknowledged that these theories developed for a homogeneous porous medium, may not be applicable for heterogeneous and anisotropic fractured rock aquifer system. It was recognition that flow through such a system is very significantly influenced by the fracture characteristics.

Investigations based on different approaches progressed in parallel, in dealing with the complexity introduced by the discontinuity in terms of fractures. One of these approaches involved adaptation of the continuum approach similar to that used in the classical study of flow through porous media. However, an attempt is made to incorporate some of the features of flow through a fractured rock, as in the case of the double porosity models. An alternative approach is based on understanding of flow fractures, based on the assumption that the rock matrix surrounding the fractures is practically impermeable. While early works in this area have focussed attention on single fracture, there have been several studies recently with complex two dimensional, and even three dimensional discrete, stochastic, fracture network models.

It might be noted that inherently the problem of understanding flow and storage in fractured rock aquifers is very complex. There is the heterogeity associated with the mechanical discontinuity resulting from the presence of fractures. The matrix or block surrounding the fracture is generally not impermeable. Though its permeability may be very small, flow transfer takes place between the fractures and the blocks. The fractured formation also generally displays anisotropic characteristics.

In the present paper, a double- porosity model based on an aquifer- water table aquitard concept is proposed for hard rock aquifers. The similarities and differences between this model and the classical double-porosity model are discussed based on typical results. The potential for using for the present model for parameter identification in an anisotropic fracture aquifer system is also illustrated.

DOUBLE POROSITY SYSTEM

Studies of fluid flow in a fractured rock mass where fissure flow is augmented by contributions from blocks, have generally adopted the "double porosity" concept proposed by Barenblatt et al(1960). In these models, the fractured rock mass is assumed to consist of two interacting overlapping continua, a continuum of low permeability, primary porosity blocks and a continuum of high permeability, secondary porosity fissures. The secondary permeability associated with the fissures is usually high, where as the secondary porosity is often low. On the other hand, the blocks are characterised by relatively high storage associated with very low permeability. The mathematical representation of the system requires two mass conservation equations, one for the fracture regime and another for the porous block regime. These two equations are coupled by a fluid transfer term which depends on the pressure difference between the block regime and fracture regime at the location.

Figure 1. presents a parallel fracture double-porosity system with a well fully penetrating the fracture aquifer over lain by a rock matrix block. The flow in the aquifer is assumed lateral and flow in the block is assumed vertical. The pumping well is cased in block zone and screened in the fracture. Such a schematization as been studied in the past by Warren and Root (1963), Kazemi et al. (1969), Gringarten(1982), and Huyakorn et al.(1983). Figure 1 also presents a fracture skin which is a thin skin of low permeability material deposited on the surface of the blocks, which impedes the free exchange of fluid between the blocks and fractures as proposed by Moench (1984). Moench has attempted to unify the prevailing theories of flow to a well in a double porosity reservoir by accounting for fracture skin.

In the following discussions, the terminology fracture or aquifer, and block or aquitard are used interchangeably. The governing equations in fracture and block in notations of the present paper, due to pumping from a well ignoring well bore storage are as follows:

Fracture:

hard_Equ1.gif (920 bytes)

Eqn.(1) represents the equation of flow in the fracture where s is the drawdown in the aquifer, s' is the drawdown in the aquiferd, z is the vertical co-ordinate , t is the time and r is the radial coordinate in the equivalent isotropic domain. Tr and St are the equivalent transmisstivity and storage coefficient of the fracture respectively and K¢ is the vertical hydraulic conductivity of the block.

The initial and boundary conditions in the fracture are:

hard_Equ2.gif (1067 bytes)

where Q is the constant pumping discharge.

Block:

hard_Equ3.gif (594 bytes)

Eqn.( 3 ) presents the equation of flow in the block where S¢ s is the specific storage of the block. The initial and boundary conditions are :

hard_Equ4.gif (1037 bytes)

Where SF is dimensionless fracture skin parameter given by ( Moench, 1984 ),

hard_Equ5.gif (525 bytes)

In eqn.(5) Ks is the permeability of fracture skin (Figure 1)

The system of eqns.(1)-(4) have been solved in the past using two approaches. One approach ignores spatial variation of hydraulic head gradients in the block, ignores the fracture skin and also assumes that the flux of fluid from the blocks to fissures occurs in response to the difference in the average hydraulic head in the fissures and the average hydraulic head in the blocks. This is the assumption of pseudo-steady state flow (Warren & Root, 1963). Models that are based on this approximation show that, well discharge consists initially of fluid derived primarily from storage in the fissures followed, at later time, by fluid derived primarily from storage in the blocks. At early and late times, draw down should therefore follow the familiar This type curve. During a sufficiently long transition from early to late time, however, draw down will approach a plateau (Gringarten, 1982; Moench, 1984 ).

Because of the simplicity of the approach and will test data that exists supporting this assumption, pseudo-steady state flow had received much attention. However, for theoretical justification a transient flow model from blocks to fissures is necessary. Kazemi et al. (1969 ) considers the second approach of transient block-fissure flow and used a numerical solution to simulate flow while ignoring the fracture skin. Type curves for this model differ from those obtained under the pseudo-steady state flow approximation during the period of transition from early to late time. Moench ( 1984 ) presents an analytical solution for the transient case incorporating the fracture skin. Typical type curves show that depending on the relative values of the hydraulic conductivity of the fracture skin, the draw down behavior resembles either the quasi-steady or unsteady state models. In both pseudo-steady and transient approaches the flow transfer term from blocks to fissures contains a convolution integral, similar to the case of "delayed yield" ( Boulton & Streltsova, 1978).

WATER TABLE IN HARD ROCK AQUIFERS

There is a widespread feeling among hydro geologists in India that ground water occurrence in hard areas is characterized by the presence of a water table in the shallow zone. It is also acknowledged that the fractures occur a little deeper, below the weathered rock zone. The water table in the weathered zone is influenced by recharge from rainfall or other sources and leakage to the fracture zone. In the double-porosity schematisation, on the block boundary which is not adjoining the fracture (upper boundary in Figure 1), a no flow condition is imposed ( eqn.4. ). If the weathered zone is looked upon as a block, There will be no fracture on the upper part of the weathered zone and hence the block schematisation as given above may not be valid. In such a case, a boundary condition corresponding to a declining water table under pumping, due to hydraulic connectivity with the fracture is required.

These concepts of an unconfined block overlying the fissure can be represented by a leaky aquifer system comprising of an aquifer-water table aquitard model. The aquifer corresponds to the fracture domain and the water table aquitard corresponds to the unconfined block domain. Cooley & Case (1973) first gave an approximate solution for the aquifer-water table aquitard problem by showing that Boulton's convolution integral can be used as the velocity of flow at the base of an incompressible water table aquitard overlying a pumped aquifer. In the present model, the effects of compressibility of the aquitard and declining water table are both considered.

AQUIFER - WATER TABLE AQUITARED SYSTEM

Conceptually, the aquifer-aquitard system is classical in ground water bydrology and so the physical meaning of the parameters associated with such models is well understood. On the contrary, this advantage does not exist for the fracture-block model. A comparison of the aquifer-water table aquitard system (Figure 2), with the fracture - block model (Figure 1), shows a difference in the domain of the upper boundary of the aquitard. The upper boundary of the aquitard is a water table instead of a no flow boundary as shown in Figure 1, and the water table is a moving boundary in response to the pumping from the aquifer.

The governing equations in aquifer and aquitard can be represented by the same fracture and block equations given by eqns. (1)-(4) earlier. However, the boundary condition at the upper boundary of the aquitard as given in eqn.(4) is no longer valid. Two boundary conditions are required on this boundary as the water table is unknown apriori, and these are given as follows.

Where ho and h are the initial saturated thickness and the transient water table elevation respectively (Figure 2); S¢ y is the specific yield of the aquitard.

An analytical solution for the governing equations incorporating eqn.(6), is not possible because of the unknown water table boundary at any instant of time. If one has to account for the compressibility of the aquitard and also satisfy the flux and varying head boundary conditions at the water table, one has to resort to numerical analysis. An efficient finite difference algorithm developed by Sridharan et al. (1900), was used here to obtain the draw downs in an aquifer-water table aquitard system. An iterative procedure is required in each time step to satisfy the continuity requirements at the aquifer-aquitard interface and the second boundary condition due to the varying water table elevation (eqn.6). The computational procedure at each iteration involves solving the aquifer draw downs at all the nodes in the radial direction and the aquitard drawdowns at all the nodes in the vertical direction, column by column (Sridharan et al.,1990). The finite difference scheme leads to a tridiagonal matrix equations for the aquifer as well as each column in the aquitard , which can be solved very conveniently by Thomas algorithm. The analysis algorithm is coupled to a parameter estimation algorithm based on weight least squares method (Sekhar et al., 1992). The portability of the software developed is enhanced by testing it on different computers ( PC, mini, mainframe).

The spatial and temporal variation of the water table aquitard and the aquifer drawdowns are governed by the following parameters in a non dimensional frame work :

In eqn.(7), W1 and W2 are the well functions; R, Z, T are the non dimensional radial distance ( r/ho) , non dimensional vertical distance in the aquitard (z/ho), and non dimensional time t/tr where tr = St h2o / Tr) respectively. S¢ t is the storage coefficient of the aquitard given by S¢ sh o. In the earlier study of Sridharan et al., ( 1990 ), the skin friction parameter SF is not considered.

A comparison of the functional relationships for the aquifer and aquitard drawdowns of eqn.(7), with the corresponding relationships for double-porosity model of Moench (1984) indicates the presence of an additional parameter, Sy/ St resulting from the decline of water table during pumping.

RESULTS AND DISCUSSION

Typical curves are presented for aquifer and water table drawdowns and comparison is made with type curves generated for fracture drawdown in an unsteady parallel fracture double-porosity model of Moench (1984). In order to generate type curves for the double-porosity model, suitable modifications in the upper boundary condition in aquitard is made in the aquifer-water table model, and drawdowns are computed numerically using the same computational procedure described in the previous section. An iterative procedure in each time step is required here for satisfying the continuity requirements at the aquifer-aquitard interface. However, the convergence was obtain more rapidly in this case in comparison to the aquifer-water table aquitard model.

Case 1 considers a variation of two orders in S¢ y / St ( = 10, 100, 1000 ) with values of S¢ t / St, R Ö K¢ ho/ Tr and SF kept constant at 100.0, 1.0 and 0.0 respectively (i.e., case 1 considers a situation without fracture skin ). Fig 3 presents the well function W1 for different S¢ y / St values. Corresponding to the negligible water table drawdown at very small time, there is no influence of specific yield of the aquitard ( S¢ y ) on the aquifer drawdown. As time increases and the water table drawdown builds up, influence of S¢ y is felt with the drawdown being lesser for larger S¢ y / St (and hence larger S¢ y ) values. It is clear that, as S¢ y / St is decreased, a matching with the double - porosity model having a no flow boundary in the block is obtained. The curves for different S¢ y / St values may not merge to a single curve even after very large times and this trend is a deviation from the double porosity model.

Case 2 considers the effect of fracture skin parameter, SF while considering the variation in S¢ y / St as in case 1. Figure 4 shows that within an increase in SF value ( low fracture skin hydraulic conductivity ), the aquifer drawdown at small time is large as expected and as time increases, the drawdown curve flattens due to small hydraulic head gradients in the block. However as time increases the influence of S¢ y is felt on the amplified drawdown, with the flat portion of the curve extending for a larger time for large S¢ y value. It is also observed that the magnitudes of the drawdown for Case 1 and Case 2 and at very large time are comparable for a particular S¢ y / St value ( Figure 3 & 4). In this case also, when S¢ y / St value is small, a good comparison could be obtained with a double porosity model considering no flow boundary in the block.

Case 3 considers both aquifer and water table drawdowns for a variation of one order in S¢ t / St ( =10,100) for a large value of S¢ y / St ( =1000). The RÖ K¢ ho/ Tr and SF kept constant at 1.0 and 0.1 respectively. Figure 5 shows a fairly significant influence of S¢ t / St on both aquifer and water table drawdowns. As S¢ t / St (and hence S¢ t) increases, both aquifer and water table drawdowns decrease. At very large times, the drawdowns are not significantly different for different S¢ t / St values, as might be expected. It is also observed that for small S¢ t / St , the shape of the aquifer drawdown curve is similar to that in pseudo-steady double- porosity model.

Case 4 considers aquifer and water table drawdowns for a variation in S¢ t / St as in the case 3, expect that the S¢ y / St parameter is considered at a low value of 10.0. It is clear from Figure 6 that the effect of low S¢ y / St is quite pronounced as the water table drawdowns, and unlike in Case 3, the drawndown curves for different S¢ t / St values are significantly separated after a very long time.

The flatness of aquifer drawdowns in well test data during a sufficiently long transition from early to late time has been a subject for theoretical justification of alternate approaches in the double porosity models. In the pseudo-steady state models, the drawdown would approach a plateau during the transition form early to late time. However, theoretical justification of this approach is not possible due to the assumption of negligible divergence of flow in the block. The transient model required the use of low permeability fracture skin to simulate such a flat behavior of the drawdown (Moench, 1994). The present aquifer-water table aquitard model shows the flat behavior of the drawdown for relatively large values of S¢ y / St even without consideration of a fracture skin. However, the flat behavior is observed at a later time due to the influence of water table, in comparison to a flat behavior in a relatively earlier time due to the effect of fracture skin. The coupled effect of both fracture skin and water table would be to have a flat behavior for a larger time.PARAMETER IDENTIFICATION

The aquifer-water table aquitard model developed in an equivalent isotropic domain can be used to estimate the anisotropic transmissivities Tx and Ty along the principal Cartesian directions x and y. In an anisotropic situation the radial co-ordinate in the equivalent isotropic domain is given by,


A weighted least squares approach based on sensitivity analysis for parameter estimation has been developed by Sekhar et al. (1992) which can estimate anisotropic transmissivities Tx and Ty if data at two observation wells are available and directions of principal axes are known.

The method is applied to a long duration pumping test in a hard rock area in the Haragonadona well field of Vedavati River Basin in India (Sridharan et al., 1990). A knowledge of hydrogeology of the region clearly indicates that there is an anisotropy, with the principal axes identified along north east-south west (x-axis) and north west-south east (y-axis) directions, with maximum transmissivity being in he y direction. Figure 7(a) presents relative location of pumping well (PW) and two observation wells (OW1 and OW2). The well is pumped at a constant discharge of 15 ls-1 for a duration of 7 days. The field observations of aquifer drawdowns and matching at both the observation wells for the estimated parameters are shown in Figure 7(b).

REFERENCES

1. Barenblatt, G.E., Zheltov, I.P., & Kochina, I.N., Basic concepts in the theory of seepage of homogeneous liquids in fissured rocks, J. Applied Mathematical Methods (USSR), 24 ,1286-1303 (1960).

2. Boulton, N.S., & Streltsova, T.D., Unsteady flow to a pumped well in an unconfined fissured aquifer, J.Hydrol., 37, 349-363. (1978).

3. Cooley, R.L., & & Case, C.M. Effect of water table aquitard on drawdwn in an underlying pumped aquifer, Water Resour. Res., 9,434-447 (1973).

4. Gringarten, A.C., Flow - test evaluation of fractured reservoirs, in Narasimhan, T.N., ed., Recent trends in hydrogeology, Special paper, Geological Society of America, 1982. pp 237-263.

Huyakorn, P.S., Lester, H., & Faust, C.R., Finite element techniques for modeling ground water flow in fractured aquifers, Water Resour. Res., 19, 1019-1035 (1983).

Kazemi, H., Leith, M.S., & Thomas, G.W., The interpretation of interference tests in naturally fractured reservoirs with uniform fracture distribution, Soc. Pet. Eng. J., 12, 463-472 (1969).

7. Moench, A.F., Double- porosity models for a fissured ground water reservoir with fracture skin, Water Resour, Res., 20, 831-846 (1984).

Sekhar, M., Mohan Kumar, M.S., & Sridharan, K., Parameter Estimation in an Aquifer-Water Table Aquitared System, J.Hydrol., 136, 177-192 (1992)

Sridharan. K., Sekhar, M., & Mohan Kumar, M.S., Analysis of Aquifer Water Table Aquitard System, J. Hydrol., 114, 175-189 (1990).

Warren, J.E., & Root, P.J., The behavior of naturally fractured reservior, Soc. Pet. Eng. J., 9 , 245-255 (1963).