Using the nonstationary spectral method to analyze
asymptotic macrodispersion in uniformly recharged
heterogeneous aquifers
Ching-Min Chang, Hund-Der Yeh
*
Institute of Environmental Engineering, National Chiao Tung University, Hsinchu, Taiwan Received 9 August 2007; received in revised form 29 November 2007; accepted 6 December 2007
KEYWORDS Nonstationary spectral approach; Nonstationary velocity fields; Asymptotic macrodispersion; Heterogeneous media
Summary This paper describes an investigation of the influence of uniformly distributed groundwater recharge on asymptotic macrodispersion in two-dimensional heterogeneous media. This is performed using a nonstationary spectral approach [Li, S.-G., McLaughlin, D., 1991. A nonstationary spectral method for solving stochastic groundwater problems: unconditional analysis. Water Resour. Res. 27 (7), 1589–1605; Li, S.-G., McLaughlin, D., 1995. Using the nonstationary spectral method to analyze flow through heterogeneous trending media. Water Resour. Res. 31 (3), 541–551] based on Fourier–Stieltjes represen-tations for the perturbed quantities. To solve the problem analytically, focus is placed on the case where the local longitudinal dispersivity aLis much smaller than the integral scale
of log transmissivity k (i.e., aL/k 1). The closed-form expressions are obtained for
describing the spectrum of flow velocity, the variability of flow velocity and asymptotic macrodispersion, in terms of the statistical properties and the integral scale of log trans-missivity, local transport parameters and a parameter b [Rubin, Y., Bellin, A., 1994. The effects of recharge on flow nonuniformity and macrodispersion. Water Resour. Res. 30 (4), 939–948] used to characterize the degree of flow nonuniformity due to the groundwater recharge. The impact of b on these results is examined.
ª 2007 Elsevier B.V. All rights reserved.
Introduction
The field-scale spreading of nonreactive solutes in porous formations is largely determined by the spatial variability
in groundwater flow velocities. From the stochastic point of view, the velocity variability is directly related to the cross-correlation between the log hydraulic conductivity perturbation and the perturbation in the hydraulic head. Therefore, the quantification of this relationship is the key in the prediction of field-scale transport processes in heter-ogeneous media.
0022-1694/$ - see front matter ª 2007 Elsevier B.V. All rights reserved. doi:10.1016/j.jhydrol.2007.12.006
* Corresponding author.
E-mail addresses:c2837@ms11.hinet.net(C.-M. Chang),hdyeh@
mail.nctu.edu.tw (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
Li and McLaughlin (1995)used the nonstationary spectral method to analyze flow in nonstationary velocity fields and concluded that the stationary spectral method (Bakr et al.,
1978) fails to capture the log conductivity-head cross covariance. Therefore, this excludes the direct applicability of the stationary spectral method to solve the problem of transport of solutes in random nonstationary velocity fields. Motivated by this, this study is devoted to the quantification of the nonstationarity in the statistics of random velocity fields using the nonstationary spectral approach (Li and
McLaughlin, 1991, 1995).
Groundwater recharge causes nonuniformity in the mean gradient of hydraulic head, and results in nonstationarity in the statistics of random velocity fields, thereby affecting the behavior of solute transport in heterogeneous aquifers (e.g.,Rubin and Bellin, 1994; Li and Graham, 1998; Butera
and Tanda, 1999; Zhu and Satish, 1999; Destouni et al., 2001; Maugis et al., 2002). This characteristic allows us to use the nonstationary spectral approach (Li and McLaughlin,
1991, 1995) to analyze the field-scale spreading process. Many studies within the framework of stochastic theory have been devoted to the investigation of field-scale solute transport in nonuniform groundwater flow (e.g.,Rubin and
Bellin, 1994; Li and Graham, 1998; Butera and Tanda, 1999; Zhu and Satish, 1999; Destouni et al., 2001; Maugis et al., 2002). The works of Rubin and Bellin (1994) and
Butera and Tanda (1999) are more directly relevant to our concern. Their investigation was carried out using the Lagrangian transport formalism in a near-source region where the effect of the pore-scale dispersion is not felt and the advective transport is strictly dominated. In other words, the effect of the pore-scale dispersion is negligible in the evaluation of their numerical results. However, in this study, focus is placed on the influence of nonuniformity in the mean flow caused by the uniform groundwater recharge on the asymptotic behavior of field-scale solute transport, which is strongly affected by the pore-scale dispersion. Chang and Yeh (2007) have shown that the prediction of the large-time behavior of macrodispersion in nonstationary velocity field made by using the advection-dominated trans-port theory, which is useful to quantify the near source transport characteristics, will not provide a good asymp-totic approximation. Therefore, the effect of the pore-scale dispersion will be included in the following analysis of the field-scale solute transport in nonstationary velocity field.
The application of the nonstationary spectral approach (Li and McLaughlin, 1991, 1995) to the investigation of the influence of uniformly distributed groundwater recharge on asymptotic macrodispersion is the task undertaken here. This will be performed by developing the closed-form expressions for the asymptotic macrodispersion coefficients within an Eulerian transport framework under consideration of the effect of the pore-scale dispersion. To the best of our knowledge, the large-time closed-form expressions have never before been presented. It is hoped that our finding will be useful for the prediction of large-time behavior of macrodispersion in nonuniform groundwater flow.
Statement of the problem
In this study we consider the problem of transport of a con-servative contaminant in a two-dimensional aquifer subject
to uniformly distributed recharge where the transmissivity is a stationary random space function. The spatial persistence of the random field of the log transmissivity can be fully characterized in terms of the covariance between two loca-tions. The theoretical analysis is developed herein for unidi-rectional mean flow with two-dimensional perturbations in transmissivity and hydraulic head. For convenience the coor-dinate axis X1is selected to be in the direction of the mean
groundwater flow so that U = < u > = (U,0), where u = (u1, u2)
is the groundwater flow velocity vector and < > stands for the expected value operator. The governing equation for hydraulic head distribution in a two-dimensional uniformly recharged aquifer can be written as (e.g., Gelhar, 1993;
Graham and Tankersley, 1994) o oXi TðXÞo/ oXi þ QR¼ 0 ð1Þ
where / is the hydraulic head, T is the transmissivity and QR
is the constant recharge rate. Assume that recharge and lnT processes are uncorrelated (e.g., Rubin and Bellin, 1994;
Butera and Tanda, 1999). In the analysis that follows, lnT(X) and /(X) are considered to be random functions.
Gelhar and Axness (1983) presented an Eulerian ap-proach to analyze the field-scale spreading of solute transport in heterogeneous media. In this approach the macrodispersion coefficients are determined by construct-ing the macroscopic dispersive flux and relatconstruct-ing it to the Fickian-type gradient transport relationship. This approach has proven useful in characterizing the larger-time behav-ior of the field-scale solute transport. We will adopt the formalism outlined by Gelhar and Axness (1983) along with (1) to investigate the influence of recharge on asymptotic macrodispersion in two-dimensional heteroge-neous media.
Head perturbation
The evaluation of asymptotic macrodispersion coefficients (Gelhar and Axness, 1983) in heterogeneous media would require relating the variability of the groundwater flow velocity to that of the local log transmissivity field. Toward the determination of the variability of the flow velocity we start with developing the hydraulic head perturbation which describes the variability of hydraulic head. We then substi-tute the head perturbation into the perturbed form of the Darcy equation in developing the spectrum of the flow velocity in the following section.
The random fields in(1), head and lnT, are decomposed into ensemble means and small perturbations around the mean, i.e.,
/ðX1; X2Þ ¼ < /ðX1; X2Þ > þhðX1; X2Þ ¼ HX1Þ þ hðX1; X2Þ
ln TðX1; X2Þ ¼ < ln TðX1; X2Þ > þfðX1; X2Þ ¼ F þ fðX1; X2Þ
ð2Þ The H in(2)is only a function of the X1direction, implying
unidirectional mean flow.
Expanding these terms and taking expectation of (1)
yields the mean head gradient equation o2H
oX21
þQR
The general solution for the ensemble mean head gradient is JðX1Þ ¼ oH oX1 ¼QR eFðx x0Þ þ J0 ð4Þ
where J0is the known value of J at X1= X0. As in the
previ-ous studies (Rubin and Bellin, 1994; Butera and Tanda,
1999), the parameter b = QRk/(eFJ0) is defined to quantify
the degree of flow nonuniformity due to the groundwater recharge, in which k is the integral scale of lnT. Thus, Eq.(4)can be rewritten as
JðX1Þ ¼ J0½1 þ bðX1 X0Þ=k ð5Þ
After subtracting the mean of the resulting Eq.(3)from(1), the result is the following first-order equation describing the hydraulic head perturbation
o2h oX2i ¼ of oX1 JðX1Þ þ QR eFf ð6Þ
Eq.(6)represents the spatial variability in head induced by aquifer heterogeneity and the groundwater recharge. It is clear from(5)that the mean hydraulic head gradient is depen-dent of X. This spatially variant mean head gradient leads the head random perturbation in(6)to be nonstationary.
The solution of(6)for the head perturbation h in terms of f and b can be developed using a nonstationary spectral approach (Li and McLaughlin, 1991, 1995) based on Fou-rier-Stieltjes representations for the perturbed quantities in wave number space. By using this approach, the random perturbations are represented by the following two-dimen-sional wave number integrals:
hðXÞ ¼ Z 1 1 UðX; KÞ dZfðKÞ ð7Þ fðXÞ ¼ Z 1 1 exp½iK X dZfðKÞ ð8Þ
where U(X, K) is a transfer function to be given, dZf(K) is
the complex Fourier amplitude of lnT, and K = K1, K2is the
wave number vector. Substituting(7) and (8) into (6)and recalling that b = QRk/(eFJ0) results in
o2U oX2i
¼ ½iK1JðX1Þ þ bJ0=k exp½iK X ð9Þ
The solution to(9)is found to be UðX; KÞ ¼ iK1K 2JðX 1Þ ðK21 K 2 1ÞbJ0=k K4 exp½iK X ð10Þ in which K¼ ðK21þ K 2 2Þ 1=2
represents the magnitude of the wave number vector K. The nonstationary representation for the head perturbation is then obtained by substituting (10) into(7) hðXÞ ¼ Z 1 1 iK1K2JðX1Þ ðK21 K 2 2ÞbJ0=k K4 exp½iK X dZfðKÞ ð11Þ where J(X1) is defined by(5).
Flow perturbation
In this section, we develop the spectrum of the flow velocity from the perturbed form of the Darcy equation, which re-lates the velocity variation to the lnT perturbations. The
spectrum of the flow velocity, the key analytical develop-ment presented in this paper, is necessary for the determina-tion of closed-form soludetermina-tions for the field-scale coefficients of transport using the formalism proposed byGelhar and
Ax-ness (1983).
Using Darcy’s equation, the first-order equation for the velocity perturbation takes the form (e.g., Gelhar, 1993;
Rubin and Bellin, 1994; Butera and Tanda, 1999) a u0 i¼ TG di1JðX1Þf oh o Xi ð12Þ where TG= exp[F], u0i= u(X) U(X), u(X) is the groundwater
flow velocity vector and U(X) = (U,0) is the mean flow veloc-ity vector. Note that the zero-order approximation for the mean velocity is in the form (Rubin and Bellin, 1994; Butera
and Tanda, 1999)
U¼ TGJ¼ TGJ0½1 þ bðX1 X0Þ=k
¼ U0½1 þ bðX1 X0Þ=k ð13Þ
where U0is the known mean velocity at X1= X0.
The last term on the right-hand side of(12) in the X1
direction is found using(11)
oh oX1 ¼ Z 1 1 K21K 2 JðX1Þ i2K1K22bJ0=k K4 exp½iK X dZfðKÞ ð14Þ
Similarly, in the X2direction
oh oX2 ¼ Z 1 1 K1K2K2JðX1Þ þ iK2ðK21 K 2 2ÞbJ0=k K4 exp½iK X dZfðKÞ ð15Þ
Substituting(8), (14), (15)and the Fourier-Stieltjes repre-sentations of velocity perturbations, i.e.,
u0 i¼
Z 1 1
exp½iK X dZuiðKÞ
into the velocity perturbation Eq. (12) gives the complex Fourier amplitudes of the longitudinal and transverse veloc-ities, respectively dZu1ðKÞ ¼ U0 ð1 K21 K2Þ 1 þ b X1 X0 k þ i2b k K1K 2 2 K4 ( ) dZfðKÞ ð16Þ dZu2ðKÞ ¼ U0 K1K2 K2 1þ b X1 X0 k þ ib k K2ðK21 K 2 2Þ K4 ( ) dZfðKÞ ð17Þ The spectra of the longitudinal and transverse velocities in terms of the lnT spectrum are obtained by multiplying each side of(16) and (17), respectively, by its complex conjugate and taking the expected value
Su1u1ðKÞ ¼ U 2 0 ð1 K21 K2Þ 2 1þ bX1 X0 k 2 þ 4b 2 k2 K21K 4 2 K8 ( ) SffðKÞ ð18Þ Su2u2ðKÞ ¼ U 2 0 K21K 2 2 K4 1þ b X1 X0 k 2 þb 2 k2 K22ðK 2 1 K 2 2Þ 2 K8 ( ) SffðKÞ ð19Þ
where Sff(K) is the spectrum of lnT. Eqs.(18) and (19)allow
for the development of macrodispersion coefficients at large times.
Velocity variances
To verify indirectly our results in(18) and (19), we compare the velocity variances, obtained using(18) and (19), with existing theoretical results ofLi and Graham (1998). In order to evaluate the velocity variation explicitly the spectrum Sff(K) in (18) and (19) must be specified. For this analysis
the random lnT perturbation field under consideration is characterized by the following spectral density function (e.g.,Mizell et al., 1982; Li and Graham, 1998)
SffðKÞ ¼ 3r2 fa 2K4 pðK2þ a2Þ4 ð20Þ where r2
f is the variance of lnT, a = 3p/16k, K is a
two-dimensional wave number vector, and K2¼ K21þ K 2 2.
With Sff(K) given in(20), integration of(18) and (19)over
the wave number domain yields the longitudinal and trans-verse velocity variances, respectively
r2 u1¼ Z 1 1 Su1u1ðKÞdK ¼ r 2 fU 2 0 3 8 1þ b X1 X0 k 2 þ32 9p2b 2 " # ð21Þ r2u 2¼ Z 1 1 Su2u2ðKÞdK ¼ rf2U 02 1 8 1þ b X1 X0 k 2 þ 32 9p2b 2 " # ð22Þ These are precisely equivalent to the results found byLi and
Graham (1998) (their(22) and (24)) using a first-order per-turbation technique under the assumption of spatially invariant recharge. This agreement between two disparate methods is important because it justifies the application of the nonstationary spectral approach to problems involved in nonstationary velocity fields.
0 1 2 3 4 σ σ 2 2 0 2 1 f u U 20 / ) (X1−X0 λ= 15 / ) (X1−X0 λ= 10 / ) (X1−X0 λ= 0 0.02 0.04 0.06 0.08 0.1 β 0 0.02 0.04 0.06 0.08 0.1 β 0 0.4 0.8 1.2 σ σ 2 2 0 2 2 f u U ( 1 0)/ =10 = −X λ X 15 / ) (X1−X0 λ 20 / ) (X1−X0 λ=
a
b
Figure 1 Dimensionless variance of the (a) longitudinal and of the (b) transverse velocity versus dimensionless recharge parameter b for various values of (X1 X0)/k.
It is clear that in the case of no recharge (b! 0), Eqs. (21) and (22) reduce to r2 u1 ¼ 3 8rf2U 2 0 ð23Þ r2 u2 ¼ 1 8rf2U 2 0 ð24Þ
which are well-known expressions for two-dimensional flow reported in the literature.
Fig. 1a and b depict the longitudinal and transverse velocity variances, respectively, as a function of b based on(21) and (22)for various values of (X1 X0)/k. These
re-sults show that an introduction of groundwater recharge leads to an increase in the velocity variance. This can be attributed to larger recharge resulting in shorter correlation distance of hydraulic head and hence larger variability of the flow velocity.
Approximate solutions for the macroscopic
dispersion coefficients
Once the relationship between the spectrum of the flow velocity and that of the local lnT field is obtained, we are in a position to develop the approximate solutions for the asymptotic macrodispersion coefficients. The goal here is to show the influence of the recharge on the asymptotic behavior of field-scale solute transport. Following Gelhar
and Axness (1983), under the unidirectional mean flow con-dition, the macrodispersion coefficient tensor, Dij, at large
times is given by: Dij¼ Z 1 1 aLK21þ aTK22 ½K2 1þ ðaLK21þ aTK22Þ 2 USujuiðKÞðKÞdK ð25Þ
where aL and aT are the local longitudinal and transverse
dispersivities and U is the mean flow velocity. The longitudi-nal and transverse macrodispersion coefficients from(25),
0 0.02 0.04 0.06 0.08 0.1 β 1 1.5 2 2.5 3 3.5 λ σ2 0 11 f U D 01 . 0 /λ= αL 1 . 0 /α = αT L 10 / ) (X1−X0 λ= 15 / ) (X1−X0 λ= 20 / ) (X1−X0 λ= 0 0.02 0.04 0.06 0.08 0.1 β 0 0.004 0.008 0.012 0.016 0.02 λ σ2 0 22 f U D 01 . 0 /λ= αL 1 . 0 /α = αT L 20 / ) (X1−X0 λ= 10 / ) (X1−X0 λ=
a
b
Figure 2 Dimensionless (a) longitudinal and (b) transverse asymptotic macrodispersion coefficients versus dimensionless recharge parameter b for various values of (X1 X0)/k and the indicated values of aL/k and aT/aL.
D11and D22, can be written by substituting(18) and (19)into (25), respectively, as D11¼ 3r2 fl 2U 0 p ( ½1 þ bðX1 X0Þ=k Z 1 1 Z kðe2m2þ gR2 2Þð1 e2m2=R 2Þ2 ½m2þ ðe2m2þ gR2 2Þ 2 R4 ðR2þ l2Þ4dm dR2 þ 4b 2 aL ½1 þ bðX1 X0Þ=k Z 1 1 Z ðR 12þ gR22Þ ½R2 1þ e2ðR 2 1þ gR 2 2Þ 2 R21R 4 2 R4ðR4þ l2Þ4dR1dR2 ) ð26Þ D22¼ 3r2 fl2U0 p ( ½1 þ bðX1 X0Þ=k Z 1 1 Z a LðR21þ gR 2 2Þ ½R21þ e2ðR 2 1þ gR 2 2Þ 2 R2 1R 2 2 ðR2þ l2Þ4dR1dR2 þ b 2 k ½1 þ bðX1 X0Þ=k Z 1 1 Z ðe2m2þ gR2 2Þ ½m2þ ðe2m2þ gR2 2Þ 2 R22ðe 2m2 R2 2Þ 2 R4ðR4þ l2Þ4 dm dR2 ) ð27Þ where Ri¼ kKi; R2¼ R21þ R 2 2¼ k 2 K2;e¼ aL=k; m¼ R1=e; g¼ aT=aL, and l = 3p/16.
The integrals of the forms of(26) and (27), for the gen-eral case, cannot be integrated analytically. Because of the form of the denominator in(26) or (27), the main con-tribution to the integral comes from mffi 0 as e ! 0 (Gelhar
and Axness, 1983). Therefore, for the physically reasonable case of relative small local dispersion (e = aL/k 1),
approximate analytical expressions can be developed by taking the limit e! 0, as illustrated byGelhar and Axness
(1983). Note that for field conditions the ratio aL/k is
typi-cally 102or small (Gelhar and Axness, 1983).
The integrals in(26) and (27)are approximated by taking the limit e! 0 in the integrands and integrating separately over each variable, to get
D11¼ U0r2fk ½1 þ bðX1 X0Þ=k þ 32 9p2 eb2ð1 þ 5gÞ ½1 þ bðX1 X0Þ=k ( ) ð28Þ D22¼ U0r2fk 1 8eð1 þ 3gÞ½1 þ bðX1 X0Þ=k þ256 9p2 b ½1 þ bðX1 X0Þ=k ð29Þ It is evident from(28) and (29)that for the case of advec-tion-dominated transport (aL! 0 and aT! 0), the
longitu-dinal and transverse macrodispersion coefficients in (28)
and (29), respectively, reduce to
D11¼ U0r2fk½1 þ bðX1 X0Þ=k ð30Þ D22¼ 256 9p2 U0r2fb 2 ½1 þ bðX1 X0Þ=k ð31Þ
Fig. 2a and b shows how the longitudinal and transverse macrodispersion coefficients, respectively, vary with b, according to(28) and (29). The increase of macrodispersion coefficient with b is related to the fact that larger recharge results in increases in the variation of flow velocity (Fig. 1a or b) and, consequently, results in more spreading of the solute plume. In the limit as b! 0 (the no-recharge case), the longitudinal and transverse macrodispersion coefficients tend to D11¼ U0r2fk ð32Þ D22¼ U0r2fke 8 ð1 þ 3gÞ ¼ U0r2faL 8 ð1 þ 3gÞ ð33Þ
which recover the results ofGelhar and Axness (1983)in the case of the two-dimensional flow.
Conclusions
The problem of solute transport in two-dimensional uni-formly recharged heterogeneous aquifers has been investi-gated from a stochastic point of view. Results here have been developed to quantify the influence of groundwater re-charge on the spectrum and the variation of flow velocity and asymptotic macrodispersion, in terms of the statistical properties and the integral scale of lnT, local transport parameters and a parameter b (Rubin and Bellin, 1994) char-acterizing the degree of flow nonuniformity due to recharge. The stochastic methodology employed to develop the results of this work is based on the nonstationary spectral approach (Li and McLaughlin, 1991, 1995). In particular, the introduc-tion of this approach allows for quantifying the nonstaintroduc-tiona- nonstationa-rity of head perturbation, and in turn developing the spectrum of the flow velocity, which is the key analytical development in the prediction of the field-scale transport coefficients in random nonstationary velocity fields. Our re-sults indicate that the increase of the variation of the groundwater flow velocity caused by the recharge leads to more spreading of the solute plume. This implies that ignor-ing the influence of the recharge in field applications leads to the erroneous conclusion in the predicted spreading of sol-ute plume. Our presented formulation for velocity variances compares well with the solutions obtained byLi and Graham
(1998) using a first-order perturbation technique.
Acknowledgements
Research leading to this work has been partially supported by the grants from Taiwan National Science Council under the contract number NSC 95-2221-E-009-017. We thank the asso-ciate editor and three anonymous reviewers for constructive comments that improved the quality of the work.
References
Bakr, A.A., Gelhar, L.W., Gutjahr, A.L., MacMillan, J.R., 1978. Stochastic analysis of spatial variability in subsurface flows: 1 comparison of one- and three-dimensional flows. Water Resour. Res. 14 (2), 263–271.
Butera, I., Tanda, M.G., 1999. Solute transport analysis through heterogeneous media in nonuniform in the average flow by a stochastic approach. Transport Porous Med. 36, 255–291. Chang, C.-M., Yeh, H.-D., 2007. Large-time behavior of
macrodi-spersion in heterogeneous trending aquifers. Water Resour. Res. 43, W11501. doi:10.1029/2007WR006017.
Destouni, G., Simic, E., Graham, W.D., 2001. On the applicability of analytical methods for estimating solute travel time statistics in nonuniform groundwater flow. Water Resour. Res. 37 (9), 2303– 2308.
Gelhar, L.W., 1993. Stochastic Subsurface Hydrology. Prentice Hall, Englewood Cliffs, New Jersey, p. 390.
Gelhar, L.W., Axness, C.L., 1983. Three-dimensional stochastic analysis of macrodispersion in aquifers. Water Resour. Res. 19 (1), 161–180.
Graham, W.D., Tankersley, C.D., 1994. Optimal estimation of spatially variable recharge and transmissivity fields under steady state groundwater flow Part 1: theory. J. Hydrol. 157, 247–266. Li, L., Graham, W.D., 1998. Stochastic analysis of solute transport in heterogeneous aquifers subject to spatially random recharge. J. Hydrol. 206, 16–38.
Li, S.-G., McLaughlin, D., 1991. A nonstationary spectral method for solving stochastic groundwater problems: unconditional analy-sis. Water Resour. Res. 27 (7), 1589–1605.
Li, S.-G., McLaughlin, D., 1995. Using the nonstationary spectral method to analyze flow through heterogeneous trending media. Water Resour. Res. 31 (3), 541–551.
Maugis, P., Mouche, E., Dewiere, L., 2002. Analysis of displacement variances in stochastic nonuniform flows by means of a first order analytical model and comparison with Monte-Carlo simu-lations. Transport Porous Med. 47, 1–27.
Mizell, S.A., Gutjhar, A.L., Gelhar, L.W., 1982. Stochastic analysis of spatial variability in two-dimensional steady groundwater flow assuming stationary and non-stationary heads. Water Resour. Res. 18 (4), 1053–1067.
Rubin, Y., Bellin, A., 1994. The effects of recharge on flow nonuniformity and macrodispersion. Water Resour. Res. 30 (4), 939–948.
Zhu, J., Satish, M.G., 1999. Stochastic analysis of macrodisper-sion in a semi-confined aquifer. Transport Porous Med. 35, 273–297.
