Determination of skin and aquifer parameters for a
slug test with wellbore-skin effect
Hund-Der Yeh
*
, Yen-Ju Chen
Institute of Environmental Engineering, National Chiao Tung University, 75 Po-Ai Street, Hsinchu 30039, Taiwan Received 15 February 2006; received in revised form 28 May 2007; accepted 29 May 2007
KEYWORDS Groundwater; Parameter estimation; Sensitivity analysis; Simulated annealing; Skin thickness; Slug test
Summary Slug test is considered to reflect the hydraulic parameters in the vicinity of the test well. The aquifer parameters are usually identified by fitting an appropriate mathe-matical solution or graphical type curves with slug test data. In this paper, we developed an approach by combining [Moench, A.F., Hsieh, P.A., 1985. Analysis of slug test data in a well with finite-thickness skin. In: Memoirs of the 17th International Congress on the Hydrogeology of Rocks of Low Permeability, U.S.A. Members of the International Associ-ation of Hydrologists, Tucson, AZ, vol. 17, pp. 17–29] and simulated annealing (SA) approach to estimate five parameters, i.e., three skin parameters and two aquifer param-eters. The three skin parameters are hydraulic conductivity, specific storage, and thick-ness of the wellbore-skin zone, while the two aquifer parameters are hydraulic conductivity and specific storage of the formation zone. It is worthy to note although the thickness of the wellbore-skin zone is usually taken as an input data for the data-ana-lyzed software, it is actually an unknown parameter that cannot be measured directly. This paper proposes a methodology for estimating the thickness of the wellbore-skin zone with other hydraulic parameters at the same time.
Eight sets of well water-level (WWL) data of aquifers with both positive and negative skins are generated by Moench and Hsieh [Moench and Hsieh, 1985] and four sets of stan-dard normally distributed noise are then added to each set of WWL data. The results indi-cate that the negative-skin cases generally give a better estimated result than that of the positive-skin cases. Sensitivity analysis is also employed to demonstrate the physical behavior when slug test was performed under positive-skin effect. For the case of an aqui-fer with a positive-skin, the use of a longer series of WWL data for analysis is strongly rec-ommended for better estimation of aquifer hydraulic conductivity. Analyzing the WWL data of the test and observation wells simultaneously could significantly improve the
esti-0022-1694/$ - see front matter ª 2007 Elsevier B.V. All rights reserved. doi:10.1016/j.jhydrol.2007.05.029
* Corresponding author. Tel.: +886 3 5731910; fax: +886 3 5726050. E-mail address:[email protected](H.-D. Yeh).
a v a i l a b l e a t w w w . s c i e n c e d i r e c t . c o m
mations on specific storages. Impetuously presuming an arbitrary value for the thickness of the wellbore-skin zone may lead to poor estimation for the other four parameters. ª 2007 Elsevier B.V. All rights reserved.
Introduction
Slug test is a quick and economical way for obtaining the changes in WWL and determining the hydraulic conductiv-ity in the vicinconductiv-ity of the test well. The test is performed by removing (adding) suddenly a small amount of water from (into) the well and measuring the change in WWL simultaneously. Subsequently, the aquifer parameters are identified by fitting an appropriate mathematical solution or graphical type curves with the slug test data. Several slug test mathematical models can be found in the ground-water literature (e.g.,Hvorslev, 1951; Cooper et al., 1967;
Bouwer and Rice, 1976; Springer and Gelhar, 1991; Hyder et al., 1994; Butler, 1998). In addition, some software, e.g., AQTESOLV (Duffield, 2002), which includes the above-mentioned models has been developed to analyze the slug test data. For a specific well, aquifer configura-tion, and slug test condiconfigura-tion, one can choose an appropri-ate mathematical model to estimappropri-ate aquifer parameters.
A wellbore skin is a zone of diverse hydraulic conductiv-ity adjacent to the wellbore face caused by well drilling or completion. A positive skin, usually arisen from the damage of well drilling, is a zone adjacent to the wellbore with a hydraulic conductivity smaller than that of the undisturbed formation. In contrast, a negative skin refers to a zone with a hydraulic conductivity larger than that of the undisturbed formation owing to excessive well completion. The well-bore-skin effect would influence the recovery rate of WWL especially in a small-scale slug test.
Regarding the investigations on wellbore-skin effect, Ramey et al. (1975) presented an analytical solution for a slug test with a skin of infinitesimal thickness. They used the type curves approach, which is developed from a cor-relating group of dimensionless storage constant and skin effect, to determine the aquifer hydraulic conductivity and the skin effect. However, there may be a large uncer-tainty in the determination of the correlating group due to the similarity in the shape of the type curves. Therefore, such a curve-matching procedure may result in an inaccu-rate estimation of aquifer hydraulic conductivity. Faust
and Mercer (1984) investigated the effect of a finite-thick-ness skin on the response of slug tests using a simple ana-lytical solution and a numerical model. They compared the normalized head responses for a slug test performed with positive-skin, negative-skin, and no-skin effects, respec-tively. Moreover, they concluded that a negative skin has little effect on the shape and does not shift the normal-ized head response curve horizontally so that it would not obviously affect the estimated result of aquifer hydraulic parameters. On the other hand, a positive skin induces a similar shape but shift the response horizontally when compared with a no-skin response curve and thus leads to an unreliable estimate of aquifer hydraulic parameters. Following the investigation of Faust and
Mer-cer (1984), Moench and Hsieh (1985) presented a
semi-analytical solution and used type curves to illustrate the influences of a finite-thickness skin on the open-well and pressurized slug tests. In addition, Yang and Gates (1997)
analyzed the effect of wellbore skin on slug test results utilizing the finite-element model and field tests. Further, Yeh and Yang (2006) proposed an analytical solution for a slug test performed in a confined aquifer with finite-thick-ness skin. Their results showed that the early- and late-time data reflect the ground-water flow within the well-bore-skin and undisturbed formation zones, respectively. The model proposed by Hyder et al. (1994), also known as the KGS model, is useful for determining the hydraulic conductivity and specific storage in a confined or uncon-fined aquifer for a partially penetrating well with skin ef-fect. The software package AQTESOLV (Duffield, 2002), which includes the KGS model, can analyze the data of slug test under skin effect. However, the AQTESOLV re-quires the ‘‘outer radius of well skin’’ or the ‘‘skin thick-ness’’ as an input data. In fact, the skin thickness is actually unknown and cannot be measured. The identifica-tion of skin thickness as well as hydraulic parameters would be more practical and useful for giving hydrologists an insight into the skin-affected slug test.
Simulated Annealing (SA) is a stochastic technique for solving optimization problems.Metropolis et al. (1953) pro-posed the SA algorithm for evaluating the properties of a material composed of interacting individual molecules. Kirkpatrick et al. (1983) were the first to present an algo-rithm based on the analogy between the annealing of solids and the problem of solving combinatorial optimization (Pham and Karaboga, 2000). Subsequently, utilization of SA in optimization problems has been widely applied in hydrological engineering. For example, the method has been employed to design the strategies of groundwater remediation (Dougherty and Marryott, 1991; Marryott
et al., 1993), identify the best parameter structure for groundwater models (Zheng and Wang, 1996), and deter-mine the optimal network design or water resource manage-ment (Pardo-Iguzquiza, 1998; Cunha and Sousa, 1999; Kuo
et al., 2001). As presented in these studies, SA is useful in solving optimization problems that consist of many degrees of freedom and several local optima.
The purpose of this paper is to develop an approach that combines the Moench and Hsieh’s slug-test-solution
(1985) and SA to analyze the slug test data for estimating three skin parameters and two aquifer parameters. The merit of this approach over existing slug test analysis is that skin thickness is treated as an unknown parameter and estimated simultaneously with the other four parame-ters. Seventeen scenarios and sensitivity analyses for se-lected scenarios of slug tests were given to demonstrate the capability of the suggested approach in identifying five unknown parameters, i.e. the hydraulic conductivities and specific storages of both the formation and skin zones as well as the skin thickness.
Methodology
In past researches, many optimization algorithms were ap-plied for finding optimum values of an objective function with several independent variables. Here, an approach that couples SA withMoench and Hsieh’s solution (1985)is devel-oped to determine three skin parameters and two aquifer parameters. This section is thus divided into two parts. We first introduce the concept, process, and parameter requirements in SA, while Moench and Hsieh’s solution
(1985), which is employed to evaluate the objective func-tion, is reviewed in the second part.
Simulated annealing
SA is known as an optimization algorithm for simulating a material crystallized in the process of annealing. The arrangement of the material molecules is initially disor-dered at high temperature. The system is gradually cooled; meanwhile, the arrangement becomes more ordered and approaches a thermodynamic equilibrium. The framework of SA is demonstrated in Fig. 1. The parameters required in SA are given in the following.
The algorithm in SA starts with an initial guess at higher temperature and regards the initial guess as the current optimal solution. The five parameters to be estimated are hydraulic conductivities k1and k2, specific storages Ss1and
Ss2 in the skin and formation zones, respectively, and the
skin thickness dsk. A lower bound (LB) as well as the upper
bound (UB) for each trial solution of the target parameter is chosen to make sure that the guessed value is physically reasonable. The hydraulic conductivity for silt sand forma-tion may range from 107to 103m/s (Freeze and Cherry,
1979). Thus, the LBs and UBs are respectively chosen as 107m/s and 103m/s for both k
1and k2. In addition, the
values of storativity generally range from 105to 103in a confined aquifer (Schwartz and Zhang, 2003). The LBs and UBs for both Ss1 and Ss2 are respectively taken to be
106m1 and 104m1 on the ground that the thickness of a confined aquifer is assumed to be 10 m. The skin thick-ness dsk, which is equal to rs rwwith rsand rwrepresenting
the outer radius of the wellbore-skin zone and the effective well radius, respectively, is calculated after rsis estimated
in SA. Typically, slug tests result in a variation in radius be-tween 1 and 2 m (Rovey, 1998). Thus, the UB for rsis set to
be 2 m; and the LB for rsis equal to rw, i.e., 0.0915 m.
Fur-thermore, the initial guessed values of k1and k2are chosen
as 103 m/s, Ss1 and Ss2 are 104m1, and rsis 0.0915 m.
The parameter estimation process beginning with a higher initial temperature (T0) may give a more complete
explora-tion of the soluexplora-tion domain but may requires much more computing time than the one with lower temperature. The value of T0is thus chosen to be 5 in this study.
The initial guessed values are first employed to generate the WWL data usingMoench and Hsieh’s slug test solution
(1985) and calculate the objective function value. Note that Nomenclature
b aquifer thickness c1 acpK0(qb) + bqK1(qb)
c2 acpI0(qb) bqI1(qb)
dsk skin thickness
ej the difference between observed and predicted
WWL well water level
h well water level of slug test in the Laplace do-main
ho(tj) observed well water level at time tj
he(tj) estimated well water level for current solution
at time tj
h0
e(tj) estimated well water level for trial solution at
time tj
H0 initial well water level
I0(Æ) modified Bessel function of the first kind of
or-der zero
I1(Æ) modified Bessel function of the first kind of
or-der one
k1 hydraulic conductivity in a skin zone
k2 hydraulic conductivity in a formation zone
K0(Æ) modified Bessel function of the second kind of
order zero
K1(Æ) modified Bessel function of the second kind of
order one
N number of parameters
NS number of searches for a parameter
NT number of times for producing NS once at a tem-perature level
p Laplace variable
P(Æ) probability q pffiffiffip
r radial distance from the centerline of well rc radius of casing
rs outer radius of wellbore-skin region
rw effective radius of test well
RD random number of a uniform (0, 1) distribution RT cooling rate
Ss1 specific storage in a skin zone
Ss2 specific storage in a formation zone
t time since the start of test T0 initial temperature
Te present temperature
T0
e new temperature V degree of freedom
VMi step length vector of the parameter i.
xi current solution of parameter i,
x0
i new trial solution of parameter i
z objective function value of current solutions z0 objective function value of trial solutions a k2/k1 b (aSs1/Ss2)1/2 c c 2pr2 wSs2b D1 aI0(qbrs)K1(qrs) + bI1(qbrs)K0(qrs) D2 aK0(qbrs)K1(qrs) bK1(qbrs)K0(qrs)
Dz difference between trial and current objective functions
the objective function is analogous to the energy in thermo-dynamics in SA. Later, SA will randomly generate a neigh-boring solution (trial solution) and evaluate the objective function. A new trial solution for parameter i namely x0
iis
produced randomly from
x0i¼ xiþ ð2 RD 1Þ VMi ð1Þ
where xiis the current solution for parameter i, RD is a
ran-dom number generated from a uniform (0, 1) distribution, and VMiis a step length vector of the parameter i. The value
of VMiis first defined as UB LB. Then VMiis automatically
adjusted according to the objective function value calcu-lated from the trial solution generated at current tempera-ture to make that approximately half of all trial solutions will be accepted at the next temperature. On the other hand, if the trial solution is out of bounds, another proce-dure for creating a trial solution within bounds is
x0
i¼ LBiþ RD ðUBi LBiÞ ð2Þ
For a given temperature, each parameter requires a number of searches (NS) and the search is repeated NT times. Thus, the total number of iterations is N· NS · NT per tempera-ture level for the target parameters, where N denotes the number of parameters. The total number of iterations at
each temperature level is taken as 500 with N = 5, NS = 20, and NT = 5. Note that a greater number of searches may give a better result but requires more computing time.
The objective functions for current and trial solutions are respectively defined as
z¼X n j¼1 ðhoðtjÞ heðtjÞÞ 2 ð3Þ and z0¼X n j¼1 ðhoðtjÞ h0eðtjÞÞ2 ð4Þ
where z and z0 are the objective functions for current and
trial solutions, respectively; ho(ti) is the observed WWL;
and he(tj) and h0e(tj) are the estimated WWL for current
solution and trial solution at time tj, respectively.
As the objective function is minimized, any downhill step is accepted or an uphill step may be accepted according to the Metropolis criterion. The Metropolis’s criterion is ex-pressed as (Pham and Karaboga, 2000)
Pðaccept trial solutionÞ ¼
1 for Dz 6 0 exp Dz Te for Dz > 0 ( ð5Þ where P is the acceptance probability of the trial solution, Dz = z0 z is the difference between the trial and current
objective functions, and Te is the present temperature.
The trial solution supplants the current optimal solution as the starting point of the next state on condition that Dz 6 0. On the other hand, if Dz > 0, the current optimal solution also has a probability to be replaced by the trial solution even though it is inferior to the current one. The replacement happens when P is larger than a random num-ber generated from a uniform (0, 1) distribution. In addition, according to the Metropolis’s criterion, the probability of accepting inferior trial solutions decreases with the cool-down process. The major advantage of SA over other optimization algorithms is its ability to accept inferior trial solution and to escape from local optimum.
After N· NS · NT iterations at a temperature level, the system decreases the temperature by a constant cooling rate (RT). Therefore, the new temperature becomes
T0
e¼ RT Te ð6Þ
The cooling rate which ranges from zero to one is chosen as 0.85 as suggested in Corana et al. (1987). The cool-down process continues until the system satisfies the termination criteria.
Two termination criteria are used in this study. The first criterion is that the difference of the objective functions between two consecutive temperatures is less than 106 for four times sequentially. The second and minor one is that the maximum number of total iterations allowed is set to be 2· 107. Once the total number of iterations
reaches this number, the identification process will be ter-minated even the solution is not optimal.
Moench and Hsieh’s slug test solution
Moench and Hsieh’s solution (1985)is chosen to couple with SA for determining the skin and aquifer parameters. Their solution assumes that a fully penetrating well is installed
in a confined, homogeneous, isotropic, and infinite lateral extent aquifer of constant thickness with a finite-thickness skin. The dimensionless form of WWL solution in the Laplace domain is expressed as h¼ac½D1K0ðqbÞ D2I0ðqbÞ c1D1 c2D2 ð7Þ with D1¼ aI0ðqbrsÞK1ðqrsÞ þ bI1ðqbrsÞK0ðqrsÞ ð8Þ D2¼ aK0ðqbrsÞK1ðqrsÞ bK1ðqbrsÞK0ðqrsÞ ð9Þ c1¼ acpK0ðqbÞ þ bqK1ðqbÞ ð10Þ c2¼ acpI0ðqbÞ bqI1ðqbÞ ð11Þ a¼ k2=k1 ð12Þ b¼ ðaSs1=Ss2Þ 1=2 ð13Þ c¼ pr 2 c 2pr2 wSs2b ð14Þ and q¼ p1=2 ð15Þ
The notation rwis the effective radius of test well; rcis the
radius of well casing; b is the aquifer thickness; and p is the Laplace variable. Moreover, I0(Æ) and I1(Æ) denote the
modi-fied Bessel functions of the first kind of order zero and one, respectively; and K0(Æ) and K1(Æ) are the modified Bessel
functions of the second kind of order zero and one, respec-tively. Eq.(7)is first employed to generate several sets of WWL data, i.e., the observed WWLs ho(tj). Then, the
esti-mated WWL data he(tj) or h0e(tj) are calculated according
to Eq.(7)to estimate the objective function, that is, either Eq.(3)or(4).
Table 1a The target values and estimated results by analyzing the WWL data within 15 s for aquifers with a positive skin
Case Estimated results
k1(m/s) k2(m/s) Ss1(1/m) Ss2(1/m) dsk(m) SEE Target value 1.00· 105 1.00· 104 1.00· 104 1.00· 104 0.9085 – 1a 1.00· 105 1.14· 104 9.94· 105 4.63· 105 0.895 3.43· 104 1b 1.03· 105 2.80· 104 8.98· 105 2.03· 105 1.051 5.97· 104 1c 1.07· 105 5.02· 104 8.17· 105 7.82· 105 1.267 8.65· 104 1d 1.04· 105 1.20· 104 8.76· 105 8.07· 105 1.024 7.41· 104 1e 1.02· 105 1.45· 104 9.51· 105 5.32· 105 0.973 5.76· 104 Mean 1.03· 105 2.32· 104 9.07· 105 5.57· 105 1.042 – RE 3.20% 132.20% 9.28% 44.26% 14.69% – Target value 1.00· 105 1.00· 104 1.00· 104 1.00· 104 0.3085 – 2a 1.53· 105 8.50· 105 3.18· 105 9.10· 105 0.856 2.82· 104 2b 1.82· 105 1.08· 104 1.68· 105 3.48· 105 1.486 6.81· 104 2c 1.55· 105 8.09· 105 3.15· 105 9.97· 105 0.907 8.56· 104 2d 1.51· 105 1.07· 104 3.17· 105 5.36· 106 0.683 7.35· 104 2e 1.58· 105 1.08· 104 2.92· 105 1.78· 105 0.869 7.78· 104 Mean 1.60· 105 9.78· 105 2.82· 105 4.97· 105 0.960 – RE 59.80% 2.22% 71.80% 50.27% 211.25% – Target value 1.00· 105 1.00· 103 1.00· 104 1.00· 104 0.9085 – 3a 9.98· 106 7.99· 104 1.00· 104 7.93· 105 0.895 3.26· 104 3b 9.98· 106 9.32· 104 9.93· 105 1.33· 105 0.891 6.01· 104 3c 1.03· 105 9.65· 104 9.20· 105 8.15· 105 1.007 8.73· 104 3d 1.05· 105 9.96· 104 8.41· 105 5.52· 105 1.039 7.45· 104 3e 1.00· 105 9.83· 104 9.98· 105 5.50· 105 0.902 7.96· 104 Mean 1.02· 105 9.35· 104 9.50· 105 5.69· 105 0.947 – RE 1.52% 6.50% 4.96% 43.14% 4.22% – Target value 1.00· 105 1.00· 103 1.00· 104 1.00· 104 0.3085 – 4a 1.07· 105 9.85· 104 8.57· 105 7.43· 106 0.346 3.15· 104 4b 1.90· 105 9.72· 104 1.53· 105 9.70· 105 1.449 7.83· 104 4c 1.14· 105 8.04· 104 7.78· 105 5.94· 105 0.400 8.17· 104 4d 1.58· 105 9.99· 104 2.82· 105 2.94· 105 0.856 7.39· 104 4e 1.22· 105 9.94· 104 6.53· 105 3.80· 105 0.462 7.91· 104 Mean 1.38· 105 9.51· 104 5.45· 105 4.62· 105 0.703 – RE 38.20% 4.92% 45.54% 53.75% 127.75% –
Results and discussion
Assume that the slug test is performed in a confined, homo-geneous, isotropic, and infinite lateral extent aquifer sys-tem. The test well is of full penetration and the radius of effective well and well casing are 0.0915 m and 0.0508 m, respectively. The initial WWL is assumed to be 1 m, while the aquifer thickness is assumed to be 10 m. The WWL data are produced according to Moench and Hsieh’s solution
(1985) and analyzed by the present approach for determin-ing the four hydraulic parameters and skin thickness for a slug test with a positive or negative skin.
Table 1presents the ‘‘true’’ aquifer conditions, i.e., the target values of the four hydraulic parameters and skin thickness, and the results estimated by the present approach for eight scenarios. In Case 1a of Scenario 1, the WWL data were directly generated usingMoench and Hsieh’s
solution (1985) and listed inTable 2. The WWL data for the other four cases in Scenario 1, i.e., Cases 1b–1e, were obtained by adding four different sets of noise to the WWL data of Case 1a. These four sets of noise are generated by the routine RNNOF ofIMSL (1997), which can generate nor-mally distributed random numbers with zero mean and unit variance. They are then divided by 1000 to represent the measurement errors of the water-level meter with the mag-nitude on the order of millimeter. Note that the standard error of the estimate (SEE) listed in the last column ofTable
1 is defined as Pnj¼1e2 j=v
1=2
, where ejrepresents the
dif-ference between the observed WWL and the predicted WWL according to the estimated parameters and v, the de-gree of freedom, equals the number of observed data point minus the number of unknowns (Yeh, 1987). In addition, the relative error (RE) is defined as the difference between the
Table 1b The target values and estimated results by analyzing the WWL data within 15 s for aquifers with a negative skin
Case Estimated results
k1(m/s) k2(m/s) Ss1(1/m) Ss2(1/m) dsk(m) SEE Target value 1.00· 104 1.00· 105 1.00· 104 1.00· 104 0.9085 – 5a 1.15· 104 9.31· 106 6.03· 105 7.35· 105 1.150 3.40· 104 5b 1.30· 104 9.94· 106 3.50· 105 3.66· 105 1.547 6.97· 104 5c 1.13· 104 9.04· 106 6.39· 105 8.06· 105 1.120 1.01· 103 5d 1.12· 104 9.74· 106 6.43· 105 6.94· 105 1.134 6.95· 104 5e 1.11· 104 1.06· 105 6.60· 105 5.88· 105 1.141 6.99· 104 Mean 1.16· 104 9.73· 106 5.79· 105 6.38· 105 1.218 – RE 16.20% 2.74% 42.10% 36.22% 34.11% – Target value 1.00· 104 1.00· 105 1.00· 104 1.00· 104 0.3085 – 6a 1.69· 104 9.85· 106 2.69· 105 3.07· 105 0.628 3.27· 104 6b 1.39· 104 9.90· 106 7.94· 105 8.95· 105 0.322 6.69· 104 6c 1.67· 104 9.42· 106 5.48· 105 7.67· 105 0.381 9.52· 104 6d 8.15· 105 1.06· 105 8.98· 105 6.77· 105 0.381 7.90· 104 6e 1.39· 104 9.94· 106 5.58· 105 6.12· 105 0.414 7.95· 104 Mean 1.39· 104 9.94· 106 6.13· 105 6.52· 105 0.425 – RE 39.10% 0.58% 38.66% 34.84% 37.83% – Target value 1.00· 103 1.00· 105 1.00· 104 1.00· 104 0.9085 – 7a 9.98· 104 1.00· 105 8.93· 105 8.88· 105 0.968 3.21· 104 7b 1.00· 103 1.09· 105 3.01· 105 2.36· 105 1.785 8.17· 104 7c 9.99· 104 1.03· 105 4.16· 105 3.73· 105 1.490 1.03· 103 7d 9.63· 104 1.03· 105 6.93· 105 6.46· 105 1.121 7.08· 104 7e 8.52· 104 1.06· 105 7.60· 105 6.41· 105 1.084 6.89· 104 Mean 9.62· 104 1.04· 105 6.13· 105 5.57· 105 1.290 – RE 3.76% 4.20% 38.74% 44.32% 41.95% – Target value 1.00· 103 1.00· 105 1.00· 104 1.00· 104 0.3085 – 8a 5.75· 104 1.00· 105 7.78· 105 7.62· 105 0.370 3.07· 104 8b 9.86· 104 1.02· 105 7.11· 105 6.61· 105 0.388 7.50· 104 8c 9.83· 104 9.64· 106 5.58· 105 6.43· 105 0.432 9.44· 104 8d 2.26· 104 1.04· 105 6.38· 105 4.96· 105 0.476 7.05· 104 8e 2.18· 104 1.07· 105 6.86· 105 4.69· 105 0.470 7.57· 104 Mean 5.98· 104 1.02· 105 6.74· 105 6.06· 105 0.427 – RE 40.24% 1.88% 32.58% 39.38% 38.48% –
mean and target value divided by the target value to reflect the magnitude of parameter uncertainty on estimation.
Positive-skin scenarios
Table 1ademonstrates the estimated results in four scenar-ios for the aquifer with a positive skin. Both Scenarscenar-ios 1 and 2 had the same target values of k1= 105m/s and
k2= 104m/s but different target values of dsk, i.e.,
0.9085 m and 0.3085 m, respectively. In Cases 1a–1e of Sce-nario 1, the estimated k1, Ss1and dskwith fairly good mean
values of 1.03· 105 m/s, 9.07· 105m1 and 1.042 m,
respectively, and RE of 3.20%, 9.28% and 14.69%. How-ever, relatively inaccurate results were obtained for k2
and Ss2 with a mean value of 2.32· 104m/s and
5.57· 105m1 and RE of 132.20% and 44.26%, respec-tively. Scenario 2 shows fairly good results on k2for a mean
value of 9.78· 105m/s and RE of2.22%. However, in this scenario, the estimated values of dsk, ranging from 0.683 to
1.486 m with a mean value of 0.960 m and RE of 211.25%, are exaggerated while the target dsk is only 0.3085 m. In
addition, slightly inaccurate results were obtained for k1,
Ss1 and Ss2 with a mean value of 1.60· 105m/s,
2.82· 105 m1 and 4.97· 105m1, and RE of 59.80%,
71.80% and 50.27%, respectively.
Scenarios 3 and 4 use the same target values of parame-ters as Scenarios 1 and 2 except that k2= 103m/s, which is
one order of magnitude higher than that of Scenarios 1 and 2. The estimated results of Scenario 3 are fairly good for skin parameters k1, Ss1 and dsk with the averages of
1.02· 105m/s, 9.50· 105 m1 and 0.947 m and RE of
1.52%, 4.96% and 4.22%, respectively. The mean value
and RE are respectively 9.35· 104m/s and6.50% for k2
and 5.69· 105 m1 and43.14% for S
s2, which are
rela-tively inaccurate but acceptable. Overall, the estimated re-sults of Scenario 3 are similar to those of Scenario 1 in that aquifer parameters, i.e., k2and Ss2were poorly estimated.
In Scenario 4, the unfavorable results of dskvary from 0.346
to 1.449 m with a mean value of 0.703 m and a RE of 127.75% while the target dskequals 0.3085 m. In addition,
relatively inaccurate estimations were also obtained for k1, Ss1, and Ss2 with a mean value of 1.38· 105 m/s,
5.45· 105 m1 and 4.62· 105 m1 and RE of 38.20%,
45.54% and 53.75%, respectively. Both Scenarios 2 and 4 represent the positive-skin scenarios that have a thinner skin zone than that simulated in Scenarios 1 and 3. More-over, relatively unfavorable estimations on k1, Ss1, Ss2,
and dskare seen in these scenarios.
Negative-skin scenarios
In contrast toTables 1a and 1blists the estimated results for four negative-skin scenarios. Scenarios 5 and 6 have the same target values of k1= 104m/s and k2= 105m/s, yet
different target values of dsk, i.e., 0.9085 m and
0.3085 m, respectively; while Scenarios 7 and 8 differ from Scenarios 5 and 6 in that one-order higher magnitude of tar-get values of k1= 103m/s. Note that the estimated k2for
the four negative-skin scenarios are accurate with each hav-ing a RE less than 5%. In addition, the estimated results of the other four parameters are acceptable. In short, the neg-ative-skin scenarios as a whole give a better estimation than that of the positive-skin scenarios. Overall, the values of the five parameters estimated using SA approach shown inTable Table 2 The WWL responses simulated byMoench and Hsieh’s solution (1985)for the hypothetical Scenarios 1–8
Time (s) WWL (m) 1a 2a 3a 4a 5a 6a 7a 8a 0.5 0.967 0.967 0.967 0.967 0.838 0.895 0.644 0.868 1.0 0.947 0.944 0.947 0.942 0.754 0.856 0.599 0.830 1.5 0.930 0.921 0.930 0.918 0.698 0.826 0.567 0.800 2.0 0.914 0.900 0.914 0.895 0.657 0.801 0.542 0.775 2.5 0.899 0.879 0.899 0.872 0.624 0.778 0.521 0.752 3.0 0.884 0.859 0.884 0.850 0.596 0.758 0.502 0.732 3.5 0.870 0.840 0.870 0.829 0.573 0.739 0.485 0.713 4.0 0.857 0.821 0.856 0.808 0.552 0.721 0.470 0.695 4.5 0.844 0.802 0.843 0.788 0.534 0.705 0.457 0.679 5.0 0.831 0.784 0.830 0.768 0.518 0.689 0.444 0.663 6.0 0.806 0.750 0.804 0.730 0.488 0.661 0.421 0.635 7.0 0.783 0.717 0.780 0.694 0.463 0.635 0.401 0.609 8.0 0.760 0.686 0.756 0.659 0.441 0.611 0.384 0.586 9.0 0.738 0.656 0.732 0.627 0.422 0.589 0.368 0.564 10.0 0.716 0.628 0.710 0.596 0.404 0.568 0.353 0.544 11.0 0.695 0.601 0.688 0.566 0.388 0.549 0.340 0.525 12.0 0.675 0.576 0.667 0.538 0.373 0.531 0.328 0.508 13.0 0.656 0.551 0.647 0.512 0.360 0.515 0.316 0.492 14.0 0.637 0.528 0.627 0.486 0.347 0.499 0.306 0.476 15.0 0.619 0.506 0.608 0.462 0.335 0.484 0.296 0.461
1 can produce a set of WWL data coinciding with those cre-ated by the target values of parameters with the SEE values being all less than 103.
Sensitivity analysis
In this section, sensitivity analysis is employed to investigate the problem of unfavorable estimated results in Scenarios 1 and 2.Fig. 2displays the normalized sensitivities, represent-ing the temporal change in WWL with respect to a relative change in one of the five parameters (Kabala, 2001). Analyzing a longer set of WWL data
Fig. 2a illustrates the normalized sensitivities of the five parameters over a 1000-s interval for Scenario 1. As can be
seen, all parameters have their own influence on the WWL. A positive change in the value of rsposes positive influence
on the WWL through the 1000-s interval, while a positive change in the values of k1, k2and Ss2has negative influence.
In addition, a positive change in the value of Ss1has negative
influence at the early 50-s interval, but exerts positive influ-ence over the 50- to 1000-s interval. The WWL is almost insensitive to the change in Ss1 and Ss2. Thus, these two
parameters allow higher degree of uncertainty in estima-tions than the other parameters. The influence of k2on the
WWL is less sensitive than that of Ss1at the prior 20-s
inter-val; afterwards, it has more sensitive characteristic than that of Ss1and Ss2. The unfavorable estimate on k2in Scenario
1 may be attributed to the short WWL period within 15 s and 20 data points used in identifying the hydraulic parameters. Scenario 9, listed inTable 3, gives the estimated results of analyzing the 180-s WWL data with the same target parameters as Scenario 1. An obvious improvement can be observed from the estimated results of k2and Ss2. It is
note-worthy that the averages of k2 and Ss2become 1.02· 104
m/s and 9.44· 105m1, which are very close to the target
values of 1.00· 104m/s and 1.00· 104m1, respec-tively. The parameter estimation took about 7.6 hours for analyzing a set of 47 WWL data of Case 9a and about 3.6 hours when analyzing a set of 20 WWL data in Case 1a when using a personal computer with 3.6 G Pentium IV CPU and 1 GB RAM. Though the use of a longer set of WWL data con-sumes more computing time, it gives a more accurate esti-mation of k2in a positive-skin case.
Fig. 2b shows the normalized sensitivities of the five parameters for Scenario 2. The influences of Ss1and Ss2on
the WWL are also insensitive over the 1000-s interval. In addition, the sensitivity of k2is apparently larger than that
of Ss1and Ss2after 4 s and the estimated result of k2in
Sce-nario 2 is fairly good. The 180-s WWL data with the same target parameters as Scenario 2 are also employed to esti-mate the five parameters. The results represented as Sce-nario 10 and listed in Table 3 reveal a more accurate estimation on k2with a mean value of 9.97· 105m/s;
nev-ertheless, there are only minor improvements in k1, Ss1, and
dskand even slight deterioration in Ss2.
Highly correlated parameters: k1and rs
The normalized sensitivities of k1and rs, shown inFig. 2, are
relatively high when compared with those of the other parameters. Nevertheless, the sensitivity curves of k1 and
rsare symmetrical in shape on the horizontal axis but have
different magnitudes, implying that the parameters k1and
rs are highly correlated. This may be the reason why SA
may sometimes misidentify the values of k1and dsk. Table 4shows the estimated results for Scenario 11 when the unknown parameters are the same as those in Scenario 2 except that dskis known. In this scenario, k1is precisely
iden-tified with a mean value of 1.06· 105m/s and RE of 6.20%; additionally, k2and Ss1are accurately estimated when
com-pared with those presented in Scenario 2. Similarly, in Sce-nario 12, k1 is regarded as a known parameter and other
parameters are the same as those given in Scenario 2. The results listed in Table 4 show a fairly good estimation on dskwith a mean value of 0.286 m and RE of7.42%. In
addi-tion, more accurate estimations on k2and Ss1are presented
when compared with those shown in Scenario 2. Scenarios 11
Figure 2 (a) The normalized sensitivities of the five parameters over a period of 1000 s for (a) Scenario 1 and (b) Scenario 2.
and 12 show a highly correlated relationship between the un-known parameters k1 and dsk in a positive-skin scenario,
which may lead to misidentification of these two parame-ters. Other positive-skin scenarios also have a normalized sensitivity plot similar to that shown inFig. 2.
Analyzing WWL data in an observation well
Assume that an observation well is located at 3 m away from the test well. In Scenario 13, we analyzed the 15-s WWL data simulated in the test and observation wells
simultaneously. The estimated results are presented in Table 5.
Fig. 3shows the normalized sensitivities of the five param-eters using the WWL data measured in the observation well of Scenario 13. The parameter k1appears to have positive
influ-ence on the WWL over the 1000-s interval, while rshas
neg-ative influence at early 100-s interval, which becomes positive afterwards. In addition, the normalized sensitivity curves of k1and rsare not highly correlated.Table 5shows
obvious improvements in the estimations of k1and dskwith Table 3 The target values and estimated results for Scenarios 1 and 2 while the analyzed WWL data is extended to 180 s
Case Estimated results
k1(m/s) k2(m/s) Ss1(1/m) Ss2(1/m) dsk(m) SEE Target value 1.00· 105 1.00· 104 1.00· 104 1.00· 104 0.9085 – 9a 1.01· 105 1.01· 104 9.67· 105 8.67· 105 0.817 3.26· 104 9b 1.10· 105 9.91· 105 7.13· 105 9.83· 105 1.206 8.99· 104 9c 1.00· 105 9.95· 105 9.95· 105 9.57· 105 0.915 7.88· 104 9d 1.00· 105 1.10· 104 9.96· 105 9.12· 105 0.936 8.48· 104 9e 1.04· 105 1.02· 104 8.68· 105 9.99· 105 1.029 1.05· 103 Mean 1.03· 105 1.02· 104 9.08· 105 9.44· 105 0.981 – RE 3.00% 2.32% 9.22% 5.64% 7.94% – Target value 1.00· 105 1.00· 104 1.00· 104 1.00· 104 0.3085 – 10a 1.51· 105 1.01· 104 3.22· 105 1.95· 105 0.756 3.08· 104 10b 1.33· 105 1.00· 104 4.69· 105 6.72· 105 0.588 9.90· 104 10c 9.99· 106 9.93· 105 9.57· 105 4.30· 106 0.243 7.92· 104 10d 1.63· 105 1.05· 104 2.64· 105 1.50· 105 0.940 8.47· 104 10e 1.63· 105 9.62· 105 2.70· 105 9.86· 105 1.104 9.87· 104 Mean 1.42· 105 1.00· 104 4.56· 105 4.09· 105 0.726 – RE 41.98% 0.30% 54.36% 59.08% 135.40% –
Table 4 The target values and estimated results of Scenario 2 on condition that dskand k1is regarded as a known parameter, respectively, in Scenarios 11 and 12
Case Estimated results
k1(m/s) k2(m/s) Ss1(1/m) Ss2(1/m) dsk(m) SEE Target value 1.00· 105 1.00· 104 1.00· 104 1.00· 104 0.3085 – 11a 1.05· 105 9.64· 105 8.44· 105 2.75· 105 – 2.85· 104 11b 1.01· 105 1.10· 104 9.50· 105 4.27· 105 – 6.70· 104 11c 1.06· 105 9.58· 105 8.82· 105 2.41· 105 – 8.41· 104 11d 1.13· 105 1.01· 104 7.14· 105 3.17· 106 – 7.51· 104 11e 1.06· 105 1.03· 104 8.73· 105 1.74· 105 – 8.12· 104 Mean 1.06· 105 1.01· 104 8.53· 105 2.30· 105 – – RE 6.20% 1.24% 14.74% 77.03% – – Target value 1.00· 105 1.00· 104 1.00· 104 1.00· 104 0.3085 – 12a – 9.55· 105 9.42· 105 3.89· 105 0.284 2.80· 104 12b – 1.09· 104 9.84· 105 7.84· 105 0.310 6.72· 104 12c – 9.08· 105 9.99· 105 8.00· 105 0.296 8.38· 104 12d – 1.00· 104 9.35· 105 1.74· 105 0.272 7.55· 104 12e – 1.01· 104 9.80· 105 1.30· 105 0.266 8.14· 104 Mean – 9.93· 105 9.68· 105 4.55· 105 0.286 – RE – 0.74% 3.20% 54.46% 7.42% –
mean values of 1.02· 105m and 0.327 m, respectively, and
RE of 2.40% and 6.00%. Moreover, the WWL data simulated in the observation well are more sensitive to the change in Ss1
and Ss2than that in rs. Particularly, the normalized
sensitiv-ities of Ss1and Ss2depicted inFig. 3become more sensitive
than those displayed inFig. 2b. The results also show signif-icant improvements in estimations of Ss1and Ss2with mean
values of 9.89· 105 m1 and 9.99· 105 m1, and RE of 1.06% and 0.08%, respectively. In this case, analyzing the composite WWL data in the test and observation wells gives better results for the five parameters.
Comparison of different positive-skin scenarios
Table 6displays the estimated results of the five parameters for four positive-skin scenarios by analyzing their 15-s WWL data. Scenario 14, whose results are shown inTable 6a, rep-resents an aquifer with a positive skin having target values of parameters same as those in Scenario 2, except that k1
equals 106. Comparing Scenario 14 with Scenario 2 shows that Scenario 14 gives a much better estimation on skin parameters k1, Ss1, and dsk, but a worse estimation on k2.
Such unfavorable estimation on k2can be improved by using
a longer set of WWL data that explicitly reflects the aquifer behavior. This phenomenon also appears in Scenario 4 in contrast to Scenario 2. In other words, it is much easier to identify skin parameters accurately from the case having a distinct difference between k1 and k2. Table 6b
demon-strates three positive-skin scenarios having different values of dskand shows that the estimations on skin parameters k1, Table 5 The target values and estimated results for Scenario 2 while the hypothetical WWL data is analyzed in the test well as well as in an observation well which is located at 3 m away from the test well
Case Estimated results
k1(m/s) k2(m/s) Ss1(1/m) Ss2(1/m) dsk(m) SEE Target value 1.00· 105 1.00· 104 1.00· 104 1.00· 104 0.3085 – 13a 1.01· 105 1.02· 104 9.98· 105 1.00· 104 0.314 3.00· 104 13b 1.02· 105 1.10· 104 9.83· 105 9.99· 105 0.329 6.89· 104 13c 1.06· 105 1.06· 104 9.95· 105 1.00· 104 0.354 7.33· 104 13d 9.70· 106 9.26· 105 1.00· 104 9.98· 105 0.282 7.91· 104 13e 1.06· 105 1.12· 104 9.71· 105 9.99· 105 0.356 7.83· 104 Mean 1.02· 105 1.05· 104 9.89· 105 9.99· 105 0.327 – RE 2.40% 4.52% 1.06% 0.08% 6.00% –
Both the analyzed WWL data in test and observation wells are within 15 s.
Figure 3 The normalized sensitivities of the five parameters over a period of 1000 s for Scenario 13.
Table 6a The target values and estimated results for a positive-skin aquifer with a distinct difference between the hydraulic conductivities of skin and formation zones
Case Estimated results
k1(m/s) k2(m/s) Ss1(1/m) Ss2(1/m) dsk(m) SEE Target value 1.00· 106 1.00· 104 1.00· 104 1.00· 104 0.3085 – 14a 1.00· 106 8.51· 104 9.70· 105 8.51· 106 0.315 3.26· 104 14b 1.00· 106 6.51· 104 9.41· 105 8.97· 105 0.315 7.38· 104 14c 1.03· 106 6.59· 104 9.94· 105 5.50· 105 0.349 8.21· 104 14d 1.38· 106 8.33· 105 4.87· 105 7.10· 105 0.623 7.49· 104 14e 1.02· 106 7.66· 104 9.97· 105 9.11· 105 0.324 6.98· 104 Mean 1.09· 106 6.02· 104 8.78· 105 6.31· 105 0.385 – RE 8.60% 502.06% 12.22% 36.94% 24.86% –
Ss1, and dsk with a larger target dskare better than those
with a smaller one. Similarly, the unfavorable estimations on the two aquifer parameters k2and Ss2can be improved
by using a longer set of WWL data. In short, it is difficult to identify the skin parameters correctly for a positive-skin aquifer having an obscure skin characteristic, e.g., a similar k1 and k2 and/or a small dsk; nevertheless, the aquifer
parameters could be identified accurately by analyzing a longer set of WWL data.
Conclusions
In this paper, we propose a methodology that couples Moench and Hsieh’s solution (1985) with SA algorithm to identify three skin parameters (k1, Ss1, and dsk) and two
aquifer parameters (k2 and Ss2) simultaneously for a slug
test performed in a skin-affected confined aquifer system. The WWL data analyzed in this paper were generated by Moench and Hsieh’s solution (1985) and with a set of stan-dard normally distributed noise added for both positive-skin and negative-skin scenarios.
From the hypothetical Scenarios 1 to 8, the negative-skin scenarios as a whole give better estimated results than those of the positive-skin scenarios. The sensitivity analysis is performed to demonstrate the WWL response to the rela-tive change in aquifer parameters when the well of slug test has a positive skin. In the sensitivity plots, the influence of
k2on the WWL is insensitive at the initial period of the test.
Accordingly, analyzing a longer set of WWL data, though consumes much computing time, would result in a better estimation on k2 for aquifers with a positive skin. In
addi-tion, the parameters k1and rsin the case of a positive-skin
aquifer are highly correlated. Hence, one may misidentify the values of k1and dskby estimating them using an
optimi-zation technique, especially for positive-skin aquifers. Note that dsk(or rs) is usually taken as an input data in some
data-analyzed software but it is actually an invisible and unknown parameter. Impetuously presuming an arbitrary value of dsk
and only estimating the other four parameters may con-found the parameter estimation for the other four parame-ters. On the other hand, the results of sensitivity analysis reveal that whichever optimization technique is applied, the insensitivities of Ss1and Ss2usually lead to inaccurately
estimated results when determining the skin and aquifer parameters simultaneously. Nevertheless, analyzing the composite WWL data, which were obtained from the test and observation wells, could significantly improve the esti-mations of Ss1and Ss2.
Acknowledgements
This study was partly supported by the Taiwan National Sci-ence Council under the grant NSC94-2211-E-009-012. The
Table 6b The target values and estimated results for positive-skin scenarios with dskrange from 0.1085 to 1.3085
Case Estimated results
k1(m/s) k2(m/s) Ss1(1/m) Ss2(1/m) dsk(m) SEE Target value 1.00· 105 1.00· 104 1.00· 104 1.00· 104 0.1085 – 15a 2.31· 105 9.97· 105 1.71· 105 1.83· 105 0.493 3.16· 104 15b 3.16· 105 1.05· 104 7.25· 106 8.91· 106 1.324 7.57· 104 15c 2.47· 105 9.34· 105 1.55· 105 2.82· 105 0.629 6.99· 104 15d 2.99· 105 9.97· 105 6.71· 106 9.38· 106 1.015 7.75· 104 15e 1.92· 105 9.66· 105 2.93· 105 5.02· 105 0.346 7.07· 104 Mean 2.57· 105 9.89· 105 1.52· 105 2.30· 105 0.761 – RE 157.00% 1.12% 84.83% 77.00% 601.75% – Target value 1.00· 105 1.00· 104 1.00· 104 1.00· 104 0.6085 – 16a 1.04· 105 1.11· 104 9.10· 105 8.79· 105 0.684 2.67· 104 16b 1.10· 105 2.23· 104 7.70· 105 1.31· 105 0.824 7.84· 104 16c 1.11· 105 1.04· 104 7.43· 105 3.72· 105 0.768 8.64· 104 16d 1.17· 105 1.16· 104 6.25· 105 4.69· 105 0.914 8.02· 104 16e 1.13· 105 1.13· 104 7.25· 105 8.08· 105 0.851 8.07· 104 Mean 1.11· 105 1.33· 104 7.55· 105 5.32· 105 0.808 – RE 11.00% 33.40% 24.54% 46.82% 32.82% – Target value 1.00· 105 1.00· 104 1.00· 104 1.00· 104 1.3085 – 17a 1.01· 105 2.40· 104 9.69· 105 3.19· 105 1.402 2.98· 104 17b 9.95· 106 1.08· 104 9.99· 105 4.30· 105 1.241 6.63· 104 17c 1.00· 105 9.02· 104 9.99· 105 2.33· 105 1.463 8.76· 104 17d 1.02· 105 7.88· 105 9.37· 105 8.36· 105 1.348 8.84· 104 17e 1.00· 105 5.12· 104 9.97· 105 8.75· 105 1.401 7.42· 104 Mean 1.01· 105 3.68· 104 9.80· 105 5.39· 105 1.371 – RE 0.50% 268.16% 1.98% 46.14% 4.78% –
authors thank three anonymous reviewers for constructive comments and suggested revisions. Section ‘‘Sensitivity analysis’’ is written based on the third reviewer’s comments and valuable suggestion.
References
Bouwer, H., Rice, R.C., 1976. A slug test method for determining hydraulic conductivity of unconfined aquifers with completely or partially penetrating wells. Water Resour. Res. 12 (3), 423–428. Butler Jr., J.J., 1998. The Design Performance and Analysis of Slug
Tests. Lewis Publishers, Boca Raton, 252pp.
Cooper, H.H., Bredehoeft, J.D., Papadopulos, I.S., 1967. Response of a finite-diameter well to an instantaneous charge of water. Water Resour. Res. 3 (1), 263–269.
Corana, A., Marchesi, M., Martini, C., Ridella, S., 1987. Minimizing multimodal functions of continuous variables with the simulated annealing algorithm. ACM Trans. Math. Soft. 13 (3), 262–280. Cunha, M.D.C., Sousa, J., 1999. Water distribution network design
optimization: simulated annealing approach. J. Water Resour. Plan. Manage. ASCE 125 (4), 215–221.
Dougherty, D.E., Marryott, R.A., 1991. Optimal groundwater man-agement 1. Simulated Annealing. Water Resour. Res. 27 (10), 2493–2508.
Duffield, G.M., 2002. AQTESOLV for Windows User’s Guide. Hydro-SOLVE, Inc., Reston, VA.
Faust, C.R., Mercer, J.W., 1984. Evaluation of slug tests in wells containing a finite-thickness skin. Water Resour. Res. 20 (4), 504–506.
Freeze, R.A., Cherry, J.A., 1979. Groundwater. Prentice-Hall, New Jersey, p. 29.
Hvorslev, M.J., 1951. Time lag and soil permeability in ground-water observations. Bull. No. 36, Waterways Experimental Station Corps of Engineers, U.S. Army, Vicksburg, MI, pp. 1–50. Hyder, Z., Butler Jr., J.J., McElwee, C.D., Liu, W.Z., 1994. Slug tests in partially penetrating wells. Water Resour. Res. 30 (11), 2945–2957.
IMSL, 1997. Stat/Library. Volume 2, Visual Numerics, Inc, Houston, TX.
Kabala, Z.J., 2001. Sensitivity analysis of a pumping test on a well with wellbore storage and skin. Adv. Water Res. 24 (5), 483–504. Kirkpatrick, S., Gelatt Jr., C.D., Vecchi, M.P., 1983. Optimization
by simulated annealing. Science 220 (4598), 671–680.
Kuo, S.F., Liu, C.W., Merkley, G.P., 2001. Application of the simulated annealing method to agricultural water resource management. J. Agric. Eng. Res. 80 (1), 109–124.
Marryott, R.A., Dougherty, D.E., Stollar, R.L., 1993. Optimal groundwater management 2. Application of simulated annealing to a field-scale contamination site. Water Resour. Res. 29 (4), 847–860.
Metropolis, N., Rosenbluth, A.W., Rosenbluth, M.N., Teller, A.H., Teller, E., 1953. Equation of state calculations by fast comput-ing machines. J. Chem. Phys. 21 (6), 1087–1092.
Moench, A.F., Hsieh, P.A., 1985. Analysis of slug test data in a well with finite-thickness skin. In: Memoirs of the 17th International Congress on the Hydrogeology of Rocks of Low Permeability, U.S.A. Members of the International Association of Hydrologists, Tucson, AZ, vol. 17, pp. 17–29.
Pardo-Iguzquiza, E., 1998. Optimal selection of number and location of rainfall gauges for areal rainfall estimation using geostatistics and simulated annealing. J. Hydrol. (210), 206– 220.
Pham, D.T., Karaboga, D., 2000. Intelligent Optimisation Tech-niques. Springer, Great Britain.
Ramey Jr., H.J., Agarwal, R.G., Martin, I., 1975. Analysis of ‘‘Slug Test’’ or DST flow period data. J. Can. Petrol Technol. 14 (3), 37–47.
Rovey II, C.W., 1998. Digital simulation of the scale effect in hydraulic conductivity. Hydrogeol. J. 6, 216–225.
Schwartz, F.W., Zhang, H., 2003. Fundamentals of Ground Water. John Wiley & Sons, Inc., New York, p. 73.
Springer, R.K., Gelhar L.W., 1991. Characterization of large-scale aquifer heterogeneity in glacial outwash by analysis of slug tests with oscillatory response. Cape Cod, MA, U.S. Geological Survey Water Reserch Investigation Report 91 (4034), pp. 36–40. Yang, Y.J., Gates, T.M., 1997. Wellbore skin effect in slug-test data
analysis for low-permeability geologic materials. Ground Water 35 (6), 931–937.
Yeh, H.D., 1987. Theis’ solution by nonlinear least-squares and finite- difference Newton’s method. Ground Water 25 (6), 710– 715.
Yeh, H.D., Yang, S.Y., 2006. A novel analytical solution for a slug test conducted in a well with a finite-thickness skin. Adv. Water Resour. 29 (10), 1479–1489.
Zheng, C., Wang, P., 1996. Parameter structure identification using tabu search and simulated annealing. Adv. Water Resour. 19 (4), 215–224.
