- Email: [email protected]
Analytic solutions for drawdowns in rectangular artesian aquifers
Journal of Hydrology, 31 (1976) 151--160
151
© Elsevier Scientific Publishing Company, Amsterdam -- Printed in The Netherlands
ANALYTIC SOLUTIONS FOR DRAWDOWNS IN RECTANGULAR ARTESIAN AQUIFERS
Y.K. CHAN', N. M U L L I N E U X ~ and J.R. REED 2
' Department of Civil Engineering, University of Birmingham, Birmingham (Great Britain) Department of Mathematics, University of Aston, Birmingham (Great Britain) (Received December 4, 1975; accepted for publication January 29, 1976)
ABSTRACT Chan, Y.K., Mullineux, N. and Reed, J.R., 1976. Analytic solutions for drawdowns in rectangular artesian aquifers. J. Hydrol., 3 1 : 1 5 1 - - 1 6 0 . Analytical solutions for the drawdowns due to a single-well pumping at a constant rate from a rectangular homogeneous artesian aquifer are obtained using the method of integral transforms for a variety of boundary conditions. The computational effort is shown to be small in comparison with that needed in employing the method of images. The method can be used in problems in which the method of images fails.
INTRODUCTION
The classical paper by Theis (1935) paved the way for the evaluation of drawdown during pumping-test analysis in a homogeneous artesian non-leaky aquifer of infinite extent. Using the method of images as, for example, in the standard texts of Walton (1970) and Bear (1972), the analysis can be extended to finite rectangular aquifers. The method of application and its limitation were pointed out by Ferris et al. (1962) and Walton (1970). The partial removal of the limitation has been discussed by Chan (1976). For a wide range of differently shaped aquifers analytical solutions are possible using the properties of the Dirac 5-function and the method of integral transforms. In this paper the possibilities are indicated using only the simple finite Fourier transforms to derive solutions for rectangular aquifers with different sets of boundary conditions. These solutions are obtained in the form of infinite series and the computational aspects of summing these are discussed. NOTATION List of symbols s s* S T
= drawdown = transform of s = storage coefficient (dimensionless) = transmissibility
L
L2/T
152
t x,y
= time
m, n
= Cartesian coordinates = coordinates of abstraction well = rate of abstraction = difference in fixed piezometric levels at x = 0,a = d r a w d o w n in case r of Fig.1 at t i m e t = steady-state d r a w d o w n in case r, r ~ 3 = pseudo-steady-state d r a w d o w n in case 3 = dimensions of aquifer = integers
OLm
= ran~a;
hra Ign o(arn,X,~ )
= (n + ~ / b ;
Q h Sr(X,y,t) Sr(X,y , ** )
s~(x,y) a, b
c(~..y,n)
L
5n = n ~ / b ; ~'m,n2 = c~n + fj~ = ( m + ~ )n/a; 2 ,n = km 2 +~ Pm vm 2 ,n = sin(amX)Sin(c~m~ )
= k m2 + I ~
= cos(~nY)COS(#n~ ) = {1 i f m = n 0 ifmen
~m,rl
+(x -
T L L L3/T L
= ( 0 ifxifX ~ =~ together with
~)
(See A p p e n d i x A)
G O V E R N I N G E Q U A T I O N AND SOLUTIONS
The differential equation governing drawdown in a homogeneous artesian aquifer due to pumping at a constant rate from an abstraction well at the point (~.~) is usually taken to be (Walton, 1960): (a2s T \ax2 - -
.$ff
a2s~ ay'----:]
~
as
S
- -
at
(1)
in regions excluding the well. The presence of the well is then represented as a b o u n d a r y condition. An alternative discussed in Appendix A is to use the equation:
[a2+
a2s I
as
T \ax:'-- +~y2J =S---Qa(x-~)a(y-~)at
(2)
which reduces to eq. 1 at all points except (~,~) but now accounts for the pumping from the well at (~,~) w i t h o u t the introduction of a b o u n d a r y condition there. For a rectangular aquifer, with 0 < x < a and 0 < y < b, eq. 2 is readily solved by the application of finite integral transforms. The m e t h o d is outlined in Appendix B for the case of fixed piezometric levels s = 0 on x = 0, s = h on x = a and impermeable boundaries at which as/ay = 0 on y = 0,b. The initial condition assumes that:
s = xh/a
(3)
153 The complete result can be written in the form:
s(x,y,t)
2Q ~ exp(--T~nt/S) s(x,y,~*)---O(am,X,~) abT m=l ~
=
4Q ~ ~ exp(--T'Y~n,nt/8) abT m'~l n=l 7~n,n a(am,X,~)C(~n,y,~?) where the steady-state solution
xh
s(x,y,~)
=--
+
a
(4)
s(x,y,~) is given by:
Q_Q_~_~, O(~rn,X,~){ cosh[otm( b aT~'= t am sinh(em b)
Irl--yl)] +
+ cosh[am(b -- r? + y)] }
(5)
The numerical evaluation of these series is discussed in a later section.
Case 1 If the piezometric levels are the same on x = 0,a it is only necessary to set h to zero in eq. 4 to obtain the solution sl(x,y,t) to case 1 of Fig.1. The transforms required for the solution of the other five cases are summarised in Appendix B, and the results listed here. Y
t................... bV case i I Jtii~L~::-"
t _:,' r
I
I
L - -
'
t
. . . . . . . . . . . . . . .
',L. . . . .2. . .
]!i 1 ,
'
'
t 4
5
......... . . . . . . . . .
J
~
I~fllilllllllllJ ........
impermeable /i//,,-//z
3
ll
l//z/.///~
!
.....
t~
./
6 III
IIII/II
I~
(i_x## . . . . .
/ z / / z / / / / / / /
Fig. 1. The six aquifers considered.
Case 2 4Q ~ ~ exp(--TT~,nt/S) s2(x,y,t) = s2(x,y,~)-- - O(am,X,~)c([Jn,Y,rl) abTm=l n=l 7~n,n where the steady-state solution is:
154
s2(x,y, oo)
~-~Q~ °(am'X'~) tcosh[am(b-- b?--yl)] - " =1 amsinh(ctmb) 't
--c°sh[am(b ~ Case 3
+ Y)] t
Qt s3(x,y,t) = s3(x,y) + abS
2Q
~
exp(--T~n
a-~ y ,
t/S)
c~n
e(am'X'~)
2Q ~ exp(--T[J~t/S) abT ~ [32n e([Jn'Y'~7) /1=1
-abT
, =
7:m,n
where the pseudo-steady-state solution is:
bQ s3(x,y) -
PaT
32
[ r t - - y [ +rt + y + n: _ _+ y2 ) b b:
+ Q~ ~_j e(am'x'~) cosh[am(b -- b? - - y ( ) ] + cosh[am(b --n +Y-~] aT D I = I ctm sinh(ct m b) Case 4 4Q o.
S4(x,y,t) = S4(X,y,~)-- - -
m=O
o~ exp(--Tp~n,nt/S)
e(km,X,~ )o(f3n,Y,q )
n=l
where the steady-state solution is:
s4(x,Y, °°) = - ~Qn ~
~ -~cosh(~na) a([Jn'Y'~7) t sinh[~n(a -- I~ -- x I)] + sinh[~n(a -- ~ + x)] 1
Case 5
ss(x,y, t) = s s ( x , y , ~ )
4Q
exp(~TV~n, n
- --
abT m = o
rl=O
t/S:e(Xm)
,x,~
122m , ri
where the steady-state solution is:
Ss(X,y,~)_,.a~ ~ C(km'X'~) tsinh[km(b_ L17_yl) ] m=o Xmeosh(Xmb) + sinh[km(b--y
+ ~)] I
)e(u.,y,~ )
155
Case 6 oo
s6(x,y,t) = s6(x,y,°°) - - 2Q ~ e x p ( - - T ~ t / S ) c(km,X,~) abT m= o X2m 4Q ~ m~ ° exp(--Tp~,nt/S) C(Xm,X,~)c(~n,y,7?) abT n= 1 = P:m,n where the steady-state solution is:
Q~ c(~rn'X'~) tcosh[~rn(b-- ,~?--y,) ] s6(x'Y'°°) = - ~ =0 kmsinh(krnb) + cosh[km(b--~? + y)] I
COMPUTATIONAL
ASPECTS
The results appear, after inversion, in the form shown in eq. B-1. Here the steady-state part may be slow to converge, whilst the time-dependent part contains negative exponential terms which ensure rapid convergence. Accordingly the series in the steady-state part are summed analytically to produce at worst a single series with an acceptable rate of convergence. These single series (e.g., eq. 5) contain hyperbolic functions and to prevent numerical overflow these must be expressed in terms of negative exponentials. The only other point warranting discussion is the number of terms to be taken in the series. In practice as each term is generated a partial sum Sn say can be formed and the process terminated when [(Sn ÷m -- Sn)/Sn +m I is less than any prescribed level of accuracy. For purposes of illustration however, the case of an abstraction well at the centre of a rectangular aquifer of sides 9 × 5 km with observation well A near the centre at the point (3.6 km, 2.5 km) and another B near the b o u n d a r y at (0.45 km, 4.75 km) is considered. The d r a w d o w n and number of terms in the series necessary to achieve 1% accuracy are given in Tables I, (A) and (B) for observation wells A and B, respectively, at very short and moderate times with S = 0.00004, T = 2,250 m2/day, Q = 7,500 m3/day. It should be noted that at moderate times very few terms of the double series are required. More terms are required for the single series giving the steady-state solution b u t being time-independent this series is only summed once for each observation well. As expected more terms are required in summing the series for observation well A than for B. This arises more because of the deliberately chosen equality of the ordinates of well A and the abstraction well than because of their proximity. Nevertheless, comparatively few terms are required and the computation time is small.
156
TABLE I C o m p u t e d results (A) For observation well A at ~ = 3.6 km and ~ = 2.5 kin:
Case
Drawdown (m) t = 0.001 (day)
1 2 3 4 5 6
0.001634 0.001634 0.001634 0.001634 0.001634 0.001634
R a t i o t i m e images t i m e fin. t r a n s
No. o f t e r m s for 1% difference t = 10.0 (day)
1.407 0.6754 42.10 0.6884 1.070 3.207
double sum t = 0.001
t = 10.0
29 30 30 28 29 30
2 2 2 2 2 2
× 29 × 30 × 30 × 28 x 29 X 30
x x x x x x
2 2 2 2 2 2
single sum
• i t = 0.001
* : t = 10.0
120 245 110 6 180 50
5:1 4:1 5:1 8:1 5:1 7:1
900:1 500:1 1,100:1 18,000:1 700:1 2,500:1
(B) F o r o b s e r v a t i o n w e l l B a t ~ = 0 . 4 5 h m a n d ~ = 4 . 7 5 k i n :
Case
Drawdown (m) t = 0.1 (day)
1 2 3 4 5 6
0.0522 0.0048 0.1887 0.01439 0.1790 0.1886
No. o f t e r m s for 1% difference t = 10.0 (day)
-0.1479 0.0055 41.42 0.02027 0.5764 2.996
double sum t = 0.1
t = 10.0
2 2 2 2 2 2
2 2 2 2 2 2
× 2 × 2 × 2 X2 × 2 × 2
x x × × x ×
2 2 2 2 2 2
R a t i o t i m e images t i m e fin. t r a n s single sum
* I t = 0.1
*2t = 10.0
5 10 4 2 6 5
700:1 500:1 800:1 1,100:1 600:1 700:1
40,000:1 26,000:1 45,000:1 60,000:1 36,000:1 40,000:1
• 5 C h a n ( 1 9 7 6 ) uses a (9 x 9) image p a t t e r n . • 2 C h a n ( 1 9 7 6 ) uses a (65 x 65) image p a t t e r n . C O M P A R I S O N WITH T H E M E T H O D O F I M A G E S
The m e t h o d of images has been successfully employed to solve all the cases discussed here except case 3 and the case with different piezometric levels on t w o boundaries. In those cases in which a comparison is possible identical results are achieved using the method of integral transforms. Nevertheless, the m e t h o d of images has its limitations, as discussed b y Chan {1976), and its application is considerably more difficult than using the m e t h o d of integral transforms for which the program is readily prepared. Moreover unlike the m e t h o d of images the same program is applicable at all observation times. As the time is increased so the efficiency of the integral transform m e t h o d increases whereas the size of the image array grows rapidly, with
157 attendant increase in c o m p u t a t i o n time. As shown in Tables I(A) and (B), the c o m p u t a t i o n times are very greatly reduced except at unrealistically short observation times when the improvement factor falls to 4:1. CONCLUSION
The use of the 6-function enables the abstraction well (or wells) to be included in the governing equation which then applies over the whole of the aquifer. The resulting domain is then simpler to handle mathematically than when the well is excluded b y the introduction of a further boundary. As illustrated this opens up the possibility of obtaining analytical solutions to problems in which the m e t h o d of images would not be directly applicable. Once the 6-function is introduced the m e t h o d of integral transforms provides a well-established procedure for solving the differential equation. Whilst the results are derived for rectangular aquifers the m e t h o d can be extended to other geometries. APPENDIX A -- DISCUSSION OF THE GOVERNING EQUATION Eq. 2 is independent of the shape of the domain and to investigate its behaviour near the well it is convenient, for this discussion, to take the origin at the well so that eq. 2 becomes: -+ T \~x2 ay2]
S--
~t = - - Q ~ ( x ) 5 ( y )
(A-I)
The 5-functions are defined as:
tO,,
in such a w a y that with e,,e2 > O: X+C
2
X-e
l
.]
f(x)~ (x -
a)dx =
f(a)
provided that f(x) is b o u n d e d and continuous at x = a. The behaviour of the -function is interpreted through this integral relation. Thus on changing to cylindrical polars it is verified that eq. A-1 becomes: T-- -r ar
-- S . . . . 6 (r) ~t 2ur
(A-2)
The transformation of the left-hand side is well known; to check the equivalence of the right-hand sides it is only necessary to integrate over a small circle including the origin, then the integral relations give from eq. A-1:
ffQ6 (x)8 (y)dx dy = Q
158
and from eq. A-2: Q ~(r) f 2~r
2~rdr
Q
completing the verification. It is now shown t h a t eq. A-2 includes the usual representation of the well of infinitesimal diameter by the b o u n d a r y condition (Walton, 1970): lira r0-~ 0
0
-
(A-3) 2~T
To this end multiply eq. A-2 t h r o u g h o u t by r and integrate over the range (0,r0) to achieve:
T
(r
° sr dr
-- S r=ro
~ t Jo
= -- -02~
As ro -~ 0 the integral vanishes giving eq. A-3. This justifies the use of eq. 2 over the whole of the aquifer w i t h o u t the introduction of a b o u n d a r y surrounding the well. Consideration of continuity conditions over the region containing the idealised well leads to the same equation. APPENDIX B -- THE APPLICATION OF INTEGRAL TRANSFORMS The problem here is to solve eq. 2 subject to the b o u n d a r y conditions: at x = 0:
s = 0;
at y = 0,b:
~s/~y = 0
a t x = a;
s=h
and the initial condition: at t = 0:
s = xh/a
The b o u n d a r y conditions suggest (MuUineux et al., 1970), applying the finite sine transform WRT x and the finite cosine transform WRT y, i.e., multiplying eq. 2 t h r o u g h o u t by sin(amX)COS(~ny) and integrating over the aquifer to obtain: ds* T 2 _, Q dt-- + -S Tm'n8 = S-- sin(am~)COS(~nW) -- (--1)m a m h b T ~ n , o Solving this equation for s* and making use of the initial condition yields:
s* = Q T
t
1 -- e x p ( - - T T ~ , n t / S )
which inverts to:
Psin(am~)COS(~n~ ) ~2m,n
( - 1 ) m hb~n~o
am
159
s(x,y,t) -
abT =1
a em
4Q ~ ,= ~= tl--exp(--TTm,nt/S) f "l-a-~-T oo
2h ~
o(C~m,X,~)C([3n,Y,~l) ~2m,n
sin(amX)
(--1) m
a m=l
(B-l)
am
Using t h e k n o w n results (Mullineux et al., 1 9 7 0 ) : cos(n0 ) n= -~ oo
.
. m=l
n cosh {p(~ -- 0 ) }
n 2 + p2
p
sinh(p~)
(_l)m Sin(mO) 0 .
.
m
2
'
'
O
IOl
t h e s t e a d y - s t a t e p a r t o f eq. B-1 can b e r e d u c e d t o t h e f o r m given in eq. 5, so t h a t eq. B-1 itself r e d u c e s t o eq. 4. F o r t h e r e m a i n i n g five cases t h e r e q u i r e d t r a n s f o r m multipliers are given in T a b l e B-I. TABLE B-I Transform multipliers Case
Transform WRT x
Transform WRT y
2
sin(amX)
sin(~ny)
3 4 5 6
cos(amX ) cos(}.mX ) cos(}..,x) cos(~.mX)
cos(~ny ) sin(~ny ) cos(u,,y) cos(~ny)
In case 5 it is necessary t o use t h e result: c o s [ ( m + ½)0]_ n s i n h [ p ( ~ - - 0 ) ] m =-~ (m + _~)2 + p2 p cosh(pn)
0<0
<2;r
in o b t a i n i n g t h e s t e a d y - s t a t e solution.
REFERENCES Bear, J., 1972. Dynamics of Fluids in Porous Media. American Elsevier, New York, N.Y., 764 pp.
160
Chan, Y.K., 1976. Improved image-well technique for aquifer analysis. J. Hydrol., 29: 149--164. Ferris, J.G., Knowles, D.B., Brown, R.H. and Stallman, R.W., 1962. Theory of aquifer tests. U.S. Geol. Surv., Water-Supply Pap., 1536-E. Mullineux, N., Reed, J.R. and Soppitt, S.A., 1970. Present mathematics syllabuses and the analytical solution of field problems. Int. J. Electr. Eng. Educ., pp. 221 and 289. Theis, C.V., 1935. The relation between the lowering of piezometric surface and the rate and duration of discharge of a well using groundwater storage. Trans. Am. Geophys. Union, 16th Annu. Meet. Pt. 2. Walton, W.C., 1970. Groundwater Resources Evaluation. McGraw-Hill, Maidenhead, 664 pp.
151
© Elsevier Scientific Publishing Company, Amsterdam -- Printed in The Netherlands
ANALYTIC SOLUTIONS FOR DRAWDOWNS IN RECTANGULAR ARTESIAN AQUIFERS
Y.K. CHAN', N. M U L L I N E U X ~ and J.R. REED 2
' Department of Civil Engineering, University of Birmingham, Birmingham (Great Britain) Department of Mathematics, University of Aston, Birmingham (Great Britain) (Received December 4, 1975; accepted for publication January 29, 1976)
ABSTRACT Chan, Y.K., Mullineux, N. and Reed, J.R., 1976. Analytic solutions for drawdowns in rectangular artesian aquifers. J. Hydrol., 3 1 : 1 5 1 - - 1 6 0 . Analytical solutions for the drawdowns due to a single-well pumping at a constant rate from a rectangular homogeneous artesian aquifer are obtained using the method of integral transforms for a variety of boundary conditions. The computational effort is shown to be small in comparison with that needed in employing the method of images. The method can be used in problems in which the method of images fails.
INTRODUCTION
The classical paper by Theis (1935) paved the way for the evaluation of drawdown during pumping-test analysis in a homogeneous artesian non-leaky aquifer of infinite extent. Using the method of images as, for example, in the standard texts of Walton (1970) and Bear (1972), the analysis can be extended to finite rectangular aquifers. The method of application and its limitation were pointed out by Ferris et al. (1962) and Walton (1970). The partial removal of the limitation has been discussed by Chan (1976). For a wide range of differently shaped aquifers analytical solutions are possible using the properties of the Dirac 5-function and the method of integral transforms. In this paper the possibilities are indicated using only the simple finite Fourier transforms to derive solutions for rectangular aquifers with different sets of boundary conditions. These solutions are obtained in the form of infinite series and the computational aspects of summing these are discussed. NOTATION List of symbols s s* S T
= drawdown = transform of s = storage coefficient (dimensionless) = transmissibility
L
L2/T
152
t x,y
= time
m, n
= Cartesian coordinates = coordinates of abstraction well = rate of abstraction = difference in fixed piezometric levels at x = 0,a = d r a w d o w n in case r of Fig.1 at t i m e t = steady-state d r a w d o w n in case r, r ~ 3 = pseudo-steady-state d r a w d o w n in case 3 = dimensions of aquifer = integers
OLm
= ran~a;
hra Ign o(arn,X,~ )
= (n + ~ / b ;
Q h Sr(X,y,t) Sr(X,y , ** )
s~(x,y) a, b
c(~..y,n)
L
5n = n ~ / b ; ~'m,n2 = c~n + fj~ = ( m + ~ )n/a; 2 ,n = km 2 +~ Pm vm 2 ,n = sin(amX)Sin(c~m~ )
= k m2 + I ~
= cos(~nY)COS(#n~ ) = {1 i f m = n 0 ifmen
~m,rl
+(x -
T L L L3/T L
= ( 0 ifxifX ~ =~ together with
~)
(See A p p e n d i x A)
G O V E R N I N G E Q U A T I O N AND SOLUTIONS
The differential equation governing drawdown in a homogeneous artesian aquifer due to pumping at a constant rate from an abstraction well at the point (~.~) is usually taken to be (Walton, 1960): (a2s T \ax2 - -
.$ff
a2s~ ay'----:]
~
as
S
- -
at
(1)
in regions excluding the well. The presence of the well is then represented as a b o u n d a r y condition. An alternative discussed in Appendix A is to use the equation:
[a2+
a2s I
as
T \ax:'-- +~y2J =S---Qa(x-~)a(y-~)at
(2)
which reduces to eq. 1 at all points except (~,~) but now accounts for the pumping from the well at (~,~) w i t h o u t the introduction of a b o u n d a r y condition there. For a rectangular aquifer, with 0 < x < a and 0 < y < b, eq. 2 is readily solved by the application of finite integral transforms. The m e t h o d is outlined in Appendix B for the case of fixed piezometric levels s = 0 on x = 0, s = h on x = a and impermeable boundaries at which as/ay = 0 on y = 0,b. The initial condition assumes that:
s = xh/a
(3)
153 The complete result can be written in the form:
s(x,y,t)
2Q ~ exp(--T~nt/S) s(x,y,~*)---O(am,X,~) abT m=l ~
=
4Q ~ ~ exp(--T'Y~n,nt/8) abT m'~l n=l 7~n,n a(am,X,~)C(~n,y,~?) where the steady-state solution
xh
s(x,y,~)
=--
+
a
(4)
s(x,y,~) is given by:
Q_Q_~_~, O(~rn,X,~){ cosh[otm( b aT~'= t am sinh(em b)
Irl--yl)] +
+ cosh[am(b -- r? + y)] }
(5)
The numerical evaluation of these series is discussed in a later section.
Case 1 If the piezometric levels are the same on x = 0,a it is only necessary to set h to zero in eq. 4 to obtain the solution sl(x,y,t) to case 1 of Fig.1. The transforms required for the solution of the other five cases are summarised in Appendix B, and the results listed here. Y
t................... bV case i I Jtii~L~::-"
t _:,' r
I
I
L - -
'
t
. . . . . . . . . . . . . . .
',L. . . . .2. . .
]!i 1 ,
'
'
t 4
5
......... . . . . . . . . .
J
~
I~fllilllllllllJ ........
impermeable /i//,,-//z
3
ll
l//z/.///~
!
.....
t~
./
6 III
IIII/II
I~
(i_x## . . . . .
/ z / / z / / / / / / /
Fig. 1. The six aquifers considered.
Case 2 4Q ~ ~ exp(--TT~,nt/S) s2(x,y,t) = s2(x,y,~)-- - O(am,X,~)c([Jn,Y,rl) abTm=l n=l 7~n,n where the steady-state solution is:
154
s2(x,y, oo)
~-~Q~ °(am'X'~) tcosh[am(b-- b?--yl)] - " =1 amsinh(ctmb) 't
--c°sh[am(b ~ Case 3
+ Y)] t
Qt s3(x,y,t) = s3(x,y) + abS
2Q
~
exp(--T~n
a-~ y ,
t/S)
c~n
e(am'X'~)
2Q ~ exp(--T[J~t/S) abT ~ [32n e([Jn'Y'~7) /1=1
-abT
, =
7:m,n
where the pseudo-steady-state solution is:
bQ s3(x,y) -
PaT
32
[ r t - - y [ +rt + y + n: _ _+ y2 ) b b:
+ Q~ ~_j e(am'x'~) cosh[am(b -- b? - - y ( ) ] + cosh[am(b --n +Y-~] aT D I = I ctm sinh(ct m b) Case 4 4Q o.
S4(x,y,t) = S4(X,y,~)-- - -
m=O
o~ exp(--Tp~n,nt/S)
e(km,X,~ )o(f3n,Y,q )
n=l
where the steady-state solution is:
s4(x,Y, °°) = - ~Qn ~
~ -~cosh(~na) a([Jn'Y'~7) t sinh[~n(a -- I~ -- x I)] + sinh[~n(a -- ~ + x)] 1
Case 5
ss(x,y, t) = s s ( x , y , ~ )
4Q
exp(~TV~n, n
- --
abT m = o
rl=O
t/S:e(Xm)
,x,~
122m , ri
where the steady-state solution is:
Ss(X,y,~)_,.a~ ~ C(km'X'~) tsinh[km(b_ L17_yl) ] m=o Xmeosh(Xmb) + sinh[km(b--y
+ ~)] I
)e(u.,y,~ )
155
Case 6 oo
s6(x,y,t) = s6(x,y,°°) - - 2Q ~ e x p ( - - T ~ t / S ) c(km,X,~) abT m= o X2m 4Q ~ m~ ° exp(--Tp~,nt/S) C(Xm,X,~)c(~n,y,7?) abT n= 1 = P:m,n where the steady-state solution is:
Q~ c(~rn'X'~) tcosh[~rn(b-- ,~?--y,) ] s6(x'Y'°°) = - ~ =0 kmsinh(krnb) + cosh[km(b--~? + y)] I
COMPUTATIONAL
ASPECTS
The results appear, after inversion, in the form shown in eq. B-1. Here the steady-state part may be slow to converge, whilst the time-dependent part contains negative exponential terms which ensure rapid convergence. Accordingly the series in the steady-state part are summed analytically to produce at worst a single series with an acceptable rate of convergence. These single series (e.g., eq. 5) contain hyperbolic functions and to prevent numerical overflow these must be expressed in terms of negative exponentials. The only other point warranting discussion is the number of terms to be taken in the series. In practice as each term is generated a partial sum Sn say can be formed and the process terminated when [(Sn ÷m -- Sn)/Sn +m I is less than any prescribed level of accuracy. For purposes of illustration however, the case of an abstraction well at the centre of a rectangular aquifer of sides 9 × 5 km with observation well A near the centre at the point (3.6 km, 2.5 km) and another B near the b o u n d a r y at (0.45 km, 4.75 km) is considered. The d r a w d o w n and number of terms in the series necessary to achieve 1% accuracy are given in Tables I, (A) and (B) for observation wells A and B, respectively, at very short and moderate times with S = 0.00004, T = 2,250 m2/day, Q = 7,500 m3/day. It should be noted that at moderate times very few terms of the double series are required. More terms are required for the single series giving the steady-state solution b u t being time-independent this series is only summed once for each observation well. As expected more terms are required in summing the series for observation well A than for B. This arises more because of the deliberately chosen equality of the ordinates of well A and the abstraction well than because of their proximity. Nevertheless, comparatively few terms are required and the computation time is small.
156
TABLE I C o m p u t e d results (A) For observation well A at ~ = 3.6 km and ~ = 2.5 kin:
Case
Drawdown (m) t = 0.001 (day)
1 2 3 4 5 6
0.001634 0.001634 0.001634 0.001634 0.001634 0.001634
R a t i o t i m e images t i m e fin. t r a n s
No. o f t e r m s for 1% difference t = 10.0 (day)
1.407 0.6754 42.10 0.6884 1.070 3.207
double sum t = 0.001
t = 10.0
29 30 30 28 29 30
2 2 2 2 2 2
× 29 × 30 × 30 × 28 x 29 X 30
x x x x x x
2 2 2 2 2 2
single sum
• i t = 0.001
* : t = 10.0
120 245 110 6 180 50
5:1 4:1 5:1 8:1 5:1 7:1
900:1 500:1 1,100:1 18,000:1 700:1 2,500:1
(B) F o r o b s e r v a t i o n w e l l B a t ~ = 0 . 4 5 h m a n d ~ = 4 . 7 5 k i n :
Case
Drawdown (m) t = 0.1 (day)
1 2 3 4 5 6
0.0522 0.0048 0.1887 0.01439 0.1790 0.1886
No. o f t e r m s for 1% difference t = 10.0 (day)
-0.1479 0.0055 41.42 0.02027 0.5764 2.996
double sum t = 0.1
t = 10.0
2 2 2 2 2 2
2 2 2 2 2 2
× 2 × 2 × 2 X2 × 2 × 2
x x × × x ×
2 2 2 2 2 2
R a t i o t i m e images t i m e fin. t r a n s single sum
* I t = 0.1
*2t = 10.0
5 10 4 2 6 5
700:1 500:1 800:1 1,100:1 600:1 700:1
40,000:1 26,000:1 45,000:1 60,000:1 36,000:1 40,000:1
• 5 C h a n ( 1 9 7 6 ) uses a (9 x 9) image p a t t e r n . • 2 C h a n ( 1 9 7 6 ) uses a (65 x 65) image p a t t e r n . C O M P A R I S O N WITH T H E M E T H O D O F I M A G E S
The m e t h o d of images has been successfully employed to solve all the cases discussed here except case 3 and the case with different piezometric levels on t w o boundaries. In those cases in which a comparison is possible identical results are achieved using the method of integral transforms. Nevertheless, the m e t h o d of images has its limitations, as discussed b y Chan {1976), and its application is considerably more difficult than using the m e t h o d of integral transforms for which the program is readily prepared. Moreover unlike the m e t h o d of images the same program is applicable at all observation times. As the time is increased so the efficiency of the integral transform m e t h o d increases whereas the size of the image array grows rapidly, with
157 attendant increase in c o m p u t a t i o n time. As shown in Tables I(A) and (B), the c o m p u t a t i o n times are very greatly reduced except at unrealistically short observation times when the improvement factor falls to 4:1. CONCLUSION
The use of the 6-function enables the abstraction well (or wells) to be included in the governing equation which then applies over the whole of the aquifer. The resulting domain is then simpler to handle mathematically than when the well is excluded b y the introduction of a further boundary. As illustrated this opens up the possibility of obtaining analytical solutions to problems in which the m e t h o d of images would not be directly applicable. Once the 6-function is introduced the m e t h o d of integral transforms provides a well-established procedure for solving the differential equation. Whilst the results are derived for rectangular aquifers the m e t h o d can be extended to other geometries. APPENDIX A -- DISCUSSION OF THE GOVERNING EQUATION Eq. 2 is independent of the shape of the domain and to investigate its behaviour near the well it is convenient, for this discussion, to take the origin at the well so that eq. 2 becomes: -+ T \~x2 ay2]
S--
~t = - - Q ~ ( x ) 5 ( y )
(A-I)
The 5-functions are defined as:
tO,,
in such a w a y that with e,,e2 > O: X+C
2
X-e
l
.]
f(x)~ (x -
a)dx =
f(a)
provided that f(x) is b o u n d e d and continuous at x = a. The behaviour of the -function is interpreted through this integral relation. Thus on changing to cylindrical polars it is verified that eq. A-1 becomes: T-- -r ar
-- S . . . . 6 (r) ~t 2ur
(A-2)
The transformation of the left-hand side is well known; to check the equivalence of the right-hand sides it is only necessary to integrate over a small circle including the origin, then the integral relations give from eq. A-1:
ffQ6 (x)8 (y)dx dy = Q
158
and from eq. A-2: Q ~(r) f 2~r
2~rdr
Q
completing the verification. It is now shown t h a t eq. A-2 includes the usual representation of the well of infinitesimal diameter by the b o u n d a r y condition (Walton, 1970): lira r0-~ 0
0
-
(A-3) 2~T
To this end multiply eq. A-2 t h r o u g h o u t by r and integrate over the range (0,r0) to achieve:
T
(r
° sr dr
-- S r=ro
~ t Jo
= -- -02~
As ro -~ 0 the integral vanishes giving eq. A-3. This justifies the use of eq. 2 over the whole of the aquifer w i t h o u t the introduction of a b o u n d a r y surrounding the well. Consideration of continuity conditions over the region containing the idealised well leads to the same equation. APPENDIX B -- THE APPLICATION OF INTEGRAL TRANSFORMS The problem here is to solve eq. 2 subject to the b o u n d a r y conditions: at x = 0:
s = 0;
at y = 0,b:
~s/~y = 0
a t x = a;
s=h
and the initial condition: at t = 0:
s = xh/a
The b o u n d a r y conditions suggest (MuUineux et al., 1970), applying the finite sine transform WRT x and the finite cosine transform WRT y, i.e., multiplying eq. 2 t h r o u g h o u t by sin(amX)COS(~ny) and integrating over the aquifer to obtain: ds* T 2 _, Q dt-- + -S Tm'n8 = S-- sin(am~)COS(~nW) -- (--1)m a m h b T ~ n , o Solving this equation for s* and making use of the initial condition yields:
s* = Q T
t
1 -- e x p ( - - T T ~ , n t / S )
which inverts to:
Psin(am~)COS(~n~ ) ~2m,n
( - 1 ) m hb~n~o
am
159
s(x,y,t) -
abT =1
a em
4Q ~ ,= ~= tl--exp(--TTm,nt/S) f "l-a-~-T oo
2h ~
o(C~m,X,~)C([3n,Y,~l) ~2m,n
sin(amX)
(--1) m
a m=l
(B-l)
am
Using t h e k n o w n results (Mullineux et al., 1 9 7 0 ) : cos(n0 ) n= -~ oo
.
. m=l
n cosh {p(~ -- 0 ) }
n 2 + p2
p
sinh(p~)
(_l)m Sin(mO) 0 .
.
m
2
'
'
O
IOl
t h e s t e a d y - s t a t e p a r t o f eq. B-1 can b e r e d u c e d t o t h e f o r m given in eq. 5, so t h a t eq. B-1 itself r e d u c e s t o eq. 4. F o r t h e r e m a i n i n g five cases t h e r e q u i r e d t r a n s f o r m multipliers are given in T a b l e B-I. TABLE B-I Transform multipliers Case
Transform WRT x
Transform WRT y
2
sin(amX)
sin(~ny)
3 4 5 6
cos(amX ) cos(}.mX ) cos(}..,x) cos(~.mX)
cos(~ny ) sin(~ny ) cos(u,,y) cos(~ny)
In case 5 it is necessary t o use t h e result: c o s [ ( m + ½)0]_ n s i n h [ p ( ~ - - 0 ) ] m =-~ (m + _~)2 + p2 p cosh(pn)
0<0
<2;r
in o b t a i n i n g t h e s t e a d y - s t a t e solution.
REFERENCES Bear, J., 1972. Dynamics of Fluids in Porous Media. American Elsevier, New York, N.Y., 764 pp.
160
Chan, Y.K., 1976. Improved image-well technique for aquifer analysis. J. Hydrol., 29: 149--164. Ferris, J.G., Knowles, D.B., Brown, R.H. and Stallman, R.W., 1962. Theory of aquifer tests. U.S. Geol. Surv., Water-Supply Pap., 1536-E. Mullineux, N., Reed, J.R. and Soppitt, S.A., 1970. Present mathematics syllabuses and the analytical solution of field problems. Int. J. Electr. Eng. Educ., pp. 221 and 289. Theis, C.V., 1935. The relation between the lowering of piezometric surface and the rate and duration of discharge of a well using groundwater storage. Trans. Am. Geophys. Union, 16th Annu. Meet. Pt. 2. Walton, W.C., 1970. Groundwater Resources Evaluation. McGraw-Hill, Maidenhead, 664 pp.