New Research In
Physical Sciences
Social Sciences
Featured Portals
Articles by Topic
Biological Sciences
Featured Portals
Articles by Topic
 Agricultural Sciences
 Anthropology
 Applied Biological Sciences
 Biochemistry
 Biophysics and Computational Biology
 Cell Biology
 Developmental Biology
 Ecology
 Environmental Sciences
 Evolution
 Genetics
 Immunology and Inflammation
 Medical Sciences
 Microbiology
 Neuroscience
 Pharmacology
 Physiology
 Plant Biology
 Population Biology
 Psychological and Cognitive Sciences
 Sustainability Science
 Systems Biology
Local molecular field theory for effective attractions between like charged objects in systems with strong Coulomb interactions

Edited by Benjamin Widom, Cornell University, Ithaca, NY, and approved March 27, 2006 (received for review January 11, 2006)
Abstract
Strong, shortranged positional correlations involving counterions can induce a net attractive force between negatively charged strands of DNA and lead to the formation of ion pairs in dilute ionic solutions. However, the long range of the Coulomb interactions impedes the development of a simple local picture. We address this general problem by mapping the properties of a nonuniform system with Coulomb interactions onto those of a simpler system with shortranged intermolecular interactions in an effective external field that accounts for the averaged effects of appropriately chosen longranged and slowly varying components of the Coulomb interactions. The remaining shortranged components combine with the other molecular core interactions and strongly affect pair correlations in dense or strongly coupled systems. We show that pair correlation functions in the effective shortranged system closely resemble those in the uniform primitive model of ionic solutions and illustrate the formation of ion pairs and clusters at low densities. The theory accurately describes detailed features of the effective attraction between two equally charged walls at strong coupling and intermediate separations of the walls. Analytical results for the minimal coupling strength needed to get any attraction and for the separation at which the attractive force is a maximum are presented.
Strong Coulomb interactions in crowded, nonuniform environments have important experimental consequences in a wide variety of biophysical applications ranging from DNA packaging in viruses to transport in ion channels (1⇓⇓–4). These interactions present major challenges to theory and computer simulations not only because of their characteristic long range but also because they can be very strong at short distances. Here, we present a local molecular field (LMF) theory (5) that averages over particular longranged and slowly varying components of the Coulomb interactions (6) while still maintaining an accurate description of the shortranged components. Our model provides a general and physically suggestive theory for strongly coupled Coulomb systems and reduces exactly to the classical Poisson–Boltzmann (PB) approximation for dilute, weakly coupled systems.
We consider a general starting point where a molecule of species i, described by a rigid body frame with center at r_{i}, interacts with an external field, φ_{fi}(r_{i}), that comes from fixed charged solutes, or walls, or particular fixed molecules of a mobile species, as illustrated below. The subscript f indicates the source of the field, which we treat as a special fixed species f. The interaction between a pair of molecules of species i and j is assumed to have the general form w_{ij}(r_{ij}) = w_{s}_{,ij}(r_{ij}) + w_{q}_{,ij}(r_{ij}), where r_{ij} ≡ r_{j} − r_{i}. The w_{s}_{,ij}(r_{ij}) denote general (repulsive core and other), shortranged intermolecular interactions. There are angular coordinates expressing orientations of the body frames that we do not denote explicitly. The w_{q}_{,ij}(r_{ij}) arises from Coulomb interactions between rigid charge distributions q_{i}(r − r_{i}) in the body frame of each molecule, so that where the caret denotes a Fourier transform and we assume that there is a uniform dielectric constant ε everywhere.
To generate uniformly slowly varying components u_{1,ij} of the full w_{q}_{,ij} ≡ w_{q}_{0,ij} + u_{1,ij} that are well suited for LMF averaging, we limit the magnitude of wave vectors making significant contributions to the integration in Eq. 2. To that end, we introduce a Gaussian function parameterized by an important length scale σ that provides a smooth cutoff in kspace and write
The first term on the right has all of the characteristic longranged Coulomb divergences as k → 0 but decays very rapidly to zero for kσ ≳ 2. The desired slowly varying components arise when only this term is used in Eq. 2 with an appropriate choice of σ. For localized charge distributions q̂_{i}(k), we expand in a Taylor series about k = 0 and take the lowest order multipole moment (7). This simplified expression defines the u_{1,ij} we consider and thus prescribes a σdependent separation of the full intermolecular potentials w_{ij}(r_{ij}) = w_{s}_{,ij}(r_{ij}) + w_{q}_{,ij}(r_{ij}) = w_{s}_{,ij}(r_{ij}) + w_{q0,ij}(r_{ij}) + u_{1,ij}(r_{ij}) ≡ u_{0,ij}(r_{ij}) + u_{1,ij}(r_{ij}) into short and longranged parts.
In rspace, Eq. 3 becomes 1/r = erf(r/σ)/r +erfc(r/σ)/r. Here, erf and erfc are the usual error and complementary error functions. The erf(r/σ)/r term is the electrostatic potential from a normalized Gaussian charge distribution with width σ. As shown in Fig. 1, this electrostatic potential remains smooth and slowly varying on the scale of σ while decaying as 1/r at large r. This use of a Gaussian charge distribution is related to the Ewald sum method, which considers periodic images of ion configurations with embedded screening and compensating Gaussian charge distributions. However, our focus is on the separation of the potential itself and not the effects of periodic boundary conditions, and our choice of σ is usually much smaller than that used in Ewald sum methods, which typically is proportional to the width of the simulation cell (8).
The shortranged components u_{0,ij}(r_{ij}) define the intermolecular interactions in the special shortranged “mimic system.” The mimic interactions are composed of the shortranged parts of the Coulomb interactions w_{q}_{0,ij} ≡ w_{q}_{,ij} − u_{1,ij} and the other shortranged core interactions w_{s}_{,ij}:
As suggested by Fig. 1, σ sets the range of w_{q}_{0,ij} and thus determines an effective Coulomb core size (6). The external potential from fixed charged solutes or walls φ_{fi}(r_{i}) ≡ φ_{0,fi}(r_{i}) + φ_{1,fi}(r_{i}) can similarly be separated into short and longranged parts, as illustrated below.
It is straightforward to arrive at explicit results for u_{1,ij}(r) and φ_{1,fi}(r_{i}). Here, we consider localized charge distributions with a net charge q̄_{i} ≡ ∫drq_{i}(r) or a net dipole moment p_{i} ≡ ∫drrq_{i}(r). If both molecules carry a net charge, we find u_{1,ij}(r) = q̄_{i}q̄_{j}erf(r/σ)/εr. The associated Coulomb core component is w_{q}_{0,ij} = q̄_{i}q̄_{j}erfc(r/σ)/εr. Thus, the results of Fig. 1, scaled by q̄_{i}q̄_{j}/ε, give examples of possible u_{1,ij} and w_{q}_{0,ij} for ionic solution models (6). For a monopole and a dipole, we find u_{1,ij}(r) = q̄_{i}(p_{j}·∇)[erf(r/σ)/εr], and for dipoles we find u_{1,ij}(r) = −(p_{i}·∇)(p_{j}·∇)[erf(r/σ)/εr]. The latter will lead to dipolar mimic systems with shortranged angulardependent interactions.
Local Molecular Field Approximation
LMF theory introduces renormalized effective fields φ_{R}_{,fi}(r_{i}) that induce nonuniform singlet densities in the mimic system (denoted by the subscript R) that are supposed to equal those induced by φ_{fi}(r_{i}) in the full system of interest: This equation defines a general mapping relating structure in the mimic and full systems. Thermodynamic properties can be determined by integration of these structural relations.
We represent the effective field φ_{R}_{,fi}(r_{i}) ≡ φ_{0,fi}(r_{i}) + φ_{R}_{1,fi}(r_{i}) as the sum of the known shortranged part φ_{0,fi}(r_{i}) of the external field in the full system and a renormalized “perturbation component” φ_{R}_{1,fi}(r_{i}) that accounts for the averaged effects of the slowly varying interactions, u_{1,ij}. As discussed in detail in refs. 5 and 6, by considering the balance of forces in the full and mimic systems when Eq. 5 holds and by making some physically motivated approximations, we find that φ_{R}_{1,fi} is determined up to a constant by the local molecular field equations: Here, the prime on the integral indicates an implicit summation over all species j and an integration over the angles of the body frames. Longranged interactions from the fixed species f are accounted for by the δ(f, j) term, which denotes products of δfunctions describing the fixed location and orientation of f.
Note that the average over the slowly varying u_{1,ji}(r_{ji}) in Eq. 6 is weighted by ρ_{R}_{,fj}(r_{j}), the singlet density for species j (in the effective field of fixed species f but with no explicit reference to species i). This neglect of pair correlations between molecules at r_{j} and r_{i} characterizes a mean field approximation (2), which would in most contexts represent a major source of error. However, a general feature illustrated by Fig. 1 is that, as σ increases, u_{1,ji} becomes progressively more slowly varying at short distances. Thus, we can ensure that all of the u_{1,ji} will slowly vary over the length scales of relevant local pair correlations in the system by choosing σ larger than some statedependent minimum value σ_{min}. This choice permits a consistent and controlled use of the mean field approximation in computing the average, and we anticipate very accurate results from the LMF theory for any σ ≥ σ_{min} if we properly describe the resulting density in the mimic system (6).
At strong coupling, we argue that σ_{min} should be on the order of a characteristic nearestneighbor spacing ā. The strong, shortranged parts of the Coulomb interactions on the scale of ā and below directly affect pair correlations between nearestneighbor molecules. These correlations will be consistently described in the mimic system if we choose σ = σ_{min} on the order of ā so that there are essentially nearestneighbor interactions between the effective Coulomb cores, and the averaged effects of the u_{1,ji} from the furthest neighbors are slowly varying on this scale. We illustrate below the utility of these ideas for simple models of ionic solutions at strong coupling.
SizeAsymmetric Primitive Model
A model of great current interest is the sizeasymmetric primitive model (SAPM) of ionic solutions, a fluid of oppositely charged hard spheres of different sizes in a uniform dielectric continuum. The different hard sphere diameters crudely account for the different core sizes of real cations and anions, and there is an interesting and not well understood dependence of the critical temperature and critical density in this model as the size or charge ratio is varied (9). We consider in particular the uniform equimolar system studied with Monte Carlo (MC) simulations by Weis and Levesque (WL) (10), with symmetric unit charges e_{0} = q̄_{1} = −q̄_{2} and a diameter ratio d_{1}/d_{2} = 0.4. Thus, w_{s}_{,ij}(r) = ∞ for r ≤ d_{ij} ≡ (d_{i} + d_{j})/2 and is zero otherwise, and w_{q}_{,ij}(r) = q̄_{i}q̄_{j}/εr.
WL characterize the states by two dimensionless parameters: a reduced density ρ* ≡ (N_{1} + N_{2})d_{2}^{3}V and an effective coupling strength Γ* ≡ l_{B}/d_{2} (called q*^{2} in the notation of WL). Here, l_{B} ≡ βe_{0}^{2}/ε is the Bjerrum length, the distance at which the interaction energy between two unit charges e_{0} equals k_{B}T and β ≡ (k_{B}T)^{−1}. The simulations of WL indicate that the critical point occurs at ρ* = 0.195 and Γ* = 15.15. We report results here for two strong coupling states with very different properties: a highdensity subcritical liquid state with ρ* = 1.4 and Γ* = 16 and a lowdensity supercritical vapor state with ρ* = 0.04 and Γ* = 9.
The competition between the Coulomb interactions and the packing arrangements of the embedded hard cores in the SAPM produces elaborate local structures in these strong coupling states as exhibited in the pair correlation functions g_{ij}(r), proportional to the density response to an external field φ_{ij}(r) = w_{ij}(r) arising from a fixed ion of type i at the origin. Thus, in LMF theory, even uniform fluid correlation functions are described from a nonuniform point of view. These characteristic features can be very accurately reproduced in the mimic system by using the strong coupling approximation (SCA). The SCA replaces φ_{R}_{,ij}(r) by the known strong shortranged component φ_{0,ij}(r) = u_{0,ij}(r) of the field from fixed ion i. This field corresponds to fixing a mimic particle at the origin, or equivalently, approximating the g_{ij}(r) in the uniform ionic system by the g_{0,ij}(r) in the uniform mimic system (6).
In Fig. 2, we compare correlation functions determined by WL for the highdensity state, with ρ* = 1.4 and Γ* = 16, to MC simulations we carried out in the uniform mimic system by using the SCA, with a “molecularsized” choice of σ = 1.2d_{2}. Simulations of the longranged system required careful and costly treatment of periodic boundary conditions using the Ewald sum method, which was not needed for the shortranged mimic system. Despite the very different range and magnitude of the mimic system interactions, all of the pair correlation functions are strikingly similar to those of the full SAPM. These functions are very different from the profiles of the associated hard sphere mixture with the charges set equal to zero, indicating the crucial importance of including the strong shortranged parts of the Coulomb interactions w_{q}_{0,ij} in defining the mimic interactions in Eq. 4. Equally good results are found for larger values of σ, as illustrated in Fig. 2, but the good agreement fails for much smaller σ, indicating that σ_{min} is ≈1.2d_{2} for this state.
Qualitatively different structures are seen in the lowdensity vapor state, with ρ* = 0.04 and Γ* = 9, as illustrated in Fig. 3. The simulations of WL show that oppositely charged ions pair together with a typical spacing close to the minimum permitted by the hard core diameters, along with some transient formation of longer chainlike structures. The correlation functions between likecharged pairs exhibit pronounced peaks of essentially the same magnitude at a separation of r = 1.4d_{2} = d_{1} + d_{2}. These results indicate the existence of small clusters of ion pairs, with the same peak position and amplitude for both “+−+” or “−+−” configurations at the minimum distance permitted by a linear arrangement of the embedded hard cores. These peaks also illustrate how counterions can induce an effective attraction between like charged objects, as discussed in detail in Charged Walls with Point Counterions. Very good agreement between full and mimic system correlation functions is achieved with a choice of σ_{min} = 3.0d_{2}, consistent with the larger average spacing between dilute ion pairs.
The clustering of the ions has probably presented the most severe challenges to theories of ionic systems. It is particularly crucial for the study of critical phenomena and vapor–liquid coexistence (9). The PB approximation and the frequently used hypernetted chain integral equation fail to predict ion clustering; indeed, the hypernetted chain equation has no solution in most of the ion pairing regime (10). In contrast, the mimic system as described by the simple SCA already builds in the most important local features of ion aggregation. This very good agreement strongly suggests that the LMF theory can accurately represent the Coulomb cores that contribute to local correlation functions in more realistic models of ionic systems. Any remaining errors can be attributed mainly to deficiencies in the description of the other shortranged core interactions, thus permitting the efficient development of more accurate models.
Charged Walls with Point Counterions
The suspension and selfassembly of highly charged polyelectrolytes (macroions) in the presence of mobile counterions (microions) is of great interest in biological systems (1). These systems usually involve charge and size asymmetries much greater than that of the SAPM and are often studied by fixing a certain macroion configuration and computing the microion distribution and resulting forces on the macroions. We discuss here the simplest such model system (11) consisting of uniformly charged infinite hard walls with neutralizing point counterions (and no coions) in a uniform dielectric environment. This model is simple enough that exact results in certain limits are known (12), but it still illustrates many fundamental issues that arise from the interplay between long and shortranged forces in an explicitly nonuniform environment. It is clear from SizeAsymmetric Primitive Model that the LMF theory can deal with more realistic models for the walls and counterions.
One Charged Wall.
We first consider the case of a single hard wall with a uniform negative charge density q_{w} at the z = 0 plane, where we take the zero of electric potential energy. Without loss of generality we can assume that the counterions have a (unit) charge e_{0} and express the wall charge density q_{w} ≡ −e_{0}/l_{w}^{2} in terms of the length l_{w} of the side of a square enclosing that amount of charge. The potential energy φ^{1w}(z) of a counterion at a distance z from the wall is 2πe_{0}^{2}z/(εl_{w}^{2}). The Gouy–Chapman length l_{G} is defined as the distance at which this potential equals k_{B}T, i.e., l_{G} ≡ k_{B}Tεl_{w}^{2}/(2πe_{0}^{2}) = l_{w}^{2}/(2πl_{B}). l_{G} characterizes the effective strength of the attractive wall–counterion interaction, and most counterions will be found near the wall in an effective slit whose width is proportional to l_{G}. Dimensionless combinations of thermodynamic variables in this simple system depend only on a single control parameter ξ ≡ l_{B}/l_{G} = l_{w}^{2}/(2πl_{G}^{2}) (13).
As ξ increases (e.g., by decreasing T at a fixed wall and counterion charge), counterions are driven increasingly close to the wall by the decreasing l_{G}. At strong coupling with ξ ≫ 1 or l_{G} ≪ l_{w}, most counterions are next to the wall and form a (“strongly correlated”) two dimensional (2D) liquid layer (14, 15) with average lateral spacing ā ≃ a_{2D} ≡ l_{w} fixed by local neutrality. There are indeed strong lateral correlations between the counterions in the 2D layer: the coupling strength Γ_{a2D} ≡ l_{B}/a_{2D} = ξ^{1/2}/(2π)^{1/2} ≫ 1 for large ξ. As discussed above, we then expect the effective Coulomb core size σ_{min} to be on the order of a_{2D} = l_{w}. However, because of these repulsive cores, particles cannot stack perpendicular to the wall and still remain near the narrow slit. Thus, the density outside the slit is very low, and there are only weak correlations normal to the wall.
In the opposite weak coupling limit, with ξ ≪ 1 or l_{G} ≫ l_{w}, counterions can take advantage of the larger effective volume of the slit and adopt a more diffuse 3D packing to reduce their repulsive interactions. Crudely assuming all counterions are found within l_{G} of the wall and using a simple cubic lattice to estimate the characteristic counterion spacing in this volume, we now have ā ≃ a_{3D} ≡ (l_{w}^{2}l_{G})^{1/3} = l_{w}/(2πξ)^{1/6}. There is weak coupling between the counterions, with Γ_{a3D} ≡ l_{B}/a_{3D} = ξ^{2/3}/(2π)^{1/2} ≪ 1, and here it is natural to take σ_{min} ≃ l_{B} = l_{w}(ξ/2π)^{1/2} as an estimate for the effective Coulomb core size (6). The crossover to strong coupling with essentially 2D packing and σ_{min} ≃ a_{2D} = l_{w} occurs for ξ on the order of unity, and the 2D packing indeed provides a larger average spacing at large ξ.
Quantitative results take an especially simple form (13) if we introduce a dimensionless rescaled density n(z/l_{G}) ≡ l_{G}l_{w}^{2}ρ(z) that incorporates the anticipated (ξdependent) scaling of the profile with l_{G}. Local neutrality requires that ∫_{0}^{∞}dz̃n(z̃) = 1, where z̃ ≡ z/l_{G}. Lengths scaled by l_{G} will generally be indicated by a tilde. Moreover, because of the impulsive δfunction force at a hard wall, there is an exact relation between the pressure and the contact density. This relation yields the well known contact theorem, which implies n(0) = 1 for the contact value of the rescaled density at a single charged hard wall (11).
Exact results (12) for n(z̃) are known in the limit ξ → 0 from a rigorous weak coupling expansion, which gives results agreeing with the PB approximation, n_{PB}(z̃) = 1/(z̃ + 1)^{2}. A different strong coupling expansion gives exact results as ξ → ∞: n_{SC}(z̃) = e^{−z̃}. However, attempts to connect these limits by analyzing higherorder terms in each expansion have had only limited success (13). We now show that the LMF theory provides a simple, accurate, and unified approach for general ξ.
LMF Equation for One Charged Wall.
We can take advantage of planar symmetry and integrate exactly over the lateral degrees of freedom in the longranged parts u_{1,ji}(r_{ji}) of the counterion–counterion and wall–counterion interactions in Eq. 6. The resulting LMF equation can be written in dimensionless form for z̃_{1}, z̃_{2} ≥ 0 as Here, φ̃_{R}_{1}(z̃_{1}) ≡ βφ_{R}_{1}(z̃_{1}l_{G}), and G(z̃_{2}, z̃_{1}) ≡ −z̃_{1} − z̃_{2}erf(z̃_{1} − z̃_{2}/σ̃) − π^{−1/2}σ̃e^{−}[(z̃_{1} − z̃_{2})/σ̃]^{2} + z̃_{2}erf(z̃_{2}/σ̃) + π^{−¿}σ̃e^{−}[z̃_{2}/σ̃]^{2} is the Green’s function associated with a normalized planar Gaussian charge distribution centered at z̃_{2}, with the zero of potential energy at z̃_{1} = 0.
The −δ(z̃_{2}) term in Eq. 7 accounts for the longranged component φ̃_{1}(z̃_{1}) = −G(0, z̃_{1}) of the dimensionless attractive potential φ̃^{1w}(z̃_{1}) = z̃_{1} between a counterion at z̃_{1} and the negatively charged wall at z̃_{2} = 0. The remaining shortranged part [φ̃_{0}(z̃_{1}) = z̃_{1} − φ̃_{1}(z̃_{1})] of the wall potential is The effective field is then φ̃_{R}(z̃_{1}) = σ̃_{0}(z̃_{1}) + φ̃_{R}_{1}(z̃_{1}), with φ̃_{R}_{1} given by Eq. 7 for z̃_{1} ≥ 0 and infinity otherwise.
mPB Approximation.
To solve Eq. 7 selfconsistently, we must accurately determine the density n_{R}(z̃) induced by φ̃_{R}(z̃). At weak coupling, neighboring ions interact weakly and the density response to the effective field is proportional to the ideal gas Boltzmann factor exp[−φ̃_{R}(z̃)]. By using this approximation in Eq. 7, we have a closed equation, which we refer to as the mimic PB (mPB) equation. Moreover, we can show (7) for all ξ that a selfconsistent solution of the mPB equation will exactly satisfy both neutrality and the contact theorem. The contact theorem implies that the density response takes the simple form n_{R}(z̃) = exp[−φ̃_{R}(z̃)] with our choice of the zero of energy.
Remarkably, however, the mPB approximation also gives accurate results at strong coupling with ξ ≫ 1, where there is an essentially 2D arrangement of the mimic particles in an effective narrow slit. Correlations normal to the wall are very weak, and the Boltzmann factor again can accurately describe the density response to the zdependent field φ̃_{R}(z̃), as can be verified by more formal arguments (13, 16).
These limits motivate our use of the mPB approximation with n_{R}(z̃) = exp[−φ̃_{R}(z̃)] for all ξ in Eq. 7. The mPB approximation is least justified at intermediate values of ξ, and it breaks down if σ is chosen much larger than σ_{min} so that there would be strong direct interactions between further neighbors in the mimic system. Thus, we provide a smooth interpolation between the known limiting values of σ_{min} in the weak and strong coupling regimes by choosing σ = σ^{1w} ≡ C min(l_{B}, l_{w}), and fix C = 0.60 by finding the best fit to simulations (13) at a moderately strong coupling state with ξ = 40. In this example, it is numerically more convenient to differentiate the resulting mPB equation and solve for the effective force, which vanishes far from the wall, and then get the effective field by integration (J. Rodgers, C. Kaur, Y.G.C., and J.D.W., unpublished data). An iterative solution is straightforward, and no other simulation data are required.
Results for One Charged Wall.
Fig. 4 gives results for n_{R}(z̃) at strong coupling, with ξ = 100. There is excellent agreement between the results of the mPB theory and MC simulations of the longranged system carried out by Moreira and Netz (13). The log–log plot emphasizes that n_{R}(z̃) has two characteristic regions. Near the wall, there is an initial exponential decay arising mainly from particles in the 2D layer, which continues until approximately z̃ ≃ σ̃_{min}/2, where the density is very low and there is a crossover to algebraic decay as in the PB solution but with a much larger effective l_{G}. A natural physical interpretation is that the small fraction of counterions outside the 2D layer adopt a diffuse PBlike profile generated by an effective wall whose charge density has been greatly reduced by the charge of the tightly bound counterions. This idea has been suggested before (15), but LMF theory provides a unified description of both limiting regions and the crossover region.
The density near the wall is very accurately described by the even simpler SCA, where φ̃_{R}(z̃) is approximated by φ̃_{0}(z̃). The resulting density n_{0}(z̃) ≡ exp[−φ̃_{0}(z̃)] can be written down immediately from Eq. 8. As shown in Fig. 4, both φ̃_{0}(z̃)] and φ̃_{R}(z̃) closely resemble the full potential z̃_{1} near the wall for z̃_{1} ≲ σ̃_{min}/2. But n_{0}(z̃) cannot describe the PBlike region at large z̃ as does the full mPB theory, and it does not obey the neutrality condition. This example highlights both the strengths and weaknesses of the SCA. When properly used to describe only shortranged correlations at strong coupling, very accurate results can be found.
Two Charged Walls.
We now briefly consider two equally charged hard walls forming a real slit with width d, with neutralizing point counterions in between. At strong coupling and intermediate widths, the counterions can induce an effective attractive force between the walls. Such effective attraction between like charged objects may have important experimental consequences, and it has generated a great deal of theoretical interest (1⇓⇓–4).
Ref. 11 gives an exact expression for the dimensionless (osmotic) pressure P̃ ≡ βl_{G}l_{w}^{2}P arising from neutralizing point counterions confined between charged hard walls at z = 0 and z = d: Thus, if the rescaled contact density n(0) is less than (greater than) 1, there is an effective attractive (repulsive) force on the walls. As d → ∞, we recover the onewall results discussed earlier, where P = 0 and n(0) = 1.
Because the total force on a counterion from equally charged walls exactly cancels for all z and all d, we now have φ̃^{2w}(z̃) = 0 for 0 ≤ z̃ ≤ d̃. As in the onewall case, it is useful to divide φ̃^{2w} into a shortranged part, given by a sum of shortranged, singlewall terms (indicated by the superscript 1w) defined in Eq. 8, and the remainder. We take the zero of energy on the left wall at z̃ = 0. The effective field φ̃_{R}(z̃) = φ̃_{0}^{2w}(z̃; d̃) + φ̃_{R}_{1}(z̃) is determined from the twowall LMF equation. This equation closely resembles Eq. 7, except that the integration is from 0 to d̃ and there is an additional −δ(d̃ − z̃_{2}) term in the integrand, accounting for interactions with the second wall at z = d. Again we use the mPB approximation n_{R}(z̃) = A exp[−φ̃_{R}(z̃)] and fix the constant A (which equals the contact density with our choice of the zero of energy) by using the neutrality condition ∫_{0}^{d̃} dz̃n_{R}(z̃) = 2. The pressure is then given by Eq. 9.
The resulting twowall mPB equation reduces exactly to an integrated form of the PB equation if σ = 0. The latter has an analytic solution and predicts a repulsive force for all d and ξ (2, 13). The mPB theory also predicts a weak repulsive force at strong coupling for sufficiently large d, arising from weak repulsions between counterions in the dilute PBlike tails of the onewall profiles discussed above. In contrast, at strong coupling and sufficiently small d, core repulsions make it unfavorable for particles in the narrow slit to stack perpendicular to the walls, and counterions will be forced into a single 2D layer with characteristic lateral spacing ã ≃ l_{w}/ fixed by neutrality. To interpolate between this limit and weak coupling, we choose σ = σ^{2w} ≡ C min(l_{B}, l_{w}/) and take the same value of C = 0.60 that we used for the onewall theory. Because σ̃_{min} = Cξ at small ξ, the PB approximation is consistent only as ξ → 0. The mPB theory naturally introduces a crucial new length scale σ_{min} that allows for a change in the functional form of n_{R}(z̃) as ξ increases.
Fig. 5 compares numerical results of mPB theory (J. Rodgers, C. Kaur, Y.G.C., and J.D.W., unpublished data) to simulation data (13) for strong coupling states. The left graph shows that, for ξ = 100, there is very good agreement between the mPB theory and computer simulations for all widths at which simulations can be performed. As shown in Fig. 5 Left Inset, the mPB theory predicts that, at still larger widths, there is a weak repulsive force between the walls. This reentrant behavior is illustrated more generally in Fig. 5 Right, which also shows that a minimal coupling strength of ξ ≥ ξ_{c} ≃ 12 is needed to get any attraction (13).
Fig. 5 Left also shows results from the analytic SCA, where n(z̃) is approximated by n_{0}(z̃) ≡ A_{0} exp[−φ̃_{0}^{2w}(z̃; d̃)], with A_{0} similarly determined by neutrality. At strong coupling and large d, φ̃_{0}^{2w}(z̃; d̃) has a deep attractive well near each wall essentially identical to that found near a single wall. However, at small enough d, the wells from the individual φ̃_{01}^{w} terms in Eq. 10 begin to overlap, and their depth decreases. As d → 0, φ̃_{0}^{2w}(z̃; d̃) vanishes for all z̃ in the slit. The SCA fails at large widths, just as it did far from the wall in the onewall case, and the reentrant behavior is completely missed. However, the SCA is very accurate at smaller widths, and it gives a good description of the location and magnitude of the maximum attractive force.
The formation of the single 2D layer at sufficiently small widths d̃ < σ̃^{2w} plays a key role in producing a strong attractive force. Because of the absence of correlations normal to the walls, the density profile will be relatively constant across the slit. We can use the SCA to describe several features analytically in this regime. The minimum pressure P̃* should occur near the largest width d̃* ≡ 2z̃*_{m} for which the single 2D layer remains stable, defined by σ̃_{0}^{2w}(z̃*_{m}; d̃*) = 1. At larger widths, there will be higher contact densities as separate 2D layers at each wall begin to form, and the nearly constant profiles at smaller widths will have higher contact densities by normalization. Expanding φ̃_{0}^{2w}(z̃; d̃) in a Taylor series, we have φ̃_{0}^{2w}(z̃; d̃) = 2z̃(d̃ − z̃)/(π^{1/2}σ̃^{2w}) to lowest order; higher order terms are negligible for d̃ ≪ σ̃^{2w}. This expansion implies d̃* = (2π^{1/2}σ̃^{2w})^{1/2} ≃ 1.94ξ^{1/4} using our expression for σ̃^{2w}. Similarly evaluating A_{0}, the minimum pressure is P̃* ≃ 3.717/d̃* − 1. When P̃* = 0, there can be no attractive forces, which determines the minimal coupling strength ξ_{c} ≃ 13.48 and associated critical spacing d̃*_{c} ≃ 3.717 needed to get any attraction. Finally, for ξ ≫ ξ_{c} and d̃ = 2, we see φ̃_{0}^{2w}(z̃; d̃) ≪ 1; therefore, the constant profile n_{0}(z̃) = 2/d̃ is very accurate. Eq. 9 then implies that, at very strong coupling, the transition from strong repulsive forces to strong attractive forces occurs near d̃ = 2, as shown in refs. 12 and 13. All these predictions are in very good agreement with numerical solutions of the mPB theory (and with simulations of the full and mimic systems, where available) for all strong coupling states, as illustrated in Fig. 5.
Discussion
In strong coupling regimes the shortranged parts of the Coulomb interactions efficiently compete with the other shortranged molecular core interactions and strongly influence pair correlations between neighboring molecules. In LMF theory, the choice of σ_{min} determines the strength and range of these important Coulomb core interactions. These core interactions play a key role in inducing an effective attraction between likecharged objects at strong coupling, as illustrated in Fig. 3 for the correlation functions between likecharged ions in the SAPM and in Fig. 5 for the osmotic pressure on two likecharged walls. In both cases, strong, shortranged forces mainly involving single counterions or a single counterion layer can mediate a strong effective attraction. Such phenomena have traditionally been interpreted as illustrating the “breakdown of mean field theory” and the need for new and more sophisticated approaches. But LMF theory using the simple SCA, in which only the shortranged parts of the Coulomb interactions are taken into account, provides very accurate results at small and moderate separations. At weak coupling and for longwavelength correlations, the averaged effects of the longranged interactions as determined by the LMF equation are needed as well (6).
LMF theory provides a general conceptual framework that connects and clarifies previous work on systems with both short and longranged interactions. It has been used to describe liquid–vapor interfaces and wetting and drying transitions for simple fluids (5) and hydrophobic interactions in water (17). LMF theory suggests simplified simulation models for general Coulomb systems based on the mimic system that do not require special treatment of periodic boundary conditions. The shortranged intermolecular interactions in the Coulomb mimic system are reminiscent of the truncated interactions used in reaction field methods (18), but it was not very clear in those approaches how to choose an appropriate cutoff and how to treat nonuniform environments. The determination of σ_{min} and the effective field in LMF theory provides a way to deal with both problems. We believe that the LMF picture will prove useful not only in formal theory but also for qualitative reasoning and in detailed simulations of biophysical systems.
Acknowledgments
We thank Michael Fisher, Charanbir Kaur, and Jocelyn Rodgers for many helpful remarks and J. J. Weis for sending us simulation data for the SAPM correlation functions. This work was supported by National Science Foundation Grants CHE0111104 and CHE0517818.
Footnotes
 ↵^{‖}To whom correspondence should be addressed at: Institute for Physical Science and Technology, University of Maryland, College Park, MD 20742. Email: jdw{at}ipst.umd.edu

Author contributions: Y.G.C. and J.D.W. designed research, performed research, and wrote the paper.

Conflict of interest statement: No conflicts declared.

This paper was submitted directly (Track II) to the PNAS office.
 Abbreviations:
 LMF,
 local molecular field;
 PB,
 Poisson–Boltzmann;
 SAPM,
 sizeasymmetric primitive model;
 WL,
 Weis and Levesque;
 MC,
 Monte Carlo;
 SCA,
 strong coupling approximation;
 mPB,
 mimic PB.
 Received January 11, 2006.
 © 2006 by The National Academy of Sciences of the USA
References
 ↵
 Gelbart W. M.,
 Bruinsma R. F.,
 Pincus P. A.,
 Parsegian V. A.
 ↵
 ↵
 ↵
 Boroudjerdi H.,
 Kim Y.W.,
 Naji A.,
 Netz R. R.,
 Schlagberger X.,
 Serr A.
 ↵
 ↵
 Chen Y.G.,
 Kaur C.,
 Weeks J. D.
 ↵
 Chen Y.G.
 ↵
 Frenkel D.,
 Smit B.
 ↵
 ↵
 Weis J. J.,
 Levesque D.
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵
 ↵