Two-dimensional Laplace-transformed power series solution
for solute transport in a radially convergent flow field
Jui-Sheng Chen
a, Chen-Wuing Liu
b,*, Chung-Min Liao
baDepartment of Environmental Engineering and Sanitation, Foo-Yin University, Kaohsiung 831, Taiwan, ROC bDepartment of Bioenvironmental Systems Engineering, National Taiwan University, Taipei 10617, Taiwan, ROC
Received 16 July 2002; received in revised form 3 April 2003; accepted 21 May 2003
Abstract
This paper presents an analytical solution for two-dimensional non-axisymmetric solute transport in a radially convergent flow field. We applied a Laplace-transformed power series (LTPS) technique to solve the two-dimensional advection-dispersion equation in cylindrical coordinates. The solution is compared with a numerical solution to evaluate its robustness and accuracy. The ap-plicable Peeclet number range of the developed power series solution is also examined. Results show that the LTPS technique can effectively and accurately handle the two-dimensional radial advection-dispersion equation for a Peeclet number up to 60. The two-dimensional power series solution is appropriate for hydrogeologic circumstances where temporally and spatially continuous so-lutions are demanded.
Ó 2003 Elsevier Ltd. All rights reserved.
Keywords: Radially convergent flow field; Laplace-transformed power series; Solute transport
1. Introduction
Tracer tests attempt to determine solute transport parameters, such as aquifer porosity, dispersion tensor and hydrogeological properties. A radially convergent tracer test facilitates the recovery of the injected mass, reduces the effect of apparent dispersion due to the flow field, and minimizes the influence of the natural hy-draulic gradient. Thus, radially convergent tracer tests are particularly useful where transport, rather than hy-draulic properties, is desired. Predictive and interpretive models are available for radially convergent tracer tests [1,2,12–15,19,20]. These studies are commonly limited to the analysis of breakthrough curves in the extraction well, although new sampling technologies are available for concentration measurements at arbitrary points in the field [11]. Because of the model complexity, the analysis is often reduced to an adjusted one-dimensional solution [17] instead of a more accurate two-dimensional transport model [22].
The one-dimensional model generally considers longitudinal dispersion only. Solutions that consider transverse dispersion of points between extraction and injection wells are not available for tracer test analysis. Unlike a divergent flow tracer test, the convergent flow one does not have axial symmetry and the transverse dispersion could be important. The magnitude of the transverse dispersion coefficient influences both the re-gion to which the pollution is extended, and the intensity of the pollution. In view of the importance of deter-mining transverse dispersion, Chen et al. [3] presented a two-dimensional mathematical model in cylindrical co-ordinates that was solved via Laplace transform finite-difference method to illustrate non-axisymmetric tracer transport in a radially convergent flow field. In addition, a curve-fitting method involving a theoretical break-through curve at an intermediate point was proposed to evaluate transverse dispersivity, a quantity could not be determined from a one-dimensional model. The nu-merical solution developed by Chen et al. [3] has the propensity to be computationally cumbersome and ex-pensive and, in the absence of interpolation, can be used to yield breakthrough curves only at discrete nodal points. An analytical solution can inherently provide a spatially continuous concentration distribution, thus is
*
Corresponding author. Tel.: +886-2-2362-6480; fax: +886-2-2363-9557.
E-mail address:[email protected](C.-W. Liu).
0309-1708/$ - see front matter Ó 2003 Elsevier Ltd. All rights reserved. doi:10.1016/S0309-1708(03)00090-3
extremely useful in generating theoretical breakthrough curves at any observation point. Beyond its obvious importance in the studying solute transport for a tracer test in a radially convergent flow field, the radially convergent dispersion problem is distinguished by its being probably the simplest case for which the disper-sion coefficient is a function of spatially varying velocity field. Accordingly, the two-dimensional analytical solu-tion to the radially convergent dispersion problem can be a valuable means of evaluating the accuracy of computer codes that simulate two-dimensional solute transport in porous medium.
It is difficult to derive the analytical solution of the two-dimensional advection-dispersion equation in cy-lindrical coordinates due to the dependence of both the longitudinal and transverse dispersion coefficients on the spatially varying velocities. To our knowledge, the two-dimensional analytical solution to the radial advection-dispersion equation in the real or Laplace domains is not currently available.
In this paper the problem is approached analytically to yield temporally and spatially continuous solutions. We used a Laplace-transformed power series (LTPS) technique [4] to solve the two-dimensional radial ad-vection-dispersion equation in cylindrical coordinates. The analytical power series solution will be compared to the numerical solution of Chen et al. [3] to examine its robustness and accuracy and to determine the applicable ranges of the Peeclet number. Moreover, the mathemat-ical characteristics of the solutions are analyzed to il-lustrate the mathematical and convergence behavior of the power series functions.
2. Problem formulation
In this paper the problem of two-dimensional tracer transport in a radially convergent flow field is consid-ered. Fig. 1 presents the conceptual configuration. For
the sake of simplicity, the following assumptions are made:
1. A vertically oriented extraction well of finite diameter is located along the vertical axis and fully penetrates a homogeneous and isotropic aquifer of constant thick-ness.
2. A steady state, horizontal flow field, radially conver-gent and axially symmetric with respect to extraction well, is established prior to the start of tracer injec-tion.
3. The tracer injection has no influence on the flow field. Based on above assumptions the seepage velocity V caused by extraction is described by
V ¼ A
r ð1Þ
where A¼ Q=2pb/ and Q is volumetric rate of water withdrawal from the extraction well. The symbols b and /represent the aquifer thickness and effective porosity, respectively, while r is the radial distance from the ex-traction well.
The two-dimensional advection-dispersion equation in a cylindrical coordinates is [22] aLA r o2C o2rþ A r oC or þ aTA r3 o2C o2h ¼ R oC ot ð2Þ
where aL and aT denote the longitudinal and transverse dispersivities, respectively, R is a retardation factor, t is time and C is the solute concentration.
The aquiferÕs initial tracer concentration is assumed to be zero before starting the test:
Cðr; h; 0Þ ¼ 0 ð3Þ
The outlet boundary condition describes the solute transport between the extraction well and aquifer. The condition which takes into account the finite mixing volume effect in the well bore is
pr2 WhW oCWðtÞ ot ¼ 2prW/baL A rW oCðr; h; tÞor at r¼ rW ð4Þ where rW is the radius of the extraction well, hW is its mixing length and CWðtÞ denotes the tracer concentra-tion in the well [9,12].
The initial condition for (4) states that the well bore contains no contaminant before pumping:
CWðt ¼ 0Þ ¼ 0: ð5Þ
Additionally, it is reasonable to assume here that CWðtÞ ¼ CðrW;h; tÞ.
Zlotnik and Logan [22] considered flow and transport in a ring-shaped domain centered at the extraction well and bounded by circles of radii r¼ rI and r¼ rL l ðl rLÞ in deriving the boundary condition at
W
L
the injection well. The physical assumption is that ad-vective transport dominates dispersive transport at a small distance l 5rIdownstream in the discharge zone of the injection well. Therefore, the boundary condition including solute transport and the well bore effect at the injection well was formulated as [22]
aL A r oCðr; h; tÞ or þ A rCðr; h; tÞ ¼ ArCIðtÞ p d < h < p 0 0 < h < p d r¼ r r L ð6Þ
where CIðtÞ is the concentration generated in the injec-tion well and transported downstream through the narrow and short (a few well diameters) discharge zone by advection. This small zone with advection-dominated flow has an aperture angle of 2d. The aperture angle of this narrow zone at this distance r r
L l from the center of extraction well is
2d¼2arI rL
ð7Þ where a is a factor that defines the distortion of distances between the two most separated stream lines that enter (or leave) the injection well. This parameter can also depend on the skin effect for an injection well (Zlotnik and Logan [22]; Eq. (3)). For a uniform isotropic aquifer with a well without a skin, a¼ 2. For skins with high conductivity, 2 6 a 6 4, whereas for skins with low conductivity, 0 < a 6 2 [7].
It remains to determine the effluent concentration from the injection well in an ambient horizontal flow. Effluent concentration from the well with initial dis-solved tracer mass M0satisfies a mass balance equation for the tracer in the borehole, namely
2arI/bjV ðrLÞjCI¼ pr2IhI dCI dt ð8Þ CIð0Þ ¼ M0 pr2 IhI ¼ C0 ð9Þ
where hI is the mixing length of injection well (Zlotnik and Logan [22]; Eq. (6)).
After integration, the known effluent concentration CIðtÞ can be substituted into boundary condition (6). The physics of the problem stipulates that C is a single-valued function in r and h coordinates. In addition, C is a continuous and symmetrical function across h¼ 0 and h¼ p. Thus, boundary conditions in the transverse directions are as follows:
oCðr; 0; tÞ
oh ¼ 0 ð10Þ
oCðr; p; tÞ
oh ¼ 0 ð11Þ
Dimensionless variables are defined in a manner similar to that used by Chen et al. [3]. Following Mo-ench [12] and substituting the definition given in Table 1 into Eq. (2), the dimensionless transport equation is presented in the following form:
1 Pe 1 rD o2C or2 D þ 1 rD oC orD þ 1 Pe aD r3 D o2C oh2 ¼ 2R 1 r2 WD oC otD ð12Þ Consequently, the initial and boundary conditions (3)– (11) become CðrD;h;0Þ ¼ 0 ð13Þ lW oCWðtDÞ otD ¼ 1 Pe oCðrD;h; tDÞ orD at r¼ rW ð14Þ CWð0Þ ¼ 0 ð15Þ 1 Pe oCðrD;h; tDÞ orD þ CðrD;h; tDÞ ¼ CIðtDÞ p d < h < p 0 0 < h < p d ð16Þ CI¼ lI dCI dtD ð17Þ CIð0Þ ¼ C0 ð18Þ oCðrD;0; tDÞ oh ¼ 0 ð19Þ oCðrD;p; tDÞ oh ¼ 0 ð20Þ
This study adopts the LTPS technique to solve the governing equation and boundary conditions (12)–(20). The LTPS has been successfully applied to solve the one-dimensional radial partial differential equation with variable-dependent coefficients [4]. The LTPS technique avoids the need to derive the full analytical solution. Table 1
Dimensionless parameters used in this study
Dimensionless quantity Expression Time tD¼ Qt ph/ðr2 L r 2 WÞ Distance rD¼ r rL
Extraction well radius rWD¼
rW
rL
Injection well radius rID¼
rI rL Peeclet number Pe¼rL aL Transverse dispersivity X¼aT aL
Extraction well mixing factor lW¼
r2 WhW
/hðr2 L r2WÞ
Injection well mixing factor lI¼
rIrLhI
/hðr2 L rW2Þ
This technique is parsimonious and easy to code into a program. The detailed application of the LTPS tech-nique is presented in Appendix A.
3. Results and discussion
3.1. Verification of power series solution
The obtained solution is compared with the numeri-cal solution from the Laplace-transformed finite differ-ence (LTFD) method [5] to demonstrate the accuracy of the LPTS solution. A hypothetical tracer test is defined for the purpose of the comparison. Table 2 provides a
summary of simulation conditions and transport pa-rameters for the hypothetical tracer test. Figs. 2–4 plot breakthrough curves at an observation well located at r¼ 5 m, h ¼ p using various Peeclet numbers, and compares them to the numerical solution of Chen et al. Table 2
Descriptive simulation conditions and transport parameters of the hypothetical tracer test
Parameter Test 1
Pumping rate (Q), m3min1 2
Aquifer thickness (h), m 10
Effective porosity (/), dimensionless 0.2 Radius of extraction well (rW), m 0.1
Extraction well mixing length (hW), m 10
Radius of injection well (rI), m 0.1
Injection well mixing length (hI), m 10
Distance to the injection well (rL), m 25
Injected mass (M), Kg 40
Longitudinal dispersivity (aL), m 25, 2.5, 0.42
Peeclet Number Pe, dimensionless 1, 10, 60 Transverse dispersivity (aT), m 5, 0.5, 0.084
Dimensionless ratio of dispersivity (X ), dimensionless 0.2 0 1000 2000 3000 4000 t [min] 0.00 0.02 0.04 0.06 0.08 Concentration [kg /m 3 ] this study Chen et al. [5] Pe=1 2 5 10
Fig. 2. Comparison of breakthrough curves at an observation well located at r¼ 5 [m], h ¼ p between the power series solution and numerical solution for Peeclet numbers of 1, 2, 5 and 10 in a radially convergent tracer test.
0 1000 2000 3000 4000 t [min] 0.00 0.10 0.20 0.30 0.40 0.50 Concentration [kg /m 3 ] 80 60 40 Pe=20 this study Chen et al. [5]
Fig. 3. Comparison of breakthrough curves at an observation well located at r¼ 5 [m], h ¼ p between the power series solution and numerical solution for Peeclet numbers of 20, 40, 60 and 80 in a radially convergent tracer test.
0 1000 2000 3000 4000 t [min] 0.00 0.40 0.80 1.20 Concentration [kg /m 3 ] Pe=100 150 200 this study Chen et al. [5]
Fig. 4. Comparison of breakthrough curves at an observation well located at r¼ 5 [m], h ¼ p between the power series solution and numerical solution for Peeclet numbers of 100, 150 and 200 in a radially convergent tracer test.
[5]. The extraction well mixing and injection well mixing factors are set to zero so that well bore mixing exerts no influence on the breakthrough curves. Concentra-tion breakthrough curves obtained from the power series solution agree well with those obtained from the numerical solution for Peeclet numbers smaller than 60. Additionally, the comparison reveals that the break-through curves obtained with LTPS are shaped roughly the same as those computed with LTFD for rising limbs and spreading tails for Peeclet numbers greater than 60, but a noticeable discrepancy occurs between the two solutions at the peak concentrations of the break-through curves for Peeclet numbers greater than 60.
The comparison of the two solutions reveals that the LPTS technique can effectively and accurately solve the two-dimensional, radial advection-dispersion equation. However, it is found that the power series method did not match the numerical solution for Peeclet numbers greater than 60. Such conditions correspond to an ad-vective-dominated solute transport and yield break-through curves of steep fronts. Transport in a single fracture or in a particularly homogeneous granular aquifer may involve large Peeclet numbers [18]. Efforts to determine solution correctness for Peeclet numbers greater than 60 will be addressed in a future study.
3.2. Mathematical behavior of power series functions The LPTS solution derived from the Appendix A includes two new functions. These two functions are in the form of infinite series, terms of which can be straightforwardly evaluated. We have to consider, however, the number of terms needed to produce ac-curate results. Chen et al. [4] illustrated that, for a fixed tolerance, the required numbers to be summed generally increase with the Peeclet numbers.
We have performed the computation for various rD and n at fixed s¼ 3 of the two new functions to examine the mathematical characteristics of the two new infinite series (some part of computational results are provided in the Appendix B). The values of the function Z1ðrD; n; sÞ increase with increasing rDand decrease with increasing n. The function values increase as Peeclet number increases. The values of oZ1ðrD; n; sÞ=orD in-crease as rD increases and decrease as n increases. For a Peeclet number of unity, however, the values of oZ1ðrD; n; sÞ=orD increase with n at rD¼ 1. The plot of Z2ðrD; n; sÞ versus rD behaves differently from that of Z1ðrD; n; sÞ : Z2ðrD; n; sÞ decreases with increasing rD and behaves irregularly with respect to n. For example, Z2ðrD; n; sÞ increases with n at rD¼ 0:2, however, no
Table 3
The values of Z2ðrD; n; sÞ for various rDand n at s¼ 3 (Pe ¼ 1)
n rD¼ 0:2 rD¼ 0:4 rD¼ 0:6 rD¼ 0:8 rD¼ 1:0 0 2.011 1012 4.156 1012 6.765 1012 1.045 1012 1.633 1012 1 3.298 1012 5.563 1012 8.165 1012 1.195 1012 1.819 1012 2 8.155 1012 9.896 1012 1.286 1012 1.837 1012 2.852 1012 3 2.431 1012 3.482 1012 7.013 1012 1.412 1012 2.725 1012 4 5.608 1012 2.439 1012 1.162 1012 1.383 1012 )1.138 1012 5 1.537 1012 4.535 1012 1.964 1012 6.649 1012 )5.521 1012 6 4.218 1012 8.213 1012 3.065 1012 1.597 1012 1.308 1012 7 1.160 1012 1.468 1012 4.248 1012 1.749 1012 1.043 1012 8 3.193 1012 2.623 1012 5.892 1012 1.895 1012 5.321 1012 9 8.802 1012 4.685 1012 8.218 1012 2.301 1012 8.057 1012 10 2.428 1012 8.363 1012 1.141 1012 2.689 1012 8.318 1012 11 6.702 1012 1.493 1012 1.583 1012 3.132 1012 8.608 1012 12 1.850 1012 2.664 1012 2.193 1012 3.636 1012 8.733 1012 13 5.111 1012 4.756 1012 3.037 1012 4.216 1012 8.844 1012 14 1.412 1012 8.488 1012 4.204 1012 4.882 1012 8.940 1012 15 3.901 1012 )1.006 1012 )5.975 1012 )1.092 1012 )1.065 1012 16 1.079 1012 2.704 1012 8.045 1012 6.529 1012 9.091 1012 17 2.982 1012 4.826 1012 1.113 1012 7.544 1012 9.152 1012 18 8.243 1012 8.613 1012 1.539 1012 8.713 1012 9.205 1012 19 2.279 1012 1.537 1012 2.127 1012 1.006 1012 9.252 1012 20 6.302 1012 2.744 1012 2.941 1012 1.161 1012 9.294 1012 21 1.743 1012 4.898 1012 4.065 1012 1.340 1012 9.331 1012 22 4.820 1012 8.742 1012 5.618 1012 1.546 1012 9.365 1012 23 1.333 1012 1.560 1012 7.765 1012 1.783 1012 9.396 1012 24 3.687 1012 2.785 1012 1.073 1012 2.056 1012 9.423 1012 25 1.020 1012 4.972 1012 1.483 1012 2.371 1012 9.448 1012 26 2.821 1012 8.874 1012 2.050 1012 2.733 1012 9.472 1012 27 7.803 1012 1.584 1012 2.832 1012 3.151 1012 9.493 1012 28 2.159 1012 2.828 1012 3.914 1012 3.632 1012 9.513 1012 29 5.971 1012 5.047 1012 5.408 1012 4.187 1012 9.531 1012 30 1.652 1012 9.009 1012 7.472 1012 4.826 1012 9.548 1012
definite rule can be concluded for Z2ðrD; n; sÞ versus n at other values of rD. Notably, an abnormally large value appears at n¼ 15 for Z2ðrD; n; sÞ and its first derivative, oZ2ðrD; n; sÞ=orD. We have scrutinized thoroughly the computer program and confirm that this extraordinary value dose not result from any computation errors but from the inherent mathematical characteristic of the new power series function. Understanding the characteristic of abnormality involves a complicated mathematical analysis, which we do not pursue here. Tables 3–8 also show that the behavior of the developed two functions are fundamentally different from that of the special Airy function solution for the one-dimensional radial advec-tion-dispersion equation.
4. Conclusions
This work presents a novel Laplace transform power series method to solve the two-dimensional variable-dependent advection-dispersion differential equation in cylindrical coordinates, accounting for non-axisymmet-rical transport during a convergent radial tracer test. The exact analytical solution in Laplace transform
do-main is derived. The new solution is compared with the numerical solution to verify its accuracy. Results show that the two solutions produce the same results for Peeclet numbers smaller than 60. Notably, it is found that the two solutions yield a noticeable difference at an in-termediate time for Peeclet numbers greater than 60. We suggest further mathematical analysis and verification of the new solution for Peeclet numbers greater than 60 are needed. The novel two-dimensional power series solu-tion is useful for quantitative hydrogeologic problems where temporally and spatially continuous solutions are demanded.
Acknowledgements
The authors would like to acknowledge the National Science Council of the Republic of China for financially supporting this work under Grant no. NSC 90-2313-B-242-004. Dr. G. Cassiani and the other anonymous reviewer are also appreciated for their constructive suggestions to improve this paper.
Table 4
The values of Z2ðrD; n; sÞ for various rDand n at s¼ 3 (Pe ¼ 10)
n rD¼ 0:2 rD¼ 0:4 rD¼ 0:6 rD¼ 0:8 rD¼ 1:0 0 2.440 1012 9.286 1012 3.875 1012 1.897 1012 1.077 1012 1 4.296 1012 1.542 1012 6.503 1012 3.219 1012 1.840 1012 2 1.941 1012 8.793 1012 4.112 1012 2.122 1012 1.239 1012 3 4.855 1012 2.732 1012 1.468 1012 8.077 1012 4.876 1012 4 5.203 1012 7.234 1012 3.858 1012 2.283 1012 1.442 1012 5 1.405 1012 1.103 1012 6.495 1012 4.234 1012 2.827 1012 6 3.871 1012 1.379 1012 8.004 1012 5.801 1012 4.135 1012 7 1.079 1012 2.285 1012 1.463 1012 1.198 1012 9.194 1012 8 3.001 1012 1.814 1012 )1.902 1012 )2.072 1012 )1.729 1012 9 8.350 1012 3.677 1012 9.640 1012 5.640 1012 5.049 1012 10 2.319 1012 6.750 1012 1.534 1012 1.093 1012 1.079 1012 11 6.435 1012 1.237 1012 3.938 1012 4.288 1012 4.693 1012 12 1.784 1012 2.232 1012 1.105 1012 )4.840 1012 )6.070 1012 13 4.946 1012 4.049 1012 1.970 1012 )6.879 1012 )1.188 1012 14 1.370 1012 7.324 1012 2.880 1012 2.646 1012 4.413 1012 15 )8.670 1012 )1.007 1012 )7.456 1012 )1.996 1012 )3.356 1012 16 1.051 1012 2.383 1012 5.792 1012 1.844 1012 )3.015 1012 17 2.911 1012 4.290 1012 8.185 1012 3.655 1012 )1.010 1012 18 8.061 1012 7.713 1012 1.153 1012 4.720 1012 )5.154 1012 19 2.232 1012 1.386 1012 1.620 1012 6.011 1012 3.659 1012 20 6.179 1012 2.488 1012 2.273 1012 6.976 1012 2.411 1012 21 1.710 1012 4.464 1012 3.184 1012 8.250 1012 2.502 1012 22 4.735 1012 8.005 1012 4.454 1012 9.733 1012 2.652 1012 23 1.311 1012 1.435 1012 6.223 1012 1.146 1012 5.241 1012 24 3.628 1012 2.571 1012 8.686 1012 1.347 1012 4.681 1012 25 1.004 1012 4.606 1012 1.211 1012 1.580 1012 4.761 1012 26 2.780 1012 8.248 1012 1.688 1012 1.851 1012 4.878 1012 27 7.693 1012 1.477 1012 2.350 1012 2.166 1012 5.018 1012 28 2.129 1012 2.643 1012 3.271 1012 2.531 1012 5.136 1012 29 5.894 1012 4.730 1012 4.550 1012 2.955 1012 5.252 1012 30 1.631 1012 8.463 1012 6.325 1012 3.447 1012 5.363 1012
Appendix A
First, proceeding with the Laplace transform of Eq. (12) and its associated boundary conditions (13)–(20), with respect to tD, we obtain
1 Pe 1 rD o2G or2 D þ 1 rD oG orD þ 1 Pe aD r3 D o2G oh2 ¼ 2R 1 r2 WD sG ðA:1Þ lWsG¼ 1 Pe oGðrD;h; sÞ orD at r¼ rW ðA:2Þ 1 Pe oGðrD;h; sÞ orD þ GðrD;h; sÞ ¼ GIðsÞ p d < h < p 0 0 < h < p d at r¼ rI ðA:3Þ GIðsÞ ¼ lI½sGIðsÞ C0 ðA:4Þ oGðrD;0; sÞ oh ¼ 0 ðA:5Þ oGðrD;p; sÞ oh ¼ 0 ðA:6Þ
where s denotes the Laplace transform parameter and G represents the Laplace transform of C, as defined by GðrD;h; sÞ ¼ Z 1 0 CðrD;h; tDÞ estDdtD ðA:7Þ GIðsÞ ¼ Z 1 0 CIðtDÞ estDdtD ðA:8Þ
The finite Fourier cosine transform with respect to h of (A.1)–(A.6) gives 1 Pe 1 rD d2W dr2 D þ 1 rD dW drD 1 Pe aDn2 r3 D þ 2R 1 r2 WD s W ¼ 0 ðA:9Þ lWsW ¼ 1 Pe oWðrD; n; sÞ orD at r¼ rW ðA:10Þ 1 Pe oWðrD; n; sÞ orD þ W ðrD; n; sÞ ¼ F ðnÞ ðA:11Þ where FðnÞ ¼ GGIðsÞd; n¼ 0 IðsÞ ð1Þ nsin nd n h i ; n¼ 1; 2; 3 . . . ( GIðsÞ ¼ lIC0 1þ lIs Table 5
The values of Z2ðrD; n; sÞ for various rDand n at s¼ 3 (Pe ¼ 60)
n rD¼ 0:2 rD¼ 0:4 rD¼ 0:6 rD¼ 0:8 rD¼ 1:0 0 7.424 1012 4.172 1012 2.934 1012 2.562 1012 2.760 1012 1 2.148 1012 1.229 1012 8.689 1012 7.608 1012 8.208 1012 2 4.393 1012 2.654 1012 1.908 1012 1.683 1012 1.824 1012 3 )7.808 1012 )5.157 1012 )3.809 1012 )3.405 1012 )3.717 1012 4 1.395 1012 1.040 1012 7.976 1012 7.259 1012 8.007 1012 5 1.265 1012 1.096 1012 8.823 1012 8.217 1012 9.186 1012 6 )1.725 1012 )1.786 1012 )1.524 1012 )1.460 1012 )1.659 1012 7 )2.387 1012 )3.023 1012 )2.762 1012 )2.734 1012 )3.166 1012 8 9.386 1012 1.485 1012 1.466 1012 1.508 1012 1.785 1012 9 )6.929 1012 )1.406 1012 )1.514 1012 )1.624 1012 )1.971 1012 10 )1.634 1012 )4.306 1012 )5.102 1012 )5.739 1012 )7.161 1012 11 )3.327 1012 )1.159 1012 )1.523 1012 )1.805 1012 )2.321 1012 12 )1.131 1012 )5.585 1012 )8.209 1012 )1.029 1012 )1.368 1012 13 )7.329 1012 )6.020 1012 )9.969 1012 )1.328 1012 )1.830 1012 14 3.234 1012 )2.984 1012 )5.606 1012 )7.966 1012 )1.141 1012 15 )4.565 1012 )5.393 1012 )1.157 1012 )1.762 1012 )2.630 1012 16 5.215 1012 )6.696 1012 )1.652 1012 )2.705 1012 )4.218 1012 17 1.607 1012 )1.621 1012 )4.624 1012 )8.176 1012 )1.335 1012 18 4.624 1012 2.896 1012 9.607 1012 1.841 1012 3.156 1012 19 1.313 1012 )1.323 1012 )5.129 1012 )1.069 1012 )1.929 1012 20 3.739 1012 )9.048 1012 )4.123 1012 )9.380 1012 )1.785 1012 21 1.060 1012 4.659 1012 2.506 1012 6.243 1012 1.256 1012 22 2.999 1012 1.296 1012 8.262 1012 2.261 1012 4.821 1012 23 8.470 1012 )6.716 1012 )5.101 1012 )1.539 1012 )3.484 1012 24 2.388 1012 5.900 1012 5.296 1012 1.766 1012 4.255 1012 25 6.722 1012 1.113 1012 1.199 1012 4.430 1012 1.138 1012 26 1.890 1012 1.467 1012 1.893 1012 7.775 1012 2.135 1012 27 5.307 1012 5.037 1012 1.470 1012 6.730 1012 1.979 1012 28 1.489 1012 1.808 1012 2.005 1012 1.025 1012 3.235 1012 29 4.173 1012 2.061 1012 1.515 1012 8.677 1012 2.943 1012 30 1.168 1012 3.187 1012 1.651 1012 1.062 1012 3.878 1012
where n denotes the finite Fourier cosine transform pa-rameter and WðrD; n; sÞ represents the finite Fourier cosine transform of GðrD;h; sÞ, as defined by
WðrD; n; sÞ ¼ Z p
0
GðrD;h; sÞ cosðnhÞ dh ðA:12Þ
Such a transform is advantageous in that the inversion is directly given by the following formula [16]:
GðrD;h; sÞ ¼ 1 pWðrD;0; sÞ þ 2 p X1 n¼1 WðrD; n; sÞ cosðnhÞ ðA:13Þ The transformed ordinary differential equation (A.9) can be directly solved using the power series method. The governing equation (A.9), however, involves an advection term that generally leads to a numerical error in numerical calculation for large Peeclet numbers. Chen et al. [4] suggested that the power series technique, with a change of variables to eliminate the first spatial de-rivative, is capable of solving the radial dispersion problem over a wide range of Peeclet numbers. There-fore, a change of variables is used to eliminate the first spatial derivative and convert Eq. (A.9) to
1 Pe d2Z dr2 D Pe 4 þ 1 Pe Xn2 r2 D þ 2rDR 1 r2 WD Z¼ 0 ðA:14Þ where Z ¼ exp Pe 2 ðrD 1Þ W
The boundary conditions (A.10) and (A.11) are in terms of Z, 1 Pe dZ drD 1 2 þ lWs Z¼ 0 ðA:15Þ and 1 Pe dWðrD; n; sÞ drD þ W ðrD; n; sÞ ¼ F ðnÞ; ðA:16Þ respectively.
The governing equation (A.14) is amenable to solu-tion using the power series method [8]. While using the power series method to solve the transformed differen-tial equation (A.14), one must consider that whether the variable-dependent coefficients of governing equation Table 6
The values ofoZ2ðrD; n; sÞ=orDfor various rDand n at s¼ 3 (Pe ¼ 1)
n rD¼ 0:2 rD¼ 0:4 rD¼ 0:6 rD¼ 0:8 rD¼ 1:0 0 1.021 1012 1.151 1012 1.507 1012 2.270 1012 3.766 1012 1 )2.070 1012 1.485 1012 1.706 1012 3.756 1012 7.102 1012 2 )1.487 1012 7.485 1012 7.594 1012 1.538 1012 2.790 1012 3 )7.842 1012 )1.633 1012 )8.785 1012 )8.487 1012 )1.184 1012 4 )2.959 1012 )3.581 1012 )1.120 1012 )6.144 1012 )5.904 1012 5 )1.045 1012 )7.988 1012 )1.636 1012 )2.491 1012 5.105 1012 6 )3.521 1012 )1.746 1012 )2.954 1012 )7.569 1012 )1.180 1012 7 )1.149 1012 )3.681 1012 )4.795 1012 )9.177 1012 2.268 1012 8 )3.665 1012 )7.589 1012 )7.751 1012 )1.508 1012 )4.290 1012 9 )1.148 1012 )1.537 1012 )1.216 1012 )1.969 1012 )4.728 1012 10 )3.550 1012 )3.069 1012 )1.881 1012 )2.540 1012 )5.194 1012 11 )1.085 1012 )6.062 1012 )2.878 1012 )3.247 1012 )5.798 1012 12 )3.288 1012 )1.186 1012 )4.364 1012 )4.112 1012 )6.417 1012 13 )9.888 1012 )2.304 1012 )6.567 1012 )5.167 1012 )7.015 1012 14 )2.955 1012 )4.445 1012 )9.815 1012 )6.448 1012 )7.622 1012 15 )8.781 1012 )7.276 1012 4.898 1012 6.787 1012 5.366 1012 16 )2.597 1012 )1.629 1012 )2.158 1012 )9.880 1012 )8.844 1012 17 )7.650 1012 )3.097 1012 )3.177 1012 )1.215 1012 )9.460 1012 18 )2.245 1012 )5.867 1012 )4.662 1012 )1.487 1012 )1.008 1012 19 )6.568 1012 )1.108 1012 )6.816 1012 )1.815 1012 )1.070 1012 20 )1.916 1012 )2.085 1012 )9.935 1012 )2.207 1012 )1.132 1012 21 )5.573 1012 )3.915 1012 )1.444 1012 )2.678 1012 )1.194 1012 22 )1.617 1012 )7.334 1012 )2.094 1012 )3.240 1012 )1.256 1012 23 )4.684 1012 )1.371 1012 )3.030 1012 )3.911 1012 )1.318 1012 24 )1.354 1012 )2.556 1012 )4.376 1012 )4.712 1012 )1.381 1012 25 )3.905 1012 )4.759 1012 )6.307 1012 )5.666 1012 )1.443 1012 26 )1.125 1012 )8.845 1012 )9.074 1012 )6.800 1012 )1.506 1012 27 )3.235 1012 )1.641 1012 )1.304 1012 )8.149 1012 )1.568 1012 28 )9.290 1012 )3.041 1012 )1.870 1012 )9.750 1012 )1.631 1012 29 )2.664 1012 )5.628 1012 )2.678 1012 )1.165 1012 )1.694 1012 30 )7.631 1012 )1.040 1012 )3.832 1012 )1.390 1012 )1.756 1012
are analytical or not. There are two different type power series solutions for the governing equation (A.14): Case 1, the coefficient is analytical for n¼ 0; and Case 2, the coefficient is not analytical for n > 0. Thus, the solution is obtained by direct application of the power series method for Case 1. The other solution is yielded using the extended power series method, or the Frobenius method, for Case 2 [10,21]. The procedures for the two approaches are described as follows:
Case 1: n¼ 0, the coefficient is analytical
Rewriting the governing equation (A.14) for n¼ 0 yields 1 Pe d2Z dr2 D Pe 4 þ 2rDR 1 r2 WD Z¼ 0 ðA:17Þ
First, the solution of Eq. (A.17) is assumed in the form of a power series with unknown coefficients,
ZðrD;0; sÞ ¼ X1 m¼0
amrmD ðA:18Þ
Termwise differentiation gives
dZðrD;0; sÞ drD ¼X 1 1 mamrm1D ðA:19Þ dZðrD;0; sÞ dr2 D ¼X 1 2 mðm 1Þamrm2D ðA:20Þ
These expressions are substituted into Eq. (A.17), yielding 1 Pe X1 2 mðm 1Þamrm2D Pe 4 þ 2Rs 1 r2 WD rD X1 0 amrmD¼ 0 ðA:21Þ By shifting summation indices, we obtain
1 Pe X1 0 ðm þ 2Þðm þ 1Þamþ2rDm Pe 4 X1 0 amrmD 2Rs 1 r2 WD X1 1 am1rmD¼ 0 ðA:22Þ
The coefficient of each power of rD to set to zero. For m¼ 0 this yields
Table 7
The values ofoZ2ðrD; n; sÞ=orDfor various rD and n at s¼ 3 (Pe ¼ 10)
n rD¼ 0:2 rD¼ 0:4 rD¼ 0:6 rD¼ 0:8 rD¼ 1:0 0 1.732 1012 6.280 1012 2.925 1012 1.579 1012 9.733 1012 1 2.421 1012 1.071 1012 5.074 1012 2.758 1012 1.707 1012 2 )1.469 1012 2.090 1012 1.201 1012 6.695 1012 4.198 1012 3 )2.491 1012 )9.155 1012 )4.691 1012 )2.688 1012 )1.722 1012 4 )1.993 1012 6.235 1012 3.743 1012 2.260 1012 1.492 1012 5 )9.838 1012 1.695 1012 1.521 1012 9.856 1012 6.771 1012 6 )3.370 1012 )7.290 1012 4.882 1012 3.516 1012 2.536 1012 7 )1.101 1012 )4.133 1012 )1.113 1012 )8.347 1012 )6.363 1012 8 )3.514 1012 )5.203 1012 1.026 1012 9.125 1012 7.423 1012 9 )1.103 1012 )1.234 1012 3.316 1012 3.856 1012 3.373 1012 10 )3.418 1012 )2.541 1012 1.070 1012 2.399 1012 2.273 1012 11 )1.047 1012 )5.085 1012 1.115 1012 3.610 1012 3.719 1012 12 )3.179 1012 )1.008 1012 )2.757 1012 9.687 1012 1.279 1012 13 )9.578 1012 )1.975 1012 )4.248 1012 1.849 1012 2.558 1012 14 )2.867 1012 )3.842 1012 )6.584 1012 1.629 1012 2.617 1012 15 )1.066 1012 )6.291 1012 )3.292 1012 )7.218 1012 )1.087 1012 16 )2.528 1012 )1.428 1012 )1.527 1012 )7.451 1012 )4.285 1012 17 )7.455 1012 )2.732 1012 )2.285 1012 )6.470 1012 )4.746 1012 18 )2.190 1012 )5.206 1012 )3.402 1012 )7.658 1012 3.915 1012 19 )6.415 1012 )9.881 1012 )5.043 1012 )1.048 1012 )1.488 1012 20 )1.873 1012 )1.869 1012 )7.444 1012 )1.235 1012 )1.853 1012 21 )5.454 1012 )3.525 1012 )1.095 1012 )1.529 1012 )8.763 1012 22 )1.584 1012 )6.630 1012 )1.605 1012 )1.889 1012 )6.922 1012 23 )4.591 1012 )1.244 1012 )2.346 1012 )2.325 1012 )5.172 1012 24 )1.328 1012 )2.328 1012 )3.418 1012 )2.854 1012 )5.831 1012 25 )3.834 1012 )4.349 1012 )4.970 1012 )3.491 1012 )6.242 1012 26 )1.105 1012 )8.107 1012 )7.209 1012 )4.259 1012 )6.697 1012 27 )3.179 1012 )1.509 1012 )1.043 1012 )5.181 1012 )7.090 1012 28 )9.135 1012 )2.803 1012 )1.507 1012 )6.288 1012 )7.565 1012 29 )2.621 1012 )5.201 1012 )2.174 1012 )7.615 1012 )8.040 1012 30 )7.512 1012 )9.636 1012 )3.129 1012 )9.202 1012 )8.521 1012
a2¼ Pe2
8 a0 ðA:23Þ
and in general and when m¼ 1; 2; 3; . . . amþ2 ¼ Pe ðm þ 2Þðm þ 1Þ Pe 4 am þ 2Rs 1 r2 WD am1 ðA:24Þ By inserting the value for the coefficients into (A.18), the general solution is obtained
ZðrD;0; sÞ ¼ b1Z1ðrD;0; sÞ þ b2Z2ðrD;0; sÞ ðA:25Þ where Z1ðrD;0; sÞ and Z2ðrD;0; sÞ are two linearly inde-pendent functions which are of the form of (A.18), with coefficients determined by (A.23) and (A.24) and set a0¼ 1, a1¼ 0 or a0¼ 0, a1¼ 1, respectively.
Case 2: n > 0, the coefficient is not analytical
The governing equation (A.14), except for n¼ 0, subjected to boundary conditions (A.15) and (A.16) is solved using the Frobenius power series method. First, assume the solution of equation (A.14) has the form of power series with undetermined coefficients,
Z¼X
1
m¼0
amrmþrD ðA:26Þ
and insert this series and the series obtained by termwise differentiation, dZ drD ¼X 1 m¼0 ðm þ rÞamrmþr1D ðA:27Þ d2Z drD ¼X 1 m¼0 ðm þ rÞðm þ r 1Þamrmþr2D ðA:28Þ
into Eq. (A.14). This yields: Pe r2 D X1 m¼0 ðm þ rÞðm þ r 1ÞamrDmþr2 Xn 2 Pe þPe 4 r 2 Dþ 2Rs 1 r2 WD r3D X 1 m¼0 amrDmþr¼ 0 ðA:29Þ By shifting summation indices, we obtain
PeX 1 m¼0 ðm þ rÞðm þ r 1ÞamrmþrD Xn2 Pe X1 m¼0 amrmþrD Pe 4 X1 m¼2 am2rDmþr 2Rs 1 r2 WD X1 m¼3 am3rmþrD ¼ 0 ðA:30Þ Table 8
The values ofoZ2ðrD; n; sÞ=orDfor various rDand n at s¼ 3 (Pe ¼ 60)
n rD¼ 0:2 rD¼ 0:4 rD¼ 0:6 rD¼ 0:8 rD¼ 1:0 0 2.308 1012 1.345 1012 9.778 1012 8.812 1012 9.779 1012 1 4.123 1012 2.435 1012 1.779 1012 1.607 1012 1.786 1012 2 )1.326 1012 )8.147 1012 )6.037 1012 )5.491 1012 )6.127 1012 3 )4.154 1012 )2.728 1012 )2.069 1012 )1.904 1012 )2.139 1012 4 1.915 1012 1.380 1012 1.081 1012 1.012 1012 1.148 1012 5 )5.434 1012 )4.409 1012 )3.602 1012 )3.441 1012 )3.951 1012 6 )1.438 1012 )1.347 1012 )1.157 1012 )1.134 1012 )1.322 1012 7 )5.281 1012 )5.842 1012 )5.330 1012 )5.380 1012 )6.385 1012 8 )5.551 1012 )7.421 1012 )7.248 1012 )7.573 1012 )9.172 1012 9 )3.531 1012 )5.817 1012 )6.134 1012 )6.662 1012 )8.257 1012 10 3.854 1012 8.217 1012 9.432 1012 1.070 1012 1.360 1012 11 2.577 1012 6.835 1012 8.607 1012 1.023 1012 1.339 1012 12 )5.255 1012 )1.766 1012 )2.459 1012 )3.078 1012 )4.154 1012 13 2.224 1012 1.844 1012 2.859 1012 3.785 1012 5.281 1012 14 1.772 1012 1.805 1012 3.139 1012 4.410 1012 6.379 1012 15 )1.965 1012 )1.605 1012 )3.152 1012 )4.720 1012 )7.094 1012 16 )1.151 1012 1.628 1012 3.633 1012 5.821 1012 9.114 1012 17 )3.516 1012 8.386 1012 2.140 1012 3.682 1012 6.020 1012 18 )1.081 1012 9.649 1012 2.833 1012 5.252 1012 8.989 1012 19 )3.280 1012 3.902 1012 1.325 1012 2.657 1012 4.771 1012 20 )9.854 1012 3.476 1012 1.372 1012 2.987 1012 5.640 1012 21 )2.946 1012 2.081 1012 9.602 1012 2.275 1012 4.528 1012 22 )8.769 1012 4.244 1012 2.300 1012 5.953 1012 1.251 1012 23 )2.599 1012 1.747 1012 1.116 1012 3.167 1012 7.047 1012 24 )7.677 1012 4.701 1012 3.560 1012 1.110 1012 2.620 1012 25 )2.260 1012 9.379 1012 8.490 1012 2.918 1012 7.320 1012 26 )6.633 1012 )2.294 1012 )2.454 1012 )9.324 1012 )2.491 1012 27 )1.941 1012 1.181 1012 1.547 1012 6.514 1012 1.857 1012 28 )5.668 1012 1.345 1012 2.718 1012 1.272 1012 3.879 1012 29 )1.651 1012 )1.270 1012 )8.639 1012 )4.505 1012 )1.471 1012 30 )4.798 1012 )5.915 1012 )9.794 1012 )5.705 1012 )2.000 1012
Equating the coefficient of each power of rD to zero yields a general formula for all m.
1 Perðr 1Þa0 Xn2 Pe a0¼ 0 ðm ¼ 0Þ ðA:31Þ 1 Peðr þ 1Þra1 Xn2 Pe a1¼ 0 ðm ¼ 1Þ ðA:32Þ 1 Peðr þ 2Þðr þ 1Þa2 Xn2 Pe a2 Pe 4 a0¼ 0 ðm ¼ 2Þ ðA:33Þ 1 Peðr þ mÞðr þ m 1Þam Xn2 Pe am Pe 4 am2 2Rs 1 r2 WD am3¼ 0 ðm ¼ 3; 4; 5; . . .Þ ðA:34Þ From (A.31), the indicial equation is obtained
rðr 1Þ Xn2¼ 0 ðA:35Þ
The roots are r1¼
1þpffiffiffiffiffiffiffiffiffiffiffi1þ4Xn2
2 and r2¼
1pffiffiffiffiffiffiffiffiffiffiffi1þ4Xn2
2 ,
respec-tively.
Eqs. (A.32)–(A.34) yield
a1¼ 0 ðA:36Þ a2¼ Pe2 4½ðr þ 2Þðr þ 1Þ Xn2a0¼ 0 ðA:37Þ am¼ Pe2 4 am2þ Peam3 ðr þ mÞðr þ m 1Þ Xn2 ðm ¼ 3; 4; 5; . . .Þ ðA:38Þ By inserting the value for the coefficients into (A.26) the general solution is obtained as
ZðrD; n; sÞ ¼ b1Z1ðrD; n; sÞ þ b2Z2ðrD; n; sÞ ðA:39Þ where Z1ðrD; n; sÞ and Z2ðrD; n; sÞ are two linearly inde-pendent functions that have the form given in (A.26), with coefficients determined by (A.36)–(A.38) and set r¼ r1 or r¼ r2, respectively.
Particular solutions can be obtained by straightfor-ward application of boundary conditions (A.15) and (A.16) to the general solution given by (A.25) and (A.39). Thus, the Laplace-finite Fourier cosine trans-form power series solution may be written as
WðrD;n;sÞ¼ exp Pe 2ð1 rDÞ
U ðSP21AP22ÞZ1ðrD;n;sÞþU ðSP11AP12ÞZ2ðrD;n;sÞ ðSP11AP12ÞðTQ21þAQ22ÞðSP21AP22ÞðTQ11þAQ12Þ
ðA:40Þ where S¼1 2þ lWs T ¼1 2 U ¼ F ðnÞ P11¼ Z1ðrWD; n; sÞ P21¼ Z2ðrWD; n; sÞ P12¼ oZ1ðrWD; n; sÞ orD P22¼ oZ2ðrWD; n; sÞ orD Q11¼ Z1ð1; n; sÞ Q21¼ Z2ð1; n; sÞ Q12¼ oZ1ð1; n; sÞ orD Q22¼ oZ2ð1; n; sÞ orD
Solutions in the original domain CðrD;h; tDÞ are the La-place and finite Fourier cosine inversion of WðrD; n; sÞ. The finite Fourier cosine inverse transform is performed first. In addition, the Laplace inverse of (A.40) has to be determined numerically. A FORTRAN subroutine, DINLAP/INLAP, provided by IMSL Subroutine Li-brary [18] and based on the de Hoog et al. [6] algorithm, is employed to perform the numerical Laplace inversion.
Appendix B
Tables 3–8 present the computed values of the func-tion, Z2ðrD; n; sÞ, and its first derivatives with respect to rD,oZ2ðrD; n; sÞ=orDfor various n and rDvalues when s is fixed at 3.
References
[1] Becker MW, Charbeneau RJ. First-passage-time transfer func-tions for groundwater tracer tests conducted in radially conver-gent flow. J Contam Hydrol 2000;40:299–310.
[2] Chen JS, Liu CW, Chen CS, Yeh HD. A Laplace transform solution for tracer test in a radially convergent flow with upstream dispersion. J Hydrol 1996;183:263–75.
[3] Chen JS, Chen CS, Gau HS, Liu CW. A two-well method to evaluate transverse dispersivity for tracer tests in a radially convergent flow field. J Hydrol 1999;223:175–97.
[4] Chen JS, Liu CW, Liao CM. A novel analytical power series solution for solute transport in a radially convergent flow field. J Hydrol 2002;266:120–38.
[5] Chen JS, Liu CW, Chen CS, Liao CM. Effect of well bore mixing volume on nonaxisymmetrical transport in a convergent radial tracer test. J Hydrol 2003;277:61–73.
[6] de Hoog FR, Knight JH, Stokes AN. An improved method of Laplace transforms using a Fourier series approximation. SIAM J Sci Stat Comput 1982;3(3):357–66.
[7] Drost W, Klotz D, Koch A, Moser H, Neumaier F, Rauert W. Point dilution methods of investigating ground water flow by means of radioisotopes. Water Resour Res 1968;4(1):125–46. [8] Gustafan KE. Introduction to partial differential equations and
Hilbert space methods. New York: John Wiley and Sons; 1980. [9] Hodgkinson DP, Lever DA. Interpretation of a field experiment
on transport of sorbed and non-sorbed tracers through a fracture in crystalline rock. Radioact Waste Manage Nucl Fuel Cycle 1983;4(2):129–58.
[10] Kreyszig E. Advanced engineering mathematics. New York: John Wiley and Sons; 1968.
[11] Mackay DM, Bianchi-Mosquera G, Kopania AW, Kianjah H. A forced-gradient experiment on solute transport in the Borden aquifer. 1. Experimental methods and moment analysis of results. Water Resour Res 1994;30(2):369–83.
[12] Moench AF. Convergent radial dispersion: a Laplace transform solution for aquifer tracer testing. Water Resour Res 1989; 25(3):439–47.
[13] Moench AF. Convergent radial dispersion: a note on evaluation of the Laplace transform solution. Water Resour Res 1991;27(12): 3261–4.
[14] Moench AF. Convergent radial dispersion in a double-porosity aquifer with fracture skin: analytical solution and application to a field experiment in fractured chalk. Water Resour Res 1995; 31(8):1823–35.
[15] Sauty JP. An analysis of hydrodispersive transfer in aquifers. Water Resour Res 1980;16(1):145–58.
[16] Sneddon IN. The use of integral transform. New York: McGraw-Hill; 1972.
[17] Thorbjarnson K, Mackay DM. A forced gradient experiment on solute transport in the Borden aquifer. 2. Transport and disper-sion of the conservative tracer. Water Resour Res 1994;30(2):385– 99.
[18] Visual Numerical, Inc. IMSL UserÕs Manual. Houston, TX, vol. 1, 1994, p. 159–61.
[19] Wang HQ, Crampon N. Method for interpreting tracer experi-ments in radial flow using modified analytical solutions. J Hydrol 1995;165:11–31.
[20] Welty C, Gelhar LW. Evaluation of longitudinal dispersivity from nonuniform flow tracer tests. J Hydrol 1994;153:71–102. [21] Wylie CR, Barrett LC. Advanced engineering mathematics. New
York: McGraw-Hill Higher Education; 1995.
[22] Zlotnik VA, Logan JD. Boundary conditions for convergent radial tracer tests and the effect of well bore mixing volume. Water Resour Res 1996;32(7):2323–8.
![Fig. 3. Comparison of breakthrough curves at an observation well located at r ¼ 5 [m], h ¼ p between the power series solution and numerical solution for P eeclet numbers of 20, 40, 60 and 80 in a radially convergent tracer test.](https://thumb-ap.123doks.com/thumbv2/9libinfo/8846165.240593/4.892.473.828.561.905/comparison-breakthrough-observation-solution-numerical-solution-radially-convergent.webp)
