- Email: [email protected]
Horizontal wells in shallow aquifers: Field experiment and numerical model
Journal of Hydrology (2006) 329, 98– 109
available at www.sciencedirect.com
journal homepage: www.elsevier.com/locate/jhydrol
Horizontal wells in shallow aquifers: Field experiment and numerical model Azuhan Mohamed 1, Ken Rushton
*
School of Civil Engineering, University of Birmingham, Birmingham B15 2TT, UK Received 10 October 2005; received in revised form 30 January 2006; accepted 3 February 2006
KEYWORDS
Summary This paper is concerned with horizontal wells in shallow aquifers; information about a field test is presented and the development of appropriate conceptual and computational models is discussed. Particular attention is paid to the processes of flow from the aquifer into the horizontal well and the hydraulic conditions within the well. A combined computational model representing the aquifer and the horizontal well is developed. This methodology is used to study a horizontal well in a shallow coastal aquifer; a detailed pumping test is analysed using the combined model. The combined model is then used to predict the consequences of increasing the abstraction, lengthening the well or using the horizontal well for an extended period when there is no recharge. The results from the combined model indicate that the commonly adopted assumptions of either a uniform flux into the well or a constant head along the well are not suitable for this field problem. c 2006 Elsevier B.V. All rights reserved.
Horizontal well; Pipe flow; Shallow aquifer; Conceptual models; Computational models; Finite differences
Introduction
Horizontal wells can have significant advantages over vertical wells in shallow aquifers because they can be positioned towards the bottom of the aquifer. Furthermore, the provision of one long horizontal well rather than a number of vertical well clusters reduces the number of pumps required. * Corresponding author. Address: 106 Moreton Road, MK18 1PW Buckingham, UK. Tel.: +44 1280 824 942. E-mail addresses: [email protected] (A. Mohamed), [email protected] (K. Rushton). 1 Present address: Economic Planning Unit, Prime Minister’s Department, Level 5, Block B6, 62502 Putrajaya, Malaysia.
For horizontal wells, any analysis of their time-variant response should include the following three component, the flow processes within the aquifer, flow from the aquifer into the horizontal well and thirdly the hydraulic conditions within the horizontal well pipe. A number of analytical solutions have been developed for horizontal wells but they do not include all of these components. For most analytical solutions the horizontal well is represented either as a constant head boundary or as a uniform flux boundary; no account is taken of hydraulic conditions within the well. For example Kawecki (2000), using approaches developed for the oil industry, considers a number of alternative approaches and presents approximate analytical expressions for drawdowns due to both constant head and uniform flux
0022-1694/$ - see front matter c 2006 Elsevier B.V. All rights reserved. doi:10.1016/j.jhydrol.2006.02.006
Horizontal wells in shallow aquifers: Field experiment and numerical model conditions at the well face; the influence of a skin factor around the well is also considered. Analytical solutions have also been developed by Zhan et al. (2001), Zhan and Zlotnik (2002) and Zhan and Park (2003); the horizontal well is treated as a line sink and the flux distribution along the well is assumed uniform. A variety of aquifer conditions are represented including leaky aquifers with and without aquitard storage; conditions on the upper boundary of the flow field include confining strata, water table conditions including delayed yield and an overlying reservoir. The issue of the distribution of flux on the well bore is considered by Cassiani and Kabala (1998) and by Zhan and Zlotnik (2002). Cassiani and Kabala (1998) present analytical solutions for vertical wells and explain that, although drawdowns are not sensitive to flux non-uniformity, the flux distribution along the well screen is non-uniform. Zhan and Zlotnik (2002) present analytical solutions for horizontal and slanted wells. They also refer to a numerical simulation of a hypothetical example with the length of the well 1.7 times the unconfined aquifer thickness with either uniform head or uniform flux well bore conditions. Differences in drawdowns from the two solutions are less than 10% at distances greater than five times the well diameter from the end of the well. The importance of hydraulic conditions in the horizontal well was recognised by Tarshish (1992) who considers steady-state conditions in a horizontal well in an aquifer underlain by an impermeable base and overlain by a water reservoir; boundaries are represented using image wells. The well is simulated as a linear sink with a non-uniform inflow. Flow within the well is represented taking into account energy losses and the momentum flux change along the well. To simulate losses due to the high approach velocities into the well, a well loss coefficient can be introduced (Kruseman and de Ridder, 1990). Chen et al. (2003) present a new treatment of in-well hydraulics for horizontal wells which avoids the assumption that the well is a constant head or a uniform inflow internal boundary. The horizontal well is represented as an equivalent hydraulic conductivity with head losses along the well. Five alternative flow regimes in the well bore are identified, a laminar regime with the head loss proportional to the velocity u, smooth turbulent flow with the head loss proportional to u1.75 and a rough turbulent regime with the loss proportional to u2.0; transitions between these regimes are also included. Finite difference solutions are developed; comparisons are made with a physical laboratory model which represents pumping from beneath a river. Chen et al. (2003) recognise that additional well losses occur due to radial flow through the slots of the well screen. Furthermore the radial flow direction into the well changes to axial at the inside wall of the pipe reducing the cross-section occupied by pure axial flow within the pipe. No corrections are made for these processes. A field test of a horizontal well in a shallow aquifer in Sarawak, Malaysia provides a focus for this study. Features which must be incorporated in the analysis of this problem include significant variations in the pumping rate from the well, the occurrence of recharge during the test, substantial falls in the water table elevation leading to a reduction in the saturated thickness, also the specification of maximum water table elevations equal to the ground surface eleva-
99
tion. The resistance to flow as the water enters the well (well loss) is an important consideration as is the estimation of hydraulic head losses in the horizontal well pipes. These complex features can be represented by three interconnected numerical models. The three components mechanisms of time-variant flow in a shallow aquifer to a horizontal well, flow from the aquifer into the horizontal well and flow within the horizontal well pipe are examined. For each component both conceptual and computational models are proposed; detailed information is presented about the computational models. Due to uncertainties about flow processes close to the well including the flows through the slots in the wall of the pipe and the resultant reduction in the area of axial flow within the horizontal pipe, coefficients are introduced to represent these features. These coefficients are adjusted during the model refinement to obtain an adequate match with field observations. In this paper the methodologies refer specifically to horizontal wells in shallow aquifers and to the analysis of a field pumping tests. However, they can be adapted for horizontal wells in deeper aquifers and horizontal wells beneath water bodies.
Field example of horizontal well The development of conceptual and computational models is strongly influenced by the intended application. Therefore a representative field problem of a horizontal well in a shallow aquifer is summarised in this section, conceptual and computational models are introduced in ‘‘Conceptual and computational models’’ section with a detailed study of the field problem in ‘‘Application to horizontal well at Loba’’ section. A horizontal well system has been installed and tested at Loba, Sarawak, Malaysia. Loba is a small coastal village, not far from Kuching, on the South China Sea. The site is suitable for tourist development since a shallow sand aquifer can provide good quality water. The sand aquifer extends about 700 m from the coast, Fig. 1(a); the base of the sand aquifer is about 4 m below mean sea level with a maximum
Figure 1 Representative cross-section and construction of horizontal well at Loba.
100 elevation of about 1.5 m above mean sea level. During extended periods of heavy rainfall the water table reaches the ground surface. Extensive fieldwork showed no evidence of the ingress of saline water, the low salinity water in the aquifer extends right up to the coast (Mailvaganam et al., 1993). The only access to the site during the construction and testing of the well was by water. The shallow sand aquifer lies above clay; apart from the immediate vicinity of the coast the quality of the water is good. The anticipated long-term requirement for water is 240 m3 d1 to cater for tourist development. Following a preliminary examination of the potential yield of vertical and horizontal wells it was concluded that a single horizontal well system feeding to a central caisson could provide the required yield, but the reliability of the system for long periods without recharge was unknown (Mailvaganam et al., 1993). The horizontal well, Fig. 1(b), consists of a central caisson 1.2 m in diameter; on each side of the caisson a pair of 150 mm diameter PVC pipes are provided with slots cut in the upper halves of the pipes. Datum is taken as the top of the clay. The slotted horizontal pipes were placed above the underlying clay with a filter sand pack extending for about 60–80 mm. Vertical tubes are located at 25 m intervals along the length of the horizontal well to permit cleaning of the well; they also allow measurements of hydraulic heads in the well during pumping. About 20 observation piezometers, at distances of between 2 and 100 m from the horizontal well, were also provided. Dewatering of the aquifer was necessary during construction of the horizontal well. The suction pumps required for dewatering were used for pumping tests.
A. Mohamed, K. Rushton
Regional groundwater flow Fig. 2(a) represents a shallow aquifer. The aquifer is of limited saturated thickness with an effectively impermeable base, hence the vertical flow crossing the lower boundary is negligible. The lateral extent of the study area is defined by boundaries on which either the groundwater head or the cross-boundary flow is defined. The upper boundary is the water table; recharge q occurs at this boundary, in addition water is released or taken into storage due to the specific yield SY. Water from the aquifer is discharged via the horizontal well or to the sea. Both horizontal and vertical flows occur within the aquifer. However, it is not necessary to use the differential equation in terms of x, y, z and t coordinates since the differential equation for regional groundwater flow in terms of the groundwater head h, which is a function of the lateral coordinates x and y and the time t, o oh o oh o oh o oh Tx þ Ty ¼ Kx h þ Ky h ox ox oy oy ox ox oy oy oh ð1Þ ¼ SY q ot does implicitly represent vertical flow components (Connorton, 1985; Pinder, 2002; Rushton, 2003). In this equation the transmissivities in the x and y directions equal the hydraulic conductivities Kx and Ky multiplied by the saturated thickness h (datum is taken as the impermeable base of
Conceptual and computational models The three component mechanisms which need to be considered when studying a horizontal well in a shallow aquifer are: (1) groundwater flow in the shallow aquifer with appropriate inflows, outflows, initial and boundary conditions, (2) flows from the aquifer into the horizontal well, (3) hydraulic conditions within the horizontal well. Diagrams of these three components are included as Fig. 2(a)–(c). For each component, conceptual models are developed and idealisations introduced so that the flow processes can be described mathematically. One possible approach for groundwater flow would be to represent the aquifer in three-dimensional space with the cross-sectional dimensions of the horizontal well simulated using a very fine mesh. However, this approach fails to allow for the complex flow processes in the vicinity of the well (for instance only the upper half of the horizontal well pipe contains slots) and the uncertain head losses through the wall of the well. Therefore an alternative procedure is adopted whereby the interaction between the aquifer and the well is represented using an empirical coefficient similar to the river coefficient used in MODFLOW (Anderson and Woessner, 1992). Furthermore, because the aquifer is relatively uniform and of limited saturated thickness, a single layer regional groundwater model is selected. The idealisations incorporated in the numerical models are influenced by the need to combine the computational model of each component to provide an overall model which represents the interactions within the total flow system.
Figure 2 Three component mechanisms when representing a horizontal well in a shallow aquifer.
Horizontal wells in shallow aquifers: Field experiment and numerical model the aquifer which is assumed to be horizontal). The implicit inclusions of vertical flow components in Eq. (1) can be appreciated by considering the governing equation for groundwater flow in terms of the coordinate directions x, y, z, integrating in the z direction and applying the boundary conditions of zero flow on the lower impermeable boundary and an inflow on the upper boundary (water table) of q SYoh/ot. Rigorous derivations can be found in Connorton (1985) and Pinder (2002). If there is an extensive and continuous zone within the aquifer having a hydraulic conductivity less than 1% of the aquifer hydraulic conductivity, it is necessary to use a multi-layer model (Rushton and Rathod, 1985). In the Loba study area the sand layer has been deposited under a sequence of wave dominated environments, a sequence of boreholes on a 100 m grid showed no evidence of clay layers or lenses. The successful use of Eq. (1) to include vertical flow components is illustrated in Appendix A by the example of a horizontal drain. There are numerous computational models available to solve Eq. (1); for instance MODFLOW automatically allows for the variable saturated thickness of the unconfined aquifer. In this study an orthogonal finite difference mesh is used with variable intervals in x and y directions. A schematic diagram of the aquifer mesh is shown in the upper part of Fig. 3; groundwater heads are functions of the coordinates x and y, hence subscripts m and n are used to indicate the location. It is also necessary to specify conditions on the lateral boundaries together with any inflows such as recharge, any natural outflows (to springs, rivers and the sea) and any man-made outflows such as to the horizontal well. In addition, initial conditions must be specified.
Flow from aquifer to the horizontal well Fig. 2(b) is a conceptual model indicating how water is drawn from the aquifer water table, through the aquifer to the horizontal well. One possible approach is to represent the actual cross-section of the horizontal well in the groundwater model using a very fine grid. However, there
Figure 3 Schematic diagram of the computational model showing part of the regional aquifer, the representation of flow from the aquifer to nodal points on the horizontal well using coefficient Ca and the increasing flows in the horizontal well.
101
are additional processes such as the effect of gravel packs, any disturbances to the aquifer during the construction of the well and the resistance to the flow of water through slots in the well pipes; these processes are difficult to quantify. Therefore an alternative approach is adopted representing all these features as an empirical coefficient. In many field studies, flow from an aquifer to a river, drain or vertical well is described by a relationship of the form Qa ¼ Caðh HÞ
ð2Þ
where Qa is the flow from the aquifer to the river, drain or well over a specified length; h is the groundwater head at the river, drain or well; H is the hydraulic head in the river, drain or well; Ca is a semi-empirical coefficient which depends on the dimensions of the river, drain or well, the hydraulic conductivity of the aquifer, any material in the river bed, drainage envelope or gravel pack and the resistance to flow through the perforations in the drain or slotted well casing. In MODFLOW, the river coefficient can be estimated in terms of the dimensions and properties of the deposits in the riverbed (Anderson and Woessner, 1992). However in most practical situations adjustments to this first estimate are made during sensitivity analyses. To explore the validity of the use of Eq. (1) to represent horizontal and vertical flow components in a shallow aquifer, and Eq. (2) to represent flows from the aquifer to a horizontal well, a related problem of groundwater flow to horizontal tile drains is considered in the Appendix A. Comparisons are made between a detailed vertical section two-dimensional time-variant analysis (profile model) of horizontal and vertical flow components from a moving water table to a drain and a simulation of drain–aquifer interaction using Eq. (2) with the aquifer represented using a one-dimensional form of Eq. (1). The discharge from the drain is specified. From information about the variation of the hydraulic head in the drain with time, the magnitude of the aquifer–drain coefficient Ca is adjusted during a sensitivity analysis. Groundwater heads in the aquifer and hydraulic heads in the drain determined using this approach show acceptable agreement with the full two-dimensional moving water table analysis, Fig. A.1. Two conclusions can be drawn from the acceptable agreement between these two alternative approaches. Firstly, the one-dimensional formulation of Eq. (1) does simulate the time-variant horizontal and vertical flow components in a shallow homogeneous aquifer system. Secondly, the interaction between the drain and aquifer can be represented by the coefficient Ca and Eq. (2); however, the magnitude of Ca has to be adjusted during a sensitivity analysis to match the hydraulic heads in the drain. For a horizontal well, since the hydraulic head varies along the horizontal pipe, the computational model for the interaction of a horizontal well with a shallow aquifer is represented as a series of discrete inflows at mesh points along the line of the well. The inflow associated with node m is indicated in Fig. 3. This inflow occurs from midway between nodes m 1 and m to midway between nodes m and m + 1, a distance Dx. The inflow over this length of pipe Qam ¼ Cam ðhm;n Hm Þ
ð3Þ
102
A. Mohamed, K. Rushton
where Cam is the aquifer-pipe coefficient for this length of pipe. Note that there are two suffices defining the location of the groundwater head hm,n since it is a function of x and y whereas the hydraulic head in the pipe Hm is only a function of x. Inflows at nodes m 1 and m + 1 are calculated in a similar manner.
Flow within horizontal well Fig. 2(c) illustrates inflows from the aquifer through the perforations in a horizontal well pipe and the increasing flows within the pipe to the discharge point. These conditions differ from those for conventional pipe theory since flows occur through the wall of the pipe causing velocities within the pipe to increase from zero at the blank end to a maximum at the discharge point. A computational model is required to represent the increasing flow along the pipe. Since the aquifer is modelled using a finite-difference mesh, the horizontal pipe is divided into discrete lengths Dx, Fig. 3. Turbulent flows often occur in pipes. For the field example discussed in ‘‘Application to horizontal well at Loba’’ section the maximum velocity in the pipe is about three times the critical ‘Reynolds velocity’, which is the maximum velocity for which the flow remains laminar. Velocities are lower away from the outlet but turbulent flows are likely to occur along the full length of the pipe due to the impact of the inflowing water through the pipe walls. There are several alternative mathematical expressions to describe pipe flow; for this analysis the Hazen-Williams expression is used (Daugherty et al., 1989), although other pipe-flow relationships could be selected, DH ¼ RQp1:85
ð4Þ
where DH is the head loss in the pipe (m), Qp is the flow rate in the pipe (m3 d1), and the coefficient of pipe resistance R¼
3:2 106 L 4:87 C1:85 H d
ð5Þ
Figure 4 well.
Detailed diagram of conditions in the horizontal
In Eq. (8) the coefficient Cpm1,m depends on the unknown flow Qpm1,m. Similarly Qpm;mþ1 ¼ Cpm;mþ1 ðHm Hmþ1 Þ
Continuity along the pipe, Fig. 4, requires that the pipe flow between nodes m 1 and m plus the inflow from the aquifer equals the pipe flow between nodes m and m + 1, Qpm1;m þ Qam ¼ Qpm;mþ1
Cpm1;m Hm1 ðCpm1;m þ Cam þ Cpm;mþ1 ÞHm þ Cpm;mþ1 Hmþ1 ¼ Cam hm;n
The manner in which the three component mechanisms are included in the total computational model is illustrated in Fig. 3; the flow chart of Fig. 5 indicates how the calculations proceed. For each time step the following procedure is adopted; note that the components are considered in the reverse order to ‘‘Regional groundwater flow’’, ‘‘Flow from aquifer to the horizontal well’’ and ‘‘Flow within horizontal well’’ sections so that they are in the same order as the flow chart. Input parameters and initial conditions 1. Solve simultaneous equns for flows in pipe, Eq. (10) to estimate Qp and H
Cpm1;m ¼
4
1
2. Calculate flows Qa from Eq. (3) 3. Solve finite difference form of Eq. (1) with appropriate 3 boundary conditions to calculate h(x, y)
2
iterative routine for Kh Repeat with updated h(x, y)
ð7aÞ
ð8Þ
.
new Cp values from Eq. (8)
ð6Þ
Next time step
where Hm1 and Hm are the hydraulic heads, whilst the pipe coefficient Cpm1,m is given by the equation 1 ðRQp0:85 m1;m Þ
ð10Þ
Combining the three component mechanisms
The term in square brackets is treated as a pipe coefficient Cp which is updated during iterative calculations. Cp is lower towards the outlet where the flow Qp increases. The pipe is divided into segments of length Dx with the hydraulic heads H defined at the nodal points, m 1, m and m + 1 but the coefficients Cp and the flows Qp are defined between nodal points as shown in Fig. 4. Hence Eq. (6) becomes Qpm1;m ¼ Cpm1;m ðHm1 Hm Þ
ð9Þ
From Eqs. (3), (7a), (7b) and (9)
In this study the pipe is divided into segments of length L = Dx; for each segment the Hazen-Williams discharge coefficient CH and the pipe diameter d are the same, hence R is unchanged for each segment. For computational purposes it is helpful to adopt the approach of McManus and Rushton (1964) by rewriting Eq. (4) as Qp ¼ DH½ðRQp0:85 Þ1 ¼ DHCp
ð7bÞ
Calculation completed
Figure 5
Flow chart for the computational model of Fig. 3.
Horizontal wells in shallow aquifers: Field experiment and numerical model Component 3: for each node along the horizontal well an equation of the form of Eq. (10) is written. This leads to a set of simultaneous equations which can be solved using Gaussian elimination. However, the coefficients Cp depend of the unknown flows in the pipe, Eq. (8), consequently an iterative technique is required, loop of Fig. 5. The right hand side of Eq. (10) involves the groundwater head hm,n. Since values of groundwater heads for the current time step have not yet been determined, it is necessary to return to Component 3 when new values of the groundwater head have been calculated. Component 2: Eq. (3) is used to calculate flows to the horizontal well from the aquifer. Component 1: groundwater heads within the aquifer are calculated from the finite-difference form of Eq. (1) with the appropriate boundary conditions and with outflows from the aquifer to the horizontal well as determined from Component 2 above. Since the transmissivities equal the hydraulic conductivity multiplied by the unknown saturated thickness, an iterative solution of the finite-difference equations is required as indicated by loop in Fig. 5; this is automatically performed by MODFLOW. After the first round of calculations for hydraulic heads in the well, groundwater heads in the aquifer and flows between aquifer and the horizontal well, the calculations are repeated using the improved values of the heads and flows until appropriate convergence criteria are met; this is indicated by iterative loop as shown in Fig. 5. For each time step the entire procedure is repeated as indicated by loop . Although preliminary estimates can be made of the coefficients Ca and Cp, they need to be refined to represent actual field conditions. Since the coefficient Ca, which is used to calculate flows from the aquifer to nodal points on the horizontal well, is an empirical coefficient, improved estimates are obtained during a sensitivity analysis which focuses on the difference between groundwater heads and hydraulic heads in the pipe. The pipe coefficient Cp depends on the coefficient of pipe resistance R and the current flow in the segment of pipe, Eq. (6). Initial estimates of R can be obtained from Eq. (5), but modifications to the Hazen-Williams discharge coefficient CH are necessary due to the effect of the slots in the wall of the pipe and the inflows through these slots. Therefore a further sensitivity analysis is required to estimate a suitable value for the coefficient of pipe resistance R to reproduce field readings of the difference in hydraulic head along the pipe. These procedures are discussed further in ‘‘Computational model for Loba’’ section.
103
flows estimated from the time taken to fill a container, variations in the pumping rate did occur. Pump breakdowns also occurred during the test. Estimates of the recharge intensity and pumping rate for periods of 6 hours are included in Fig. 6(a). Field readings of drawdowns in two of the aquifer observation piezometers, in the caisson and in the horizontal well are indicated by the discrete symbols in Figs. 6(b) and 8; all measurements are relative to datum which is 4.0 m below mean sea level. The estimated tolerance in the observed readings is ±1 cm due to uncertainties in the datum for measurement, the limited experience of the operators and the fact that some readings were taken at night. Measurements of the hydraulic head in the horizontal well are also influenced by adjustments to the pumping rate from the caisson.
Computational model for Loba The groundwater model for the shallow sand aquifer at Loba extends for a distance of 700 m perpendicular to the coast and 2600 m parallel to the coast, Fig. 7(a). The horizontal well is located in the centre of the modelled area and is parallel to the coast. The mesh spacing is indicated in Fig. 7(a); the minimum spacing of 5.0 m perpendicular to the well is too small to be identified on the diagram. A mesh spacing of 25 m is used along the horizontal well hence there are 13 nodes along the 300 m length of the horizontal well. Recharge and discharge data are input at intervals of 0.25 d, two time steps of 0.125 d are used for each data input. For the aquifer, the boundary conditions are a specified groundwater head of 4.0 m above datum along the
Application to horizontal well at Loba Pumping tests in the horizontal well Several pumping tests were carried out on the horizontal well at Loba; this discussion will concentrate on the first six days of a 14-day test. During the first six days there was only limited rainfall; thereafter, substantial rainfall led to total recovery in many of the observation piezometers. Since four dewatering pumps were used for the pumping test with gate valves at the discharge outlets and with
Figure 6 Pumping test at Loba: (a) assumed recharge and pumping rate, (b) field readings and modelled results, x and y axes defined in Fig. 7.
104
A. Mohamed, K. Rushton
Figure 7
Computational model for Loba.
coastline, in this shallow aquifer the presence of any saline interface can be ignored. With a transmissivity of 50 m2/d and a storage coefficient of 0.033, the tidal effect spreads no more than 70 m from the coast; since all the observation wells and the horizontal well are more than 250 m from the coastline, tidal effects are not significant. No-flow conditions are applied on the other three boundaries to represent the aquifer extending about 700 m perpendicular to the coast with the effect of pumping extending no more than 1 km on either side of the horizontal well. Initial conditions are based on the observation that maximum groundwater heads approximately coincide with the general shape of the land surface. Therefore initial heads are set equal to a smoothed land surface profile; the simulation starts with a period of 6 days without recharge or pumping to allow any inconsistencies in the initial heads to dissipate. Iterative routines are required both for the representation of hydraulic conditions in the horizontal well, loop of Fig. 5, and for the regional groundwater model, loop of Fig. 5. For the computational model of the horizontal well, improved values of the coefficients Cp are obtained from iterative loop ; convergence is achieved in four cycles. For the groundwater model the iterative calculation continues until the out of balance of flows at every node is equivalent to less than 0.005 mm d1. After three cycles of iterative loop , convergence is achieved. In the absence of information to the contrary, the hydraulic conductivity and specific yield are assumed to be uniform over the aquifer. The aquifer-pipe coefficient Ca is also assumed to take a constant value at all nodal points on the horizontal well. In addition, the coefficient
of pipe resistance R needs to be refined to allow for the effect of the inflows of water along the horizontal well. Using a sensitivity analysis, trial parameter values are explored; the criteria used to estimate suitable values of the parameters are listed in Table 1. Estimates of the aquifer hydraulic conductivity and specific yield were obtained from earlier studies (Mailvaganam et al., 1993); therefore initial trial solutions concentrated on adjustments to the aquifer-pipe coefficient Ca and the coefficient of pipe resistance R. An increase in R from the initial estimate (based on Eq. (5) for a smooth straight pipe, 25 m in length) by a factor of 10 was required to reflect the disturbances to flow cause by the slots and the water flowing through the slots into the pipe. Reductions in Ca, from a first estimate of 65 m2 d1 (deduced from Eq. (2) and field information about flows and head differences) were required to lead to greater differences between groundwater and hydraulic head. After achieving a reasonable simulation, further refinements were made to each of the parameters in turn. The final parameter values are K ¼ 12 m d1 ; SY ¼ 0:033 ðvalue is influenced by fine sand near ground surfaceÞ Ca ¼ 56 m2 d1 for a 25 m length of pipe, R ¼ 2:0 106 d1:85 m4:55 for 25 m of pipe It would be possible to use an alternative expression to Hazen-Williams with the loss proportional to Qp2; slightly different values of the pipe coefficient would be obtained but this would not have a significant effect on the overall model results.
Horizontal wells in shallow aquifers: Field experiment and numerical model Table 1
105
Parameters refined during the sensitivity analysis
Parameter
Key criterion
Hydraulic conductivity K
Groundwater head gradient towards the horizontal well depends primarily on K, reductions in K increase the gradient The overall decline in groundwater heads is strongly influenced by SY, a lower SY leads to greater drawdowns Difference between groundwater heads and hydraulic heads at the horizontal well Variation of hydraulic head along pipes of horizontal well, reduce R to increase hydraulic head difference along pipe
Specific yield SY Aquifer-pipe coefficient Ca Coefficient of pipe resistance R
Computational model results Fig. 6(b) includes comparisons between field readings (the discrete symbols) and computational model results (continuous lines) for the hydraulic head in the caisson and groundwater heads in observation piezometers at distances of 10 m and 30 m perpendicular to the horizontal well. The agreement is encouraging. Small differences, which occur in the measured groundwater heads, can be due to the tolerances in measurements discussed in ‘‘Pumping tests in the horizontal well’’ section and local variations in aquifer parameters such as varying elevations of the clay base to the aquifer. Changes in hydraulic heads in the caisson are also represented satisfactorily by the model. Since the effect of pump failure is represented as occurring throughout the time step between 3.75 and 4.0 days, there are small discrepancies in the timing of this recovery. Field and modelled hydraulic head distributions within the horizontal well at 3.5 and 6.0 days are compared in Fig. 8. The discrete symbols represent the hydraulic heads measured in the cleaning tubes of the horizontal well, the continuous line is constructed through the numerical model results for hydraulic head H at 25 m intervals along the horizontal well. Hydraulic heads in the caisson are sensitive to the pumping rate at the time when the reading was taken since the four dewatering pumps did not maintain a steady discharge. The difference in hydraulic head between the ends of the pipe and the caisson depends on the choice of the coefficient of pipe resistance R; the value selected gives an adequate simulation for each of the times included in Fig. 8. From the computational model results, Fig. 9 is prepared which shows both calculated inflows from the aquifer to the horizontal well at 25 m intervals and flows within the pipe towards the caisson. Three times are selected to illustrate
Figure 8
Hydraulic heads in horizontal well.
a range of pumping rates. Fig. 9(a) refers to a time of 1.25 d when 230 m3 d1 is pumped from the horizontal well; the vertical arrows, which indicate the flows from the aquifer into the well, are almost uniform. However when the pumping rate falls to 160 m3 d1, Fig. 9(b), inflows from the aquifer towards the centre of the horizontal well reduce significantly. This occurs due to the substantial decrease in hydraulic head drawdown at the caisson compared to the much smaller change in groundwater heads as shown in Fig. 6(b). An increase in total pumping rate to 280 m3 d1 also leads to a non-uniform inflow distribution to the horizontal well, Fig. 9(c). For these three times, substantial changes occur in flow in the horizontal pipe with flows to the caisson varying form 75 to 129 m3 d1. Results from Figs. 8 and 9 indicate that the horizontal well should not be represented either as a constant head or as a uniform inflow. Hydraulic heads at the caisson are typically 0.05 m lower than at the ends of the well; this variation is significant compared to the differences between groundwater head and hydraulic head of 0.25–0.40 m. Furthermore, inflows from the aquifer to the horizontal well are sensitive to any variations in the pumping rate; the maximum inflow is often 1.3 times and occasionally 1.5 times the minimum inflow.
Predictive simulations Since the combined numerical model provides an adequate representation of field conditions in the aquifer and in the horizontal well, it can be used for predictive simulations. Three predictive simulations are described. For the first predictive simulation the pumping rates are set at double the historical values. It is essential to consider whether modifications should be made to any of the three component models due to this increase in abstraction rate. For the groundwater component, decreases in saturated thickness are included in Eq. (1). The second model refers to flows from the aquifer into the horizontal well, Fig. 2(b); due to increased approach velocities at the horizontal well there is likely to be an increased resistance to flow which can be represented by decreasing Ca. Substantial increases in the well loss factor occur for vertical wells following increases in abstraction rates (Rushton, 2003); consequently the value of Ca is halved from the value deduced from the historical value (this is equivalent to doubling the well loss factor). For the third component model relating to flow in the horizontal pipe, Fig. 2(c), increased inflows through the wall of the pipe result in higher values of the pipe resistance R. In the absence of field information, R is doubled.
A. Mohamed, K. Rushton aquifer to pipe flow (m3 d-1 )
106 25 20
(a) 1.25 day
230 m3 d-1 on day 1.25 19.2
15
17.9
17.2
37
54
aquifer to pipe flow (m 3 d-1 ) aquifer to pipe flow (m 3 d-1 )
17.8
17.6
17.6
17.3
88
71
17
17.2
17.9
37
19
19.2
10 5 0 19
-150 25
71
-100
88
106 106
-50
0
50
54
100
150
(b) 3.875 day
20 15
160 m3 d-1 on day 3.875
16.2 13.6
10
12.1
11.3
10.8
10.6
10.7
10.6
10.8
11.3
12.1
42
53
64
75
75
64
53
42
30
16.2 13.6
5 0 16
30
16
(c) 4.50 day
25
280 m3 d-1 on day 4.50
25.1 22.6
20
21.1
20.4
20.1
68
89
109
20.3
20.9
20.3
129
129
109
20.1
20.4
21.1
89
68
48
25.1 22.6
15 10 5 0 25
-150
Figure 9
17.3
17
48
-100
-50
0
50
100
25
150
Calculated inflows to horizontal well (vertical arrows) and flow through horizontal well at specific times, units m3 d1.
When these modifications are included in the simulation the predicted maximum drawdown in the caisson is 2.76 m compared to 0.99 m for the historical simulation (see Fig. 8). From the end of the horizontal well to the caisson, the difference in hydraulic head is 0.24 m compared to 0.034 m for the historical simulation. This means that the well would be operating under a pressure head of about 2.0 m, hence this increase in pumping rate appears to be feasible. Nevertheless, since the calculated drawdowns in the caisson are based on assumed values of Ca and R, any proposed increase in the pumping rate should be confirmed by field-testing. In the second predictive simulation, the pumping rates are again double the historical values but the horizontal pipes are also doubled in length; Ca and R are returned to their original values. Lengthening the horizontal wells leads to a reduction in the water table drawdowns above the well. The hydraulic head at the caisson is also reduced to 1.79 m compared to 2.76 m for the first predictive run. However, the hydraulic head differences along each well become 0.62 m compared to 0.24 m for the first predictive run. This
results in a variation in the fluxes into the well with 140% of the average flux at the caisson and a minimum of 75% of the average flux towards the end of the well. These findings indicate that a lengthening of the well reduces the water table drawdown but the well becomes less efficient at collecting water at locations distant from the caisson. In a third predictive simulation the long-term yield of the aquifer system is examined. The pumping rate is set at the design value of 240 m3 d1, the simulation continues for 100 days without any recharge occurring (this would be an exceptional situation for a humid tropical climate). Since inflows into the horizontal well remain of similar magnitude to the historical simulation, an unchanged value of R is appropriate. It is possible that Ca will also remain unchanged but in case some increase in the resistance to flow does occur in the vicinity of the horizontal well, Ca is reduced to two-thirds of the value for the historical simulation. In Fig. 10, groundwater heads are plotted on a crosssection perpendicular to the horizontal well extending from the coast to the no-flow boundary. After 30 d the groundwa-
Horizontal wells in shallow aquifers: Field experiment and numerical model
Figure 10 Predictive response due to pumping for 100 days with no recharge.
ter head above the horizontal well falls to 4.0 m which coincides with sea level. However there is a groundwater mound between the horizontal well and the coast with a maximum head of 4.32 m; consequently fresh water continues to flow to the coast. By day 62 this mound has disappeared with a very small groundwater head gradient from the coast towards the horizontal well. At day 100 the gradient from the coast in the direction of the horizontal well becomes more significant. The groundwater head gradient from the coast at day 100 is 0.0014. When multiplied by the hydraulic conductivity of 12.0 m d1, this gives a Darcy velocity of 0.017 m d1. If the effective porosity is 0.10, the seepage velocity for the ingress of saline water equals 0.17 m d1. Consequently saline water would move a few metres into the aquifer; once substantial recharge occurs this saline water would be flushed back to the coast. This simulation suggests that in the very unlikely circumstances of no recharge occurring for 100 days, the horizontal well would continue to provide good quality water with only a small ingress of saline water at the coast.
Conclusions This paper describes a methodology for analysing the response of horizontal wells in shallow aquifers without the assumptions either of a uniform hydraulic head in the well or a uniformly distributed inflow from the aquifer along the length of the well. A computational model is introduced which represents the three components of flow in the aquifer, the hydraulic conditions in the horizontal well, and the interaction between the aquifer and the well. This model is used to study the field situation of a horizontal well in a shallow costal aquifer. From a sensitivity analysis, values are deduced for aquifer parameters and coefficients describing the well-aquifer interaction and the resistance to flow in the pipe. Satisfactory agreement is achieved between field and modelled values of both hydraulic heads and groundwater heads. The model is then used to predict the consequences of increasing the pumping rate, lengthening the well, or operating the existing well for 100 days without recharge.
107
Certain of the references quoted in the Introduction suggest that it may be acceptable to assume a uniform flux into the well. From the limited number of simulations described in this paper relating to a specific location, it appears that the uniform flux assumption may be valid for wells of limited length but for longer wells (and drains) the hydraulic head losses are significant leading to substantial variations in the flux. Since the aim of the current study is to develop numerical models which can be used to represent horizontal wells without any assumptions about the distribution of heads or fluxes, no firm conclusions can be drawn about situations where the uniform flux assumption is acceptable. The field study and subsequent mathematical modelling demonstrate the suitability of horizontal wells for abstracting substantial quantities of water from shallow aquifers. Although the cost of constructing a horizontal well is greater than for a number of vertical wells, the advantages include the need for only one pump and the ability to increase the discharge during periods of high demand. Furthermore, the horizontal well can continue to operate efficiently during periods of zero recharge. For an aquifer with a greater saturated thickness, a multi-layered regional groundwater model is required. Due to the uncertainties about the flow processes close to the horizontal well and the difficulties in representing the precise geometry of the well cross-section in a regional groundwater model, Eq. (3) should be used to represent the aquifer-well interaction with the coefficient Ca deduced during sensitivity analyses. Hydraulic conditions in the pipe should be simulated using Eq. (7). Alternative conditions on the upper and lateral boundaries can be included in the regional groundwater model. The methodology described in this paper is also useful when studying conditions in interceptor drains which are pumped horizontal drains positioned to intercept water lost from irrigation canals (Rushton, 2003).
Appendix A. Representation of drain–aquifer interaction in a shallow aquifer The purpose of this appendix is to explore the adequacy of the representation of flow in a shallow aquifer using a single layer model (‘‘Regional groundwater flow’’ section) and also representation of conditions in the vicinity of the horizontal well using an aquifer-well coefficient (‘‘Flow from aquifer to the horizontal well’’ section). This is achieved by considering the comparable problem of flow to horizontal tile drains. A study by Khan and Rushton (1996) examines time-variant flow in a vertical section to horizontal tile drains from which constant discharges are withdrawn. Fig. A.1(a) illustrates a two-dimensional vertical section from which the discharge from each drain is Qd = 0.125 m3 d1 per metre length of the drain. Since the tile drains are equally spaced at 50 m, lateral no-flow boundaries are located at 25 m from the drain centreline. The drains are 0.1 m in diameter and 2.0 m above the aquifer base. Initial conditions are that the water table and all the groundwater heads equal 4.5 m relative to the base of the aquifer. The differential equation describing the flow in a vertical section in terms of the groundwater head h(y, z, t) is
108
A. Mohamed, K. Rushton H ¼ h Qd=Ca
Figure A.1 Comparison of two-dimensional vertical section model with one-dimensional solution in which aquifer–drain interaction is represented by a drain coefficient: (a) twodimensional formulation, (b) results from two-dimensional solution and one-dimensional approximation, (c) one-dimensional formulation showing shape of water table.
o oh o oh oh Ky þ Kz ¼ SS oy oy oz oz ot
ðA:1Þ
in which SS is the specific storage coefficient. Khan and Rushton (1996) describe numerical solutions to Eq. (A.1) using time-variant finite-difference approximations in which both horizontal and vertical flow components are represented; the upper boundary is a moving water table which is a function of the specific yield. The continuous lines in Fig. A.1(b) show the calculated water table elevations midway between drains, at the quarter points and immediately above the drain; the head in the drain is also plotted. Aquifer parameter values are listed in Fig. A.1(b). The same problem is analysed using a one-dimensional approximation of Eq. (1) in which the groundwater head h(y, t), Fig. A.1(c), is a function of the horizontal coordinate y and the time t. The governing equation is o oh oh Kh ¼ SY ðA:2Þ oy oy ot where SY is the specific yield. The drain is represented as a constant discharge Qd per unit length of the drain; the hydraulic head in the drain is calculated by rearranging Eq. (2) as
ðA:3Þ
where Ca, the aquifer–drain coefficient, is adjusted during a sensitivity analysis to match the hydraulic head in the drain. A best fit is obtained with Ca = 0.32 m2 d1 m1; this coefficient applies for all times. Numerical solutions to the problem of the drain with a constant discharge are obtained using the finite-difference form of Eq. (A.2) with a mesh spacing of 0.5 m and a time step of 0.5 d. Since in Eq. (A.2) the hydraulic conductivity is multiplied by the unknown saturated thickness h, iterative solutions are required for each time step. In Fig. A.1(b), results from this numerical solution are indicated by the discrete symbols, results are plotted for alternate time steps. Water table elevations at the lateral boundaries and at the quarter points show no difference between the results for the two-dimensional and one-dimensional approximations. For the water table above the drain, the one-dimensional solution predicts elevations slightly lower than the two-dimensional solution. An examination of Fig. A.1(a) shows that for the two-dimensional solution the water table above the drain is horizontal. However, for the one-dimensional approximation, Fig. A.1(c), there is a discontinuity in the groundwater head because an outflow occurs from the drain node, consequently the water table is not horizontal. The hydraulic head in the drain calculated from the simulation based on Eq. (A.2) is similar to that determined from the two-dimensional simulation, although there are small differences for early and late times. These findings support the proposition that groundwater flow through a shallow aquifer towards a horizontal drain can be studied using a single layer approximation with the aquifer–drain interaction represented by a coefficient which is determined during a sensitivity analysis.
References Anderson, M.P., Woessner, W.W., 1992. Applied Groundwater Modeling: Simulation of Flow and Advective Transport. Academic Press, San Diego, 381pp. Cassiani, G., Kabala, Z.J., 1998. Hydraulics of a partially penetrating well: solution to a mixed-type boundary value problem via dual integration equations. J. Hydrol. 211, 100–111. Chen, C., Wan, J., Zhan, H., 2003. Theoretical and experimental studies of coupled seepage-pipe flow to a horizontal well. J. Hydrol. 281, 159–171. Connorton, B.J., 1985. Does the regional groundwater-flow equation model vertical flow. J. Hydrol. 79, 279–299. Daugherty, R.L., Franzini, J.B., Finnemore, E.J., 1989. Fluid Mechanics with Engineering Applications, SI Metric Edition. McGraw-Hill, Singapore, 596pp. Kawecki, M.W., 2000. Transient flow to a horizontal water well. Ground Water 38, 842–850. Khan, S., Rushton, K.R., 1996. Reappraisal of flow to tile drains: III Drains with limited flow capacity. J. Hydrol. 183, 383–395. Kruseman, G.P., de Ridder, N.A., 1990. Analysis and Evaluation of Pumping Test Data, second ed. International Institute for Land Reclamation and Improvement, Wageningen, The Netherlands, Publ. no. 47, 377pp. Mailvaganam, Y., Ramili, M.Z., Rushton, K.R., Ong, B.Y., 1993. Groundwater Exploitation of a Shallow Coastal Sand Aquifer in Sarawak, Malaysia. IAHS Publ. no. 216, pp. 451–462.
Horizontal wells in shallow aquifers: Field experiment and numerical model McManus, B.D., Rushton, K.R., 1964. A pure resistance electrical analogue for the analysis of pipe networks. Civ. Eng. Pub. Works Rev. London, 585–589. Pinder, G.F., 2002. Groundwater Modeling Using Geographical Information Systems. Wiley, New York, 233pp. Rushton, K.R., 2003. Groundwater Hydrology: Conceptual and Computational Models. Wiley, Chichester, 416pp. Rushton, K.R., Rathod, K.S., 1985. Horizontal and vertical flow components deduced from groundwater heads. J. Hydrol. 79, 261–278.
109
Tarshish, M., 1992. Combined mathematical model of flow in an aquifer-horizontal well system. Ground Water 30, 931–935. Zhan, H., Park, E., 2003. Horizontal well hydraulics in leaky aquifers. J. Hydrol. 281, 129–146. Zhan, H., Zlotnik, V.A., 2002. Groundwater flow to a horizontal or slanted well in an unconfined aquifer. Water Resour. Res. 38 (7), 13-1–13-11. Zhan, H., Wang, L.V., Park, E., 2001. On the horizontal-well pumping tests in anisotropic confined aquifers. J. Hydrol. 252, 37–50.
available at www.sciencedirect.com
journal homepage: www.elsevier.com/locate/jhydrol
Horizontal wells in shallow aquifers: Field experiment and numerical model Azuhan Mohamed 1, Ken Rushton
*
School of Civil Engineering, University of Birmingham, Birmingham B15 2TT, UK Received 10 October 2005; received in revised form 30 January 2006; accepted 3 February 2006
KEYWORDS
Summary This paper is concerned with horizontal wells in shallow aquifers; information about a field test is presented and the development of appropriate conceptual and computational models is discussed. Particular attention is paid to the processes of flow from the aquifer into the horizontal well and the hydraulic conditions within the well. A combined computational model representing the aquifer and the horizontal well is developed. This methodology is used to study a horizontal well in a shallow coastal aquifer; a detailed pumping test is analysed using the combined model. The combined model is then used to predict the consequences of increasing the abstraction, lengthening the well or using the horizontal well for an extended period when there is no recharge. The results from the combined model indicate that the commonly adopted assumptions of either a uniform flux into the well or a constant head along the well are not suitable for this field problem. c 2006 Elsevier B.V. All rights reserved.
Horizontal well; Pipe flow; Shallow aquifer; Conceptual models; Computational models; Finite differences
Introduction
Horizontal wells can have significant advantages over vertical wells in shallow aquifers because they can be positioned towards the bottom of the aquifer. Furthermore, the provision of one long horizontal well rather than a number of vertical well clusters reduces the number of pumps required. * Corresponding author. Address: 106 Moreton Road, MK18 1PW Buckingham, UK. Tel.: +44 1280 824 942. E-mail addresses: [email protected] (A. Mohamed), [email protected] (K. Rushton). 1 Present address: Economic Planning Unit, Prime Minister’s Department, Level 5, Block B6, 62502 Putrajaya, Malaysia.
For horizontal wells, any analysis of their time-variant response should include the following three component, the flow processes within the aquifer, flow from the aquifer into the horizontal well and thirdly the hydraulic conditions within the horizontal well pipe. A number of analytical solutions have been developed for horizontal wells but they do not include all of these components. For most analytical solutions the horizontal well is represented either as a constant head boundary or as a uniform flux boundary; no account is taken of hydraulic conditions within the well. For example Kawecki (2000), using approaches developed for the oil industry, considers a number of alternative approaches and presents approximate analytical expressions for drawdowns due to both constant head and uniform flux
0022-1694/$ - see front matter c 2006 Elsevier B.V. All rights reserved. doi:10.1016/j.jhydrol.2006.02.006
Horizontal wells in shallow aquifers: Field experiment and numerical model conditions at the well face; the influence of a skin factor around the well is also considered. Analytical solutions have also been developed by Zhan et al. (2001), Zhan and Zlotnik (2002) and Zhan and Park (2003); the horizontal well is treated as a line sink and the flux distribution along the well is assumed uniform. A variety of aquifer conditions are represented including leaky aquifers with and without aquitard storage; conditions on the upper boundary of the flow field include confining strata, water table conditions including delayed yield and an overlying reservoir. The issue of the distribution of flux on the well bore is considered by Cassiani and Kabala (1998) and by Zhan and Zlotnik (2002). Cassiani and Kabala (1998) present analytical solutions for vertical wells and explain that, although drawdowns are not sensitive to flux non-uniformity, the flux distribution along the well screen is non-uniform. Zhan and Zlotnik (2002) present analytical solutions for horizontal and slanted wells. They also refer to a numerical simulation of a hypothetical example with the length of the well 1.7 times the unconfined aquifer thickness with either uniform head or uniform flux well bore conditions. Differences in drawdowns from the two solutions are less than 10% at distances greater than five times the well diameter from the end of the well. The importance of hydraulic conditions in the horizontal well was recognised by Tarshish (1992) who considers steady-state conditions in a horizontal well in an aquifer underlain by an impermeable base and overlain by a water reservoir; boundaries are represented using image wells. The well is simulated as a linear sink with a non-uniform inflow. Flow within the well is represented taking into account energy losses and the momentum flux change along the well. To simulate losses due to the high approach velocities into the well, a well loss coefficient can be introduced (Kruseman and de Ridder, 1990). Chen et al. (2003) present a new treatment of in-well hydraulics for horizontal wells which avoids the assumption that the well is a constant head or a uniform inflow internal boundary. The horizontal well is represented as an equivalent hydraulic conductivity with head losses along the well. Five alternative flow regimes in the well bore are identified, a laminar regime with the head loss proportional to the velocity u, smooth turbulent flow with the head loss proportional to u1.75 and a rough turbulent regime with the loss proportional to u2.0; transitions between these regimes are also included. Finite difference solutions are developed; comparisons are made with a physical laboratory model which represents pumping from beneath a river. Chen et al. (2003) recognise that additional well losses occur due to radial flow through the slots of the well screen. Furthermore the radial flow direction into the well changes to axial at the inside wall of the pipe reducing the cross-section occupied by pure axial flow within the pipe. No corrections are made for these processes. A field test of a horizontal well in a shallow aquifer in Sarawak, Malaysia provides a focus for this study. Features which must be incorporated in the analysis of this problem include significant variations in the pumping rate from the well, the occurrence of recharge during the test, substantial falls in the water table elevation leading to a reduction in the saturated thickness, also the specification of maximum water table elevations equal to the ground surface eleva-
99
tion. The resistance to flow as the water enters the well (well loss) is an important consideration as is the estimation of hydraulic head losses in the horizontal well pipes. These complex features can be represented by three interconnected numerical models. The three components mechanisms of time-variant flow in a shallow aquifer to a horizontal well, flow from the aquifer into the horizontal well and flow within the horizontal well pipe are examined. For each component both conceptual and computational models are proposed; detailed information is presented about the computational models. Due to uncertainties about flow processes close to the well including the flows through the slots in the wall of the pipe and the resultant reduction in the area of axial flow within the horizontal pipe, coefficients are introduced to represent these features. These coefficients are adjusted during the model refinement to obtain an adequate match with field observations. In this paper the methodologies refer specifically to horizontal wells in shallow aquifers and to the analysis of a field pumping tests. However, they can be adapted for horizontal wells in deeper aquifers and horizontal wells beneath water bodies.
Field example of horizontal well The development of conceptual and computational models is strongly influenced by the intended application. Therefore a representative field problem of a horizontal well in a shallow aquifer is summarised in this section, conceptual and computational models are introduced in ‘‘Conceptual and computational models’’ section with a detailed study of the field problem in ‘‘Application to horizontal well at Loba’’ section. A horizontal well system has been installed and tested at Loba, Sarawak, Malaysia. Loba is a small coastal village, not far from Kuching, on the South China Sea. The site is suitable for tourist development since a shallow sand aquifer can provide good quality water. The sand aquifer extends about 700 m from the coast, Fig. 1(a); the base of the sand aquifer is about 4 m below mean sea level with a maximum
Figure 1 Representative cross-section and construction of horizontal well at Loba.
100 elevation of about 1.5 m above mean sea level. During extended periods of heavy rainfall the water table reaches the ground surface. Extensive fieldwork showed no evidence of the ingress of saline water, the low salinity water in the aquifer extends right up to the coast (Mailvaganam et al., 1993). The only access to the site during the construction and testing of the well was by water. The shallow sand aquifer lies above clay; apart from the immediate vicinity of the coast the quality of the water is good. The anticipated long-term requirement for water is 240 m3 d1 to cater for tourist development. Following a preliminary examination of the potential yield of vertical and horizontal wells it was concluded that a single horizontal well system feeding to a central caisson could provide the required yield, but the reliability of the system for long periods without recharge was unknown (Mailvaganam et al., 1993). The horizontal well, Fig. 1(b), consists of a central caisson 1.2 m in diameter; on each side of the caisson a pair of 150 mm diameter PVC pipes are provided with slots cut in the upper halves of the pipes. Datum is taken as the top of the clay. The slotted horizontal pipes were placed above the underlying clay with a filter sand pack extending for about 60–80 mm. Vertical tubes are located at 25 m intervals along the length of the horizontal well to permit cleaning of the well; they also allow measurements of hydraulic heads in the well during pumping. About 20 observation piezometers, at distances of between 2 and 100 m from the horizontal well, were also provided. Dewatering of the aquifer was necessary during construction of the horizontal well. The suction pumps required for dewatering were used for pumping tests.
A. Mohamed, K. Rushton
Regional groundwater flow Fig. 2(a) represents a shallow aquifer. The aquifer is of limited saturated thickness with an effectively impermeable base, hence the vertical flow crossing the lower boundary is negligible. The lateral extent of the study area is defined by boundaries on which either the groundwater head or the cross-boundary flow is defined. The upper boundary is the water table; recharge q occurs at this boundary, in addition water is released or taken into storage due to the specific yield SY. Water from the aquifer is discharged via the horizontal well or to the sea. Both horizontal and vertical flows occur within the aquifer. However, it is not necessary to use the differential equation in terms of x, y, z and t coordinates since the differential equation for regional groundwater flow in terms of the groundwater head h, which is a function of the lateral coordinates x and y and the time t, o oh o oh o oh o oh Tx þ Ty ¼ Kx h þ Ky h ox ox oy oy ox ox oy oy oh ð1Þ ¼ SY q ot does implicitly represent vertical flow components (Connorton, 1985; Pinder, 2002; Rushton, 2003). In this equation the transmissivities in the x and y directions equal the hydraulic conductivities Kx and Ky multiplied by the saturated thickness h (datum is taken as the impermeable base of
Conceptual and computational models The three component mechanisms which need to be considered when studying a horizontal well in a shallow aquifer are: (1) groundwater flow in the shallow aquifer with appropriate inflows, outflows, initial and boundary conditions, (2) flows from the aquifer into the horizontal well, (3) hydraulic conditions within the horizontal well. Diagrams of these three components are included as Fig. 2(a)–(c). For each component, conceptual models are developed and idealisations introduced so that the flow processes can be described mathematically. One possible approach for groundwater flow would be to represent the aquifer in three-dimensional space with the cross-sectional dimensions of the horizontal well simulated using a very fine mesh. However, this approach fails to allow for the complex flow processes in the vicinity of the well (for instance only the upper half of the horizontal well pipe contains slots) and the uncertain head losses through the wall of the well. Therefore an alternative procedure is adopted whereby the interaction between the aquifer and the well is represented using an empirical coefficient similar to the river coefficient used in MODFLOW (Anderson and Woessner, 1992). Furthermore, because the aquifer is relatively uniform and of limited saturated thickness, a single layer regional groundwater model is selected. The idealisations incorporated in the numerical models are influenced by the need to combine the computational model of each component to provide an overall model which represents the interactions within the total flow system.
Figure 2 Three component mechanisms when representing a horizontal well in a shallow aquifer.
Horizontal wells in shallow aquifers: Field experiment and numerical model the aquifer which is assumed to be horizontal). The implicit inclusions of vertical flow components in Eq. (1) can be appreciated by considering the governing equation for groundwater flow in terms of the coordinate directions x, y, z, integrating in the z direction and applying the boundary conditions of zero flow on the lower impermeable boundary and an inflow on the upper boundary (water table) of q SYoh/ot. Rigorous derivations can be found in Connorton (1985) and Pinder (2002). If there is an extensive and continuous zone within the aquifer having a hydraulic conductivity less than 1% of the aquifer hydraulic conductivity, it is necessary to use a multi-layer model (Rushton and Rathod, 1985). In the Loba study area the sand layer has been deposited under a sequence of wave dominated environments, a sequence of boreholes on a 100 m grid showed no evidence of clay layers or lenses. The successful use of Eq. (1) to include vertical flow components is illustrated in Appendix A by the example of a horizontal drain. There are numerous computational models available to solve Eq. (1); for instance MODFLOW automatically allows for the variable saturated thickness of the unconfined aquifer. In this study an orthogonal finite difference mesh is used with variable intervals in x and y directions. A schematic diagram of the aquifer mesh is shown in the upper part of Fig. 3; groundwater heads are functions of the coordinates x and y, hence subscripts m and n are used to indicate the location. It is also necessary to specify conditions on the lateral boundaries together with any inflows such as recharge, any natural outflows (to springs, rivers and the sea) and any man-made outflows such as to the horizontal well. In addition, initial conditions must be specified.
Flow from aquifer to the horizontal well Fig. 2(b) is a conceptual model indicating how water is drawn from the aquifer water table, through the aquifer to the horizontal well. One possible approach is to represent the actual cross-section of the horizontal well in the groundwater model using a very fine grid. However, there
Figure 3 Schematic diagram of the computational model showing part of the regional aquifer, the representation of flow from the aquifer to nodal points on the horizontal well using coefficient Ca and the increasing flows in the horizontal well.
101
are additional processes such as the effect of gravel packs, any disturbances to the aquifer during the construction of the well and the resistance to the flow of water through slots in the well pipes; these processes are difficult to quantify. Therefore an alternative approach is adopted representing all these features as an empirical coefficient. In many field studies, flow from an aquifer to a river, drain or vertical well is described by a relationship of the form Qa ¼ Caðh HÞ
ð2Þ
where Qa is the flow from the aquifer to the river, drain or well over a specified length; h is the groundwater head at the river, drain or well; H is the hydraulic head in the river, drain or well; Ca is a semi-empirical coefficient which depends on the dimensions of the river, drain or well, the hydraulic conductivity of the aquifer, any material in the river bed, drainage envelope or gravel pack and the resistance to flow through the perforations in the drain or slotted well casing. In MODFLOW, the river coefficient can be estimated in terms of the dimensions and properties of the deposits in the riverbed (Anderson and Woessner, 1992). However in most practical situations adjustments to this first estimate are made during sensitivity analyses. To explore the validity of the use of Eq. (1) to represent horizontal and vertical flow components in a shallow aquifer, and Eq. (2) to represent flows from the aquifer to a horizontal well, a related problem of groundwater flow to horizontal tile drains is considered in the Appendix A. Comparisons are made between a detailed vertical section two-dimensional time-variant analysis (profile model) of horizontal and vertical flow components from a moving water table to a drain and a simulation of drain–aquifer interaction using Eq. (2) with the aquifer represented using a one-dimensional form of Eq. (1). The discharge from the drain is specified. From information about the variation of the hydraulic head in the drain with time, the magnitude of the aquifer–drain coefficient Ca is adjusted during a sensitivity analysis. Groundwater heads in the aquifer and hydraulic heads in the drain determined using this approach show acceptable agreement with the full two-dimensional moving water table analysis, Fig. A.1. Two conclusions can be drawn from the acceptable agreement between these two alternative approaches. Firstly, the one-dimensional formulation of Eq. (1) does simulate the time-variant horizontal and vertical flow components in a shallow homogeneous aquifer system. Secondly, the interaction between the drain and aquifer can be represented by the coefficient Ca and Eq. (2); however, the magnitude of Ca has to be adjusted during a sensitivity analysis to match the hydraulic heads in the drain. For a horizontal well, since the hydraulic head varies along the horizontal pipe, the computational model for the interaction of a horizontal well with a shallow aquifer is represented as a series of discrete inflows at mesh points along the line of the well. The inflow associated with node m is indicated in Fig. 3. This inflow occurs from midway between nodes m 1 and m to midway between nodes m and m + 1, a distance Dx. The inflow over this length of pipe Qam ¼ Cam ðhm;n Hm Þ
ð3Þ
102
A. Mohamed, K. Rushton
where Cam is the aquifer-pipe coefficient for this length of pipe. Note that there are two suffices defining the location of the groundwater head hm,n since it is a function of x and y whereas the hydraulic head in the pipe Hm is only a function of x. Inflows at nodes m 1 and m + 1 are calculated in a similar manner.
Flow within horizontal well Fig. 2(c) illustrates inflows from the aquifer through the perforations in a horizontal well pipe and the increasing flows within the pipe to the discharge point. These conditions differ from those for conventional pipe theory since flows occur through the wall of the pipe causing velocities within the pipe to increase from zero at the blank end to a maximum at the discharge point. A computational model is required to represent the increasing flow along the pipe. Since the aquifer is modelled using a finite-difference mesh, the horizontal pipe is divided into discrete lengths Dx, Fig. 3. Turbulent flows often occur in pipes. For the field example discussed in ‘‘Application to horizontal well at Loba’’ section the maximum velocity in the pipe is about three times the critical ‘Reynolds velocity’, which is the maximum velocity for which the flow remains laminar. Velocities are lower away from the outlet but turbulent flows are likely to occur along the full length of the pipe due to the impact of the inflowing water through the pipe walls. There are several alternative mathematical expressions to describe pipe flow; for this analysis the Hazen-Williams expression is used (Daugherty et al., 1989), although other pipe-flow relationships could be selected, DH ¼ RQp1:85
ð4Þ
where DH is the head loss in the pipe (m), Qp is the flow rate in the pipe (m3 d1), and the coefficient of pipe resistance R¼
3:2 106 L 4:87 C1:85 H d
ð5Þ
Figure 4 well.
Detailed diagram of conditions in the horizontal
In Eq. (8) the coefficient Cpm1,m depends on the unknown flow Qpm1,m. Similarly Qpm;mþ1 ¼ Cpm;mþ1 ðHm Hmþ1 Þ
Continuity along the pipe, Fig. 4, requires that the pipe flow between nodes m 1 and m plus the inflow from the aquifer equals the pipe flow between nodes m and m + 1, Qpm1;m þ Qam ¼ Qpm;mþ1
Cpm1;m Hm1 ðCpm1;m þ Cam þ Cpm;mþ1 ÞHm þ Cpm;mþ1 Hmþ1 ¼ Cam hm;n
The manner in which the three component mechanisms are included in the total computational model is illustrated in Fig. 3; the flow chart of Fig. 5 indicates how the calculations proceed. For each time step the following procedure is adopted; note that the components are considered in the reverse order to ‘‘Regional groundwater flow’’, ‘‘Flow from aquifer to the horizontal well’’ and ‘‘Flow within horizontal well’’ sections so that they are in the same order as the flow chart. Input parameters and initial conditions 1. Solve simultaneous equns for flows in pipe, Eq. (10) to estimate Qp and H
Cpm1;m ¼
4
1
2. Calculate flows Qa from Eq. (3) 3. Solve finite difference form of Eq. (1) with appropriate 3 boundary conditions to calculate h(x, y)
2
iterative routine for Kh Repeat with updated h(x, y)
ð7aÞ
ð8Þ
.
new Cp values from Eq. (8)
ð6Þ
Next time step
where Hm1 and Hm are the hydraulic heads, whilst the pipe coefficient Cpm1,m is given by the equation 1 ðRQp0:85 m1;m Þ
ð10Þ
Combining the three component mechanisms
The term in square brackets is treated as a pipe coefficient Cp which is updated during iterative calculations. Cp is lower towards the outlet where the flow Qp increases. The pipe is divided into segments of length Dx with the hydraulic heads H defined at the nodal points, m 1, m and m + 1 but the coefficients Cp and the flows Qp are defined between nodal points as shown in Fig. 4. Hence Eq. (6) becomes Qpm1;m ¼ Cpm1;m ðHm1 Hm Þ
ð9Þ
From Eqs. (3), (7a), (7b) and (9)
In this study the pipe is divided into segments of length L = Dx; for each segment the Hazen-Williams discharge coefficient CH and the pipe diameter d are the same, hence R is unchanged for each segment. For computational purposes it is helpful to adopt the approach of McManus and Rushton (1964) by rewriting Eq. (4) as Qp ¼ DH½ðRQp0:85 Þ1 ¼ DHCp
ð7bÞ
Calculation completed
Figure 5
Flow chart for the computational model of Fig. 3.
Horizontal wells in shallow aquifers: Field experiment and numerical model Component 3: for each node along the horizontal well an equation of the form of Eq. (10) is written. This leads to a set of simultaneous equations which can be solved using Gaussian elimination. However, the coefficients Cp depend of the unknown flows in the pipe, Eq. (8), consequently an iterative technique is required, loop of Fig. 5. The right hand side of Eq. (10) involves the groundwater head hm,n. Since values of groundwater heads for the current time step have not yet been determined, it is necessary to return to Component 3 when new values of the groundwater head have been calculated. Component 2: Eq. (3) is used to calculate flows to the horizontal well from the aquifer. Component 1: groundwater heads within the aquifer are calculated from the finite-difference form of Eq. (1) with the appropriate boundary conditions and with outflows from the aquifer to the horizontal well as determined from Component 2 above. Since the transmissivities equal the hydraulic conductivity multiplied by the unknown saturated thickness, an iterative solution of the finite-difference equations is required as indicated by loop in Fig. 5; this is automatically performed by MODFLOW. After the first round of calculations for hydraulic heads in the well, groundwater heads in the aquifer and flows between aquifer and the horizontal well, the calculations are repeated using the improved values of the heads and flows until appropriate convergence criteria are met; this is indicated by iterative loop as shown in Fig. 5. For each time step the entire procedure is repeated as indicated by loop . Although preliminary estimates can be made of the coefficients Ca and Cp, they need to be refined to represent actual field conditions. Since the coefficient Ca, which is used to calculate flows from the aquifer to nodal points on the horizontal well, is an empirical coefficient, improved estimates are obtained during a sensitivity analysis which focuses on the difference between groundwater heads and hydraulic heads in the pipe. The pipe coefficient Cp depends on the coefficient of pipe resistance R and the current flow in the segment of pipe, Eq. (6). Initial estimates of R can be obtained from Eq. (5), but modifications to the Hazen-Williams discharge coefficient CH are necessary due to the effect of the slots in the wall of the pipe and the inflows through these slots. Therefore a further sensitivity analysis is required to estimate a suitable value for the coefficient of pipe resistance R to reproduce field readings of the difference in hydraulic head along the pipe. These procedures are discussed further in ‘‘Computational model for Loba’’ section.
103
flows estimated from the time taken to fill a container, variations in the pumping rate did occur. Pump breakdowns also occurred during the test. Estimates of the recharge intensity and pumping rate for periods of 6 hours are included in Fig. 6(a). Field readings of drawdowns in two of the aquifer observation piezometers, in the caisson and in the horizontal well are indicated by the discrete symbols in Figs. 6(b) and 8; all measurements are relative to datum which is 4.0 m below mean sea level. The estimated tolerance in the observed readings is ±1 cm due to uncertainties in the datum for measurement, the limited experience of the operators and the fact that some readings were taken at night. Measurements of the hydraulic head in the horizontal well are also influenced by adjustments to the pumping rate from the caisson.
Computational model for Loba The groundwater model for the shallow sand aquifer at Loba extends for a distance of 700 m perpendicular to the coast and 2600 m parallel to the coast, Fig. 7(a). The horizontal well is located in the centre of the modelled area and is parallel to the coast. The mesh spacing is indicated in Fig. 7(a); the minimum spacing of 5.0 m perpendicular to the well is too small to be identified on the diagram. A mesh spacing of 25 m is used along the horizontal well hence there are 13 nodes along the 300 m length of the horizontal well. Recharge and discharge data are input at intervals of 0.25 d, two time steps of 0.125 d are used for each data input. For the aquifer, the boundary conditions are a specified groundwater head of 4.0 m above datum along the
Application to horizontal well at Loba Pumping tests in the horizontal well Several pumping tests were carried out on the horizontal well at Loba; this discussion will concentrate on the first six days of a 14-day test. During the first six days there was only limited rainfall; thereafter, substantial rainfall led to total recovery in many of the observation piezometers. Since four dewatering pumps were used for the pumping test with gate valves at the discharge outlets and with
Figure 6 Pumping test at Loba: (a) assumed recharge and pumping rate, (b) field readings and modelled results, x and y axes defined in Fig. 7.
104
A. Mohamed, K. Rushton
Figure 7
Computational model for Loba.
coastline, in this shallow aquifer the presence of any saline interface can be ignored. With a transmissivity of 50 m2/d and a storage coefficient of 0.033, the tidal effect spreads no more than 70 m from the coast; since all the observation wells and the horizontal well are more than 250 m from the coastline, tidal effects are not significant. No-flow conditions are applied on the other three boundaries to represent the aquifer extending about 700 m perpendicular to the coast with the effect of pumping extending no more than 1 km on either side of the horizontal well. Initial conditions are based on the observation that maximum groundwater heads approximately coincide with the general shape of the land surface. Therefore initial heads are set equal to a smoothed land surface profile; the simulation starts with a period of 6 days without recharge or pumping to allow any inconsistencies in the initial heads to dissipate. Iterative routines are required both for the representation of hydraulic conditions in the horizontal well, loop of Fig. 5, and for the regional groundwater model, loop of Fig. 5. For the computational model of the horizontal well, improved values of the coefficients Cp are obtained from iterative loop ; convergence is achieved in four cycles. For the groundwater model the iterative calculation continues until the out of balance of flows at every node is equivalent to less than 0.005 mm d1. After three cycles of iterative loop , convergence is achieved. In the absence of information to the contrary, the hydraulic conductivity and specific yield are assumed to be uniform over the aquifer. The aquifer-pipe coefficient Ca is also assumed to take a constant value at all nodal points on the horizontal well. In addition, the coefficient
of pipe resistance R needs to be refined to allow for the effect of the inflows of water along the horizontal well. Using a sensitivity analysis, trial parameter values are explored; the criteria used to estimate suitable values of the parameters are listed in Table 1. Estimates of the aquifer hydraulic conductivity and specific yield were obtained from earlier studies (Mailvaganam et al., 1993); therefore initial trial solutions concentrated on adjustments to the aquifer-pipe coefficient Ca and the coefficient of pipe resistance R. An increase in R from the initial estimate (based on Eq. (5) for a smooth straight pipe, 25 m in length) by a factor of 10 was required to reflect the disturbances to flow cause by the slots and the water flowing through the slots into the pipe. Reductions in Ca, from a first estimate of 65 m2 d1 (deduced from Eq. (2) and field information about flows and head differences) were required to lead to greater differences between groundwater and hydraulic head. After achieving a reasonable simulation, further refinements were made to each of the parameters in turn. The final parameter values are K ¼ 12 m d1 ; SY ¼ 0:033 ðvalue is influenced by fine sand near ground surfaceÞ Ca ¼ 56 m2 d1 for a 25 m length of pipe, R ¼ 2:0 106 d1:85 m4:55 for 25 m of pipe It would be possible to use an alternative expression to Hazen-Williams with the loss proportional to Qp2; slightly different values of the pipe coefficient would be obtained but this would not have a significant effect on the overall model results.
Horizontal wells in shallow aquifers: Field experiment and numerical model Table 1
105
Parameters refined during the sensitivity analysis
Parameter
Key criterion
Hydraulic conductivity K
Groundwater head gradient towards the horizontal well depends primarily on K, reductions in K increase the gradient The overall decline in groundwater heads is strongly influenced by SY, a lower SY leads to greater drawdowns Difference between groundwater heads and hydraulic heads at the horizontal well Variation of hydraulic head along pipes of horizontal well, reduce R to increase hydraulic head difference along pipe
Specific yield SY Aquifer-pipe coefficient Ca Coefficient of pipe resistance R
Computational model results Fig. 6(b) includes comparisons between field readings (the discrete symbols) and computational model results (continuous lines) for the hydraulic head in the caisson and groundwater heads in observation piezometers at distances of 10 m and 30 m perpendicular to the horizontal well. The agreement is encouraging. Small differences, which occur in the measured groundwater heads, can be due to the tolerances in measurements discussed in ‘‘Pumping tests in the horizontal well’’ section and local variations in aquifer parameters such as varying elevations of the clay base to the aquifer. Changes in hydraulic heads in the caisson are also represented satisfactorily by the model. Since the effect of pump failure is represented as occurring throughout the time step between 3.75 and 4.0 days, there are small discrepancies in the timing of this recovery. Field and modelled hydraulic head distributions within the horizontal well at 3.5 and 6.0 days are compared in Fig. 8. The discrete symbols represent the hydraulic heads measured in the cleaning tubes of the horizontal well, the continuous line is constructed through the numerical model results for hydraulic head H at 25 m intervals along the horizontal well. Hydraulic heads in the caisson are sensitive to the pumping rate at the time when the reading was taken since the four dewatering pumps did not maintain a steady discharge. The difference in hydraulic head between the ends of the pipe and the caisson depends on the choice of the coefficient of pipe resistance R; the value selected gives an adequate simulation for each of the times included in Fig. 8. From the computational model results, Fig. 9 is prepared which shows both calculated inflows from the aquifer to the horizontal well at 25 m intervals and flows within the pipe towards the caisson. Three times are selected to illustrate
Figure 8
Hydraulic heads in horizontal well.
a range of pumping rates. Fig. 9(a) refers to a time of 1.25 d when 230 m3 d1 is pumped from the horizontal well; the vertical arrows, which indicate the flows from the aquifer into the well, are almost uniform. However when the pumping rate falls to 160 m3 d1, Fig. 9(b), inflows from the aquifer towards the centre of the horizontal well reduce significantly. This occurs due to the substantial decrease in hydraulic head drawdown at the caisson compared to the much smaller change in groundwater heads as shown in Fig. 6(b). An increase in total pumping rate to 280 m3 d1 also leads to a non-uniform inflow distribution to the horizontal well, Fig. 9(c). For these three times, substantial changes occur in flow in the horizontal pipe with flows to the caisson varying form 75 to 129 m3 d1. Results from Figs. 8 and 9 indicate that the horizontal well should not be represented either as a constant head or as a uniform inflow. Hydraulic heads at the caisson are typically 0.05 m lower than at the ends of the well; this variation is significant compared to the differences between groundwater head and hydraulic head of 0.25–0.40 m. Furthermore, inflows from the aquifer to the horizontal well are sensitive to any variations in the pumping rate; the maximum inflow is often 1.3 times and occasionally 1.5 times the minimum inflow.
Predictive simulations Since the combined numerical model provides an adequate representation of field conditions in the aquifer and in the horizontal well, it can be used for predictive simulations. Three predictive simulations are described. For the first predictive simulation the pumping rates are set at double the historical values. It is essential to consider whether modifications should be made to any of the three component models due to this increase in abstraction rate. For the groundwater component, decreases in saturated thickness are included in Eq. (1). The second model refers to flows from the aquifer into the horizontal well, Fig. 2(b); due to increased approach velocities at the horizontal well there is likely to be an increased resistance to flow which can be represented by decreasing Ca. Substantial increases in the well loss factor occur for vertical wells following increases in abstraction rates (Rushton, 2003); consequently the value of Ca is halved from the value deduced from the historical value (this is equivalent to doubling the well loss factor). For the third component model relating to flow in the horizontal pipe, Fig. 2(c), increased inflows through the wall of the pipe result in higher values of the pipe resistance R. In the absence of field information, R is doubled.
A. Mohamed, K. Rushton aquifer to pipe flow (m3 d-1 )
106 25 20
(a) 1.25 day
230 m3 d-1 on day 1.25 19.2
15
17.9
17.2
37
54
aquifer to pipe flow (m 3 d-1 ) aquifer to pipe flow (m 3 d-1 )
17.8
17.6
17.6
17.3
88
71
17
17.2
17.9
37
19
19.2
10 5 0 19
-150 25
71
-100
88
106 106
-50
0
50
54
100
150
(b) 3.875 day
20 15
160 m3 d-1 on day 3.875
16.2 13.6
10
12.1
11.3
10.8
10.6
10.7
10.6
10.8
11.3
12.1
42
53
64
75
75
64
53
42
30
16.2 13.6
5 0 16
30
16
(c) 4.50 day
25
280 m3 d-1 on day 4.50
25.1 22.6
20
21.1
20.4
20.1
68
89
109
20.3
20.9
20.3
129
129
109
20.1
20.4
21.1
89
68
48
25.1 22.6
15 10 5 0 25
-150
Figure 9
17.3
17
48
-100
-50
0
50
100
25
150
Calculated inflows to horizontal well (vertical arrows) and flow through horizontal well at specific times, units m3 d1.
When these modifications are included in the simulation the predicted maximum drawdown in the caisson is 2.76 m compared to 0.99 m for the historical simulation (see Fig. 8). From the end of the horizontal well to the caisson, the difference in hydraulic head is 0.24 m compared to 0.034 m for the historical simulation. This means that the well would be operating under a pressure head of about 2.0 m, hence this increase in pumping rate appears to be feasible. Nevertheless, since the calculated drawdowns in the caisson are based on assumed values of Ca and R, any proposed increase in the pumping rate should be confirmed by field-testing. In the second predictive simulation, the pumping rates are again double the historical values but the horizontal pipes are also doubled in length; Ca and R are returned to their original values. Lengthening the horizontal wells leads to a reduction in the water table drawdowns above the well. The hydraulic head at the caisson is also reduced to 1.79 m compared to 2.76 m for the first predictive run. However, the hydraulic head differences along each well become 0.62 m compared to 0.24 m for the first predictive run. This
results in a variation in the fluxes into the well with 140% of the average flux at the caisson and a minimum of 75% of the average flux towards the end of the well. These findings indicate that a lengthening of the well reduces the water table drawdown but the well becomes less efficient at collecting water at locations distant from the caisson. In a third predictive simulation the long-term yield of the aquifer system is examined. The pumping rate is set at the design value of 240 m3 d1, the simulation continues for 100 days without any recharge occurring (this would be an exceptional situation for a humid tropical climate). Since inflows into the horizontal well remain of similar magnitude to the historical simulation, an unchanged value of R is appropriate. It is possible that Ca will also remain unchanged but in case some increase in the resistance to flow does occur in the vicinity of the horizontal well, Ca is reduced to two-thirds of the value for the historical simulation. In Fig. 10, groundwater heads are plotted on a crosssection perpendicular to the horizontal well extending from the coast to the no-flow boundary. After 30 d the groundwa-
Horizontal wells in shallow aquifers: Field experiment and numerical model
Figure 10 Predictive response due to pumping for 100 days with no recharge.
ter head above the horizontal well falls to 4.0 m which coincides with sea level. However there is a groundwater mound between the horizontal well and the coast with a maximum head of 4.32 m; consequently fresh water continues to flow to the coast. By day 62 this mound has disappeared with a very small groundwater head gradient from the coast towards the horizontal well. At day 100 the gradient from the coast in the direction of the horizontal well becomes more significant. The groundwater head gradient from the coast at day 100 is 0.0014. When multiplied by the hydraulic conductivity of 12.0 m d1, this gives a Darcy velocity of 0.017 m d1. If the effective porosity is 0.10, the seepage velocity for the ingress of saline water equals 0.17 m d1. Consequently saline water would move a few metres into the aquifer; once substantial recharge occurs this saline water would be flushed back to the coast. This simulation suggests that in the very unlikely circumstances of no recharge occurring for 100 days, the horizontal well would continue to provide good quality water with only a small ingress of saline water at the coast.
Conclusions This paper describes a methodology for analysing the response of horizontal wells in shallow aquifers without the assumptions either of a uniform hydraulic head in the well or a uniformly distributed inflow from the aquifer along the length of the well. A computational model is introduced which represents the three components of flow in the aquifer, the hydraulic conditions in the horizontal well, and the interaction between the aquifer and the well. This model is used to study the field situation of a horizontal well in a shallow costal aquifer. From a sensitivity analysis, values are deduced for aquifer parameters and coefficients describing the well-aquifer interaction and the resistance to flow in the pipe. Satisfactory agreement is achieved between field and modelled values of both hydraulic heads and groundwater heads. The model is then used to predict the consequences of increasing the pumping rate, lengthening the well, or operating the existing well for 100 days without recharge.
107
Certain of the references quoted in the Introduction suggest that it may be acceptable to assume a uniform flux into the well. From the limited number of simulations described in this paper relating to a specific location, it appears that the uniform flux assumption may be valid for wells of limited length but for longer wells (and drains) the hydraulic head losses are significant leading to substantial variations in the flux. Since the aim of the current study is to develop numerical models which can be used to represent horizontal wells without any assumptions about the distribution of heads or fluxes, no firm conclusions can be drawn about situations where the uniform flux assumption is acceptable. The field study and subsequent mathematical modelling demonstrate the suitability of horizontal wells for abstracting substantial quantities of water from shallow aquifers. Although the cost of constructing a horizontal well is greater than for a number of vertical wells, the advantages include the need for only one pump and the ability to increase the discharge during periods of high demand. Furthermore, the horizontal well can continue to operate efficiently during periods of zero recharge. For an aquifer with a greater saturated thickness, a multi-layered regional groundwater model is required. Due to the uncertainties about the flow processes close to the horizontal well and the difficulties in representing the precise geometry of the well cross-section in a regional groundwater model, Eq. (3) should be used to represent the aquifer-well interaction with the coefficient Ca deduced during sensitivity analyses. Hydraulic conditions in the pipe should be simulated using Eq. (7). Alternative conditions on the upper and lateral boundaries can be included in the regional groundwater model. The methodology described in this paper is also useful when studying conditions in interceptor drains which are pumped horizontal drains positioned to intercept water lost from irrigation canals (Rushton, 2003).
Appendix A. Representation of drain–aquifer interaction in a shallow aquifer The purpose of this appendix is to explore the adequacy of the representation of flow in a shallow aquifer using a single layer model (‘‘Regional groundwater flow’’ section) and also representation of conditions in the vicinity of the horizontal well using an aquifer-well coefficient (‘‘Flow from aquifer to the horizontal well’’ section). This is achieved by considering the comparable problem of flow to horizontal tile drains. A study by Khan and Rushton (1996) examines time-variant flow in a vertical section to horizontal tile drains from which constant discharges are withdrawn. Fig. A.1(a) illustrates a two-dimensional vertical section from which the discharge from each drain is Qd = 0.125 m3 d1 per metre length of the drain. Since the tile drains are equally spaced at 50 m, lateral no-flow boundaries are located at 25 m from the drain centreline. The drains are 0.1 m in diameter and 2.0 m above the aquifer base. Initial conditions are that the water table and all the groundwater heads equal 4.5 m relative to the base of the aquifer. The differential equation describing the flow in a vertical section in terms of the groundwater head h(y, z, t) is
108
A. Mohamed, K. Rushton H ¼ h Qd=Ca
Figure A.1 Comparison of two-dimensional vertical section model with one-dimensional solution in which aquifer–drain interaction is represented by a drain coefficient: (a) twodimensional formulation, (b) results from two-dimensional solution and one-dimensional approximation, (c) one-dimensional formulation showing shape of water table.
o oh o oh oh Ky þ Kz ¼ SS oy oy oz oz ot
ðA:1Þ
in which SS is the specific storage coefficient. Khan and Rushton (1996) describe numerical solutions to Eq. (A.1) using time-variant finite-difference approximations in which both horizontal and vertical flow components are represented; the upper boundary is a moving water table which is a function of the specific yield. The continuous lines in Fig. A.1(b) show the calculated water table elevations midway between drains, at the quarter points and immediately above the drain; the head in the drain is also plotted. Aquifer parameter values are listed in Fig. A.1(b). The same problem is analysed using a one-dimensional approximation of Eq. (1) in which the groundwater head h(y, t), Fig. A.1(c), is a function of the horizontal coordinate y and the time t. The governing equation is o oh oh Kh ¼ SY ðA:2Þ oy oy ot where SY is the specific yield. The drain is represented as a constant discharge Qd per unit length of the drain; the hydraulic head in the drain is calculated by rearranging Eq. (2) as
ðA:3Þ
where Ca, the aquifer–drain coefficient, is adjusted during a sensitivity analysis to match the hydraulic head in the drain. A best fit is obtained with Ca = 0.32 m2 d1 m1; this coefficient applies for all times. Numerical solutions to the problem of the drain with a constant discharge are obtained using the finite-difference form of Eq. (A.2) with a mesh spacing of 0.5 m and a time step of 0.5 d. Since in Eq. (A.2) the hydraulic conductivity is multiplied by the unknown saturated thickness h, iterative solutions are required for each time step. In Fig. A.1(b), results from this numerical solution are indicated by the discrete symbols, results are plotted for alternate time steps. Water table elevations at the lateral boundaries and at the quarter points show no difference between the results for the two-dimensional and one-dimensional approximations. For the water table above the drain, the one-dimensional solution predicts elevations slightly lower than the two-dimensional solution. An examination of Fig. A.1(a) shows that for the two-dimensional solution the water table above the drain is horizontal. However, for the one-dimensional approximation, Fig. A.1(c), there is a discontinuity in the groundwater head because an outflow occurs from the drain node, consequently the water table is not horizontal. The hydraulic head in the drain calculated from the simulation based on Eq. (A.2) is similar to that determined from the two-dimensional simulation, although there are small differences for early and late times. These findings support the proposition that groundwater flow through a shallow aquifer towards a horizontal drain can be studied using a single layer approximation with the aquifer–drain interaction represented by a coefficient which is determined during a sensitivity analysis.
References Anderson, M.P., Woessner, W.W., 1992. Applied Groundwater Modeling: Simulation of Flow and Advective Transport. Academic Press, San Diego, 381pp. Cassiani, G., Kabala, Z.J., 1998. Hydraulics of a partially penetrating well: solution to a mixed-type boundary value problem via dual integration equations. J. Hydrol. 211, 100–111. Chen, C., Wan, J., Zhan, H., 2003. Theoretical and experimental studies of coupled seepage-pipe flow to a horizontal well. J. Hydrol. 281, 159–171. Connorton, B.J., 1985. Does the regional groundwater-flow equation model vertical flow. J. Hydrol. 79, 279–299. Daugherty, R.L., Franzini, J.B., Finnemore, E.J., 1989. Fluid Mechanics with Engineering Applications, SI Metric Edition. McGraw-Hill, Singapore, 596pp. Kawecki, M.W., 2000. Transient flow to a horizontal water well. Ground Water 38, 842–850. Khan, S., Rushton, K.R., 1996. Reappraisal of flow to tile drains: III Drains with limited flow capacity. J. Hydrol. 183, 383–395. Kruseman, G.P., de Ridder, N.A., 1990. Analysis and Evaluation of Pumping Test Data, second ed. International Institute for Land Reclamation and Improvement, Wageningen, The Netherlands, Publ. no. 47, 377pp. Mailvaganam, Y., Ramili, M.Z., Rushton, K.R., Ong, B.Y., 1993. Groundwater Exploitation of a Shallow Coastal Sand Aquifer in Sarawak, Malaysia. IAHS Publ. no. 216, pp. 451–462.
Horizontal wells in shallow aquifers: Field experiment and numerical model McManus, B.D., Rushton, K.R., 1964. A pure resistance electrical analogue for the analysis of pipe networks. Civ. Eng. Pub. Works Rev. London, 585–589. Pinder, G.F., 2002. Groundwater Modeling Using Geographical Information Systems. Wiley, New York, 233pp. Rushton, K.R., 2003. Groundwater Hydrology: Conceptual and Computational Models. Wiley, Chichester, 416pp. Rushton, K.R., Rathod, K.S., 1985. Horizontal and vertical flow components deduced from groundwater heads. J. Hydrol. 79, 261–278.
109
Tarshish, M., 1992. Combined mathematical model of flow in an aquifer-horizontal well system. Ground Water 30, 931–935. Zhan, H., Park, E., 2003. Horizontal well hydraulics in leaky aquifers. J. Hydrol. 281, 129–146. Zhan, H., Zlotnik, V.A., 2002. Groundwater flow to a horizontal or slanted well in an unconfined aquifer. Water Resour. Res. 38 (7), 13-1–13-11. Zhan, H., Wang, L.V., Park, E., 2001. On the horizontal-well pumping tests in anisotropic confined aquifers. J. Hydrol. 252, 37–50.