J. Phys. Chem. 1996, 100, 1847-1851
1847
Ion Adsorption in the Electrical Double Layer: Variational Solution to the Hypernetted Chain Equations Jeffery A. Greathouse and Donald A. McQuarrie* Department of Chemistry and Institute of Theoretical Dynamics, UniVersity of California, DaVis, California 95616 ReceiVed: August 4, 1995; In Final Form: October 24, 1995X
In order to investigate the effect of specific adsorption of counterions onto a charged planar surface, the HNC/MSA equations were solved using a variational method. Ion adsorption is introduced through a sticky potential. The variational method involves a simple computational algorithm, and although it does not yield exact solutions, the contact condition and electroneutrality are automatically satisfied. Results are presented for a surface charge density of -0.10 C m-2 and for a range of bulk ion concentrations of 1:1 electrolyte and counterion adsorption tendencies. While the sticky potential provides a simple model of wall-ion interactions, the effect of increasing ion adsorption on the structure of the double layers can be seen. Comparisons are also made with the modified Gouy-Chapman theory, which neglects the finite size of the ions. At high bulk ion concentrations, ion size effects result in nonmonotonic ion distributions.
Introduction The distribution of ions about a charged surface is known as the electrical double layer and is the focus of much research in physical chemistry. The standard theoretical treatment that has enjoyed great success is that of Gouy and Chapman,1 in which the ions are modeled as point charges. The modified GouyChapman (MGC) theory2,3 introduces ionic size in that the ions are excluded from a region around the surface representing one ionic radius. However, the ions still interact with each other as point charges. As a result, MGC theory fails for cases of high surface charge, multivalent ions, and large ion concentrations.4,5 Theories of the electrical double layer based on statistical mechanical theories of liquids have introduced ion size effects. These treatments begin with the Ornstein-Zernicke equation and involve modeling the charged surfaces as planes of infinite extent. The Ornstein-Zernicke equation separates the total correlation between the wall and an ion into a direct correlation and an indirect correlation. However, an approximation must be made to provide closure to the resulting integral equations. The hypernetted chain/mean spherical approximation (HNC/ MSA) has been applied to the restricted primitive model of planar double layers,6-8 showing excellent agreement with Monte Carlo computer simulations. In this theory, the hypernetted chain approximation is applied to the wall-ion direct correlation functions while the mean spherical approximation is used for the ion-ion direct correlation functions. The MSA treats the ion-ion correlations according to their bulk values, neglecting the presence of the wall. However, the MSA has received widespread use because it has an analytic solution and it yields excellent results. Improvements of the HNC/MSA theory have included wall-ion bridge functions,9,10 which are neglected in the HNC closure. These results compare well with Monte Carlo calculations presumably because the bridge functions provide a more accurate description of the ionic structure near the wall, while the HNC/MSA uses the bulk ion correlations in this area. Despite these improvements in HNC theory, very little has been done to investigate another source of interaction near the X
Abstract published in AdVance ACS Abstracts, January 1, 1996.
0022-3654/96/20100-1847$12.00/0
wall, namely, that of ion adsorption onto the charged surface. Results from clay swelling experiments11 show dramatic effects due to the type of counterion used, and these differences have been attributed to the different adsorption tendencies of the ions. Recently, Glosli and Philpott12,13 have performed molecular dynamics simulations on double layer systems containing a single ion pair in a molecular solvent. Their results show that some ions contact adsorb onto the surface while others remain fully solvated and do not adsorb. One model of ionic adsorption that has been used because of its ease of incorporation into the HNC/MSA theory is Baxter’s sticky potential,14 in which an ion experiences an additional attraction to the wall only when the two come into contact. Using this approach, Carnie and Chan15 solved the HNC/MSA equations for adsorbing ions near a charged wall. More recently, Marcelja16 used the anisotropic HNC theory for interacting double layers to study counterion adsorption effects. The anisotropic HNC theory uses the hypernetted chain approximation for both wall-ion and ionion correlations, and counterion adsorption is included in the form of a negative adsorption energy for counterions that come into contact with the surface. In a comparison of this theory to atomic force microscope data on mica surfaces, it was found that an adsorption energy of 2.5kT (k is the Boltzmann constant, T is the kelvin temperature) provided the best fit to the experimental data.17 The aim of this work is to solve the HNC/MSA equations for a planar double layer with counterion adsorption using a variational method of solution, which employs a simple computational scheme involving numerical integration and multivariate optimization. Variational methods have been used previously to solve integral equations for bulk electrolytes18-20 and for planar double layers involving constant surface charge8 and constant surface potential.21 The variational solutions are expressed in terms of analytic trial functions, allowing trends in the solutions to be seen as double layer conditions change. This technique avoids the problems of convergence present in numerical solutions to these equations,22 and boundary conditions such as electroneutrality can be built into the algorithm so that they are always satisfied. However, the variational solutions are not exact, achieving only the accuracy of the trial functions. The model consists of a planar sheet of infinite extent © 1996 American Chemical Society
1848 J. Phys. Chem., Vol. 100, No. 5, 1996
Greathouse and McQuarrie second and third terms on the right side of eq 1 vanish, and the HNC/MSA equations reduce to the modified Gouy-Chapman theory. Using the electroneutrality condition for this system
σ0 + σs + ∑zjeρj∫0 dy hj(y) ) 0 ∞
(3)
j
where σs ) ∑jzjeFjSj[1 + hj(0)] is the adsorbed surface charge density, the mean electrostatic potential can be expressed in terms of the charge densities σ0 and σs and the functions hi(x),
ψ(x) )
Figure 1. Geometry and coordinate system used in the HNC/MSA equations.
with a constant surface charge density σ0 in contact with an infinite reservoir of a 1:1 electrolyte solution. The solvated ions are represented as hard spheres of diameter D with charges +1 or -1 embedded in their centers, and the solvent is taken to be a continuum of dielectric constant . The effect of images and the molecular nature of the solvent are neglected. Counterion adsorption onto the surface results in an adsorbed charge density σs forming at the plane of closest approach, located at a distance D/2 from the wall. Theory In order to describe the structure of a liquid about a planar interface, the Ornstein-Zernicke equation is used. Using the procedure of Henderson and Blum,2 the charged wall is represented as a species at infinite dilution whose diameter and charge approach infinity while the surface charge density σ0 remains constant. The interaction potential between the wall and ions includes a Coulombic term as well as an adsorption term. This adsorption term is controlled by a parameter, S, which has units of length and is related to the adsorption potential by a Dirac delta function. A relationship between the adsorption potential and S can be derived,15 and for the case of a square well a plot of the well depth versus S is found in a previous study.23 Using the hypernetted chain approximation for the wall-ion closure and the mean spherical approximation for the ion-ion direct correlation functions, one obtains the HNC/MSA equations for this system15
ln[1 + hi(x)] ) -zieψ(x) kT
-x 1 ∞ (σ0 + σs) - ∑zjeρj∫0 dy hj(y) Max(x,y) (4) 0 0 j
where 0 is the permittivity of free space and Max(x,y) ) (x + y + |x - y|)/2. Substitution of eq 4 into eq 1 results in a set of coupled, nonlinear, integral equations for h+(x) and h-(x). Following the method of Feller and McQuarrie,8 the first step in the variational method of solution is to derive a functional F[h(x)] such that the functional derivative δF[h(x)]/δh(x) ) 0 for all values of x. If δF[h(x)]/δhi(x) ) B[hi(x)], then the functional can be written18
F[h(x)] ) F0 + ∑∫0 dx hi(x)∫0 dt B[thi(x)] ∞
where F0 is a constant. In this case, B[hi(x)] is found by rearranging eq 1:
B[hi(x)] ) ln[1 + hi(x)] zieψ(x) ∞ - 2π∑ρj∫-∞dy hj(y) CijMSA(|x - y|) kT j
2π∑ρjSj[1 + h+(0)]CijMSA(x) (6) j
Written in this manner, B[hi(x)] ) δF[h(x)]/δh(x) ) 0, and eq 5 can be used to obtain the functional F[h(x)]. The functional thus obtained is
∞
∞
i
x > 0 (1)
j
hi(x) ) -1
x 0, where x is the distance from the plane of closest approach to the centers of the ions. Figure 1 shows the coordinate system. In the case where D ) 0, the
2πhi(x)∑ρjSj[1 + hj(0)]CijMSA(x) j
}
(7)
After constructing a trial function with variational parameters hTi(x, ai1, ..., ain), the functional is minimized when
( ) δF
∑i δa
dai ) 0
(8)
i
Equation 8 is satisfied at the minimum of F with respect to the variational parameters {an}, and the functions hTi represent variational solutions to eq 1. In constructing a suitable trial function, several factors must be considered. First, the electroneutrality condition must be automatically satisfied. Second, the contact condition15 which relates the values of h+(0) and h-(0) to the charge densities σ0
Ion Adsorption in Electrical Double Layers
J. Phys. Chem., Vol. 100, No. 5, 1996 1849
and σs, should also be satisfied. Without adsorption a reasonable approximation is to set hi(0) equal to the MGC values,6,8 but in the present case the contact condition involves an integral of the ion-ion direct correlation functions which extends to x ) D. Therefore the contact condition must be included in the trial function to control the choices for h+(0) and h-(0). Finally, because the functional depends on the hi at all values of x, the terms containing hj(0) in eq 7 complicate the selection of a trial function. To rectify this, eq 1 can be used to write h+(0) as integrals over h+(x) and h-(x), and this can be substituted into eq 7. While a trial function consisting of the solution to the Poisson-Boltzmann equation plus an expansion in Laguerre polynomials was successfully used for the case in which there is no adsorption,8 this type of trial function was useful in the present work only at bulk ion concentrations of 0.010 and 0.050 M, where the HNC/MSA results and MGC results are nearly identical. For concentrations greater than 0.050 M, the following trial functions are used:
h-(x) )
{
a1-ex/2D + a2-ex/D + a3-e-κx, x < D e-κx(a4- + a5-κx), x>D
h+(x) )
{
(9)
a1+e-x/D + a2+e-3x/D, x < D e-κx(a3+ + a4+κx), x>D
where the ai are the variational parameters and κ ) (2Fe2/ 0kT)1/2 is the Debye-Huckel screening parameter. The use of the exponential function simplifies the integral evaluations needed to enforce the restrictions on the trial functions discussed previously. The presence of separate solutions for x < D and x > D is a consequence of the mean spherical approximation. The last term in eq 1 is due to ion-ion interactions between adsorbed ions and diffuse ions. This term is nonzero only when x < D, however, because CijMSA(x) is discontinuous at x ) D. Additionally, the parameters are chosen so that the hi are continuous across x ) D. As shown in the Results section, the trial functions given by eqs 9 are flexible enough to describe wall-ion distributions for a wide range of double-layer conditions, especially when the differences between the MGC results and the HNC/MSA results are large. Finally, it is important to note that any choice of the nine parameters in eqs (9) will not automatically satisfy the constraints mentioned above. In implementing the optimization algorithm, four of the parameters in eqs 9 are independent variables: a3-, a1+, a2+, and a3+. The parameter a4+ is found by enforcing the constraint on h+(D). The remaining parameters (a1-, a2-, a4-, and a5-) are found by simultaneously solving the four equations produced by electroneutrality (eq 3), the contact condition, eq 1, and the constraint on h-(D). Choosing the independent parameters in this manner is effective because the value of h+(0) ) a1+ + a2+ is chosen first. This choice then determines the value of the adsorbed surface charge, which affects the net surface charge acting on the diffuse ions. Results Calculations have been performed on a 1:1 electrolyte with bulk electrolyte concentrations from 0.010 to 1.0 M at a temperature of 25 °C and a surface charge density of -0.10 C m-2. The ionic diameter D was assigned a value of 4.25 Å. The co-ion adsorption parameter, S-, was zero for all cases, and the counterion adsorption parameter, S+ ) S, ranged from 0 to 10D. Comparisons are made with the MGC theory, which neglects the finite size of the ions but maintains a plane of closest approach for the ions. The MGC results involve solving
Figure 2. Wall-ion total correlation functions vs the distance from the plane of closest approach for σ0 ) -0.10 C m-2 and F ) 0.10 M. The labeled curves correspond to adsorption parameters of S/D ) 0.0, 5.0, and 10.0.
the nonlinear Poisson-Boltzmann equation with adsorption using a procedure similar to that of Andrietti et al.24 The ion distributions depend solely on the potential at the plane of closest approach, ψs, which is related to the surface charge density and adsorbed charged density through a cubic equation. All HNC/MSA results presented here for bulk ion concentrations greater than 0.050 M were obtained using the trial functions from eq 9. The functional was minimized with respect to the variational parameters using a multivariate optimization algorithm from the IMSL library.25 The integrations performed in the functional evaluations were carried out with a Gauss quadrature numerical integration subroutine, also from the IMSL library.26 All work was done on a Vax station 3100 M38. Figure 2 shows the dependence of the ionic distributions h+ and h- on the reduced counterion adsorption parameter S/D for a 0.10 M electrolyte solution. The distance from the plane of closest approach x is given in units of the ionic diameter D, and hi ) 0 corresponds to the bulk ion concentration. For a comparison of S to the adsorption energies used by Marcelja,16,17 S/D values of 5.0 and 10.0 correspond to square well depths of 1.8kT and 2.4kT, respectively. For the S ) 0 case, h+(0) ) 30.0, although it has been truncated. As noted previously, the cusps in h+ and h- at x/D ) 1 for S/D values of 5.0 and 10.0 are a consequence of the mean spherical approximation. As S increases, the value of h+(0) decreases, while the value of h-(0) increases. The adsorbed charge densities are 0.082 and 0.094 C m-2 for S/D ) 5.0 and 10.0, respectively. This near cancellation of the surface charge density σ0 results in counterion concentrations almost equal to the bulk values even at close distances to the wall. Even for the case of S/D ) 10.0 with a net surface charge density of -0.0060 C m-2 acting on the diffuse ions, the co-ions are almost completely excluded from the wall due to ion-ion interactions with adsorbed counterions. Another feature of these results is that, as S increases, the ion distributions reach their bulk values at closer distances to the
1850 J. Phys. Chem., Vol. 100, No. 5, 1996
Figure 3. Ratio of adsorbed charge density σs to surface charge density σo vs the bulk ion concentration. The solid and dashed curves represent S/D ) 5.0 and 10.0, respectively.
Figure 4. Wall-ion total correlation functions vs the distance from the plane of closest approach for σo ) -0.107 C m-2, F ) 0.50 M, and S/D ) 5.0. The solid and dashed curves represent the HNC/MSA and MGC results, respectively.
wall, indicating that the wall-counterion adsorption strength influences the thickness of the double layer. Figure 3 shows the ratio of adsorbed charge density to surface charge density (σ0 ) -0.10 C m-2) versus the bulk ion concentration for S/D ) 5.0 and 10.0. A value of 1.0 corresponds to a cancellation of surface charge, and values greater than 1.0 represent a charge reversal where a net positive charge acts on the diffuse layer. This reversal occurs around 0.35 M for S/D ) 5.0 and 0.15 M for S/D ) 10.0. In both cases the curves become steep as F f 0 since σs f 0 also. The lowest bulk ion concentration at which calculations were performed was 0.010 M. In this case, the trial function used consisted of the MGC solution plus Laguerre polynomials, which was used for calculations without adsorption.8 The results presented in Figures 4 and 5 show a comparison between the HNC/MSA results and the MGC results. Figure 4
Greathouse and McQuarrie
Figure 5. Reduced mean electrostatic potential eψ/kT vs the distance from the plane of closest approach. The conditions and meaning of the curves are the same as in Figure 4.
shows the ion distributions for σo ) -0.107 C m-2, F ) 0.50 M, and S/D ) 5.0. The HNC/MSA results are solid lines, and the MGC results are represented by dashed lines. For the MGC results, the adsorbed charge density is 0.106 C m-2, which produces a very small but net negative charge acting on the diffuse layer. As a result, the ion concentrations are nearly zero even at contact. For the HNC/MSA results, the adsorbed charge density is 0.135 C m-2, indicating that a net positive charge is acting on the diffuse layer. This reversal of charge is evident in the ion distributions at very close distances from the wall. In this case ion size effects are responsible for the nonmonotonic behavior in the HNC/MSA theory, a feature which is not possible in MGC theory. For both h+ and h- in the region x > D, there are no longer any correlations between adsorbed counterions and diffuse ions, so the distributions resemble those of nonadsorbing ions near a positive wall. Figure 5 shows the mean electrostatic potential versus distance from the plane of closest approach for the same case. Again, the solid curve corresponds to the HNC/MSA results while the dashed curve corresponds to the MGC results. The values of eψs/kT are -0.030 (MGC) and 0.51 (HNC/MSA). While the potential in both cases approaches zero at nearly the same distance from the wall, the opposite signs of ψ clearly demonstrate the differences between the two theories. As a check of the accuracy of the variational solution, the numerical results obtained by Carnie and Chan15 for this case are σs ) 0.131 C m-2 and eψs/kT ) 0.52, in excellent agreement with the variational results presented here. Discussion We have used a variational method to solve the HNC/MSA equations for a planar electrical double layer, including the effects of counterion adsorption onto the charged wall. The variational method is a useful tool in solving integral equations of this type because they avoid many of the problems associated with a numerical solution and they do not require large amounts of computer resources necessary for simulations. While numerical solutions do exist for this system,15 to our knowledge no computer simulations have been performed on a similar system for comparison. Additionally, the contact condition, which was satisfied to 1% at 0.10 M and 5% at 0.50 M in the numerical solutions, is automatically satisfied by the trial functions. Results were presented for a surface charge density of -0.10 C m-2 and a bulk ion concentration range of 0.0101.0 M. These are double-layer conditions where the HNC/MSA has showed agreement with Monte Carlo results.
Ion Adsorption in Electrical Double Layers While no attempt has been made to associate a particular value of the adsorption parameter S with a particular cation, an increase in S in our calculations corresponds to an increasing tendency for an ion to adsorb onto the surface. The result of this was an increase in the adsorbed charge density σs as well as a decrease in the thickness of the double layer. We have also shown that, at higher bulk ion concentrations, a charge reversal occurs in which a net positive charge acts on the diffuse ions. However, results of Feller and McQuarrie8 and Kjellander et al.17 without adsorption show that charge reversal occurs in systems with high concentration and divalent ions. Finally, a comparison with the MGC calculations which neglect the size of the ions shows that ion-ion interactions near the adsorbed region of the double layer play a significant role in the doublelayer structure. Acknowledgment. This work has been supported by National Science Foundation Grant NSF EAR-8910530. References and Notes (1) Verwey, E. J. W.; Overbeek, J. T. G. The Theory of the Stability of Lyophobic Colloids; Elsevier: Amsterdam, 1948. (2) Henderson, D.; Blum, L. J. Chem. Phys. 1978, 69, 5441. (3) Torrie, G. M.; Valleau, J. P. Chem. Phys. Lett. 1979, 65, 343. (4) Valleau, J. P.; Ivkov, R.; Torrie, G. M. J. Chem. Phys. 1991, 95, 520.
J. Phys. Chem., Vol. 100, No. 5, 1996 1851 (5) Torrie, G. M.; Valleau, J. P. J. Chem. Phys. 1980, 73, 5807. (6) Lozada-Casssou, M.; Saavedra-Barrera, R.; Henderson, D. J. Chem. Phys. 1982, 77, 5150. (7) Lozada-Cassou, M.; Henderson, D. J. Phys. Chem. 1983, 87, 2821. (8) Feller, S. E.; McQuarrie, D. A. J. Phys. Chem. 1992, 96, 3454. (9) Ballone, P.; Pastore, G.; Tosi, M. P. J. Chem. Phys. 1986, 85, 2943. (10) Attard, P.; Miklavic, S. J. J. Chem. Phys. 1993, 99, 6078. (11) Lubetkin, S. D.; Middleton, S. R.; Ottewill, R. H. Philos. Trans. R. Soc. London, A 1984, 311, 353. (12) Glosli, J. N.; Philpott, M. R. J. Chem. Phys. 1992, 96, 6962. (13) Glosli, J. N.; Philpott, M. R. J. Chem. Phys. 1993, 98, 9995. (14) Baxter, R. J. J. Chem. Phys. 1968, 49, 2770. (15) Carnie, S. L.; Chan, D. Y. C. J. Chem. Phys. 1981, 75, 3485. (16) Marcelja, S. Biophys. J. 1992, 61, 1117. (17) Ke´kicheff, P.; Marcelja, S.; Senden, T. J.; Shubin, V. E. J. Chem. Phys. 1993, 99, 6098. (18) Olivares, W.; McQuarrie, D. A. J. Chem. Phys. 1976, 65, 3604. (19) Baxter, R. J. J. Chem. Phys. 1970, 52, 4559. (20) Andersen, H.; Chandler, D. J. Chem. Phys. 1972, 57, 1918. (21) Feller, S. E.; McQuarrie, D. A. J. Colloid Interface Sci. 1994, 162, 208. (22) Carnie, S. L.; Torrie, G. M. AdV. Chem. Phys. 1984, 56, 141. (23) Greathouse, J. A.; McQuarrie, D. A. J. Colloid Interface Sci., 1995, 175, 219. (24) Andrietti, F.; Peres, A.; Pezzotta, R. Biophys. J. 1976, 16, 1121. (25) Gill, P. E.; Murray, W. Minimization Subject to Simple Bounds on the Variables; NPL Report NAC 72; National Physics Laboratory, England, 1976. (26) Piessens, R.; deDoncker-Kapenga, E.; Uberhuber, C. W.; Kahaner, D. K. QUADPACK; Springer-Verlag: New York, 1983.
JP952238T
