7870
Ind. Eng. Chem. Res. 2004, 43, 7870-7879
Single-Parameter Model for Binary Ion-Exchange Equilibria John L. Provis, Grant C. Lukey, and David C. Shallcross* Department of Chemical and Biomolecular Engineering, University of Melbourne, Victoria 3010, Australia
A model is presented to describe binary ion exchange equilibria with the use of a single fitted parameter per binary system. The proposed model incorporates the Pitzer model for calculation of ion activity coefficients in solution, a description of ion association in a concentrated solution, and a statistical thermodynamic treatment of exchange reactions. The behavior of all binary combinations of Ca2+, Mg2+, Na+, K+, and H+ ions in the presence of the Cl- ion as a nonexchanging anion is modeled, and the model fit to all data sets is found to be acceptable. Fitted parameters for each binary system are independent of the solution concentration, and the model is internally consistent. Model results compare favorably with those of other available models, particularly when taking into consideration the fact that only a single fitted parameter per binary system is used. Accurate prediction of the behavior of binary systems by use of data from other binary systems containing each component in the presence of a common “reference component” is also shown to be possible. Introduction A wide variety of methods have been proposed over the past 80 years for the description and analysis of multicomponent ion-exchange equilibria. The primary purpose of many of these models is to enable the prediction of the behavior of multicomponent ionexchange systems by analysis of their constitutive binary systems. These models were reviewed recently by Vo and Shallcross1 and generally rely upon the use of thermodynamic equilibrium constants2,3 in combination with other empirically determined parameters to describe exchange reactions. Models for the prediction of ion-exchange equilibria generally combine one of a selection of existing methods for calculating solutionphase activity coefficients with an empirical description of exchanger-phase nonidealities. Thus, the proportion of each exchanging ion present that is sorbed onto the exchanger when the system reaches equilibrium may be predicted. The model presented in the current work also follows this approach. The models used to describe the solution-phase activity range from a simple assumption of unity activity coefficients4,5 to the use of the extended Debye-Hu¨ckel equation3,6 or the model of Meissner and Kusik.7 The use of the Pitzer model8,9 for the prediction of solutionphase ion activities in ion-exchange modeling was first suggested by Elprince and Babcock,10 but it has only been widely implemented more recently11,12 because increases in computing power and parameter availability permitted its application to complex multicomponent systems. The equilibrium-based approach of Kester and Pytkowicz13 is often used in addition to the activity coefficient models to account for ion association in solution.1,12 The description of exchanger-phase nonidealities has most commonly utilized a Wilson-type model, following the approach of Elprince and Babcock,10 to calculate exchanger-phase activity coefficients from empirically determined “binary interaction parameters”. This model, coupled with each of the possible treatments of solutionphase activity coefficients, has been applied by a number of workers in several slightly different forms to provide a range of models of varying complexity and accuracy.14 * To whom correspondence should be addressed. Fax: +61 3 8344 4153. E-mail:
[email protected].
The Wilson equation15 uses two fitted interaction parameters per binary system to describe the deviation from ideality of the exchanger phase. Attempts to reduce these to a single parameter by use of symmetry relations16 provided fitted parameters that were dependent on the solution concentration. However, from the definitions of these parameters as equilibrium constants, they must necessarily be independent of the concentration. This then provides serious difficulties in the application of a model to be used as a predictive tool over a range of solution compositions and ionic strengths, as well as calling into question the fundamental validity of the model.17 The method of Melis et al.18 utilizes two fitted parameters to describe adsorption onto a heterogeneous surface, describing the equilibrium constants for exchange at each of two types of active sites. In the original model formulation, the two types of sites are assumed to be present in equal proportion on the exchanger phase. An extension of this model to include Pitzer’s activity coefficient model8,9 and the ion association treatment of Kester and Pytkowicz13 was found to provide an accurate description of ternary ion-exchange equilibria, based on the use of two fitted parameters from each component binary system.1 Myers and Byington19 also developed a three-parameter model from fundamental thermodynamic quantities, with two of these parameters used to describe resin heterogeneity and ionic strength effects incorporated explicitly. However, this model was based on the Langmuir isotherm, which has been shown not to be the most appropriate description of adsorption on heterogeneous resins.20 It is considered desirable to develop a model without such prior assumptions, which will then be applicable in a more general sense. An alternative model for exchanger-phase nonidealities was presented by de Kock and van Deventer21 and refined and extended by Lukey et al.5,22 This model was based on the Metropolis Monte Carlo method for statistical thermodynamics and was initially applied to complex multicomponent sorption in hydrometallurgical processes. As was noted by de Kock and van Deventer,21 its applicability is not limited to these complex systems, and a slightly simplified version of their model was found to model the ternary H+/Na+/K+ system with a high degree of accuracy at low ionic strengths. However, this model has to date never been combined with a
10.1021/ie049581h CCC: $27.50 © 2004 American Chemical Society Published on Web 10/27/2004
Ind. Eng. Chem. Res., Vol. 43, No. 24, 2004 7871
rigorous description of the solution-phase behavior because all previous authors implementing this model to describe the exchanger phase have assumed unity activity coefficients in solution. Shallcross17 recently presented a set of criteria that must be met by any empirical or semiempirical model describing multicomponent ion-exchange equilibria. Those criteria of interest in the description of binary equilibria may be summarized as follows: (i) Model parameters must be independent of the solution-phase concentration. (ii) Model parameters must remain constant in the presence of different nonexchanging anions. (iii) Internal consistency of the model is essential. Most ion-exchange models are overspecified, with more equations to be solved than unknowns present. The parameter values must not be significantly dependent on which of these equations are used in the solution procedure. Further criteria were also presented for the extension of a model from binary to ternary and higher systems, which will be addressed in a subsequent paper detailing the application of this model to multicomponent ionexchange equilibria.23 The aim of this study is to develop and present a model utilizing a single empirically determined but physically meaningful parameter for each binary system, satisfying these three criteria. It must be acknowledged that a single-parameter model cannot possibly provide as good a fit to experimental data as the more commonly used models with two or more parameters per binary system, particularly in highly nonideal systems or where multiple exchange mechanisms are taking place simultaneously. However, the fundamental value of the model developed in this work is centered on its ability to demonstrate that a single parameter is sufficient to obtain at least a qualitatively accurate description, and in some cases a quantitatively accurate prediction, of the equilibria present in ion-exchange processes. Model Development The model presented in the current work consists essentially of three components: (i) a simplified form of the model of de Kock and van Deventer21 for a description of the behavior of the exchanger phase; (ii) an implementation of the Pitzer model8,9 for the prediction of ion activities in the solution phase; (iii) use of an approach similar to that of Kester and Pytkowicz13 for the prediction of ion association in solution. This model resembles those developed by Mehablia et al.12 and Vo and Shallcross1 in the approach used to predict the solution-phase behavior. However, the novelty of this work lies in the coupling of this proven method of describing solution properties to a very different means of predicting the exchanger-phase behavior. In particular, the progressive reduction of the number of empirically determined parameters per binary system from the three used by Mehablia et al.12 to the two used by Vo and Shallcross1 continues with the presentation of this model. The prediction of binary ion-exchange equilibria is undertaken with the use of only a single fitted parameter per system. This parameter represents the mean sorption energy of each component. Parameter values therefore have physical significance rather than simply representing a multiplying factor applied to provide some degree of realism to a model otherwise assuming ideal behavior. The fundamental basis of the proposed exchangerphase model is the observation that any ion-exchange
process will involve a change in the total system energy. The position of the equilibrium in any system will therefore be determined by the relative sorption energies of all components participating in ion-exchange reactions.24 The difference in the relative sorption energies between the two exchanging components of a binary system is used as the single fitted parameter in the formulation of the model, giving physical significance to any observed variations in parameter values between different systems. The nature of the exchanging species and of the exchanger phase used will determine the sorption energies, which are currently not able to be predicted from first principles. Formulation of the Proposed Exchanger-Phase Model The assumptions made in formulating the proposed model for the exchanger-phase behavior are based on the refinements of Lukey et al.22 to those originally presented by de Kock and van Deventer.21 However, small but highly significant modifications have been made to some of the assumptions previously outlined to allow for nonideal behavior in the solution phase. Therefore, the full set of assumptions used in the current model formulation as well as the derivation of all equations used will be presented here. In the formulation of the proposed model, the solvent is treated as a species that is capable of exchange, meaning that a binary system actually involves exchange between three species. This assumption is made primarily as a mathematical convenience because the final loading of the solvent on the resin is negligible. The solvent is denoted species 0, while the exchanging species in a binary system are denoted A and B. The mean sorption energy of the solvent is set to be zero, with the sorption energy of each exchanging species being defined as negative. Consequently, the total system energy decreases with the replacement of a sorbed solvent molecule by an exchanging species on the exchanger phase. Without taking into account nonideal (selective) sorption, in its most simple form the exchanger-phase model assumes that all active sites in the exchanger phase are available for exchange by all exchanging species. Therefore, the probability Pr(Iload) that a particular active site is loaded with species I is given by eq 1, where yI
Pr(Iload) ) yI
(1)
represents the fraction of all active sites occupied by species I. This assumption means implicitly that variation in the selectivity between exchanging species at specific active sites is not incorporated into this model and will be revisited in the analysis of the model results. In the previous application of this exchanger-phase model to the highly nonideal selective sorption of metal cyanide complexes, it was necessary to include selectivity effects,21,22,25 effectively modifying the form of eq 1. This has not been incorporated into the current model, which is applied to ions exhibiting close to ideal sorption behavior. Ion exchange is assumed to occur by collision of dissolved exchanging species, with active sites occupied either by exchanging species or by solvent molecules. Species dissolved in the solution phase will participate in collisions with a probability proportional to their solution-phase activities. This may be represented symbolically by eq 2, where Pr(Jcol) represents the probability that a particular collision will involve dissolved species J and aJ is the activity of species J in
7872
Ind. Eng. Chem. Res., Vol. 43, No. 24, 2004
Pr(Jcol) )
aJ aQ ∑ Q
(2)
the solution phase as calculated from the Pitzer and ion association models. “Species Q” will be used throughout the model formulation as a dummy summation index, where all summations using dummy index Q include all exchanging species present in the system as well as species 0, the solvent. The probability that a change in the system configuration via replacement of species I with species J on the exchanger phase will occur is calculated by use of the Boltzmann criterion as expressed in terms of dimensionless sorption energies by eq 3. This equation is expressed in terms of Pr(IfJ)col, which is defined as the probability that, if it is known that a collision occurs between sorbed species I and dissolved species J, replacement of I by J will occur.
Pr(IfJ)col ) min [exp(-∆EIfJ), 1]
(3)
∆EIfJ is defined as the difference between the dimensionless sorption energies of species I and J. It should be noted that ∆EIfJ ) -∆EJfI. All sorption energies are expressed in units of kBT, where kB is Boltzmann’s constant and T is the absolute temperature, fixed at 294 K throughout this work.1 Following the usage of work by Lukey et al.,22 all sorption energies will be divided by the constant value kBT ()4.06 × 10-21 J) and therefore treated as dimensionless numbers throughout the development and analysis of this model. For the purposes of this study, it has been determined that a value of ∆EIf0 in the range of 10-20 is sufficient to match experimental observations of negligible solvent loading on the exchanger at equilibrium. Sorption energy reference points have therefore been set so as to force all calculated energies to fall approximately in this range. The probability of a particular collision between a solution-phase species and an occupied active site resulting in the replacement of sorbed species I with solution species J at that active site is then given by eq 4. Here Pr(IfJ) is the probability that a collision between any dissolved species and any sorbed species (i.e., any attempt at a change in the system configuration) will lead to the replacement of sorbed species I by dissolved species J.
Pr(IfJ) ) Pr(Jcol) Pr(Iload) Pr(IfJ)col
(4)
This relationship arises because three conditions must all be met for a given collision to result in the replacement of species I by species J on a particular active site. First, the dissolved species involved in the collision must be species J. Second, the sorbed species involved must be species I. Finally, given that both of these conditions have been met, exchange must occur. The method used to solve this model resembles that used in Monte Carlo simulations, with proposed changes in the system configuration accepted or otherwise according to a particular rule. In fact, if this model were to be used in the complete kinetic description of an ionexchange process, initialization of the model with conditions matching those used experimentally with incorporation of mass balance equations for both phases is possible, and the progression of the system toward equilibrium is monitored using a Metropolis Monte Carlo formulation. However, in this investigation focus-
ing on the description of ion-exchange equilibrium, kinetic effects are not modeled and the complete Metropolis Monte Carlo method formalism is not required. Rather, the model equations are solved analytically, starting from the requirement that, at equilibrium, the probability of an increase or a decrease in the loading of any sorbed species must be equal. This requirement is expressed by eq 5. Substituting eq 4 into both sides
Pr(IfQ) ) ∑Pr(QfI) ∑ Q Q
for each species I (5)
of eq 5:
Pr(Iload)
[Pr(Qcol) Pr(IfQ)col] ) ∑ Q Pr(Icol)∑[Pr(Qload) Pr(QfI)col] Q
(6)
Using eqs 1 and 2, simplifying, rearranging, and using the fact that ∑QyQ ) 1 yield eq 7. When eq 7 is expressed
yI
∑ Q
[
]
aQ Pr(IfQ)col + aI
{yQ[Pr(0fI)col ∑ Q Pr(QfI)col]} ) Pr(0fI)col (7)
in terms of each component of a multicomponent system, of which the binary systems of interest in this investigation are simply special cases, the entire system may be conveniently presented in a matrix formulation (eq 822),
h +N h )S (M h )D h
(8a)
h is a square matrix with elements defined by where M eq 8b.
MIJ ) Pr(0fI)col - Pr(JfI)col
(8b)
h is a diagonal matrix defined by eq 8c. N
NII )
aQ
Pr(IfQ)col; ∑ Q a
NIJ ) 0 for I * J (8c)
I
S h is a column vector containing the exchanger-phase loadings
SI ) yI
(8d)
and D h is a column vector containing the probabilities that a sorbed molecule of the solvent, species 0, will be replaced by species I in the event of a collision with species I.
DI ) Pr(0fI)col
(8e)
From the Boltzmann criterion and the observation made previously that ∆EIf0 values in the range of 1020 correspond to experimental observations of negligible solvent loadings on the exchanger phase, all entries in the column vector D h will have the value of unity. This means that the exchanging species will always replace the solvent on the resin if a collision occurs. It must also be noted that Pr(IfI) ) 1. This is a direct consequence of the Boltzmann criterion because if one molecule of species I is replaced by another, there is no
Ind. Eng. Chem. Res., Vol. 43, No. 24, 2004 7873
change in the system energy. Hence, this change in the system configuration will always be accepted.
D(n) )
Description of Exchanger-Phase Heterogeneity
∫-∞ψ exp n
[
]
(ξ - 1)2 2σ2
dξ
[ ]
(9)
MX(aq) T M+(aq) + X-(aq) KMX S )
[MX] [M+]f[X-]f
(10a)
(10b)
min [exp(-ψn∆E0fA), 1] - 1 min [exp(-ψn∆E0fB), 1] - min [exp(-ψn∆EAfB), 1]
N(n) )
[
0
(11b)
]
min [exp(-ψn∆E0fA), 1] - min [exp(-ψn∆EBfA), 1] min [exp(-ψn∆E0fB), 1] - 1 (10c)
aQ min [exp(-ψn∆EAfQ), 1]
∑ Q a
(11a)
defined by eq 11b is known as the stability constant of the compound MX, and tabulated values of this parameter for a wide range of compounds are available in the literature. Subscript f in eq 11b denotes free ion concentrations, taking into account only those ions that are dissociated according to eq 11a. Dissociation of compounds containing divalent cations and monovalent anions is treated as a sequential process, with an equilibrium established at each step of the dissociation. No assumptions are made regarding completeness or otherwise of the first dissociation of these compounds, but rather all equilibria are calculated explicitly. Given that this model is being used to describe exchange processes involving monovalent and divalent cations in the presence of a monovalent nonexchanging anion, mass balance equations for each type of cation may then be written in terms of free ion concentrations. Substitution of eq 11b into the mass balance equations yields a relationship for the free concentration of each ion in terms of all other free ion concentrations, total ion concentrations (denoted by subscript t, including both associated and dissociated ions), and stability constants. Equation 12 displays this for a monovalent cation M+.
M(n) )
[
(10e)
The method used to describe nonideal behavior in the solution phase is similar to that utilized by Mehablia et al.12 and by Vo and Shallcross.1 In solutions as concentrated as those used in some of the experimental work to which this model is applied,1 dissociation of ionic compounds will not be complete. Rather, following the approach of Kester and Pytkowicz,13 an equilibrium will be established according to eq 11. The value KMX S
with
yA(n) S(n) ) y B(n)
]
Description of Solution-Phase Nonideality
A preliminary analysis of the model output showed that N ) 30 and σ ) 0.25 provide an adequate description of the heterogeneity of the ion-exchange resin (Amberjet 1200H, Rohm and Haas, Philadelphia, PA) used throughout the experiments on which this work is based.1 The value of N describes the closeness of the approximation to a true Gaussian distribution, and increasing this above approximately 30 requires additional computing time for very little change in the model results. The value of σ is a property of the resin and so is held constant throughout this work. The model sensitivity to σ is also relatively low, so it has not been the subject of detailed investigation in this paper. The proposed model requires that eq 8 must therefore be solved separately for each of the N patches, multiplying the values of ∆EIfJ in eq 3 by the appropriate ψn for each patch when substituting the Boltzmann criterion into the matrix formulation of eq 8. For a binary system with exchanging ions A and B, this results in eq 10, the final formulation of the exchanger-phase section of the proposed model. Equation 10 is presented here with all matrix entries explicitly displayed for the sake of completeness.
S(n) ) (M(n) + N(n))-1D(n)
min [exp(-ψn∆E0fA), 1] min [exp(-ψn∆E0fB), 1]
This matrix equation is defined and solved for each of the N patches. The overall results are obtained by summing the calculated loadings on all patches. This model may then be used to predict the equilibrium exchanger-phase behavior of systems containing exchanging species of approximately equal size undergoing a reversible sorption process with negligible lateral interaction between sorbed species.22
Heterogeneity of the active sites in the exchanger phase is incorporated into the model by dividing the resin into N uniformly sized “patches”, with the sorption energy of the exchanging species on each patch following a normal (Gaussian) distribution.24 That is, the energy of sorption of each exchanging species onto the active sites of the nth patch of the resin is multiplied by a factor ψn calculated by the numerical solution of eq 9, the right-hand side of which is a Gaussian cumulative distribution function (CDF) with a mean of 1 and a standard deviation of σ.
1 2n - 1 ) 2N σx2π
[
0
A
∑ Q
aQ min [exp(-ψn∆EBfQ), 1] aB
]
(10d)
7874
Ind. Eng. Chem. Res., Vol. 43, No. 24, 2004 [M]t ) [M+]f(1 + KMX S [X ]f)
(12)
Equation 13 is the corresponding relationship for a divalent cation L2+, with both first and second dissocations described explicitly. +
+
LX2 - 2 [L]t ) [L2+]f(1 + KLX KLX S [X ]f + KS S [X ]f ) (13)
A similar mass balance may be carried out to determine [X-]f, the free concentration of the nonexchanging anion. The exact form of this equation depends on the valencies of the cations present. Equation 14 shows the equation used for a mixture of M+ and L2+, referred to as a 1-2 system, with analogous formulations possible for monovalent/monovalent (1-1) and divalent/divalent (2-2) systems. +
+ 2+ LX [X]t ) [X-]f(1 + KMX S [M ]f + KS [L ]f + +
2+ + 2 2KLX KLX S S [L ]f[X ]f) (14)
For each binary system, an equation of the form of eq 12 or eq 13 may be formulated for each exchanging cation and solved simultaneously with the appropriate version of eq 14 according to the valencies of the cations present. Numerical solution of these three equations yields the free ion concentrations present in the solution phase of the ion-exchange system of interest.13 The cation activities may then be calculated from these values according to eq 15. The activity coefficients γM
aMm+ ) γM[Mm+]f
(15)
are calculated by use of the Pitzer model8,9 in the form described by Harvie and Weare.26 The Pitzer model, as applied to a binary cation-exchange system between monovalent and/or divalent cations in the presence of a single monovalent anion X-, is presented as eq 16.
|zK|mK)CMX] + ∑ K ∑c mc(2ΦMc + mXψMcX) + |zM|∑c mcmXCcX 2
ln γM ) zM F + mX[2BMX + (
[
xI
1 + 1.2xI
+
ln(1 + 1.2xI) 0.6
]
(16)
+
∑c mcmXB′cX + mAmBΦ′AB
(17)
where Aφ is a lumped combination of physical parameters equal to 0.392 for water at 25 °C. The coefficients BMX follow the functional form given by Pitzer8
BMX ) β(0) MX +
2β(1) MX 4I -1 + 1 + 2xI + exp(-2xI) 4I 2
[
(
)
]
(18b)
and the terms CMX and ΦAB, representing the contribution of the third virial coefficients to single- and mixedelectrolyte effects, respectively, are defined by
CMX )
CφMX
(19)
2x|zMzX|
ΦAB ) θAB + EθAB(I)
(20a)
Φ′AB ) Eθ′AB(I)
(20b)
The values β(0), β(1), Cφ, θ, and ψ are empirical parameters characteristic of each single- or mixed-electrolyte system. The term EθAB accounts for the electrostatic effects of mixing electrolytes containing ions whose charge is of the same sign but different magnitude and is zero for 1-1 and 2-2 systems. For 1-2 systems, these “unsymmetrical mixing terms” are given by eq 21,27 expressed in terms of the variable xIJ ) 6zIzJAφxI E
E
θAB )
z AzB 1 1 Ω0(xAB) - Ω0(xAA) - Ω0(xBB) 4I 2 2
[
]
(21a)
θ′AB ) E θAB 1 1 Ω Ω (21b) Ω (x ) (x ) (x ) 1 AB 1 AA 1 BB 2 2 2 I 8I
z AzB
[
]
where Ω0 and Ω1 are defined by
∫0∞ξ2[1 - exp(- ξxe-ξ)] dξ
Ω0(x) )
x 1 -1+ 4 x
Ω1(x) )
x 1 4 x
(22a)
∫0∞ξ2[1 - (1 + ξxe-ξ) exp(- ξxe-ξ)] dξ (22b)
The symbol zJ represents the charge of species J, mJ is the molality of species J, and I is the solution-phase ionic strength. The dummy summation index c represents a sum over all cations present. Where terms relating to particular cations are required, those cations are denoted as A and B. The term F is a combination of a Debye-Hu¨ckel term and terms incorporating the dependence of the second virial coefficients on ionic strength
F ) -Aφ
B′MX )
2β(1) MX [1 - (1 + 2xI) exp(-2xI)] (18a) 4I
The integrals in eq 22 were evaluated using the QAGS algorithm as implemented for semiinfinite integrals in the GNU Scientific Library (GSL) version 1.4 (Free Software Federation, Boston, MA).28 Model Application The proposed model has been tested by comparison with the experimental data of Vo14 for a range of binary cation-exchange systems on the commercially available resin Amberjet 1200H (Rohm and Haas, Philadelphia, PA), with Cl- as a nonexchanging anion. Full experimental details are presented by Vo and Shallcross.1 Data were available for the 1-1 systems Na+-K+-Cl-, Na+-H+-Cl-, and K+-H+-Cl-, the 2-1 systems Ca2+Na+-Cl-, Ca2+-K+-Cl-, Ca2+-H+-Cl-, Mg2+-Na+Cl-, Mg2+-K+-Cl-, and Mg2+-H+-Cl-, and the 2-2 system Ca2+-Mg2+-Cl-. For each system, a minimum of 12 data points (compositions) were available at each of four concentrations between 0.10 and 1.00 N. Appropriate stability constant and Pitzer model parameters for these systems were obtained from the literature and are tabulated in Tables 1-3. As was noted previously, the sorption energy of the solvent was set to zero in the application of this model. The requirement for ∆EIf0 to fall approximately in the range of 10-20 for each species I was met in the
Ind. Eng. Chem. Res., Vol. 43, No. 24, 2004 7875 Table 1. Stability Constant Data Used in the Proposed Model MX
KMX S
ref
CaCl2 CaCl+ MgCl2
0.60256 exp(1.073 - 0.464I) (2(1 - R)/IR(1 + R)), where R ) 1 - 0.294xI - 0.406I + 0.222I3/2 0.77 exp(-0.537 - 1.002I) exp(-0.491 - 0.464I) exp(-1.179 - 0.982I)
33 34 35
MgCl+ NaCl KCl HCl
36 34 34 34
Table 2. Single-Salt Pitzer Model Parameters9 MX
β(0) MX
β(1) MX
CφMX
CaCl2 MgCl2 NaCl KCl HCl
0.3159 0.3524 0.0765 0.04835 0.1775
1.614 1.6815 0.2664 0.2122 0.2945
-0.000 17 0.002 596 0.001 27 -0.000 84 0.000 80
Figure 1. Comparison of model predictions to experimental data for the Na+-H+-Cl- binary system, ∆ENafH ) 1.00.
Table 3. Mixed-Electrolyte Pitzer Model Parameters9 A
B
X
θAB
ψABX
Ca2+
Na+
Cl-
0.07 0.092 0.032 0.007 0.07 0 0.10 0.036 -0.0012 0.005
-0.007 -0.015 -0.025 -0.012 -0.012 -0.022 -0.011 -0.004 -0.0018 -0.011
Ca2+ Ca2+ Ca2+ Mg2+ Mg2+ Mg2+ Na+ Na+ H+
H+ K+ Mg2+ Na+ K+ H+ H+ K+ K+
ClClClClClClClClCl-
investigation of each reaction system by setting ∆E0fA to a constant value in the appropriate range and allowing ∆EAfB and therefore ∆E0fB to vary freely. In fitting the model to each system, the optimal value of ∆EAfB for each system was determined by use of a simple descent direction search algorithm. The objective function minimized was the weighted sum of squared errors given in eq 23, where the sum is taken over all points in the relevant data set. Subscript mod refers to model predictions, and exp denotes experimental values.
S)
Figure 2. Comparison of model predictions to experimental data for the Ca2+-K+-Cl- binary system, ∆ECafK ) 3.77.
∑j W(xj,yj) [y0,j,mod2 + (yA,j,mod - yA,j,exp)2 + (yB,j,mod - yB,j,exp)2] (23a)
where
W(x,y) ) exp(-25xAxByAyB)
(23b)
The use of this weighting function minimizes artifacts in the fitting procedure due to inaccuracies in the determination of the composition in binary systems containing much more of one component than the other.14 Figures 1-3 show the results obtained by fitting the proposed model by optimization of ∆EAfB to experimental data for representative 1-1, 1-2, and 2-2 ionexchange systems. Fitted energy parameter ∆EAfB values for all systems are presented in Table 4. It is observed that the model predicts these equilibria to a remarkable degree of accuracy considering that it is based on the use of only a single fitted parameter per binary system. In each case, data are obtained at solution concentrations from 0.10 to 1.00 N. It may be noted that the model tends to slightly underestimate the sorption of components with yM < 0.5 and overestimate those with yM > 0.5. This trend
Figure 3. Comparison of model predictions to experimental data for the Ca2+-Mg2+-Cl- binary system, ∆ECafMg ) 0.85. Table 4. Fitted Values of ∆EAfB for Each Binary System Investigated A-B
∆EAfB
A-B
∆EAfB
Na-H Na-K Ca-Na Mg-Na K-H
1.00 0.10 4.51 3.81 1.29
Ca-H Mg-H Ca-K Ca-Mg Mg-K
5.85 4.89 3.77 0.85 3.16
suggests that this discrepancy may be due to variations in selectivity between components at different active sites on the resin. A more pronounced form of this type of behavior was predicted by de Kock and van Deven-
7876
Ind. Eng. Chem. Res., Vol. 43, No. 24, 2004
ter21 when some active sites were not available for sorption of the most strongly sorbed species present in a binary system. The isotherms predicted by de Kock and van Deventer to occur in those situations resemble to some extent those observed by Vo14 in the obtention of the experimental data used in this work. However, it appears that the variations in selectivity observed experimentally are less pronounced than the complete exclusion from a proportion of active sites proposed by de Kock and van Deventer. Rather, sorption of the most strongly sorbed species on some proportion of the active sites in the resin is somewhat less favorable than that on the remainder of the active sites. Explicit incorporation of this effect into the proposed model on an empirical level, as was done by Valverde et al.29 in development of their many-parameter model, would improve the fit of predictions to the experimental results. However, the aim of this investigation is to develop a model with only a single fitted parameter per binary system. Therefore, this extension to the model has not been undertaken at this timebecause it is beyond the scope of the current investigation. A comparison of model output for different model formulations in the same binary system is presented in Figure 4. The results of fitting the full model to these data, for Ca2+-K+ exchange, are displayed in Figure 2. Figure 4a shows the model fit to the experimental data where ideal solution behavior (unity activity coefficients and no ion association) is assumed. Figure 4b shows the results of including ion association but not the Pitzer activity coefficient model, and Figure 4c includes the Pitzer model but no ion association. Sorption energies were fitted to the experimental data in each case, using the appropriate form of the model. It is clear from these plots that both components of the solution-phase model are required for an accurate description of exchange equilibria in this system. Harvie et al.30 investigated the requirement for explicit incorporation of ion association in association with the Pitzer model and found that such a treatment is desirable in the case of moderately or strongly associating ions. The stability constants used in this work generally fall slightly below the regime in which Harvie et al. state that an explicit treatment of ion association is required. However, it is seen from Figure 4 that inclusion of these effects improves the model fit significantly in this system. Such effects become much more significant in more strongly associating systems, to which this model has also been shown to fit accurately.23,31 The ability of the proposed model to accurately predict the effect of the solution concentration on the position of equilibrium is displayed in Figure 5. Sorption isotherms for Mg2+ in the presence of different monovalent cations at different normalities are presented as an example of this, with the model accuracy seen to be unaffected by changes in the ionic strength. The accuracy of the model predictions compares very favorably with those of other models. The behavior of the binary system Mg2+-Na+-Cl- is actually predicted significantly more accurately by the proposed model than by the three-parameter model of Shallcross et al.11 Other binary systems are predicted slightly less accurately, as would be expected in the comparison of a singleparameter model with a three-parameter model. Analysis and Validation of the Proposed Model From the isotherms for Mg2+-containing systems as presented in Figure 5, it may be seen that the proposed model fits experimental data accurately over a range
Figure 4. Comparison of model predictions to experimental data for the Ca2+-K+-Cl- system for different model formulations: (a) assuming an ideal solution phase; (b) incorporating ion association but not the Pitzer model; (c) incorporating the Pitzer model but not ion association. Parameters are refitted using each model formulation. The results for a full model formulation are shown in Figure 2.
of compositions with a single value of ∆EAfB for each pair of cations. This invariance of fitted parameters with changes in the composition is a fundamental requirement17 for any model of ion-exchange equilibria. Similar behavior is observed for every cation pair investigated, with the accuracy of the model fit to the experimental data never observed to be significantly affected by the solution concentration.
Ind. Eng. Chem. Res., Vol. 43, No. 24, 2004 7877
Figure 5. Experimental and predicted sorption isotherms for the system Mg2+-M+-Cl-, where M+ represents H+, Na+, or K+, at504 concentrations of (a) 0.10 N, (b) 0.20 N, (c) 0.50 N, and (d) 1.00 N.
The second previously identified requirement that must be met by an ion-exchange model is that the fitted parameters are independent of the nonexchanging anion present. Unfortunately, the data set being utilized in this investigation contains only data for systems with Cl- as a nonexchanging anion, and so a full examination of the effect of different nonexchanging anions is not possible in the current work. However, this model has been applied to a limited range of data sets for Na+H+ exchange in the presence of Cl-, F-, and NO3-, and the model parameters were found to be independent of the anion present.31 Therefore, it may be considered that this criterion is met by the proposed model. The final requirement that must be met by a model for binary ion-exchange equilibria is that the model must be internally consistent. Internal consistency is examined using a modified form of the test of Sircar and Myers.32 This test states that when the model is fitted to the three binary systems A-B, A-R, and B-R, where R is a third, “reference”, component, eq 24 must hold.
∆EAfB ) ∆EAfR - ∆EBfR
(24)
Fundamentally, if this criterion holds, the predictive ability of the model is such that the ion-exchange equilibrium between species A and B may be predicted accurately using data obtained from mixtures of each of A and B with a third species R. This has important implications for the application of the model in “real” contexts, where equilibria must often be predicted for systems with little or no experimental data available.
Figure 6 compares the fitted values of ∆EAfB to those calculated using each of the various cations as a reference. This graph shows that there is little difference in each system between the fitted value of ∆EAfB and those calculated using each of the other cations as a reference. For example, in the Na+-H+-Cl- system, the fitted value of ∆ENafH is 1.00, compared to 1.39 calculated using K+ as a reference, 1.34 using Ca2+, or 1.08 using Mg2+. The largest discrepancies appear to be in the prediction of ∆ENafK using a divalent cation (Ca2+ or Mg2+) as a reference. These difficulties can largely be attributed to the fact that the fitted value of ∆ENafK is very small, meaning that the selectivity between these two ions is low when they are present together in a binary system. However, in combination with a more strongly sorbed divalent cation L2+, a significant difference between ∆ELfNa and ∆ELfK is observed. Similar effects are observed in the analysis of the model parameters of Valverde et al.,29 where the behavior of the Na-K system is significantly different from that which would be predicted from the Ca-Na and Na-K systems. This effect is attributable to the difference in hydrated cation size between Na+ and K+, which leads to variations in the accessibility of different active sites within the resin by these cations when in competition with a divalent cation. However, it may be concluded from Figure 6 that the model is generally internally consistent and therefore meets the requirements of this criterion. Figure 7 shows the effect of the choice of the reference cation on the overall weighted error of the model
7878
Ind. Eng. Chem. Res., Vol. 43, No. 24, 2004
Figure 6. Comparison of fitted values of ∆EAfB with those calculated using each component as a reference, plotted by system (in the form A-B) and reference components. The values for each system shown with either component A or B as the reference are the fitted values.
predictions when compared to experimental data. The error values are calculated according to eq 25. This equation differs from the objective function used in the model fitting (eq 23) only in the use of the absolute value rather than the square of the difference between model predictions and in the normalization by the sum of the weighting function used. This provides a result viewable as an actual percentage error. % Err ) 100 ×
∑W(x ,y ) [y j
j
0,j,mod
+ |yA,j,mod - yA,j,exp| + |yB,j,mod - yB,j,exp|]
j
∑W(x ,y ) j
j
j
(25)
However, the use of the absolute value rather than the square of the error in reporting of the errors as displayed in Figure 7 leads to the slightly unusual situation that, in some cases, the parameters fitted to a particular system provide a comparison with experimental data that appears to be less accurate than the parameters obtained by use of certain reference components. This is solely due to the different measures of error used and has no real physical significance. The results presented in Figure 7 further support those in Figure 6, showing that predicting the behavior of Na+-containing systems using K+ as a reference component, and vice versa, leads to results that are relatively inaccurate. This effect is particularly pronounced when the second component of the system being modeled is divalent. The effect of the choice of the reference cation on the Ca2+-Mg2+-Cl- system is seen to be minimal. The prediction of the 1-1 system behavior using monovalent reference components gen-
Figure 7. Errors in model application to each data set for (a) 2-2 and 2-1 systems and (b) 1-1 systems. Points where a component of the system is shown as the reference component represent results of the parameter fitting.
erally produces results that are better than those with divalent reference components, although this effect is less pronounced in the Na+-H+-Cl- system than in the others investigated. It must also be noted that the system Ca2+-Na+Cl- appears from Figure 7a to be described somewhat poorly by this model. Difficulties in modeling this system occur primarily in the intermediate composition region, where the model predicts a greater selectivity for Ca2+ than is observed experimentally. However, the experimentally observed variation in selectivity with ionic strength is reproduced accurately by the model, a result possibly almost as important as the actual fit to the experimental data sets. Conclusion A model is presented that provides a generally accurate description of binary ion-exchange equilibria between alkali and alkaline earth cations on a commercially available strong acid cation-exchange resin with the use of a single fitted parameter per binary system. The proposed model incorporates the Pitzer model for calculation of ion activity coefficients in solution, a description of ion association in a concen-
Ind. Eng. Chem. Res., Vol. 43, No. 24, 2004 7879
trated solution, and a statistical thermodynamic treatment of the exchange reactions occurring on a heterogeneous resin phase. The model is able to predict the behavior of all binary combinations of Ca2+, Mg2+, Na+, K+, and H+ in the presence of Cl- as a nonexchanging anion. Fitted parameters are independent of the solution concentration, and the model is seen to be internally consistent. Model results compare favorably with those of other available models, particularly when taking into consideration the fact that only a single fitted parameter per binary system is used. The prediction of the behavior of binary systems by use of data from other binary systems containing each component in the presence of a common reference component is also possible, although care is required in some cases in the selection of an appropriate reference component. Extension of the model to multicomponent systems, other nonexchanging anions, and different exchanger phases is now possible. Literature Cited (1) Vo, B. S.; Shallcross, D. C. Multi-component ion exchange equilibria prediction. Chem. Eng. Res. Des. 2003, 81, 1311. (2) Gaines, G. L.; Thomas, H. C. Adsorption studies on clay minerals. 2. A formulation of the thermodynamics of exchange adsorption. J. Chem. Phys. 1953, 21, 714. (3) Smith, R. P.; Woodburn, E. T. Prediction of multicomponent ion exchange equilibria for the ternary system SO42--NO3--Clfrom data of binary systems. AIChE J. 1978, 24, 577. (4) Dranoff, J. S.; Lapidus, L. Equilibrium in ternary ion exchange systems. Ind. Eng. Chem. 1957, 49, 1297. (5) Lukey, G. C.; van Deventer, J. S. J.; Shallcross, D. C. Equilibrium model for the sorption of gold cyanide and copper cyanide on trimethylamine ion-exchange resin in saline solutions. Hydrometallurgy 2001, 59, 101. (6) Shehata, F. A.; El-Kamash, A. M.; El-Sorougy, M. R.; Aly, H. F. Prediction of multicomponent ion-exchange equilibria for a ternary system from data of binary systems. Sep. Sci. Technol. 2000, 35, 1887. (7) de Lucas, A.; Zarca, J.; Can˜izares, P. Ion-exchange equilibrium of Ca2+, Mg2+, K+, Na+, and H+ ions on Amberlite IR-120: Experimental determination and theoretical prediction of the ternary and quaternary equilibrium data. Sep. Sci. Technol. 1992, 27, 823. (8) Pitzer, K. S. Thermodynamics of electrolytes. I. Theoretical basis and general equations. J. Phys. Chem. 1973, 77, 268. (9) Pitzer, K. S. In Activity Coefficients in Electrolyte Solutions; Pitzer, K. S., Ed.; CRC Press: Boca Raton, FL, 1991; pp 75-153. (10) Elprince, A. M.; Babcock, K. L. Prediction of ion-exchange equilibria in aqueous systems with more than two counter-ions. Soil Sci. 1975, 120, 332. (11) Shallcross, D. C.; Herrmann, C. C.; McCoy, B. J. An improved model for the prediction of multicomponent ion exchange equilibria. Chem. Eng. Sci. 1988, 43, 279. (12) Mehablia, M. A.; Shallcross, D. C.; Stevens, G. W. Prediction of multicomponent ion exchange equilibria. Chem. Eng. Sci. 1994, 49, 2277. (13) Kester, D. R.; Pytkowicz, R. M. Sodium, magnesium, and calcium sulfate ion-pairs in seawater at 25 °C. Limnol. Oceanogr. 1969, 14, 686. (14) Vo, B. S. Ph.D. Thesis, The University of Melbourne, Melbourne, Australia, 2004. (15) Wilson, G. M. Vapor-liquid equilibrium. XI. A new expression for the excess free energy of mixing. J. Am. Chem. Soc. 1964, 86, 127. (16) Allen, R. M.; Addison, P. A.; Dechapunya, A. H. The characterization of binary and ternary ion-exchange equilibria. Chem. Eng. J. Biochem. Eng. J. 1989, 40, 151.
(17) Shallcross, D. C. Modelling multi-component ion exchange equilibrium behaviour. J. Ion Exch. 2003, 14 (Suppl), 5. (18) Melis, S.; Cao, G.; Morbidelli, M. A new model for the simulation of ion exchange equilibria. Ind. Eng. Chem. Res. 1995, 34, 3916. (19) Myers, A. L.; Byington, S. In Ion Exchange: Science and Technology; Rodrigues, A. E., Ed.; Nijhoff: Dordrecht, The Netherlands, 1986; pp 119-145. (20) de Kock, F. P.; van Deventer, J. S. J. An evaluation of isotherms in the description of competitive adsorption equilibria. Chem. Eng. Commun. 1997, 160, 35. (21) de Kock, F. P.; van Deventer, J. S. J. Statistical thermodynamic model for competitive ion exchange. Chem. Eng. Commun. 1995, 135, 21. (22) Lukey, G. C.; van Deventer, J. S. J.; Shallcross, D. C. Equilibrium model for the selective sorption of gold cyanide on different ion-exchange functional groups. Miner. Eng. 2000, 13, 1243. (23) Provis, J. L.; Lukey, G. C.; Shallcross, D. C. Singleparameter model for ion exchange equilibria applied to multicomponent systems. Ind. Eng. Chem. Res. 2004, in press. (24) Ross, S.; Olivier, J. P. On Physical Adsorption; Interscience: New York, 1964. (25) Lea˜o, V. A.; Lukey, G. C.; van Deventer, J. S. J.; Ciminelli, V. S. T. The dependence of sorbed copper and nickel cyanide speciation on ion-exchange resin type. Hydrometallurgy 2001, 61, 105. (26) Harvie, C. E.; Weare, J. H. The prediction of mineral solubilities in natural waters: The Na-K-Mg-Ca-Cl-SO4-H2O system from zero to high concentration at 25 °C. Geochim. Cosmochim. Acta 1980, 44, 981. (27) Pitzer, K. S. Thermodynamics of electrolytes. V. Effects of higher order electrostatic terms. J. Solution Chem. 1975, 4, 249. (28) Galassi, M.; Davies, J.; Theiler, J.; Gough, B.; Jungman, G.; Booth, M.; Rossi, F. GNU Scientific Library Reference Manual, 2nd ed.; Network Theory Ltd.: Bristol, U.K., 2003. (29) Valverde, J. L.; de Lucas, A.; Rodriguez, J. F. Comparison between heterogeneous and homogeneous MASS action models in the prediction of ternary ion exchange equilibria. Ind. Eng. Chem. Res. 1999, 38, 251. (30) Harvie, C. E.; Møller, N.; Weare, J. H. The prediction of mineral solubilities in natural waters: The Na-K-Mg-Ca-HCl-SO4-OH-HCO3-CO3-H2O system to high ionic strength at 25 °C. Geochim. Cosmochim. Acta 1984, 48, 723. (31) Provis, J. L.; Lukey, G. C.; Shallcross, D. C. An improved model for binary system ion exchange equilibria. IEX 2004sIon Exchange Technology for Today and Tomorrow, Cambridge, U.K., 2004. (32) Sircar, S.; Myers, A. L. A thermodynamic consistency test for adsorption from binary liquid mixtures on solids. AIChE J. 1971, 17, 186. (33) Roine, A. HSC Chemistry 4.1; Outokumpu Research Oy: Pori, Finland, 1999. (34) Mehablia, M. A. Ph.D. Thesis, The University of Melbourne, Melbourne, Australia, 1994. (35) de Robertis, A.; Rigano, C.; Sammartano, S.; Zerbinati, O. Ion association of Cl- with Na+, K+, Mg2+ and Ca2+ in aqueous solution at 10 < T < 45 °C and 0 < I < 1 mol L-1. A literature data analysis. Thermochim. Acta 1987, 115, 241. (36) Majer, V.; Sˇ tulı´k, K. A study of alkaline-earth metal complexes with fluoride and chloride ions at various temperatures by potentiometry with ion-selective electrodes. Talanta 1982, 29, 145.
Received for review May 17, 2004 Revised manuscript received September 6, 2004 Accepted September 11, 2004 IE049581H
