1840
The Journal of Physical Chemistry, Vol. 83, No. 14, 1979
T. L. Croxton and D. A. McQuarrie
Numerical Solution of the Born-Green-Yvon Equation for the Restricted Primitive Model of Ionic Solutions Thomas L. Croxtont and Donald A. McQuarrIe" Departments of Chemistry, Indiana University, Bloomington, Indiana 4740 1, and Universlty of Callfornia, Davis, California 956 16 (Received February 13, 1979) Publication costs assisted by the National Institutes of Health
The Born-Green-Yvon integral equation under the Kirkwood superposition approximation is formulated for the restricted primitive model of electrolyte solutions. An iterative solution of this equation is discussed and results are presented for 1:l and 2:2 electrolytes at concentrationsup to 2 M. Various thermodynamic quantities are calculated and compared to hypernetted chain and Monte Carlo data. Agreement is found to be very good for both 1:l electrolyte and 2:2 electrolyte solutions,except at the higher concentrations of 2:2 electrolyte solutions. It is suggested that the BGY formalism can provide numerically accessible theories for more complex Coulombic systems.
Introduction Coulombic chemical systems involve an intermolecular potential which may be written in the form ucj(R) = u,*(R)
+ u,le1(R)
(1)
where u,,*(R) is a short-ranged potential, ubjel(R)is the familiar electrostatic term, and the subscripts index the species. The particular difficulty in the theoretical treatment of such systems stems from the fact that the two contributions are of comparable importance but have quite dissimilar character. For the case of electrolyte solutions, consideration of the short-ranged potential would suggest a virial expansion, but the virial coefficients diverge because of the electrostatic term. On the other hand, a simple electrostatic approach yields the Debye-Huckel theory which is numerically inaccurate at other than very low concentrations partly because of the importance of the core potential. In practice a rigorous theory is required to treat both parts of the potential in a satisfactory manner. Many theories have been tested for a special model of ionic solutions called the restricted primitive model. This model consists of hard spheres of equal diameters with point charges at their centers immersed in a continuum solvent of uniform dielectric constant E . This model is important since Monte Carlo (MC) data are available for it and since it is believed to represent the significant interaction in real electrolytes. Indeed, it has been shown that the simple addition of a square well potential and different ionic radii provides good agreement with experimental data for 1:l electrolytes up to 1.0 M.1,2 The quality of a given theory is indicated not only by its agreement with MC calculations but also by its internal consistency, indicated by the consistency of thermodynamic quantities when calculated by different thermodynamic routes. By these criteria most theories, such as the mean spherical approximation, show serious deficiencies as the ionic charge and concentration increase. The hypernetted chain (HNC) integral equation, however, has proven to be an excellent theory. A disadvantage of the HNC equation is the difficulty in obtaining numerical solutions. Also, it +Supportedin part by the Medical Science Department, Indiana University, Bloomington, IN. *Address correspondence to this author at the University of California. 0022-3654/79/2083-1840$0 1.OO/O
is not always obvious how to extend the HNC equation to systems of more complex geometry such as that involved in a calculation of surface tension of ionic solutions. We feel it worthwhile, then, to investigate a somewhat simpler theory, hoping that it will provide thermodynamic calculations comparable to HNC while maintaining conceptual and numerical simplicity. To that end, we have effected a numerical solution of the Born-Green-Yvon (BGY) integral equation for the restricted primitive model since its derivation suggests that equivalent theories may be more easily developed for systems of less symmetric geometry. For example, a theory of the surface tension of electrolyte solutions3 and a theory of the electrical double layer4 have been developed previously through the BGY and Kirkwood formalism and the superposition approximation. The BGY equation has been solved for a variety of shorbranged potentials. For example, hard sphere: square and Lennard-Jones' fluids have been investigated into the region of very high density. Generally, problems of numerical convergence develop as the density is increased. The study of ionic solutions requires two extensions: (1)generalization to multiple components and (2) treatment of the long-ranged potential. The first is well known and, by using the notation of the BGY equation for w components is
uij(R) In gij(R)= -- kT
+
where
where g,(R) is the radial distribution function of type j particles about a central particle of type i, h,j(R) = g@) - 1, P k is the bulk density of component k, and other variables have their usual meaning. Ono9 has obtained Debye-Huckel results from an equivalent equation in the limit of infinite dilution and MollerlO has solved a linearized form for ionic solutions. The form of eq 2 is most commonly used for computations which involve only a short-ranged p ~ t e n t i a l but ~ - ~ is unsuitable for ionic solutions since the kernals, K,(t),are individually divergent. It is possible, however, to algebraically cancel infinite 0 1979 American Chemical Society
The Journal of Physical Chemistty, Vol. 83, No. 14, 1979
Restricted Primitive Model of Ionic Solutions portions of the kernal difference and obtain a finite expression. The BGY equation for the restricted primitive model with ~ij*(R)= R