Table Of Content

Computing Counterion Densities at Intermediate Coupling Christian D. Santangelo1,2,∗ 1Department of Physics and Astronomy, University of Pennsylvania, Philadelphia, 19104 2Department of Physics, University of California, Santa Barbara 93106 (Dated: February 2, 2008) By decomposing the Coulomb interaction into a long distance component appropriate for mean- 5 fieldtheory,andanonmean-fieldshortdistancecomponent,wecomputethecounteriondensitynear 0 a charged surface for all values of the counterion coupling parameter. A modified strong-coupling 0 expansion that is manifestly finite at all coupling strengths is used to treat the short distance 2 component. Wefindanonperturbativecorrectionrelatedtothelateralcounterioncorrelationsthat n modifies thedensity at intermediate coupling. a J PACSnumbers: 82.70.-y,82.45.-h,87.15.Kg 0 2 The rise of biological physics has rekindled the long- exponentialforminthe strongcouplinglimitfailstoelu- ] standing interest in aqueous electrostatics [1]. Poisson- cidate the physics behind the corrections to that limit t f Boltzmann mean-field theory fails to describe a number in a clear and satisfactory manner [16]. Furthermore, I o ofstrikingphenomena,suchaschargeinversion[2,3]and find a nonperturbative correction to the density related s . counterion-mediated attraction [4, 5, 6, 7, 8, 9, 10], that to the lateral counterion correlations that becomes im- t a occur when strong correlations develop between multi- portant at intermediate coupling. Unlike in SC theory, m valent counterions. Although there has been some suc- I can unambiguously compute an expression for the free - cess understanding counterion correlations using both a energy and show that these lateral correlations play a d phenomenological Wigner crystal theory [3], and sys- role. n o tematic weak-coupling (WC) [11] and strong-coupling Considertheprimitivemodelforachargedsurfaceneu- c (SC) [12, 13, 14] expansions, a complete quantitative tralized by pointlike counterions of the opposite charge [ theory spanning the entire range of counterion behav- in a dielectric medium with dielectric constant, ǫ. To 1 ior is still lacking. This letter introduces a method to proceed, introduce a length scale, ℓ, and define Vs(r) = v computethe counteriondensityatintermediatecoupling l Q2e−r/ℓ/r and V (r)=l Q2(1 e−r/ℓ)/r, where l = B l B B − 8 by decomposing the Coulomb interaction into long and e2/(ǫk T)istheBjerrumlengthandQisthecounterion B 9 short distance components in the spirit of the Weeks- valence. The length ℓ is currently arbitrary and will be 4 Chandler-Andersentheory of simple fluids [15]. This de- chosen later to optimize the calculation. The Hamilto- 1 0 composition not only gives good quantitative agreement nian for ions of charge Qe centered at positions Rα in- 5 withsimulations,italsoprovidesanaturalframeworkto teracting with a surface of charge density nf(r)=σδ(z) 0 understand both the success of SC expansion, as well as is given by / t the role that lateral correlations play in the counterion 1 ma density. H= 2Z d3rd3r′ρ(r)[Vs(r−r′)+Vl(r−r′)]ρ(r′), (1) Here, I use mean-field theory for the long distance d- interaction, and introduce a modified SC expansion at where ρ(r)=nf(r)/Q− αδ(r−Rα). It is understood n short distances. The traditional SC expansion [12] is in this expression that ioPn self-interactions, which arise o problematicbecauseitisformallyavirialexpansion,and only from Vs(r), are to be neglected. The long range c one naively expects it to be invalid precisely when the interactioncanbedecoupledbyintroducingacontinuous : v counterions are strongly interacting. Nonetheless, nu- field φ through a Hubbard-Stratonovich transformation, Xi merical simulations have demonstrated that it not only resulting in the action S′ =Ss+Sl, where ar cstorrornecgt-lcyoupprleidnigctlsimtihte, bauvteraalgseoccooumnptuertieosnthdeenfosirtmy ionf tthhee Ss = 1 d3rd3r′ nf(r)Vs(r r′ )nf(r) 2Z Q | − | Q corrections [14]. In contrast, the modified SC expansion introducedhereismanifestlyfiniteinthelimitofinfinite d3r n (r)V (r R )/Q (2) f s α − Z | − | counterion coupling, but recovers the SC corrections at Xα large, finite coupling. + V (R R ), s α β | − | In addition, this decomposition correctly reproduces αX<β the counterion distribution around a charged surface at ℓ =4πl Q2, and both strong- and weak- coupling, and agrees very well B B withsimulationsatintermediatecoupling. Atest-charge 1 ρ(r)φ(r) S = d3r[( φ)2+ℓ2( 2φ)2]+i d3r . theory (TCT) which also computes approximatecounte- l 2ℓ Z ∇ ∇ Z Q B rion densities at intermediate coupling and explains the (3) 2 The counterion positions, R , are restricted to be over reproducedby expanding lnZ in powers of ρ (z). How- α s 0 the volume of space that canbe occupied by the counte- ever,forcounterionsin the presenceofa chargedsurface rions. with charge density σ, each term of this expansion di- In the Grand Canonical ensemble for the counterions, vergesas the coupling constantΓ=2πl σQ3 indi- B →∞ the partition function is catingthatthecounterioninteractionsarenotsmall;this divergence can be absorbed by shifting κ2 and utilizing Z 1 κ2 N φ N d3R e−S′, (4) overallcharge neutrality [12]. α ∝ N!(cid:18)ℓ (cid:19) Z D Instead, expand in powers of δρ (z) = ρ (z) XN B αY=1 2/(ℓ λ)δ(z), where λ = (2πσl Q)−01 is the 0Gouy−- B B where the length κ−1 is defined by κ2/ℓB = eβµ/a3, a Chapman length for a surface of charge density σ. This is the counterion radius, and µ is the chemical poten- has the property that dz δρ (z) = 0 due to overall 0 tial [11]. Now define the partial partition function with charge neutrality, and yRields integrals only over the counterion positions, n 1 ∞ N Zs = d2rαdzα δρ˜0(rα,zα) (8) 1 n!Z Z = d3R ρ (R ) Xn αY=1 Yα s α 0 α N!Z NX=0 αY=1(cid:2) (cid:3) exp Vs(rα rβ;zα zβ)Zp, × − − − ×exp− Vs(|Rα−Rβ|) (5)  αX<β   αX<β  whererα indicatesthepositionofacounterionprojected where ρ = (κ2/ℓ )eF+iφ, and F(r) = d3r′ V (r′ to the surface, zα its distance from the surface, and r)nf(r′)0/Q. ThisBleaves Z φe−S, wRhere s | − δρ˜0(rα,zα) = hexp[− iVs(rα − ri;zα)]ipδρ0(rα,zα). | ∝R D Here,h···ip is the averaPgetakenwith respectto the par- S = 1 d3r 1( φ)2+ ℓ2 2φ 2 lnZ [ρ (r)]. tition function s 0 ℓ Z (cid:20)2 ∇ 2 ∇ (cid:21)− B (cid:0) (cid:1) 1 2 m m (6) Zp = d2riexp Vs(ri rj;0). The average counterion density can be formally com- Xm m!(cid:18)ℓBλ(cid:19) Z iY=1 −Xi<j − puted as ρ(r) =δlnZ /δF(r).   s (9) h i Since this formulation is exact, the partition function Theexpression exp[ V (r r ;z )] representsthe is independent of the choice of ℓ. To proceed, make the interactionofahcharg−eaPticoosrdαin−ateis(αrα,ipzα)withalayer followingapproximations: (1)the mean-fieldapproxima- of counterions at positions (r ,z = 0). In deriving Eq. i i tion for the long-distance interaction (saddle point in (8), I have assumed exp[ V (r r ;z )] iφn)g, aanmdo(d2i)fieedxpSaCnd-litkheeceluffsetcetriveexppaontesniotnialdlenscZrsib[ρe0d] bues-- αhexp[− iVs(rα−hQri;αzα)]ip−,Pwhiichs isαt−rueiasαlonigp a≈s Qthe counterPions at z >0 are far enough apart compared low. Making these approximations, the theory will lose to ℓ. its independence on the choice of ℓ, and there will be a Performing a cluster expansion with respect to δρ˜ 0 “best”ℓwhosevaluegivestheclosestagreementwiththe yields lnZ [φ] = d3r δρ˜ (r)[1 + 1 d3r′ v (r full theory. In principle, its value should be determined r′)δρ˜0(r′)+s ], wheRre v2(r) =0 1 exp[2RVs(r)]. T2erm−s by optimizing the errorbetween a loopexpansiononthe ··· − − higher than zeroth order vanish in the limit Γ as long-distanceinteractionandperturbativecorrectionsto → ∞ the density becomes delta function-like and δρ˜ 0. As 0 theshort-distanceexpansion. Onphysicalgrounds,how- Γ 0, these corrections also vanish because ℓ3δ→ρ˜ 0, 0 ever,we arguethatℓ l Q2, or,inotherwords,ℓ is the → → B andthe interactionsbecome predominantlylong-ranged. ∼ distance that fixed counterions interact with an energy Itisusefultodefineρ˜ =eζ(z)−φ(z) withexp[ζ(r,z)]= 0 of kBT. Counterions at separations larger than this in- exp[F(z) V (r r ;0)] . The function ζ(z) can be teractweakly,andmean-fieldtheoryis likelyvalidabove h − i s − i i interpreted aPs the short-distance interaction potential of this length scale. In addition, the distinction between a chargeatheight z with the chargedsurface andwith a short and long length scales should not depend on the layer of counterions at z = 0. Therefore, it encodes the geometryofthe fixedchargedistribution, andcanthere- response of the z = 0 layer to the presence of a charge fore only depend on the Bjerrum length and counterion at some z > 0, and is reminiscent of theTCT [16]. One valence. For concreteness, I will choose ℓ=l Q2. B difference between ζ(z) and the TCT, however, is that The mean-field approximation, given by δS/δφ = 0, ζ(z) also depends, at least in principle, on the short- results in the equation range structure of the counterions induced by the short- 2φ ℓ2 4φ+ℓ ρ(r) =n (r)ℓ /Q, (7) range interaction. To develop a simple approximation B f B ∇ − ∇ h i for ζ(z) which I will use throughout the remainder of where I have performedthe Wick rotationφ iφ in the the letter, assume that each counterion at z > 0 in- → complexplaneforconvenience. TheSCexpansioncanbe teracts with a uniform distribution of charge at z = 0 3 containing an induced circular correlation hole of radius 0.15 ar0v=acapncℓyBλin/2a=lopcaQlly/σo.rdTerheids alaptptricoexiomfactoeusnttheerisoinzse aotf 0.15 (z) B 0.1 tζzh(≪ze)su=ℓr,f(ζaℓ(c/zeλ,))w[≈eh−iΓzc/(hℓ1−s−heo√eu−lrdr020+b/zℓe2)/v−ℓa]l.zidI/tλwishisienndtoeΓrmeissintliaantrgegdeth.baTyt,htufhoser, n(z) - n(z)PB 0.00.51 n(z) - nP-00..00550 1 2 3 4 5 6 7 interaction of the counterions with the bare surface and z / λ not the z = 0 layer of counterions whose contribution is 0 of order z2/(λℓ) z/λ. ∼ ≪ 1 2 3 4 5 6 7 Tolowestorder, ρ =ρ˜ , andthe mean-fieldequation 0 for a charged surfacheinow reads z / λ ∂2φ ℓ2∂4φ+κ2Θ(z)eζ(z)−φ(z) = ℓBσδ(z), (10) FIG. 1: Normalized density difference, n(z)−nPB(z), as a z − z Q function of z/λ for Γ = 100 (Γ = 10 in inset). Solutions to Eq. (10) (solid line) are compared to numerical simulations whereσ isthesurfacechargedensity. Forsmallz,ζ(z) fromRef.[14](diamonds)andtheTCTfromRef.[16](dashed − F(z) is analytic and can be expanded as a power series line). in (z/ℓ)2. On the other hand, F(z) is nonanalytic at z = 0 and contributes to the boundary conditions. This additional contribution can be disentangled by defining Φ=φ F. IntermsofΦ,ζ(z)isreplacedwithζ(z) F(z) can be neglected in this limit because κ is exponentially − − andthereisanadditionalsourcetermontherightofEq. suppressedbyζ(0)beinglarge. Therefore,thedensityat (10)oftheform ℓ σℓ2∂2δ(z)/Q. Thisequationencodes largedistancesisPB-likeandiscontrolledbyarenormal- two boundary co−ndBitionsz: ∂zΦ|00+−−ℓ2∂z3Φ|00+− = ℓBσ/Q, ized Gouy-Chapman length, λren =√2/κ=λexp[ζ(0)]. and ∂ Φ0+ = ℓ σ/Q. In terms of the original φ, the This is in agreement with the arguments of Burak et boundzary|0−conditBions are ∂ φ0+ = 0, and ℓ2∂3φ0+ = al. [16], which also exhibits a crossover to a slow decay σℓ /Q. Charge neutralityz, |0−dz κ2exp[ζ(−z) zφ(|z0)−] = far from the surface. Using ζ(z), I obtain the estimate B ℓBσ/Q, is ensured for any solRution to Eq. (10)−. ln(λren/λ)=Γ[1−e−√ℓB/(2ℓΓ)]. Analytical approximations to Eq. (10) can be found Eq. (10) has been solved numerically for Γ = 100 in the WC andSC limits. Inthe WC limit, I assume the and Γ = 10 (inset). Fig. 1 plots n(z) n (z), where PB − solution will decay with characteristic length λ. Thus, n(z) = ρ(z)ℓ λ2/2 and n (z) = 1/(1+z/λ)2 is the B PB the fourth order derivative is negligible. Since ζ(r) 1, normalizedPoisson-Boltzmanndensity. Thesenumerical ≪ φ has the PB equation form, given by solutions are compared to actual simulation data from Ref. [14] (courtesy of A. Moreira) and show quite good φ(z >0)=2ln 1+κz/√2 . (11) agreement. Furthermore,Eq. (10)outperformstheTCT (cid:16) (cid:17) (shown as dashed lines using data providedby Y. Burak The boundaryconditions aresatisfiedby choosingφ(z < from Ref. [16]). The nonperturbative function ζ(z) is 0) = 2(ℓ/λ)A(ez/ℓ 1), where ℓ2A3/λ2 A = 1. This an important component of this numerical agreement; − − requires κλ/2 = A. As Γ 0, φ(z < 0) 0, and Eq. when a virial expansion in ρ is used to compute lnZ → → 0 s (11) becomes exact. (equivalent to setting ζ(z) = F(z) at lowest order), the In the SC limit, the fourth order term dominates over agreementwiththesimulationdataisonlyslightlybetter the second order term near the surface. I make the ad- than the TCT and not nearly as good as Fig. 1. A more ditional assumption that φ 1, and solve carefulevaluationofζ(z)islikelytoimprovetheseresults ≪ further. To be clear, the value for r used in Fig. 1 is ℓ2∂4φ=κ2eζ(z) κ2eζ(0)−z/λ, (12) 0 z ≈ simplyanestimateoftherealcorrelationholesize,which may differ up to a factor of order one depending on the which has the solution model used. Though there is still very good agreement φ(z >0)=κ2eζ(0) e−z/λ 1 . (13) for other reasonable values of r0, r0 = ℓBλ/2 seems to (cid:16) − (cid:17) give the best agreement with simulatiopns. In contrastto Applying the boundary conditions, φ(z < 0) = the Γ = 10 and 100 results, the density at Γ = 1 always (2/ℓ2)/(1 λ2/ℓ2)[exp(z/ℓ) 1] and κ2 =(2/λ2)/(1 decays faster than the Poisson-Boltzmann density (not λ−2/ℓ2)exp[−ζ(0)]. The expon−ential SC density is recov−- shown). This is preciselywhere boththe shortandlong- ered for la−rge Γ. This solution also becomes exact as distance expansions in the interpolation scheme attain Γ . their maximum error, and it is possible that including →For∞z ℓ, ζ(z) 1 and Eq. (10) has an approximate higher order terms may improve this. ≫ ≪ solution given by Eq. (11). The fourth order derivative The SC expansion can be reconstructed in this frame- 4 workasanasymptoticexpansionaroundtheΓ= limit well with the simulation data, and depends crucially on ∞ by considering the higher order terms in δρ˜ for lnZ . a nonperturbative correlation correction whose form we 0 s Here we sketch this result: substituting the conjectured haveestimated. Thesecorrelationsalsoplayaroleinde- asymptoticallyexactresultρ˜ =2e−z/λ/(ℓ λ2) into ρ , terminingtherenormalizedGouy-Chapmanlengthwhen 0 B h i thefirstordercorrectioncanbecomputed. Itisstraight- thedensitiesrecoversitsPoisson-Boltzmanformfarfrom forwardtoshowthatthisgivesonlythefinite partofthe the surface. first order SC correction as Γ , and therefore that → ∞ this correction vanishes in this limit [12]. This occurs IwouldliketothankA.Moreiraforprovidingthesim- because the delta function in δρ˜ exactly cancels the di- 0 ulationdataandY.Burakforprovidingsolutionsforthe vergencefromρ˜ asΓ . Higherordertermswillalso 0 →∞ TCT. I also acknowledge many helpful discussions with vanishinthis limit, andanasymptotic expansioncanbe P.Pincus,M.Henle,D.Needleman,A.W.C.Lau,A.Liu, constructed in powers of 1/Γ which agrees exactly with and R. Kamien. This work was supported by the NSF SC up to first order and which I conjecture also agrees under Award No. DMR02-03755, Award No. DMR01- at all orders. It is interesting that, although no explicit 29804, and by the Materials Research Laboratory at renormalization is necessary in the modified expansion, UCSB under Award No. DMR00-80034. ζ(z)arisestoshiftthe fugacityκ2 inasimilarmanneras the SC renormalization. A more systematic accounting of these corrections is left for future work. Despite the quantitative agreement, the modified SC expansion suggests a different physical picture of the ∗ Electronic address: [email protected] strong-coupling limit: the corrections to the SC limit [1] A.Yu.Grosberg,T.T.Nguyen,andB.I.Shklovskii,Rev. arise from the interactions between only those counte- Mod. Phys. 74, 329 (2002); Y. Levin, Rep. Prog. Phys. rions that have made excursions away from the wall, as 65, 1577 (2002); A.G. Moriera, R.R. Netz in Electro- measured by δρ˜ . It is clear that the density of these static Effects in Soft Matter and Biophysics, edited by 0 C. Holm, P.Kekicheff,and R. Podgornik (KluwerAcad. excited counterions becomes small for large Γ, and that Pub., Boston, 2001). the SC limit becomes exact. The function ζ(z) encodes [2] V.I.PerelandB.I.Shklovskii,PhysicaA274,446(1999). the interaction of these excitations with their z =0 cor- [3] B.I. Shklovskii,Phys. Rev.E 60, 5802 (1999). relation holes, and becomes important at intermediate [4] P.Kékicheff,S.Marôelja, T.J.Senden,andV.E.Shubin, coupling, once counterions get far enough from the sur- J. Chem. Phys. 99, 6098 (1993). face. Interestingly, this correlation hole picture has al- [5] B.-Y.HaandA.J.Liu,Phys.Rev.Lett.79,1289 (1997) readybeendescribedinthecontextofWignercrystal-like [6] V.A. Bloomfield, Biopolymers 31, 1471 (1991) [7] R. Podgornik, D. Rau, and V.A. Parsegian, Biophys. J. structural correlations for multivalent ions [3, 18] and in 66, 962 (1994) the interpretation of the TCT [16]. [8] A.E. Larson and D.G. Grier, Nature (London) 385, 230 Thisschemealsoprovidesanaturalframeworktocom- (1997) pute the counterion free energy, which is obscured by [9] A.W.C. Lau and P. Pincus, Phys. Rev. E 66, 041501 the traditional SC expansion [12]. This is given by F = (2002) S(iφ) µN/(k T),whereN isthenumberofcounterions, [10] A.W.C. Lau, D.Levine,and P.Pincus, Phys.Rev.Lett. B φ is t−he mean-field long range potential, and the chem- 84, 4116 (2000) [11] R.R.NetzandH.Orland,Eur.Phys.J.E1,203(2000). ical potential is related to κ2 by µ = k T ln(κ2a3/ℓ ). B B [12] R.R. Netz,Eur. Phys. J. E 5, 557 (2001). Since κ2 depends exponentiallyonζ(0),nonperturbative [13] A.G.MoreiraandR.R.Netz,Phys.Rev.Lett.87,078301 correlationsalsoplayaroleinthefreeenergy,andsubse- (2001). quently in the interaction between two surfaces at sepa- [14] A.G. Moreira and R.R. Netz, Eur. Phys. J. E 8, 33 rations where the SC expansion cannot be applied. This (2002). will be explored in a future publication. [15] J.D.Weeks,D.Chandler,H.C.Andersen,J.Chem.Phys. 54, 5237 (1971); D. Chandler, J.D. Weeks, H.C. Ander- To summarize, I have computed the counterion den- sen, Science 220, 787 (1983). sity around a charged surface using a scheme to decom- [16] Y. Burak, D. Andelman, and H. Orland, Phys. Rev. E, pose the Coulomb interaction into short and long dis- 70, 016102 (2004). tance components. Each is treated with different ap- [17] For a simple and elegant proof, see reference [11]. proximations. ForlargeΓ,we recoverthe SC resultsand [18] T.T. Nguyen, A. Yu. Grosberg, B.I. Shklovskii, Phys. for small Γ we recover the WC Poisson-Boltzmann den- Rev. Lett. 85, 1568 (2000); T.T. Nguyen, A. Yu. Gros- sity. At intermediate coupling, the model agrees very berg, B.I. Shklovskii, J. Chem. Phys. 113, 3 (2000).