Geosphere
Abstract
Most stochastic analyses of reactive transport in physically and geochemically heterogeneous aquifers have focused on the analysis of a single reactive species. Here we conduct the stochastic analysis of multicomponent competitive monovalent cation exchange. Transport equations for dissolved cations are coupled with nonlinear cation exchange terms, which, for chemical equilibrium, are described by massaction law expressions. These equations can be effectively decoupled by assuming that the weighted sum of cation concentrations is constant. The weight of each cation is equal to the reciprocal of its selectivity. Randomness of cation exchange capacity (CEC) leads to random retardation factors. Analytical expressions for effective retardation factors, longitudinal macrodispersivities, and concentration spatial moments are derived for a chemical system made of three monovalent cations (Na^{+}, K^{+}, and Cs^{+}) using the stochastic analytical solution of MirallesWilhelm and Gelhar (1996). Our results indicate that effective retardation factor, R_{C,i}, spatial moments, and macrodispersivities of K^{+} are significantly different from those of Na^{+}. Effective retardation factors asymptotically attain their mean values after a transient phase of cationdependent duration. They strongly depend on the correlation between logpermeability (log K) and CEC. Preasymptotic effective retardation factor values for a negative correlation are smaller than the mean value, regardless of the value of the coefficient of variation of CEC (CV_{CEC}). The smaller (larger) the variance of log K, , the greater (smaller) the effective retardation factor for a negative (positive) correlation. Cation macrodispersivities for a positive correlation structure are smaller than that of a nonreactive species and increase with increasing . On the other hand, for a negative correlation structure they are larger than that of a nonreactive species. Macrodispersivity of Na^{+} is smaller (larger) than that of K^{+} for a negative (positive) correlation structure. The macrodispersivity of K^{+} increases with k_{1} and decreases with the Na^{+}/Cs^{+} selectivity, k_{2}, for a positive correlation. The firstorder spatial moment of Na^{+} is greater than that of K^{+}, but is smaller than that of a nonreactive species at the asymptotic phase. Secondorder spatial moments of cations are smaller than that of a nonreactive species for a positive correlation structure, but are larger for a negative correlation.
 stochastic analysis
 cation exchange
 physical
 geochemical heterogeneities
 spatial moment
 macrodispersivity
 effective retardation factor
INTRODUCTION
A vast amount of research has been conducted in the last two decades to analyze the impact of geological heterogeneity on contaminant transport (Bosma et al., 1993; Brusseau, 1998; Dentz et al., 2000a, 2000b; Espinoza and Valocchi, 1997; Gelhar, 1993; Ginn et al., 1995; Kapoor et al., 1997; MirallesWilhelm and Gelhar, 1996; Mojid and Vereecken, 2005; Pang et al., 2003; Simmons et al., 1995; Srivastava et al., 2002, 2004). Experimental evidence and theoretical results indicate that the transport of reactive solutes in porous media can be significantly influenced by the spatial variability in physical and geochemical parameters, such as hydraulic conductivity, porosity, and sorption coefficient. Ptak et al. (2004) summarized recent developments on field tracer tests in heterogeneous porous media and stochastic modeling of flow and transport.
Contaminant mixtures are usually considered to contain two or more distinct compounds that chemically interact with one another as well as with solids in a highly nonlinear manner. One of the most challenging aspects is the derivation of effective parameters at field scale for reactive species that migrate through heterogeneous porous formations. MirallesWilhelm et al. (1997) argued that the assumption of equal dispersivities for substrate and oxygen could be inadequate because they found longitudinal macrodispersivities for the contaminant and dissolved oxygen to be significantly affected by biodegradation. In comparison with the large amount of stochastic studies for a single reactive solute transport with adsorption, much less attention has been paid to multicomponent competitive cation exchange reactive transport in physically and geochemically heterogeneous aquifers. Yang and Samper (2004) and Samper and Yang (2004) analyzed the behavior of multicomponent cation exchange in a heterogeneous column by MonteCarlo simulations. They pointed out that apparent dispersivities are different for each reactive cation.
In this paper, we evaluate the effects of hydrodynamic and geochemical heterogeneities on the transport of a multicomponent set of chemical species undergoing competitive monovalent cation exchange. A synthetic case with three monovalent cation exchange reactions is considered. K^{+} and Na^{+} are injected into a physically (hydraulic conductivity) and geochemically (cation exchange capacity [CEC]) heterogeneous aquifer in which Cs^{+} is the initially dominant cation occupying the exchange positions of solid phase. Na^{+} and K^{+} competitively exchange with Cs^{+} and enter into the solid phase. Cation exchange in groundwater is a surface reaction in which mass transfer between solution and solid commonly proceeds according to nonlinear massaction equations. Cation transport equations, which are generally coupled due to nonlinear cation exchange sink/source terms, can be effectively decoupled when the total weighted dissolved cation concentration (Q = C_{1} + C_{2}/k_{1} + C_{3}/k_{2}) is constant, where C_{i} is dissolved cation concentration, and k_{i} is selectivity coefficient. Retardation factors for each cation are derived in terms of CEC and selectivity coefficients. Randomness of retardation factors stems from the randomness of CEC. Mean and variance of retardation factors are calculated from those of CEC. The original problem is decoupled into three transport equations with random groundwater velocity and different random retardation factors. Analytical expressions of effective retardation factor, longitudinal macrodispersivity, and concentration spatial moments derived by MirallesWilhelm and Gelhar (1996) are employed to analyze the stochastic behavior of Na^{+} and K^{+}.
We describe here the mathematical formulation of the problem, and then we present the analytical expressions for effective retardation factors, macrodispersivities, and concentration spatial moments for Na^{+} and K^{+}. Effective retardation factors, macrodispersivities, and concentration spatial moments of Na^{+} and K^{+} are evaluated in terms of correlation structure, variance, and correlation length of log K and selectivity coefficients.
MATHEMATICAL STATEMENT OF STOCHASTIC CATION TRANSPORT AND EXCHANGE
Cation Transport Equations
Transport equations for N dissolved cations undergoing cation exchange in porous media can be expressed as: where C_{i} is concentration of the ith dissolved cation (mol/L); W_{i} is concentration of the ith exchanged cation in the solid phase (mol/L); D is the local dispersion tensor (L^{2}/T); and v is the water velocity vector (L/T).
Cation Exchange Reactions
Cation exchange takes place when free dissolved cations exchange with interlayer cations of clay minerals. This process can be described as an equilibrium reaction between dissolved and exchanged cations. According to the GainesThomas convention (Appelo and Postma, 1993), cation exchange reactions of Na^{+}, K^{+}, and Cs^{+} can be expressed as where (XNa), (XK), and (XCs) denote exchanged cations. It should be noted that these three reactions may take place simultaneously; however, they are linearly dependent. In fact, the third one can be obtained by adding the first two equations. Therefore, only two out of three possible exchange reactions are sufficient to describe all exchange reactions. Application of massaction laws to reactions 2 and 3 leads to: where k_{1} and k_{2} are selectivity coefficients of Na^{+} to K^{+} and Na^{+} to Cs^{+}, respectively, and β_{Na}, β_{K}, and β_{Cs} are equivalent fractions of exchanged Na^{+}, K^{+}, and Cs^{+} in the solid phase, respectively, which add to one:
Equivalent fractions are related to exchanged cation concentrations in equation 1, W_{i}, through where CEC is cation exchange capacity (meq/100 g of solids), 𝛉 is porosity (−), z_{i} is charge (in this case, z_{i} = 1 ∀ i), ρ_{d} is dry density of medium (kg solid/dm^{3}), and i = Na^{+}, K^{+}, and Cs^{+}.
Random Fields
Concentrations of dissolved cations, C_{i}, and exchanged cations, W_{i}, as well as groundwater velocity, v, in equation 1 are modeled as random fields. These random fields are driven by the inherent physical and geochemical heterogeneities of aquifers. Hydraulic conductivity, K, is typically found to be lognormally distributed, so that in order to describe its spatial variability, it is convenient to work with its logarithm:
The log conductivity field f(x) can be expressed as the sum of its ensemble mean f̄ and a random perturbation f′: Here f′ is assumed to be a statistically homogeneous random field, described by a threedimensional anisotropic exponential autocovariance function with correlation lengths λ_{i} (i = 1, 2, 3) and a variance (Gelhar, 1993).
Several researchers (Cassel et al., 2000; Jacques et al., 1999; Kuzyakova et al., 2001; Matschonat and Vogt, 1996; Pucci et al., 1997; Wirth, 2001) have reported that CEC of soils can be regarded as a stochastic field showing distinct spatial correlation as well as crosscorrelation with log K according to where a and b are the intersect and the slope of the linear correlation of CEC and logconductivity, f, and η is a zeromean random residual accounting for the fact that some of the variance of CEC is not related to the variability of f. It should be noticed that b > 0 indicates a positive correlation, b < 0 indicates a negative correlation, and b = 0 means no correlation between CEC and f.
The mean,, and fluctuation, CEC′, of CEC are given by
In our analysis, CEC and f are assumed to be perfectly correlated, and therefore η is equal to 0.
SOLUTION METHODOLOGY
Decoupled Cation Transport Equations
Contrary to the linear sorption problems usually discussed in the literature, the problem described in the previous section has nonlinear chemical source/sink terms involving concentrations of several components. Therefore, the first step is to linearize source/sink terms so that transport equations can be decoupled.
Concentrations of exchanged species (W_{i}) in equation 1 can be expressed in terms of solute concentrations by solving equations 5, 6, 7, and taking into account equation 8. According to the derivation given in Appendix A, exchanged concentrations, W_{Na}, W_{K}, and W_{Cs} can be written as where
In the case study analyzed here, it is assumed that Cs^{+} occupies initially all exchange sites. After the injection of Na^{+} and K^{+} into the system, these cations are sorbed to the solid phase, while Cs^{+} is released to the aqueous solution. Dissolved concentrations of Na^{+} and K^{+} are much smaller than that of Cs^{+}, so Na^{+} and K^{+} are trace cations compared to Cs^{+}. The focus of this paper is the stochastic analysis of Na^{+} and K^{+}.
By substituting equations 14 and 15 back into equation 1, and assuming that Q is constant, one finds where L(C_{i}) = ∇ · (D∇C_{i}) – v∇C_{i}, and
It should be noted that the two equations in equation 17 are uncoupled stochastic partial differential equations similar to those usually found in the literature. R_{Na} and R_{K} are retardation factors for Na^{+} and K^{+}, respectively, while CEC/100Q and CEC/100Qk_{1} are the distribution coefficients of Na^{+} and K^{+}, respectively.
Validity of Constant Q
Our solution assumes that the weighted sum of cation concentrations, Q, is constant. The weight of each cation is equal to the reciprocal of its selectivity. Jin and Ye (1999) used a similar assumption when deriving an analytical solution for binary cation exchange.
The validity of constant Q is discussed by Appelo and Postma (1993). They present calculations corresponding to a set of divalent cations including Ca^{+2}, Mg^{+2}, and Sr^{+2}. The latter is taken as a trace cation element, and its distribution coefficient is calculated according to: It should be noted that the denominator in this expression is Q. They compared the K_{d} values computed with this expression with values measured by Johnston et al. (1985) at an experimental site. Both measured and computed values agreed well and showed that K_{d} (and therefore Q) remained constant for a wide interval of Sr^{+2} concentrations, ranging from 10^{−6} to 1 mg/L. These results provide additional evidence for the validity of our assumption of constant Q for trace cations, satisfying the condition of C_{i} /k_{i} ≪ Q.
The validity of our assumption has also been tested numerically by solving numerically for nonlinear reactive cation exchange with the reactive transport code CORE^{2D} (Samper et al., 2003; Molinero et al., 2004; Dai and Samper, 2004). The case used to test the assumption corresponds to uncorrelated random conductivity and CEC in a rectangular 40 × 10 m domain. A random realization of conductivity and CEC was generated using MonteCarlo simulation. Both hydraulic conductivity and CEC had a correlation length of 2.8 m. Variances of log K and log CEC were both equal to 1. Na^{+}, K^{+}, and Cl^{−} were instantaneously injected at the center of the domain (see Fig. 1, Animation 1). Initially, the concentration of Cs^{+} was 0.01 mol/L. Other parameters of this numerical simulation can be found in Table 1. Figure 1 shows contour plots of the spatial distributions of Na^{+} and K^{+} in the aquifer 45 d after the pulse injection of Na^{+} and K^{+} into the system. For comparison purposes, the plot also shows the spatial distribution of a conservative species, Cl^{−}. It can be seen clearly that the displacement of the center of mass for each species is different. Reactive species are retarded compared to the conservative species, Cl^{−}. K^{+} suffers more retardation than Na^{+} because its selectivity, k_{1}, is less than 1.
Temporal fluctuations of Q around its initial value, Q_{0}, were evaluated at a selected node located downstream from the injection point. Relative fluctuations, computed as 100 (Q − Q_{0} / Q_{0})%, were < 0.5% (see Fig. 2).
Effective Retardation Coefficients, Macrodispersivities, and Spatial Moments of Na^{+} and K^{+}
Taking Na^{+} as the reference cation, retardation factors of other cations, such as K^{+}, can be written in terms of those of Na^{+}. From equations 18 and 19, it follows that:
According to the definition of cation retardation factors given in equations 18 and 19, randomness of CEC causes retardation factors to become random. Therefore, the mean and variance of cation retardation factors can be computed from those of CEC and where and are the mean and variance of CEC, respectively; and
From equations 21 and 22, it follows that the coefficient of variation of R̄_{i} is related to statistics of CEC through CV_{Ri} = [σ_{CEC}(P/k_{i})]/[1 + (P/k_{i})] with k_{i} = 1 for i = 1, which corresponds to Na^{+}. The righthand side of this expression can be rewritten in terms of the coefficient of variation of CEC, CV_{CEC}, and R̄_{i}:
It can be seen that in general CV_{Ri} depends on both CV_{CEC} and R̄_{i}. For a strongly sorbing cation (R̄_{i} ≫ 1), both coefficients of variation coincide. However, a cation having a retardation factor close to 1 has CV_{Ri} much smaller than that of CEC.
According to the results of MirallesWilhelm and Gelhar (1996), the macrodispersivities of Na^{+} or K^{+} are given by where γ is the flow factor, while the socalled effective retardation factors, R_{C,i}, are given by where R̄_{i} and σ_{R,i} are given by equations 21 and 22; t_{D,i} = vt/λ_{1}R̄_{i} is dimensionless time; and b_{i} is a parameter that defines the correlation structure between the retardation factor of each cation and log K. The values of b_{Na} and b_{K} can be calculated from equations 11, 21, and 22 by taking into account equation 23, such that
As stated by MirallesWilhelm and Gelhar (1996), the field scale retardation factor has a macrokinetic term that is caused by the crosscorrelation of retardation and the log hydraulic conductivity field and by the local retardation factor field. It should be noted that for a large value of time, effective retardation factors are equal to the mean retardation factors. Macrodispersivities in equation 25 also include a transient term.
Firstorder spatial moments of Na^{+} and K^{+}, X̄, are given by while longitudinal, secondorder spatial moments can be computed from:
(Gelhar, 1993) and
The ratio of the firstorder spatial moment of Na^{+} to that of K^{+} at the asymptotic phase is given by:
It can be seen that at the asymptotic phase, this ratio depends on the selectivity coefficient (k_{1}), P, and the mean CEC.
The ratio of the secondorder spatial moment of Na^{+} to that of K^{+} at the asymptotic phase is given by
Similarly, the ratio of the macrodispersivity of Na^{+} to that of K^{+} at the asymptotic phase can be written as:
One can see that this ratio depends on k_{1}, the mean CEC, and the correlation structure between log K and CEC.
RESULTS
A synthetic case inspired from the flow configuration and heterogeneity conditions of the Borden site (MirallesWilhelm and Gelhar, 1996; Sudicky, 1986; Roberts et al., 1986) is employed here to evaluate the behavior of multicomponent competitive monovalent cation exchange in heterogeneous aquifers by using the stochastic approach outlined in the previous section. Both hydraulic conductivity and CEC are assumed to be random processes having correlation lengths of 2.8 m along the mean direction of groundwater flow. The variance of log K is 0.29, and the average velocity of groundwater is 0.1 m/d. CEC has a mean value of 0.1 meq/100 g (typical of a sandy aquifer) and a coefficient of variation of 70%. Selectivity coefficients k_{1} and k_{2} are equal to 0.2 and 0.08, respectively, which correspond to published values for selectivities of Na^{+} to K^{+} and Na^{+} to Cs^{+} (Appelo and Postma, 1993). Correlation between log K and CEC is defined by parameter b in equation 12, which takes values of −0.1 (negative correlation), 0 (uncorrelated), and 0.1 (positive correlation). Table 1 summarizes the main parameters for this synthetic case.
It is worth noting that retardation factors, spatial moments, and macrodispersivities of each cation depend on the variance and correlation length of log K, the correlation structure between log K and CEC, and selectivity coefficients. Effective retardation factors, macrodispersivities, and spatial moments of each cation during both preasymptotic and asymptotic phases are analyzed next.
Effective Retardation Factors
According to equation 26, the effective retardation factor, R_{C,i} for the ith cation shows a transient trend that depends on the coefficient of variation of CEC (CV_{CEC}), the variance of f (), the correlation structure (b), and the mean retardation factor (R̄_{i}). R_{C,i} depends on a cationspecific dimensionless time, t_{D,i}, which measures the time needed for a cation to migrate a distance equal to the correlation length of log K.
Figure 3 shows the effective retardation factors of Na^{+} and K^{+} as a function of dimensionless time for two values of CV_{CEC} (30% and 70%). One can see that R_{C,i} depends on time and reaches the asymptotic mean value, which equals 1.125 for Na^{+} and 1.61 for K^{+}. Dimensionless times needed to reach the mean values are different for each cation. Such times are approximately equal to 3 for Na^{+} and 4 for K^{+}. As stated by MirallesWilhelm and Gelhar (1996), the effective retardation factor is significantly influenced by the correlation structure of log K and retardation. One can see in Figure 3 that for a negative correlation structure (b < 0), preasymptotic values of R_{C,i} are smaller than mean values, regardless of the value of CV_{CEC}. The reverse is true for a positive correlation structure (b > 0).
Curves in Figure 3 illustrate also that preasymptotic retardation factors, R_{C,i}, for CV_{CEC} = 30% are always greater than those corresponding to CV_{CEC} = 70%. Therefore, the smaller the CV_{CEC}, the larger the effective retardation coefficient. Figure 4 illustrates the sensitivity of effective retardation factors to changes in . During the preasymptotic phase, retardation factors, R_{C,i}, for a negative correlation structure and = 0.29 are greater than those corresponding to = 1. Therefore, the smaller the value of , the greater the effective retardation factor. However, the opposite holds true for a positive correlation structure. For a positive correlation structure and a large , all R_{C,i} values reach large values at t = 0.
Macrodispersivities
MirallesWilhelm et al. (1997) claimed dispersivities of reactive species in heterogeneous systems may be different. Next we show that our results confirm their statement. According to equation 25, cation macrodispersivities show a transient preasymptotic phase and depend on CV_{CEC}, , b, R_{C}, and λ_{1}.
Figure 5 shows the macrodispersivities of Na^{+} and K^{+} together with that of a nonreactive species as a function of dimensionless time (t_{D,K}) for positive and negative correlation structures. One can see that macrodispersivities for Na^{+} and K^{+} are different from each other and depend on the correlation structure between log K and CEC. Macrodispersivities of all species increase with increasing dimensionless time and asymptotically attain constant values, which are equal to 0.77 m for Na^{+} and 1.26 m for K^{+} for a negative correlation structure. Macrodispersivities are smaller for a positive correlation structure: 0.46 m for Na^{+} and 0.19 m for K^{+}. It should be noted that the asymptotic macrodispersivity of a conservative species is equal to 0.6 m. Macrodispersivity of the nonreactive species is smaller (greater) than that of reactive species for a negative (positive) correlation structure. The macrodispersivity of Na^{+} for a negative correlation structure is ∼1.5 times that for a positive correlation structure. A negative correlation structure results in a greater macrodispersivity of reactive species than a positive correlation structure.
Figure 5 also shows the ratio of the macrodispersivities of Na^{+} to K^{+} as a function of dimensionless time for three correlation structures. Such ratios slightly decrease with time for all three correlation structures and reach constant values. The ratio is greater than two for a positive correlation structure but less than one (∼0.7) for a negative correlation structure (see Fig. 5). At the asymptotic phase, the ratio is 0.7 for a negative correlation, around 1 for an uncorrelated structure, and ∼ 2.4 for a positive correlation structure.
Stochastic studies of heterogeneous porous media have focused on the major role played by the variance of log K. This parameter measures the degree of physical heterogeneity. To evaluate the effect of this parameter on the behavior of reactive species involved in competitive cation exchange, macrodispersivities at the asymptotic phase were computed for a wide range of values of the variance of log K from 0.001 to 1. Figure 6 shows macrodispersivities of Na^{+}, K^{+}, and a nonreactive species as a function of the variance of log K at the asymptotic phase for negative and positive correlation structures. Macrodispersivities of all species increase with increasing variance of log K. The macrodispersivity of the nonreactive species is greater (smaller) than that of Na^{+} and K^{+} for a positive (negative) correlation structure. For = 1, macrodispersivities are equal to 0.35 m for Na^{+} and 0.12 m for K^{+} for a positive correlation structure. They increase to 1.4 m for Na^{+} and 2.7 m for K^{+} for a negative correlation structure.
Figure 6 also shows the ratio of the macrodispersivities of Na^{+} to K^{+} as a function of the variance of log K at the asymptotic phase. One can see that the ratio is nearly constant and does not depend on the variance of log K for uncorrelated and negative correlation structures. From equation 32, one can see that for b = 0, the macrodispersivity ratio is equal to 1. For a positive correlation structure, the ratio of the macrodispersivity of Na^{+} to K^{+} increases from 2 to 4.5 when increases from 0.001 to 1.
The coefficient of variation of CEC is the parameter that measures the degree of geochemical heterogeneity. Macrodispersivities are not sensitive to CV_{CEC} because the case considered here assumes a perfect relationship between f and CEC without a contribution from η, the random residual accounting for the variability of CEC not related to the variability of log K. In the case of perfect correlation between f and CEC, the coefficient of variation of CEC has almost no effect on the macrodispersivities of Na^{+} and K^{+}. The contribution of η to the macrodispersivity could be significant depending on the effective retardation factor and the correlation length. Future studies should address the case of imperfect correlation between CEC and f.
According to equation 25, cation macrodispersivities at the asymptotic phase increase linearly with increasing correlation length.
Cation selectivity coefficients are parameters that control the distribution of mass between solution and solid in cation exchange reactions. From equation 32, it follows that macrodispersivities of Na^{+} and K^{+} are equal when k_{1} is 1. Figure 7 shows the macrodispersivities of reactive and nonreactive species as a function of selectivity coefficient, k_{1}, at the asymptotic phase. Macrodispersivity of the nonreactive species is 0.6 m. Macrodispersivity of Na^{+} does not depend on k_{1} and is equal to 0.46 m for a positive correlation and 0.77 m for a negative correlation. For a negative correlation structure, the macrodispersivity of K^{+} decreases with increasing k_{1} and asymptotically approaches the macrodispersivity of a nonreactive species. It is greater than that of Na^{+} if k_{1} is less than 1 and smaller otherwise. For a positive correlation structure, it increases with increasing k_{1} and asymptotically reaches the macrodispersivity of a nonreactive species.
The effect of k_{2} on macrodispersivities is illustrated in Figure 8. One can see that cation macrodispersivities increase with increasing selectivity coefficient, k_{2}, for a negative correlation structure. The opposite holds true for a positive correlation structure. Macrodispersivities of Na^{+} and K^{+} reach asymptotically identical values for a negative correlation structure and for a large value of k_{2}, but the difference in macrodispersivities of Na^{+} and K^{+} becomes more significant for a positive correlation structure.
Spatial Moments
According to equations 28 and 29, spatial moments for the ith cation depend on CV_{CEC}, , b, R_{C}, and λ_{1}. The results of our analysis of spatial moments are consistent with published results, which indicate that the firstorder spatial moment of a single reactive species does not depend on the correlation structure, while the secondorder spatial moment does. A positive correlation structure leads to a decrease in the secondorder spatial moment, whereas a negative correlation structure causes an increase in the secondorder moment (Bellin et al., 1993; Bosma et al., 1993; Hu et al., 1995).
Reactive species are retarded and move more slowly than a nonreactive species because Na^{+} and K^{+} are exchanged with Cs^{+} and get sorbed onto the solid. For this reason the firstorder spatial moment of K^{+} is smaller than that of Na^{+}, which in turn is smaller than that of the conservative species. Figure 9 shows the secondorder spatial moments of Na^{+}, K^{+}, and a nonreactive species as a function of dimensionless time. One can see that the secondorder spatial moment of a nonreactive species is smaller (greater) than that of reactive species for a negative (positive) correlation structure. Differences in secondorder spatial moments of Na^{+} and K^{+} are more significant for a positive correlation structure than for a negative correlation structure.
Figure 10 shows asymptotic first and secondorder spatial moments as a function of selectivity coefficient, k_{1}, for both positive and negative correlation structures. Firstorder spatial moments of reactive species are smaller than those of nonreactive species. However, one can see that the firstorder spatial moments of Na^{+} and the nonreactive species do not depend on the selectivity coefficient. The firstorder moment of K^{+}, however, increases with increasing selectivity coefficient and reaches a value equal to that of Na^{+} when k_{1} = 1. Therefore, the ratio of the firstorder spatial moments of Na^{+} and K^{+} increases with increasing k_{1} and is equal to 1 when k_{1} = 1. It can be seen that the differences in secondorder spatial moments between Na^{+} and K^{+} are greater for a positive correlation structure than they are for negative correlation.
CONCLUSIONS
We have presented here the stochastic analysis of multicomponent competitive monovalent cation exchange. Coupled transport equations for dissolved cations have been decoupled by assuming that the weighted sum of cation concentrations, Q, is constant. Derived expressions for cation retardation factors depend on CEC and selectivity coefficients. Randomness of CEC leads to random retardation factors. We have derived analytical expressions for effective retardation factors, longitudinal macrodispersivities, and concentration spatial moments for a chemical system made of Na^{+}, K^{+}, and Cs^{+} using the stochastic analytical solution of MirallesWilhelm and Gelhar (1996). The validity of the assumption of constant Q has been tested in a synthetic case by numerically solving the reactive transport of the three cations. The maximum relative change in Q was found to be <0.5%. Effective retardation factors, spatial moments, and macrodispersivities of Na^{+} and K^{+} and their ratios have been analyzed in terms of the correlation structure of logconductivity (log K) and CEC, the variance and correlation length of log K and selectivity coefficients. The analysis of the results leads to the following major conclusions:

Effective retardation factors, spatial moments, and macrodispersivities of K^{+} differ significantly from those of Na^{+}.

Cation effective retardation factors are time dependent and asymptotically attain their mean values. Times needed to reach mean values are different for each cation. Effective retardation factors are influenced by the correlation structure between log K and CEC: the smaller the CV_{CEC} value, the larger the effective retardation factor. In the case of perfect correlation between log K and CEC, the coefficient of variation of CEC has a minor effect on the macrodispersivities of Na^{+} and K^{+}. The smaller (larger) the variance of log K, , the greater (smaller) the effective retardation factor for a negative (positive) correlation.

Macrodispersivities of reactive species at the asymptotic phase increase with increasing variance of log K and depend on the correlation structure. The macrodispersivity of a nonreactive species is greater than that of reactive species for positive correlation, while it is smaller than that of reactive species for negative correlation. The macrodispersivity of Na^{+} is smaller than that of K^{+} for a negative correlation structure. The opposite holds true for the positive correlation structure. Macrodispersivity of K^{+} increases with k_{1} and decreases with k_{2} for the positive correlation structure, while it decreases with k_{1} and increases with k_{2} for the negative correlation structure.

The firstorder spatial moment of Na^{+} is greater than that of K^{+} but smaller than that of the nonreactive species at the asymptotic phase.

The secondorder spatial moment of nonreactive species is greater than that of a reactive species for a positive correlation structure, but smaller than that for a negative correlation structure. At the asymptotic phase, the secondorder spatial moment of Na^{+} is greater than that of K^{+} for a positive correlation structure. The opposite holds true for a negative correlation structure. Differences in secondorder spatial moments of Na^{+} and K^{+} are more significant for a positive correlation structure than for a negative correlation structure.
The method outlined in this paper can be applied to any number of monovalent cation exchange reactions as long as the weighted sum of dissolved cation concentrations, Q, is constant. Conditions for the validity of this assumption have been provided. Future work is needed to test the results of the stochastic analysis presented in this paper with both MonteCarlo simulations and field data.
APPENDIX A. DERIVATION OF EQUATIONS 14, 15, AND 16
Equation 5 can be rewritten so that β_{k} is on the lefthand side:
Similarly, from equation 6, one can get the expression for β_{Cs}:
By substituting A1 and A2 back into equation 7, one obtains:
This firstorder equation can be solved for β_{Na} in terms of cation concentrations and selectivity coefficients, such that
Substitution of A4 into A1 and A2 provides explicit expressions for β_{k} and β_{Cs}, respectively:
Substituting A4, A5, and A6 back into equation 8, and taking into account that we are dealing with monovalent cations (z_{i} = 1), gives equations 14, 15, and 16.
Acknowledgments
This work was supported by the Spanish Nuclear Waste Company (ENRESA) through contracts 770045 and 770125, the European Union through projects NFPRO (FI6WCT200302389) and FUNMIG (FP6516514), the Spanish Ministry of Science and Technology (CICYT Project HID98282), Xunta de Galicia (Incentive PGIDT04PXIC11801PM), and the University of La Coruña through a research scholarship awarded to the second author. We thank the two reviewers for their constructive comments and suggestions, which resulted in an improved version of the paper.
Footnotes

↵*jsamper{at}udc.es
 Accepted 1 January 2006.
 Received 25 July 2005.
 Revision received 15 December 2005.
Attribution: You must attribute the work in the manner specified by the author or licensor (but no in any way that suggests that they endorse you or your use of the work).
Noncommercial ‒ you may not use this work for commercial purpose.
No Derivative works ‒ You may not alter, transform, or build upon this work.
Sharing ‒ Individual scientists are hereby granted permission, without fees or further requests to GSA, to use a single figure, a single table, and/or a brief paragraph of text in other subsequent works and to make unlimited photocopies of items in this journal for noncommercial use in classrooms to further education and science.
 Geological Society of America