Semi-analytical and approximate solutions for contaminant transport
from an injection well in a two-zone confined aquifer system
Ping-Feng Hsieh, Hund-Der Yeh
⇑Institute of Environmental Engineering, National Chiao Tung University, Hsinchu, Taiwan
a r t i c l e
i n f o
Article history:Received 23 February 2014
Received in revised form 20 August 2014 Accepted 25 August 2014
Available online 4 September 2014 This manuscript was handled by Peter K. Kitanidis, Editor-in-Chief, with the assistance of Adrian Deane Werner, Associate Editor
Keywords: Radial transport Robin boundary condition Skin zone
Laplace transform Groundwater pollution
s u m m a r y
This study develops a mathematical model for contaminant transport due to well injection in a radial two-zone confined aquifer system, which is composed of a wellbore skin zone and a formation zone. The model contains two transient equations describing the contaminant concentration distributions; one is for contaminant transport in the skin zone while the other is for transport in the formation zone. The contaminants are injected into the well with given dispersive and advective fluxes; therefore, the well boundary is treated as a third-type (Robin) condition. The solution of the model derived by the method of Laplace transforms can reduce to a single-zone solution in the absence of the skin zone. In addition, an approximate solution in the time domain is also developed by neglecting dispersion for the case that the contaminants move away from the injection well. Analysis of the semi-analytical solution showed that the influence of the skin zone on the concentration distribution decreases as time elapses. The distribution will be over-estimated near the wellbore if the constant concentration (Dirichlet) condition is adopted at the well boundary. The approximate solution has advantages of easy computing and yield reasonable predictions for Peclet numbers larger than 50, and thus is a practical extension to existing methods for designing aquifer remediation systems or performing risk assessments. Ó 2014 Elsevier B.V. All rights reserved.
1. Introduction
Radial contaminant transport problems have been intensively studied.Ogata (1958)was the first to develop an analytical solution using the complex integral method for radial transport problems with an injection of constant contaminant concentration at the well-bore, yet his solution was in terms of an integral form and cannot be evaluated numerically.Tang and Babu (1979)presented a complete solution in terms of Bessel functions and modified Bessel functions for radial transport; however, their solution was in a very compli-cated form and difficult to evaluate.Moench and Ogata (1981) solved the radial transport equation by Laplace transforms. They obtained the Laplace-domain solution involving Airy functions and then inverted the solution numerically using the Stehfest algorithm to the time domain. Later,Hsieh (1986)gave an analytical solution for the radial transport problem; this solution consisted of an inte-gral form with Airy functions, and was computed by the 20-point Gaussian quadrature.Chen et al. (2002)applied a Laplace transform power series (LTPS) technique to solve the radial solute transport
equation with a spatially variable coefficient. The aquifer pumping induced a convergent flow field and the groundwater velocity there-fore varied with radial distance. Their analytical results indicated that the LTPS technique can effectively and accurately handle the radial transport equation under the condition of high Peclet number. Later,Chen (2010)presented a mathematical model for describing three-dimensional transport of a contaminant originating from an area centered within a radial, non-uniform flow field. The solution of the model was developed by coupling the power series technique, the Laplace transform and the finite Fourier cosine transform. The comparison between this solution and the Laplace transform solu-tion showed excellent agreement.Chen et al. (2012a)considered radial contaminant transport problems in a two-zone confined aqui-fer system and presented a semi-analytical solution for describing the concentration distribution in the aquifer system, which consists of a formation zone and a skin zone resulted from well drilling and/ or well completion.Liu et al. (2013)derived a semi-analytical solu-tion to the problem of groundwater contaminasolu-tion in an aquifer– aquitard–aquifer system, considering both advective transport and diffusive transport of contaminants in the aquifers and the interven-ing aquitard. All of the solutions mentioned above adopt a constant concentration condition at the inlet boundary, implying that the contaminants are well mixed and continuously enter the aquifer systems. In other words, those solutions were derived under the
http://dx.doi.org/10.1016/j.jhydrol.2014.08.046
0022-1694/Ó 2014 Elsevier B.V. All rights reserved.
⇑ Corresponding author. Address: 300 Institute of Environmental Engineering, National Chiao Tung University, 1001 University Road, Hsinchu, Taiwan. Tel.: +886 3 5731910; fax: +886 3 5725958.
E-mail address:[email protected](H.-D. Yeh).
Contents lists available atScienceDirect
Journal of Hydrology
first-type boundary (or Dirichlet boundary) condition at the well boundary.
The third-type boundary condition, also called Robin boundary condition, which considers the effects of both the dispersive and advective fluxes, may also be adopted at the rim of the wellbore. This condition leads to conservation of mass inside the formation, while the wellbore has a well-mixed concentration at a constant flow rate entering the formation (Bear, 1972). Chen (1987) pre-sented an analytical solution for radial dispersion problems with Robin conditions at the rim of the wellbore. The analysis ofYeh and Yeh (2007)indicates that the solution obtained from the con-taminant transport equation with the Dirichlet boundary condi-tion over-estimates the concentracondi-tion near the wellbore if the flow regime is dispersion dominant. Pérez Guerrero and Skaggs (2010) presented a general analytical solution depicting solute transport with a distance-dependent dispersivity in a heteroge-neous medium subject to a general boundary condition, which can be a first-, second-, or third-type.Veling (2012)presented a mathematical model describing the solute concentration distribu-tion in a radial groundwater velocity field due to well extracdistribu-tion or injection. The model was composed of a radial transport equation with a Dirichlet, Neumann, or inhomogeneous mixed boundary condition at the well boundary. The solution of the model was solved using the methods of Laplace transform and generalized Hankel transform. Chen et al. (2011) developed an analytical model depicting two-dimensional radial transport in a finite-domain medium subject to either the first- or the third-type boundary condition at the well boundary. Recently, Chen et al. (2012b) derived a generalized analytical solution for the problem of coupled multi-species contaminant transport in a finite-domain medium under an arbitrary time-dependent third-type boundary condition. Wang and Zhan (2013) developed a mathematical model for describing radial reactive solute transport due to well injection in an aquifer-aquitard system con-sisting of a main aquifer and overlying and underlying aquitards. The well boundary was specified as either Dirichlet or Robin type. In fact, the Dirichlet boundary can be considered as a special case of the Robin boundary because the solution developed with the Robin boundary can reduce to the one with the Dirichlet boundary if the dispersion mechanism is negligible.
In the past, many studies have been devoted to the develop-ment of approximate solutions for radial dispersion contaminant transport problems.Raimondi et al. (1959)derived an approximate solution based on two assumptions. One was that the total deriva-tive of the contaminant concentration with respect to time is equal to zero when the solute is far away from the well. The other was to neglect the radius of the injection well. Hoopes and Harleman (1967) presented a summary of earlier works and provided an approximate solution by neglecting the effect of a finite well radius. Later, both Dagan (1971) and Gelhar and Collins (1971) obtained approximate solutions by employing the perturbation method. Tang and Babu (1979) also presented an approximate solution based on the work ofRaimondi et al. (1959)with the con-sideration of the well radius.
The objective of this study is to develop a mathematical model for describing radial contaminant transport in a two-zone confined aquifer with a Robin boundary condition specified at the injection well. The solution (i.e., in Laplace domain) of the model is derived by the method of Laplace transform, and the time-domain results (hereinafter referred to as ‘‘semi-analytical solution’’) are obtained by the Crump algorithm (1976). In addition, an approximate solution in the time domain is also developed in terms of error and complementary error functions. The impacts of the skin zone and the use of different boundary conditions on the contaminant concentration distribution in the aquifer system are investigated based on the developed solution.
2. Methodology 2.1. Analytical solution
Some assumptions are made for the mathematical model describing radial transport of the injected contaminant in a two-zone confined aquifer system. They are: (1) the aquifer is homogeneous, isotropic and of uniform thickness, (2) the injection well has a finite radius and a finite thickness skin, and therefore the aquifer can be considered as a two-zone system, (3) the well fully penetrates the aquifer, and (4) the effect of molecular diffusion is negligible. The groundwater velocity in a steady-state radial flow system can be written as:
v
¼ Q =2p
rbn ð1Þwhere Q is a constant injection rate (L3/T), r is a radial distance from
the center of the wellbore (L), b is the aquifer thickness (L), and n is the aquifer porosity (–).
The governing equation describing the concentration distribu-tions in the wellbore skin zone and formation zone are expressed, respectively, as: @C1 @t þ
v
@C1 @r ¼ D1 @2C1 @r2 for rw<r 6 r1 and t > 0 ð2Þ and @C2 @t þv
@C2 @r ¼ D2 @2C2 @r2 for r1<r < 1 and t > 0 ð3Þ where C1 and C2 are the contaminant concentrations in the skinzone (or called first zone) and formation zone (or called second zone) [M/L3], respectively; D
1and D2are the dispersion coefficients
in the first zone and second zone defined by D1=
a
1v
and D2=a
2v
,respectively [L2/T];
a
1anda
2are the radial dispersivities in the firstzone and second zone, respectively [L]; t is the time since injection [T]; rwis the well radius [L]; r1is the radial distance from the central
of the well to the outer radius of the skin zone [L].
For the sake of convenience, the dimensionless forms for Eqs. (2) and (3)can be formulated, respectively, as:
j
@ 2G 1 @q
2 @G1 @q
¼q
@G1 @s
forq
w<q
6q
1 ands
>0 ð4Þ and @2G2 @q
2 @G2 @q
¼q
@G2 @s
forq
1<q
<1 ands
>0 ð5Þ where G1= C1/C0 and G2= C2/C0are the dimensionlessconcentra-tions in the first and second zones, respectively;
q
is a dimension-less radial distance defined asq
= r/a
2 and the other twodimensionless radial distances are
q
w= rw/a
2andq
1= r1/a
2;s
is adimensionless time defined as
s
¼ Qt=ð2p
bna
22Þ;
j
=a
1/a
2 is theratio of the skin-zone dispersivity to the formation-zone dispersiv-ity. Initially, the aquifer system is considered to be free from con-tamination; i.e., the contaminant concentrations in both the skin and formation zones are equal to zero and expressed as:
G1ð
q
;0Þ ¼ G2ðq
;0Þ ¼ 0 forq
Pq
w ð6Þ The Robin condition specified at the wellbore boundary is expressed asG1ð
q
w;s
Þj
@G1ð
q
w;s
Þ@
q
¼ 1 fors
>0 ð7ÞThe condition at the remote boundary is considered to be free from contamination and thus described as:
The continuity requirements for the contaminant concentration and flux at the interface of the skin zone and formation zone are, respectively: G1ð
q
1;s
Þ ¼ G2ðq
1;s
Þ fors
>0 ð9Þ andj
@G1ðq
1;s
Þ @q
¼ @G2ðq
1;s
Þ @q
fors
>0 ð10ÞApplying the Laplace transform to Eqs.4–10results in:
j
d 2 G1 dq
2 dG1 dq
¼q
s G1 ð11Þ d2G2 dq
2 dG2 dq
¼q
s G2 ð12Þ G1ðq
w;sÞj
dG1ðq
w;sÞ dq
¼ 1 s ð13Þ G2ð1; sÞ ¼ 0 ð14Þ G1ðq
1;sÞ ¼ G2ðq
1;sÞ ð15Þj
dG1ðq
1;sÞ dq
¼ dG2ðq
1;sÞ dq
ð16Þwhere G is the dimensionless Laplace-domain concentration and s is the transform parameter. Eqs.11–16can be solved as:
G1ð
q
;sÞ ¼ 1 sexpq
q
w 2j
2f ðq
;q
1Þ 2j
2=3gðq
;q
1Þ f ðq
w;q
1Þj
2=3gðq
w;q
1Þ þ 2s1=3j
2=3hðq
1;q
wÞ 2s1=3j
4=3iðq
1;q
wÞ ð17Þ and G2ðq
;sÞ ¼ 1 sexpq
1q
w 2j
þq
q
1 2 2j
2=3jðq
1;q
1Þ f ðq
w;q
1Þj
2=3gðq
w;q
1Þ þ 2s1=3j
2=3hðq
1;q
wÞ 2s1=3j
4=3iðq
1;q
wÞ ð18Þwhere f(x,y), g(x,y), h(x,y), i(x,y) and j(x,y) are functions composed of the Airy functions Ai(z), Bi(z), and the derivatives of the Airy func-tions Ai’(z) and Bi’(z). Detailed derivation for Eqs.(17) and (18)is given inAppendix A.
For the absence of the wellbore skin (i.e.,
q
1=q
w,a
1=a
2,and
j
= 1), both Eqs. (17) and (18) reduce to the same result, expressed as: G1ðq
;sÞ ¼ G2ðq
;sÞ ¼1 sexpq
q
w 2 2Aiðz1;q
Þ Aiðz1;q
wÞ 2s1=3Ai 0 ðz1;q
wÞ ð19Þwhich is identical to Chen’s (1987) solution for a single-zone aquifer. Obviously,Chen’s (1987)solution can be considered as a special case of the present solution.
Eqs.(17) and (18)are in the Laplace domain and expressed in terms of the Airy functions. The inversion of those two equations to the time domain may not be tractable due to the complexity of the Airy functions. TheCrump algorithm (1976)is therefore adopted to obtain the time-domain solution. Based on Abramowitz and Stegun (1972), the Airy functions are associated with the modified Bessel functions for positive arguments. These functions and their derivatives can be written as:
AiðzÞ ¼1
p
ffiffiffi z 3 r K1=3ðnÞ ð20Þ BiðzÞ ¼ ffiffiffi z 3 r ½I1=3ðnÞ þ I1=3ðnÞ ð21Þ Ai0 ðzÞ ¼p
1 zffiffiffi 3 p K2=3ðnÞ ð22Þ Bi0 ðzÞ ¼ zffiffiffi 3 p ½I2=3ðnÞ þ I2=3ðnÞ ð23Þwhere I and K are the first kind and second kind modified Bessel functions for n = 2/3 z3/2.
2.2. Approximate solution
Consider that the contaminant concentration in the region far away from the injection well does not change with time. That is to say dC/dt = 0 and:
@C @t ¼
v
@C
@r ð24Þ
Accordingly, the right-hand side terms in both Eqs.(2) and (3) are the same and can be written, respectively, as:
D1 @2C @r2 ¼
a
1v
@2C @t2 ð25Þ and D2 @2C @r2 ¼a
2v
@2C @t2 ð26ÞBased on Eqs.(25) and (26), Eqs.(4) and (5)can be transformed, respectively, to: @G1 @
s
þ 1q
@G1 @q
¼j q
@2G1 @s
2 ð27Þ and @G2 @s
þ 1q
@G2 @q
¼q
@2G2 @s
2 ð28ÞAlso, Eqs.(27) and (28)can, respectively, be reduced to:
d2G1 dW21 þ 2 W1 dG1 dW1 ¼ 0 ð29Þ and d2G2 dW22 þ 2 W2 dG2 dW2¼ 0 ð30Þ
if introducing the new variables W1and W2, respectively, defined
as: W1ð
q
;s
Þ ¼q
2 2s
ffiffiffiffiffiffiffiffiffiffiffiffiffiffi4 3j q
3 r ð31Þ and W2ðq
;s
Þ ¼q
2 2s
ffiffiffiffiffiffiffiffiffiffiffi4 3q
3 r ð32ÞWith Eqs.6–10, Eqs.(29) and (30)can be found, respectively, as:
G1ð
q
;s
Þ ¼g
erfcðW2;q1Þ þ f ½erf ðW1;q1Þ erf ðW1Þg
erfcðW2;q1Þ þ f ½erf ðW1;q1Þ erf ðW1;qwÞ þ hð33Þ
and
G2ð
q
;s
Þ ¼g
erfcðW2Þ
g
erfcðW2;q1Þ þ f ½erf ðW1;q1Þ erf ðW1;qwÞ þ hwith
g
¼ 4pffiffiffiffiffiffiffijp
q
4 w expðW 2 1;wþ W 2 2;1Þ ð35Þ f¼ 4pffiffiffiffip
q
4 wexpðW 2 1;wþ W 2 1;1Þ ð36Þ h¼ ffiffiffiffiffiffiffiffiffiffiffiffiffi3jq
3 w q ð6s
þq
2 wÞ expðW 2 1;1Þ ð37Þwhere erf() and erfc() are, respectively, the error function and the complementary error function with the arguments W1;q1,
W1;qw and W2;q1, respectively, representing W1(
q
1,s
), W1(q
w,s
)and W2(
q
1,s
).3. Results and discussion
The temporal distribution curves of the dimensionless contaminant concentration predicted by the present solution with dispersivity ratios
j
= 0.5, 1, and 2 are shown inFig. 1aforq
w= 1and
q
1= 4 andFig. 1bforq
w= 10 andq
1= 40. These two plotsindicate that the skin zone with a smaller dispersivity has a lower concentration at the early injection period but a higher concentration at the late period. In addition, the effect of the skin-zone dispersivity on the concentration distribution is more significant in the skin zone than in the formation zone.
The spatial distribution curves of the dimensionless concentra-tion predicted by the present soluconcentra-tion at small dimensionless times
Fig. 1a. Temporal distributions of dimensionless concentration atq= 2, 4 and 6 for qw= 1,q1= 4 andj= 0.5, 1 and 2.
Fig. 1b. Temporal distributions of dimensionless concentration atq= 20, 40 and 60 forqw= 10,q1= 40 andj= 0.5, 1 and 2.
Fig. 2a. Spatial distributions of dimensionless concentration whens= 0.5, 4.5, and 18 forqw= 1,q1= 2 andj= 0.5, 1 and 2.
Fig. 2b. Spatial distributions of dimensionless concentration whens= 0.5, 4.5, and 18 forqw= 1,q1= 4 andj= 0.5, 1 and 2.
are shown inFig. 2.Fig. 2aindicates that the influence of the skin-zone dispersivity on the concentration decreases quickly with increasing time. The concentration curves for different values of
j
tend to merge into one line as the time elapses. A smaller disper-sivity ratio has a higher concentration near the well but a lower concentration away from the well. Such a phenomenon can be attributed to the use of the Robin condition at the well boundary. Fig. 2bshows the concentration distributions for the aquifer with a larger skin thickness. Compared toFig. 2a, the effect of the skin thickness on the concentration distribution becomes large at large times.Fig. 3 illustrates the comparison of the spatial dimensionless concentration distributions predicted by the present semi-analytical solution and approximate solution. The figure indicates that the
approximate solution predicts poorly in the regions of high and low concentrations, but accurately in the region of intermediate concentrations (i.e., dimensionless concentrations in the range 0.2–0.8) as compared with those from the semi-analytical solution. The difference between these two solutions arises because the effect of dispersion is neglected in the development of the approximate solution. From a remediation perspective, the approximate solution is a convenient tool for providing useful information in designing aquifer clean-up systems or performing risk assessments.
The accuracy of the approximate solution depends on the mag-nitude of the Peclet number, defined as Pe= vL/D, where
v
isdefined in Eq.(1), and L is a characteristic length chosen as the distance between the injection well and the observation well. Pe
reduces to r/
a
, the dimensionless radial distance (q
), because D ¼a
v
and L = r.Fig. 4 shows the temporal distributions of the dimensionless concentrations at r = 20 m predicted by the semi-analytical solution and the approximate solution, for rw= 0.1 m,r1= 1 m, and Pe= 10, 20 and 50 (i.e.,
a
= 2, 1 and 0.4 m). WhenPe= 50, both solutions agree well, and the largest difference in
the predicted concentration is less than 0.05, indicating that the approximate solution gives good predictions when PeP50.
4. Conclusions
A mathematical model is presented for describing the concen-tration distribution in a radial two-zone aquifer system due to well injection at a constant rate and well mixed contaminant concentra-tion. The solution of the model is derived based on the methods of the Laplace transform and the Crump algorithm. The present solution reduces toChen’s (1987)solution in the absence of the wellbore skin. In addition, the present solution allows for the investigation of the influences of wellbore skin and different boundary conditions on the spatiotemporal dimensionless concentration distributions.
It was found that the dimensionless concentration distributions in the skin and formation zones differ from those in the homoge-neous (single-zone) system. For the skin zone with a small disper-sivity ratio, the concentration is lower at early injection periods but higher at late injection periods. In contrast, the concentration will be higher at the early period but lower at the late period for the case that the two-zone aquifer system has a large dispersivity ratio. The effect of skin thickness on the concentration distribution is large if the skin zone is thick and/or the time is large. In addition, the influence of skin zone on the dimensionless concentration decreases with increasing dimensionless time.
An approximate solution is also developed by considering that the contaminant concentration remains constant at significant distances from the well. When PeP50, the concentrations
predicted by the approximate solution have good agreement with those of the semi-analytical solution. The approximate solution has a much simpler form and therefore more easily evaluates the numerical value than the semi-analytical solution.
The predicted results from the semi-analytical solution demonstrate that the contaminant concentration at the wellbore during the early period of injection will be less than the injected concentration if the Robin boundary condition is adopted in the radial transport model. This is an important deviation from the models that adopt the Dirichlet condition, which causes the wellbore rim concentration to be equal to the injected concentration. Acknowledgements
Research leading to this paper has been partially supported by the grants from the Taiwan National Science Council under the contract numbers NSC 101-2221-E-009-105-MY2 and Fig. 3. Spatial distributions of dimensionless concentration predicted by the
semi-analytical solution and the approximate solution whens= 50, 225, 450, 900 and 1800.
Fig. 4. Comparison between the semi-analytical solution and the approximate solution for Pe= 10, 20, and 50.
102-2221-E-009-072-MY2. We are grateful to the editor, the associate editor Dr. Adrian Werner, and three anonymous reviewers for constructive comments that improved the quality of the work. Appendix A. Derivation of Eqs.(17) and (18)
Assume that:
G1¼ U1expðm
q
Þ ðA:1Þwhere m = 1/2
j
. Substituting Eq.(A.1)into Eq.(11)results in:d2U1 d
q
2 1 4j
2þq
sj
U1¼ 0 ðA:2ÞDefining Z1(
q
, s) = (s/j
)1/3(q
+ 1/4j
s), Eq. (A.2) can betransformed to the Airy equation expressed as:
d2U1
dZ21
Z1U1¼ 0 for
q
w<q
6q
1 ðA:3Þ Also, let:G2¼ U2expðn
q
Þ ðA:4Þwhere n = 1/2. With Eq.(A.4), Eq.(12)leads to:
d2U2 d
q
2 1 4þq
s U2¼ 0 ðA:5ÞSetting Z2(
q
, s) = s1/3(q
+ 1/4s), Eq.(A.5)becomes: d2U2dZ22
Z2U2¼ 0 for
q
1<q
<1 ðA:6Þ To solve Eqs.(A.3) and (A.6), we assume:U1ð
q
;sÞ ¼ aAiðZ1Þ þ bBiðZ1Þ ðA:7Þ andU2ð
q
;sÞ ¼ cAiðZ1Þ þ dBiðZ1Þ ðA:8Þ Based on the boundary conditions (Eqs.13–16), the coefficients a and b in Eq.(7)as well as c and d in Eq.(8)can be simultaneously determined as: a ¼1 sexp qw 2j
½2Ai0ðZ2;q1ÞBiðZ1;q1Þ 2j
2=3AiðZ 2;q1ÞBi 0 ðZ1;q1Þ 1 W ðA:9Þ b ¼1 sexp qw 2j
½2j
2=3Ai0ðZ1;q1ÞAiðZ2;q1Þ 2AiðZ1;q1ÞAi
0 ðZ2;q1Þ 1 W ðA:10Þ c ¼1 sexp q1qw 2
j
q1 2j
2=3½2Ai0ðZ1;q1ÞBiðZ1;q1Þ 2AiðZ1;q1ÞBi
0 ðZ1;q1Þ 1 W ðA:11Þ d ¼ 0 ðA:12Þ where
W
¼ f ðq
w;q
1Þj
2=3gðq
w;q
1Þ þ 2s 1=3j
2=3hðq
1;q
wÞ 2s1=3j
4=3iðq
1;q
wÞ ðA:13Þand f(x,y), g(x,y), h(x,y), i(x,y) and j(x,y) are functions composed of the Airy functions and expressed as:
f ðx; yÞ ¼ Ai0ðZ2;q1Þ½AiðZ1;xÞBiðZ1;yÞ AiðZ1;yÞBiðZ1;xÞ ðA:14Þ
gðx; yÞ ¼ AiðZ2;q1Þ½AiðZ1;xÞBi
0
ðZ1;yÞ Ai0ðZ1;yÞBiðZ1;xÞ ðA:15Þ
hðx; yÞ ¼ Ai0ðZ2;q1Þ½AiðZ1;xÞBi
0
ðZ1;yÞ Ai0ðZ1;yÞBiðZ1;xÞ ðA:16Þ
iðx; yÞ ¼ AiðZ2;q1Þ½Ai
0
ðZ1;xÞBi0ðZ1;yÞ Ai0ðZ1;yÞBi0ðZ1;xÞ ðA:17Þ
jðx; yÞ ¼ AiðZ2;qÞ½Ai0ðZ1;xÞBiðZ1;yÞ AiðZ1;yÞBi0ðZ1;xÞ ðA:18Þ where the arguments Z1;q1, Z1;qw and Z2;q1 represent Z1(
q
1, s),Z1(
q
w, s) and Z2(q
1, s), respectively.References
Abramowitz, M., Stegun, I.A., 1972. Handbook of mathematical functions: with formulas, graphs, and mathematical tables. Dover, New York.
Bear, J., 1972. Dynamics of fluids in porous media. American Elsevier, New York.
Chen, C.S., 1987. Analytical solutions for radial dispersion with Cauchy boundary at injection well. Water Resour. Res. 23 (7), 1217–1224.
Chen, J.S., 2010. Analytical model for fully three-dimensional radial dispersion in a finite-thickness aquifer. Hydrol. Process. 24 (7), 934–945.
Chen, J.-S., Liu, C.-W., Liao, C.-M., 2002. A novel analytical power series solution for solute transport in a radially convergent flow field. J. Hydrol. 266 (1), 120–138.
Chen, J.-S., Chen, J.-T., Liu, C.-W., Liang, C.-P., Lin, C.-W., 2011. Analytical solutions to two-dimensional advection–dispersion equation in cylindrical coordinates in finite domain subject to first-and third-type inlet boundary conditions. J. Hydrol. 405 (3), 522–531.
Chen, J.-S., Liu, C.-W., Liang, C.-P., Lai, K.-H., 2012. Generalized analytical solutions to sequentially coupled multi-species advective–dispersive transport equations in a finite domain subject to an arbitrary time-dependent source boundary condition. J. Hydrol. 456, 101–109.
Chen, Y.-J., Yeh, H.-D., Chang, K.-J., 2012. A mathematical solution and analysis of contaminant transport in a radial two-zone confined aquifer. J. Contam. Hydrol. 138–139, 75–82.
Crump, K.S., 1976. Numerical inversion of Laplace transforms using a Fourier series approximation. J. Assoc. Comput. Mach. 23 (1), 89–96.
Dagan, G., 1971. Perturbation solutions of the dispersion equation in porous mediums. Water Resour. Res. 7 (1), 135–142.
Gelhar, L., Collins, M., 1971. General analysis of longitudinal dispersion in nonuniform flow. Water Resour. Res. 7 (6), 1511–1521.
Hoopes, J.A., Harleman, D.R., 1967. Dispersion in radial flow from a recharge well. J. Geophys. Res. 72 (14), 3595–3607.
Hsieh, P.A., 1986. A new formula for the analytical solution of the radial dispersion problem. Water Resour. Res. 22 (11), 1597–1605.
Liu, C.-T., Yeh, H.-D., Yeh, L.-M., 2013. Modeling contaminant transport in a two-aquifer system with an intervening aquitard. J. Hydrol. 499, 200–209.
Moench, A., Ogata, A., 1981. A numerical inversion of the Laplace transform solution to radial dispersion in a porous medium. Water Resour. Res. 17 (1), 250– 252.
Ogata, A., 1958. Dispersion in Porous Media. Northwestern University.
Pérez Guerrero, J., Skaggs, T., 2010. Analytical solution for one-dimensional advection–dispersion transport equation with distance-dependent coefficients. J. Hydrol. 390 (1), 57–65.
Raimondi, P., Gardner, G., Petrick, C., 1959. Effect of pore structure and molecular diffusion on the mixing of miscible liquids flowing in porous media. In: American Institute of Chemical Engineers-Society of Petroleum Engineers Conference.
Tang, D., Babu, D., 1979. Analytical solution of a velocity dependent dispersion problem. Water Resour. Res. 15 (6), 1471–1478.
Veling, E., 2012. Radial transport in a porous medium with Dirichlet, Neumann and Robin-type inhomogeneous boundary values and general initial data: analytical solution and evaluation. J. Eng. Math. 75 (1), 173–189.
Wang, Q., Zhan, H., 2013. Radial reactive solute transport in an aquifer–aquitard system. Adv. Water Resour. 61, 51–61.
Yeh, H.-D., Yeh, G.-T., 2007. Analysis of point-source and boundary-source solutions of one-dimensional groundwater transport equation. J. Environ. Eng. 133 (11), 1032–1041.