Slug Test in a Large Dip Angle Fracture Zone : Model and Field Experiment

A single-porosity model is developed to deal with a fracture zone slug test in a large dip angle by assuming the fracture zone causes a downward regional flow. For the oscillatory response, a larger dip angle causes larger amplitude while introduces little impact on period. The effective water length, an important parameter necessary for analyzing the oscillatory response, is proven to be independent of the dip angle and can be evaluated using the available horizontal formation methods. The dip angle effect is more pronounced for a larger storage coefficient. An empirical relationship is developed to evaluate the limiting dip angle, below which the dip angle effect is negligible. Field data analysis of a slug test in a 47° dip angle fracture zone indicates that neglecting the dip angle can result in a 27% transmissivity over estimation and a 53% storage coefficient under estimation.


INTRODUTION
The slug test is a relatively simple field method for estimating pertinent hydrogeologic parameters.Slug testing has been widely used in either granular aquifers or fractured formations.There are many studies for the slug test under various hydrogeological test conditions.It is well recognized (e.g., Van der Kamp 1976;Kipp 1985;Butler 1998;Zurbuchen et al. 2002; and so on) that the slug test response is oscillatory under a relatively high permeable condition while non-oscillatory under a relatively low permeable condition.Most of the works prior to 1998 can be found in Butler (1998), and those after 1998 include; e.g., McElwee and Zenner (1998); Zlotnik and McGuire (1998); Zurbuchen et al. (2002); Butler et al. (2003); Butler and Zhan (2004); Ostendorf et al. (2005); Chen (2006); Ross and McElwee (2007); Alexander et al. (2011); Rozos et al. (2015); and so on.Cooper et al. (1967) conducted one pioneering work in slug testing.There are double-porosity models, specifically for fractured formations, that take into account the interactive flow between the fracture and adjacent porous matrices based on various flow mechanism assumptions (e.g., Barker and Black 1983;Dougherty and Babu 1984;Mateen and Ramey 1984;Moench 1984;Black 1985;Sageev and Ramey 1986;Grader and Ramey 1988;Barker et al. 2000), and the fractal geometry models of a fractional flow dimensionality between 1 and 3 (e.g., Barker 1988;Novakowski and Bickerton 1997;Audouin and Bodin 2008).All of these models assume that the aquifer or formation being investigated is horizontal with little or no dip angle effect, in which the flow field caused by the slug test is assumed to be radially symmetric with respect to the test well.
However, field investigations conducted in a Cenozoic folded sandstone formation revealed that the dip angle of a fracture zone is as large as 47°.For a large dip angle, a uniform regional flow prevails through the sloping fracture zone and its presence will cause the flow field due to the slug test to be asymmetrical with respect to the test well.This asymmetric flow is similar to that of a regional flow superimposed onto a radial flow; e.g., see McWhorter and Sunada (1981, Figs. 4 -12) for the steady-state condition.As shown by these figures, the asymmertic flow is characteristic of a "groundwater divide" that consists of all the points of a zero hydraulic gradient because of the balance of the radial flow caused by the pumping well and the uniform regional flow.This divide extends to infinity up gradient of the pumping well while terminates at a stagnant point somewhere down gradient of the pumping well.All the groundwater within the divide can be withdrawn by the pumping well.All the groundwater outside the divide; however, bypasses the pumping well.Clearly, all of the currently available slug test models are not suitable for such a flow field, and thus how the large dip angle will influence the slug test response warrants study.It is important to have appropriate data analysis methods for both the oscillatory and non-oscillatory test responses with the dip angle effect, and to know when the dip angle effect can be neglected (more specifically, the limiting dip angle smaller than which the dip angle effect can be neglected.)Therefore, the purpose of this paper is to develop a slug test model that takes into account the dip angle effects and investigates the dip angle influence on the oscillatory and non-oscillatory slug test response.The test response solution and analysis in nearby observation wells will be presented elsewhere.

Model Development
The problem of interest is schematically illustrated in Fig. 1a, where a double-packer slug test is conducted in a sloping fracture zone with a dip angle a and aperture bcosa , b being the vertical uniform thickness of the fracture zone.Groundwater in the fracture zone can enter or leave the borehole only through the test section between the packers.The borehole radius is r w .The test response H(t) is measured using a pressure transducer submerged in the riser pipe of radius r c .The model is developed based on the single-porosity approach, considering that a double-porosity characteristic is usually evident only at large times where the test response is too small to be practically useful (Butler 1998).As a result, the fracture zone is assumed to be a homogeneous and isotropic continuous porous medium under the confined condition, where Darcy's law is applicable.It is understood that the flow in close vicinity of the test well can be larger than non-Darcian flow (e.g., Wu 2002).This possible condition is not taken into account in the current model.
Owing to the dip angle, the elevation head Z(x, y) of the fracture zone varies as tan , where cos x l a = is the horizontal projection of l.Accordingly the piezometric head can be expressed as ( , , ) ( , , ) tan h x y t P x y t x a = + (1) where P(x, y, t) represents the pressure head.It is assumed that the dipping angle induces a downward regional flow of uniform hydraulic gradient tan Z x 2 2 a = in the fracture zone.The problem associated with an upward regional flow in the fracture zone will be presented elsewhere.
Note that the current model is unavailable for a vertical fracture zone because 2 a r = leads to an infinity hydraulic gradient for the downward regional flow [i.e., ( ) tan 2 " 3 r ].As shown in Fig. 1b, it is further assumed that the streamlines in the fracture zone are nearly parallel to the sloping bed (Boussinesq 1877).Applying the mass balance principle to the control volume gives where S is the storage coefficient of the fracture zone.The symbol q l [L 2 /T] represents the Darcy flow rate per unit fracture zone width along the l-direction, as defined by The symbol q y [L 2 /T] represents the Darcy flow rate per unit fracture zone width along the y-direction, as defined by The flow equation in the fracture zone is formulated on the horizontal xy plane.This is done by first introducing cos to Eq. (3).Substituting the resultant equation for q l along with Eq. ( 4 where v = S/T with T = Kb.Since the regional flow in the fracture zone has a uniform hydraulic gradient tan , it is known from Eq. (1) that the pressure head before the test does not change in space ( P x P y 0 2 2 2 2 = = ); that is P(x, y, 0) = constant.Without loss of generality, the initial condition for Eq. ( 5) can thus be appropriately prescribed as ( , , ) P x y 0 0 = (6) The fracture zone is modeled as having an infinite extent.This assumption is invoked based on the consideration that the slug test is initiated by limiting the initial water level drop H 0 less than 0.5 m in order to secure a hydrostatic water pressure response of H(t); e.g., see Butler et al. (2003), Chen et al. (2012).For H 0 < 0.5 m, the radius of influence of the slug testing flow may not reach any of the hydrogeological boundaries of the finite fracture zone, and thus the fracture zone can be assumed to have an infinite extent.However, if H(t) measured exhibits the interference from a hydrogeological boundary, the well-known image well method can be applied to deal with the boundary effect.The outer boundary conditions for Eq. ( 5) are set forth as ( ,) , Due to the dipping angle, both h(x, y, t) and its hydraulic gradient are asymmetric with respect to the test well during the slug test.Note that Eq. ( 5) is mathematically equivalent to the governing equation of groundwater flow in a homogeneous anisotropic aquifer under the confined condition when Tcos 2 a is visualized as T x and T as T y .Hantush (1966) investigated the constant-rate pumping problem in a homogeneous and anisotropic leaky aquifer where the pumping well is assumed to be a line sink (i.e., r 0 w " ), and indicated that the flux around the well varies with the flow directions because of aquifer anisotropy [e.g., see his Eq. ( 11 , indicating that the net flow rate of the regional flow entering and leaving the test section is zero.As a result, the boundary condition of the testing well is Equations ( 5) -( 7) and (9b) form the pressure head variation model in the fracture zone due to the slug test on the horizontal xy plane.The model solution determination is made easier using an appropriate polar coordinate system, which is established through the following two steps.First, letting cos x x a = l enables us to express Eq. ( 5) in its standard Laplacian form P x P y v P t . Second, this Lapacian equation for ( , , ) P x y t l is transformed to the typical radial flow equation as  6) and ( 7) are ( , ) While the well geometry on the xy plane is a circle ( x y r w ) on the x y -l plane.As derived in the Appendix, this ellipse on the r i l l plane is described using the equation of ( ) r r w m i = l l with ( ) cos cos sin Correspondingly the boundary condition of Eq. ( 9b) is modified to (see the appendix for a more detailed derivation) As stated earlier, the dipping angle cannot be 2 r , so ( ) . As a result, the integration for Eq. ( 13) exists for 0 2 < # a r .For a horizontal fracture ( 0 a = ), ( ) 1 m i = l and the test well recovers its circular geometry.Eqs. ( 10) -( 13) are used to determine the solution for the pressure head distribution.Now we proceed to deal with the slug test model for H(t).It is recognized that for a relatively low-K condition, the water level in the riser pipe recovers slowly and H(t) is considered to be instantaneously equilibrated to h w (t), the averaged value of the pressure head in the aquifer around the wellbore of the test section.The solution obtained for h w (t) is therefore set equal to H(t).For a relatively high-K condition, however, the water level in the riser pipe recovers relatively fast and H(t) is oscillatory due to the significant inertial force.In this event, H(t) is not the same as h w (t), and they can be related to each other through the following linearized momentum equation (e.g., Van in which the averaged elevation head around the test section wellbore vanishes because of 0 tan cos r d . The initial conditions for Eq. ( 14) are prescribed as where H 0 is the slug test initial head drop.The mathematical model of interest consists of Eqs. ( 10) -( 17).As Eq. ( 14) is typical of free oscillation, it is valid for either the underdamped condition where the solution is oscillatory or the over damped condition where the solution is non-oscillatory (e.g., see Wylie and Barrett 1982).That is, we use the above model for both the relatively high-and low-K conditions.

Model Solution
The model solution will be determined using appropriate dimensionless parameters.The dimensionless mathematical model is ( ) 1 < where ( ) is the dimensionless pressure head in the fracture zone, 2r S r 2 z = is the dimensionless transmissivity, and L e is the effective water column length in the test well.
Applying the Laplace transform with respect to x to Eqs. ( 18) -(21) yields 0 s 1 where s is the Laplace transform variable and ( , ) represents the dimensionless pressure head ( , ) after the Laplace transform is defined as # (e.g., see Kreyszig 1998).The solution to Eqs. ( 26) -( 28) is where and K 0 (x), K 1 (x) are the modified Bessel functions of the second order zero kind and one, respectively.Note that ( ) , and the evaluation of f 1 (s) can be numerically carried out without difficulty.
( ) w x is determined using Eqs.( 22) -( 24).The ( ) Laplace transform in Eq. ( 25) is obtained by substituting Eq. ( 29) into Eq.( 25).In a straightforward manner the Laplace-domain solution for ( ) where , and the f 2 (s) evaluation can be numerically carried out without difficulty.When 0 a = , Eq. ( 31) becomes identical to its counterpart in a horizontal fracture zone [e.g., see Eq. ( 16) from Kipp (1985) while setting both the skin effect and water surface initial velocity to zero].The Laplace inverse of Eq. ( 31) gives the solution for ( ) w x .The Laplace inverse is numerically calculated using the De Hoog et al. ( 1982) method.

THEORETICAL ANALYSIS
The theoretical analysis focuses on how the dip angle a influences the test response ( ) w x under various v and z. Figure 2 pertains to the relatively low-K non-oscillatory response conditions and Fig. 3 to the relatively high-K oscillatory response conditions.For the set v and z, the symbol * a denotes the upper limit below which the dip angle effect can be neglected.For the same set v and z, the test response recovers faster when the dip angle becomes larger, as shown in Fig. 2a.This is observed in the solid curves group associated with v = 0.01 and z = 0.05, where ( ) w x of a = 85° recovers faster than does that of a = 70° and ( ) w x of a = 70° recovers faster than does that of a = 34°.This faster recovery is attributed to the fact that the wellbore flow rate ( ) q * x associated with the slug test increases as a increases, as shown in Fig. 2b.The determination of ( ) = is made by numerically inverting its Laplace-domain counterpart, Eq. ( 28).When z remains unchanged, * a increases as v decreases, as indicated by the group of dashed curves of z = 0.05, which show * a increases from 34 -45° as v = 0.01 decreases to 10 -4 .When v remains constant, * a is not significantly influenced by the variation of z.This is observed by noting * a remains as 45° for both the dashed curve (z = 0.05) and the broken curve (z = 0.5), when v is kept to be 10 -4 .
For the relatively high-K conditions, the oscillation amplitude increases as a increases while the period shows little change, as shown in Fig. 3, where the solid curve amplitude of a = 85° is greater than that of a = 70° and the solid curve amplitude of a = 70° is greater than that of a = 34°.This is because a larger dip angle causes smaller ( ) w h x .Using an electric analog model, Bredehoeft et al. (1966) demonstrated that h w (t) functions as a frictional force to the oscillatory H(t).Below it further proves that h w (t) indeed introduces a damping force to H(t).As far as * a is concerned, it decreases as v increases and has little influence from the variation of z; the same as for the relatively low-K conditions.
We determine * a for various values of v , and the results are plotted in Fig. 4. Using regression analysis an empirical relation for * a (in radian) as a function of v is determined as 1.95 2.36 As S is unknown a priori while a is usually measured before the test, Eq. ( 34) is more useful to "predict" whether or not the test response would be influenced by the dip angle.Take the well field site in the Cenozoic folded sandstone formation for example.It is known that S of sandstone rock is approximately 10 -4 -10 -2 (Singhal and Gupta 2010).As r w = 0.0508 m and r c = 0.0135 m, the values of ( ) r S r 2 w c 2 2 v = range from 2.83 × 10 -3 to 0.283.Equation ( 34) then gives 26° ≤ * a ≤ 37°.That is, the dip angle effect can be quite safely neglected in the data analysis if a measures less than 26°, and should be included if a measures greater than 37°.When a is measured between 37° and 26°, the dip angle effect may or may not be negligible.For practical purposes the dip angle effect should be considered in the data analysis when a > 26°.
The above theoretical analysis focuses primarily on the test response in the test well.For a cross-borehole slug test, the test response can be transmitted from the test well to the nearby observation wells.The test response analysis in an observation well requires the response solution in an observation well and the momentum balance in the observation well (e.g., see Zhan and Butler 2003;Audouin and Bodin 2008).The flow field complexity is primarily due to the fact that the test response in the fracture zone changes with distance as well as direction.Observation well response modeling and data analysis will be presented elsewhere.

DATA ANALYSIS
In data analysis, the effective water length L e is recommended to be treated as a fitting parameter (e.g., Kabala et al. 1985;Butler et al. 2003;Chen 2006).For a horizontal formation L e is defined using the following equation where T 1 ~-6 @ and T 1 b - 6 @ are the frequency and the damping coefficient of oscillation, respectively (Van der Kamp 1976).Both ~ and b can be estimated from the measured test response (Chen 2006;Chen and Wu 2006), and thus L e can be uniquely evaluated using the field data.Whether or not L e is so determined can be used for a sloping fracture zone depending on whether Eq. ( 35) is valid for a sloping formation.To this end, we start with the oscillatory response solution in a horizontal formation (e.g., Wylie and Barrett 1982;Springer and Gelhar 1991) In terms of the dimensionless parameters used in this paper, the Laplace transform of Eq. ( 36) is  31) to (37) discloses the following two facts.First, the ( ) g s 2 z term in Eq. ( 31) is equivalent to the damping coefficient * b in Eq. ( 37).Given that ( ) g s 2 z is from the Laplace transform of h w (t), the electric analog observation by Bredehoeft et al. (1966), as discussed above, is confirmed.Second, the constant term of 1 2 z in Eq. ( 31) is equal to the constant term of 0.25 * * 2 2 b + in Eq. ( 37).The dimensional analysis of this equality proves Eq. ( 35) is also valid for a non-zero dip angle.Therefore, L e is evaluated using Eq. ( 35) with ~ and b determined using the field data is applicable to the sloping fracture zone case.After L e is known the two parameters T and S can be determined using the least root mean square error fitting method.
There is a well field (Fig. 5) in a Cenozoic folded sandstone formation overlain by a weathered zone about 20 m in thickness in northern Taiwan.This site was about 15 m above a nearby stream and had seven wells (C, N, E, S, W, W1, and W2).The layout of the wells is shown in Fig. 5.All of the wells were 4-inch in diameter (r w = 0.0508 m) and cased in the weathered zone and uncased as an open borehole from 20 -50 m in sandstone.Wells W1 and W2 had not been in use due to certain construction problems.The image taken by a borehole optical televiewer in a well reveals that at 31 m depth a fracture zone exists consisting of a set of fractures distributed over the 0.1 m packed interval and dipping to the northwest (344°) with a dip angle of 47°.The actual extent of this fracture zone is unknown.A doublepacker slug test was conducted at well C in the a = 47° dip angle fracture zone.No response was observed in the four surrounding wells during the test.Therefore, it is appropriate to assume the sloping fracture zone to be infinite in extent for the test.In general, the current model assumptions adequately satisfy the field conditions.
The slug test was initiated pneumatically with an initial head drop of H 0 = 0.48 m.The pressure transducer for H(t) measurement is placed 0.5 m below the initial water level.The measured H(t) is nearly associated with hydrostatic water pressure and the linear model of Eq. ( 14) is adequate (e.g., Butler et al. 2003;Chen et al. 2012).Since the dip angle is greater than the upper limit of * a = 37°, as determined earlier, the data analysis needs to take into account the dip angle effect.The fracture zone thickness is b = 0.1 m.The effective water column length L e is determined using the method given by (Chen 2006), as shown in Fig. 6.The frequency is determined as 2 ( ) , where t k is the time corresponding to the occurrence of the k th extremity (peak or valley of the oscillation), where the subscript k =1, 2, ….As indicated in Fig. 6, t k + 2t k = 8.74 sec gives ~ = 0.72 sec -1 .The damping coefficient is calculated by 4 , where H k denotes the k th extremity displacement.As H k /H k + 1 = 1.97, b is 0.31 sec -1 .As a result, L e = 18.2 m, as determined using Eq. ( 35).Using Eq. ( 32) and L e = 18.2 m, the best oscillatory test response fit, as shown by the solid curve, gives T = 4.12 × 10 -4 m 2 s -1 and S = 5.34 × 10 -4 .Now we proceed to investigate the effect of neglecting the dip angle in data analysis.This is done by setting a = 0° in Eq. ( 31) when analyzing the field data.As a result, the best fit represented by the dashed curve in Fig. 6 is associated with T = 5.24 × 10 -4 m 2 s -1 and S = 2.53 × 10 -4 .Although the solid and dashed curves are nearly coincident, the T and S estimates are quite different.The calculated dip angle negligence causes a 27% transmissivity over estimation and a 53% storage coefficient under estimation.Furthermore, we used S = 5.34 × 10 -4 and Eq. ( 34) to obtain * a = 34°, which is indeed less than 47°, supporting using the current model to analyze the field data.
The sensitivity of the current model to T (for a given a and a given S) is demonstrated in Fig. 7a, where three T values, 5.24 × 10 -4 m 2 s -1 , 4.1 × 10 -4 and 3.0 × 10 -4 m 2 s -1 are used to represent a ±27% change in transmissivity from 4.1 × 10 -4 m 2 s -1 .A ±27% change in T causes a significant variation in the calculated test response.Figure 7b shows the model solution sensitivity to S, where three values for S = 8.15 × 10 -4 , 5.3 × 10 -4 , and 2.53 × 10 -4 are used to represent a ±53% variation in the storage coefficient from the one of 5.3 × 10 -4 .In comparison to Fig. 7a it is clear that the current model is less sensitive to S.

CONCLUSIONS
Slug test models in the literature assume a horizontal formation.This research developed a new slug test model that takes the dip angle effects into account.The dip angle is assumed to cause a downward regional flow in the fracture zone.The problems associated with an upward regional flow will be presented elsewhere.The Laplace-domain solution was determined for both the oscillatory and non-oscillatory test response.When the fracture zone transmissivity is relatively low (i.e., T < 1.0 × 10 -4 m 2 s -1 ), a larger dip angle causes a faster non-oscillatory test response recovery.When the fracture zone transmissivity is relatively high (i.e., T > 1.0 × 10 -4 m 2 s -1 ), a larger dip angle causes an increase in oscillatory test response amplitude.It is proven that the effective water length, an important parameter necessary for oscillatory test response analysis, is independent of the dip angle and can be evaluated using the methods currently available for horizontal formations.
In general, neglecting the dip angle may lead to hydraulic conductivity over estimation and storage coefficient under estimation.The dip angle effect is more pronounced for a larger storage coefficient, being less sensitive to the change in transmissivity.An empirical relationship as a function of the dimensionless storage coefficient is developed for the limiting dip angle, below which the dip angle effect is negligible.A number of cross-borehole slug tests were conducted in a Cenozoic folded sandstone formation, where a dip angle fracture zone as large as 47° was measured.The slug test data in the 47° dip angle fracture zone was analyzed using the current model.Neglecting the dip angle can result in a 27% transmissivity over estimation and

APPENDIX A
Derivation of Eq. ( 13) from Eq. (9b).The well circular geometry on the xy plane needs to be expressed on the r i l l plane.By substituting cos x r i = and sin r y which after rearranging terms can be expressed as Using the definitions of i and il, it is known that tan tan cos which is the integral appeared in Eq. ( 13).
)].Similarly, the flux around the well varies with the flow directions in the current problem and the mass flow rate continuity across the test well test section can be formulated as Fig. 1.(a) Schematics of a double-packer slug test in a dip angle a fracture zone, which induces a uniform regional flow in a constant hydraulic gradient of tana .(b) The control volume diagram for deriving the flow equation.

Fig. 2 .
Fig. 2. Influence of a on the test response for relatively low-K conditions: (a) the head recovers faster for a larger a ; (b) the wellbore flow rate increases with increasing a .The limiting angle * a increases as v decreases, while it remains relatively insensitive to z for the same v .

Fig. 5 .
Fig. 5. Schematics of the well field in northern Taiwan, which has seven wells (C, N, E, S, W, W1, and W2).

Fig. 6 .
Fig. 6.Analysis of the slug test data in a dip angle equal to 47° with a fracture zone, where t k and H k are used for the effective water length determination L e .The solid ( a = 47°) and dashed curves ( a = 0°) are nearly coincident while the T and S estimates are different.
Fig. 7. (a) Sensitivity of the model solution to a ±27% change in transmissivity.(b) Sensitivity of the current model to a variation of ±53% in the storage coefficient.
shown earlier in the text.Now the relationship for d d i il is derived by taking the first derivative with respect to i (or il for the same result) of both of the second kind of order 0 K 1 (x) modified Bessel function of the second kind of order 1 l distance along the direction paralleling to the sloping bed [L] L e effective water length of the slug test [L] P(x, y, t) pressure head response in the fracture [L]q l (t) Darcy flow rate per unit fracture zone width along the l-direction [L 2 /T] q y (t) the Darcy flow rate per unit fracture zone width along the y-direction [