Comparison of density functional and modified Poisson-Boltzmann structural properties for a spherical double layer

The density functional and modified Poisson-Boltzmann descriptions of a spherical (electric) double layer are compared and contrasted vis-a-vis existing Monte Carlo simulation data (for small ion diameter 4.25 · 10 m) from the literature for a range of physical parameters such as macroion surface charge, macroion radius, valencies of the small ions, and electrolyte concentration. Overall, the theoretical predictions are seen to be remarkably consistent between themselves, being also in very good agreement with the simulations. Some modified Poisson-Boltzmann results for the zeta potential at small ion diameters of 3 and 2 · 10 m are also reported.


Introduction
The accumulation of an ionic charge cloud in the vicinity of an electrode in an electrolyte is commonly referred to as an electric double layer with the geometry of the electrode determining the geometry of the double layer.The microscopic characterization of the charge cloud has wide ranging significance in important biological, industrial, and technological processes.Such relevance has made the double layer phenomenon a major research area in colloid and polyelectrolyte science for over three decades, see the recent reviews [1][2][3][4][5][6].
Recently there has been a resurgence of interest in the use of the density functional theory (DFT) as a theoretical probe to study the properties of electric double layers in all geometries.For example, Boda et al. [7,8] have applied the theory to a planar double layer (PDL) containing a restricted primitive model (RPM) electrolyte (equisized rigid ions moving in a continuum dielectric) and have obtained structural and thermodynamic results for different ion sizes.The application to the cylindrical (electric) double layer (CDL) has been considered by Patra and Yethiraj [9,10].Yu et al. [11] have employed the technique to spherical double layers (SDL).In another study Valisko et al. [12] have further treated the PDL containing a primitive model electrolyte but now with unequal ionic radii.In many cases the DFT results from these studies have been compared with the corresponding Monte Carlo (MC) simulation data.A consensus that emerges from these comparisons is that the DFT is a viable approach to the electric double layer theory.It ought to be mentioned here that the DFT was originally applied to the PDL in the early 1990s by various groups [13][14][15][16].
It is natural to wonder how the DFT compares with other more established theories of the electric double layer, viz., integral equation methods and potential based theories.To this end some work has already been reported.In an earlier paper [17] we have compared the DFT results of Boda et al. [7,8] for a RPM planar double layer at ionic diameter of a = 4.25 • 10 −10 m with the corresponding modified Poisson-Boltzmann (MPB) results.With few exceptions the two theories were seen to be on par with each other in reproducing the MC simulation data of Boda et al. [7,8] for the model systems.More recently, Patra and Bhuiyan [18] have compared the DFT and MPB results for a CDL where the calculations were done at physical parameters for a double stranded DNA.Notwithstanding the fact that MC results for the CDL are presently rather sketchy in the literature, the two theories again showed notable consistency in their predictions.
The above findings have encouraged us to extend such comparisons to a double layer involving spherical symmetry.In this paper we will compare the structural descriptions of a SDL arising out of the DFT and the MPB theory.As indicated previously, the DFT for the SDL has recently been solved by Yu et al. [11], the calculations having been made at the ionic diameter of a = 4.25•10 −10 m.These DFT results will be utilized in the present study, and in the rest of the paper "DFT" will refer to this version.The MC simulations at this ionic diameter were obtained earlier by Degreve et al. [19] and Degreve and Lozada-Cassou [20].González-Továr et al. [21] have reported simulation and hypernetted chain/mean spherical approximation (HNC/MSA) integral equation results at higher ionic diameters.
As the name suggests the modified Poisson-Boltzmann formalism arose from modifications to the classical Poisson-Boltzmann (PB) theory to account for (i) the interionic correlations, and (ii) the ionic exclusion volume effects, both of which are ordinarily neglected in the classical mean field theory.The MPB theory was first applied to the SDL by Outhwaite and Bhuiyan [22,23].Some zeta potential results from this study were subsequently used by Degreve and co-workers [19,20] in their MC studies for comparison purposes.Indeed the MPB zeta potentials compared rather well with the simulations and the HNC/MSA results of González-Továr and Lozada-Cassou [24].

Model and theory
The SDL model utilized in this study consists of an isolated spherical macroion of radius R and uniform surface charge density σ, bathed in a RPM electrolyte, the model being equivalent to a macroion at infinite dilution.The common diameter of the small, simple ions is a, while the mean number density and charge of an ion of species s are ρ s and e s = z s |e| (z s is the valency and e the electron charge), respectively.
The development of the MPB theory for the SDL has been described in detail in reference [23].We will therefore outline here only the principal equations.
The mean electrostatic potential ψ(r) satisfies Poisson's equation where r is the radial distance from the macroion centre, and 0 , r are the vacuum and relative permittivities, respectively.In spherical symmetry, the macroion-small ion radial distribution function g s (r) reads where with u(r) = rψ(r) and β = 1/(k B T ) (k B the Boltzmann's constant and T the absolute temperature).The exclusion volume term ξ s (r) is approximated using the Bogoliubov-Born-Green-Yvon hierarchy DFT results from reference [11] and MC data from references [19,20]. and is the fluctuation potential evaluated on the surface of the exclusion volume V of the discharged ion.The equations ( 1)-( 9) are collectively known as the MPB equation for the SDL.The classical PB equation is obtained upon taking ξ(r) = H(r − R), F = F 0 , and a → 0. The MPB and the PB equations were solved numerically using a quasilinearization iterative procedure [25,26] for a range of R, σ, z s , a, and electrolyte concentration c.

Results and discussion
The physical parameters used in the present work follow closely those used in the earlier MC [19,20] and DFT [11] calculations.Thus unless otherwise mentioned  The meaning of the symbol and curves as in figure 1. DFT results from reference [11] and MC data from reference [20].
The singlet distribution functions of the simple ions around the spherical colloid of radius R = 15 • 10 −10 m are shown in figure 1 for 1:1, 1:2 and 2:2 electrolytes.The 1:1 electrolyte has σ = 0.102 C/m 2 with c = 1 mol/dm 3 while the divalent electrolytes have σ = 0.204 C/m 2 with c = 0.5 mol/dm 3 .Comparison with the MC results indicates that both the MPB and DFT theories accurately reproduce the simulation results.For the 1:2 and 2:2 cases both the theories predict a nonmonotonic behaviour in the distribution functions indicative of overscreening.This non-monotonic behaviour cannot be produced by the standard PB theory.Figures 2  and 3 illustrate the ionic singlet distribution functions and the mean electrostatic potential around the colloid for 1:1 and 2:1 electrolytes at c = 1 and 0.5 mol/dm 3 , respectively, with σ = 0.102 and 0.204 C/m 2 .Note that there are misprints in the captions of figures 4 and 6 in reference [11] where it should read 0.204 C/m 2 rather than 0.306 C/m 2 [27].Again the MPB and DFT theories are in close agreement, with both theories having a shallow minimum in the mean electrostatic potential (see insets).As the surface charge is increased we expect the steric effect of the counterions to become important and a second layer of counterions to form at r ∼ R+3a/2.This is seen in the planar electric double layer in both theory and simulation [28][29][30].Provided R is not too small this layering effect is predicted by the DFT theory (figures 4-6 of Yu et al. [11]) but not by the MPB theory due to the breakdown of the numerical procedure at high σ.It is expected that the MPB theory will predict counterion layering on using an alternative numerical technique.In practice, however, such high surface charges are not observed in typical colloidal suspensions.DFT results from reference [11] and MC data from reference [31].
In figure 4 we plot the zeta potential ζ = ψ(R + a/2) as a function of surface charge σ for 1:1, 2:1/1:2, and 2:2 electrolytes at various electrolyte concentrations.Again we use the colloid radius R = 15 • 10 −10 m with also R = 5 • 10 −10 m for the 1:1 case.There is good agreement between the MPB, DFT, and simulation results with all the three electrolyte cases, although the MPB is closer to the MC values for divalent counterions at c = 0.005 mol/dm 3 .The PB theory gives the correct qualitative response for 1:1 electrolytes but fails to predict the maxima and minima for the 2:2 and 2:1 electrolytes, respectively.Figure 5 illustrates the variation of ζ with a/R for 1:2 and 2:1 electrolytes at σ = 0.204 C/m 2 with c = 0.005 and 0.5 mol/dm 3 .The MPB results for 1:1 and 2:2 electrolytes with the same parameters are very similar to the 2:1 and 1:2 graphs, respectively, in figure 5 and hence are not shown.Increasing R for divalent counterions leads to a maximum in the ζ plot which cannot be predicted by the PB theory.The recent HNC/MSA calculations of González-Továr et al. [21] further illustrate the severe limitations of using the PB zeta potential to analyse electrophoretic mobility experiments.These authors studied the behaviour of ζ at varying ionic radii and surface charge for fixed R and c and found that for both 1:1 and 2:2 electrolytes there is a change in curvature sign as the ion size is increased.The non-monotonic behaviour of g i (r) and ψ(r) seen in figures 1-3 indicate that the local density of the counterions is such that the colloidal charge is overcompensated giving rise to overcharging or overscreening.To illustrate this overscreening we consider the integrated charge distribution [11]   At the lowest concentration of 0.05 mol/dm 3 the integrated charge is monotonic, but when the concentration is increased to the higher concentrations of 1 and 1.5 mol/dm 3 a damped oscillatory behaviour sets in.As expected the divalent counterion charge inversion response is more pronounced than that of the univalent counterion at the same concentration.
Finally, in figures 8-10 we present some MPB and PB zeta potential results when the small ion diameter is 3 • 10 −10 and 2 • 10 −10 m, respectively.(The macroion radius has been kept at R = 15 • 10 −10 m.)A motivation for these calculations is the fact that in their DFT and MC work on the PDL at these diameters, Boda et al. [7,8] had noted an increasing discrepancy between the classical PB and DFT (or MC) diffuse layer potentials with decreasing ion size.The behaviour was confirmed by the MPB calculations in the PDL [17], and by the DFT and MPB calculations in the CDL [18].In the present case the trend is clearly visible at the higher absolute values of the surface charge for each of the 1:1, 2:1/1:2, and 2:2 situations with the deviation between the PB and the MPB ζ increasing as the ion size decreases from 3 • 10 −10 m to 2 • 10 −10 m.For divalent counterions the MPB ζ continues to show a maximum (or a minimum), while the PB ζ remains monotonic as in the planar or cylindrical geometries [17,18].

Concluding remarks
This paper has focussed on a comparative study of the structure of an electric double layer in spherical symmetry for a range of physical parameters using the density functional and modified Poisson-Boltzmann theories.Although quite distinct in their origins, the theories show a remarkable consistency overall for the cases studied, and their further predictions closely follow the simulation data.
The DFT and MPB approaches to the electric double layer have now been compared for planar [17], cylindrical [18], and spherical geometries.The accumulated evidence clearly suggests that the aforementioned consistency is global, extending to all the three symmetries, and that under ordinary laboratory conditions, viz., σ 0.3 C/m 2 and electrolyte solution regimes, the two theories are comparable in their predictive properties.

Figure 5 .
Figure5.Variation of the zeta potential for 2:1 and 1:2 electrolytes as a function of inverse macroion radius at a fixed macroion surface charge density of σ = 0.204 C/m 2 .The meaning of the symbol and curves as in figure1.DFT results from reference[11] and MC data from reference[20].

Figure 6 .
Figure 6.Charge inversion in a spherical double layer with a 2:2 RPM electrolyte at c = 1.25 mol/dm 3 and a = 4 • 10 −10 m.The macroion radius is R = 10 • 10 −10 m and contains 20 elementary charges.The meaning of the symbol and curves as in figure 1.DFT results from reference[11] and MC data from reference[31].

Figure 8 .
Figure 8. Zeta potential in a spherical double layer for a 1:1 RPM electrolyte at different concentrations with small ion diameter (a) a = 3 • 10 −10 m, and (b) a = 2 • 10 −10 m.The macroion radius is R = 15 • 10 −10 m.The meaning of the curves as in figure 1.

Figure 6
Figure 6 compares the theoretical predictions with the MC results of Terao and Nakayama [31] for a 2:2 electrolyte with R = 10 • 10 −10 m, a = 4 • 10 −10 m and macroion charge Z|e| = −20|e|.The damped oscillatory behaviour of the charge density profile is clearly seen with the MPB and DFT giving excellent agreement with the simulations, while the PB is qualitatively different, showing a monotonic behaviour.The graphs in figure 7 give P (r) for 1:1 and 2:1 electrolytes respectively at 3 different concentrations with Z = −18 and R = 15 • 10 −10 m.

Figure 10 .
Figure 10.Zeta potential in a spherical double layer for a 2:2 RPM electrolyte at different concentrations with small ion diameter (a) a = 3 • 10 −10 m, and (b) a = 2 • 10 −10 m.The macroion radius is R = 15 • 10 −10 m.The meaning of the curves as in figure 1.