The Electrostatic Repulsion between Charged Spheres from Exact Solutions to the Linearized Poisson-Boltzmann Equation A. B. GLENDINNING 1 AND W. B. RUSSEL Department of Chemical Engineering, Princeton University, Princeton, New Jersey 08544 Received June 25, 1982; accepted August 27, 1982 An exact solution to the linearized Poisson-Boltzmann equation for the electrostatic potential around two equal, charged spheres is constructed via a multipole expansion. Results for the repulsive forces acting on the spheres are presented for 0.1 ~< aK ~ 20 over the full range of separations. Comparison with the linear superposition and Derjaguin approximations defines their ranges of validity and highlights the failure of the latter at small separations for the case of fixed surface charge densities. INTRODUCTION
surface potentials, or charges, over a broad range of separations and Debye lengths. An exact solution to the linearized PoissonBoltzmann equation is constructed in terms of multipole expansions centered on the two spheres and the boundary conditions are applied through use of an addition theorem for Bessel functions. Solutions for the interparticle force are presented for approach at constant surface potential and at constant charge density.
In aqueous colloidal dispersions electrostatic interactions play a major role in determining the stability of the dispersion and the magnitudes of various transport properties. Unfortunately, exact results for the repulsive force, even for the case of two equal spheres, are available only from computations at a limited number of conditions (1, 2) and asymptotic results for weak interactions or thin double layers. The linear superposition approximation (3) is valid in the far-field limit, i.e., for separations of several double layer thicknesses, while the Derjaguin approximation (4) is exact when both the Debye length and the separation are small relative to the sphere radius. Both approaches can accommodate the nonlinear PoissonBoltzmann equation. The regime of intermediate double layer thicknesses with arbitrary separations and surface potentials remains despite attempts (5) to patch together the Derjaguin approximation, as an inner solution, with the far-field superposition approximation. Here we address the interaction between equal spheres of radius a with small, equal
MULTIPOLE EXPANSIONS
For equilibrium double layers with small surface potentials the electrostatic potential 4/is determined by V2~ = (aK)Z~p
[ 1]
with the potential scaled with 60 and lengths with the radius; K-1 is the Debye length. In spherical coordinates with axial symmetry [1] has the general solution for 6 ~ 0 as r ----* oO
~p(r, O) = ~ a,K,(aKr)Pn(cos O)
[2]
n=0
with Kn the modified spherical Bessel function of the second kind and Pn the Legendre polynomial.
Present address: Raychem Corporation, Menlo Park, California. 95
0021-9797/83 $3.00 Journal of Colloid and Interface Science, Vol. 93, No. 1, May 1983
Copyright © 1983 by Academic Press, lnc, All rights of reproduction in any form reserved.
96
GLENDINNING
For two equal spheres the appropriate solution to [2] consists of two infinite sums of the form shown, one based on the (rl, 01) coordinate system and the other on (r2, 02) as defined in Fig. 1. The coefficients are equal due to the symmetry of the problem. To apply the boundary condition at the surface of either sphere we transform one of the multipole expansions into the other set of spherical coordinates through an addition theorem
AND
RUSSEL
for Bessel functions (6) as described by Marcelja et al. (7) so that oo
~b = ~ a~[Kn(axrOPn(cOs 01) n=O
+ X (2m + 1)B,mlm(arrOPm(cOs 01)], [3] m=O
where
c~3
Brim = ~ AnmKn+m-z~(arR) v=0
1
Anm =
+m-v)l(n+rn-2v+
~rr(m+n-v+~)(n-v)!(m-v)lv!
Im is the modified spherical Bessel function 0fthe first kind, r the Gamma function, and R the center-to-center separation. Note that the expression for A~n in (7) has several typographical errors. (a) Constant surface potential. For a fixed surface potential 6o on both spheres the dimensionless potential satisfies ~b= 1 on
1)
r l = 1.
With [3] this determines the linear system of equations (L + 6).A = e, where
Aj = ajKj(aK) Ljn = (2j + l)Bnslj(aK)/Ko(aK)
ej=
1, O,
j =0 j>0
and ~ is the unit tensor. Except at touching where the multipole expansion fails to converge, the aj can be determined by truncating the infinite system and then inverting numerically as discussed later. (b) Constant charge density. For fixed charge densities ¢ on each sphere, setting 6o = a¢l~ with ~ the dielectric permittivity provides the boundary condition
06 --1 Or1
on
rl-- 1.
The resulting linear system of equations has the same form but with
70 3_ FIG. 1. C o o r d i n a t e s y s t e m for t w o - s p h e r e interactions. Journal of Colloid and Interface Science, VoL 93, No, 1, May 1983
[4]
POISSON-BOLTZMANN
Aj = --(aK)ajK~(aK) Lj, = (2j + 1)B,,jI'j(a~c)/K'(aK).
[5]
The prime on the Bessel functions indicates differentiation with respect to the argument.
with the argument of all the Bessel functions being aK and CI (n, m) 2(m + 1) (2m+3)(2m+
FORCE CALCULATION
+ 1 .E2)n-(E.n)E]
1)' n = m +
1
( 2 m + l ) ( 2 m - 1)' n = m -
1
n#m+
1.
2m
Both the pressure and the Maxwell stresses contribute to the force (3) as
_
97
EQUATION
dA. [61
0,
For constant surface charge
The electric field E is scaled with ~o/a and the force with ~b2. n is the unit normal on sphere 1 and 1
H = ~(aK)
22
¢
[lO]
+ (,- .)4°ffL \ OM _1 +2(1-Fz 2)~
[7]
du
oo
is the local excess osmotic pressure due to ions in the double layer. The sole nonzero component of the force, of course, acts along the line of centers with magnitude
= 7r ~ an[(a~c)2Kn ~ aiKiCl(n, i) n=0
i=0 oo
+ 2(aK)ZKn Z ai E (2j + 1)BijljCl(rt, j) i=0
F = ~
{[(a~)2¢ ~ + E0~ - E l ] . [81
with t~ = cos 0t. For constant potential [8] reduces to F = -~r
f_
j =0
oo
I
+ 2ERE0(1 - u2) 1/2} & ,
co
oo
+ (aK)2 ~, (2m + 1) BnmIm ~, m=0
i=0
× ai ~ (2j + 1)BqI~G(m, j) j=0 of)
E~ u d .
o(3
+ K. ~ aiKiC3(n, i) + 21£.
1
i=0
oo
i=0
oo
= (a~) 2 ~, a, [K',~ ~ aiK}Cffi, n) n=0
i=0
co
× ai ~ (2j + 1)BijIj C3(n, j) j=0
oo
co
+ 2K" ~ ai ~ (2j + 1)BqI~Cl(j, n) i=0
j =0
o~
m=0
o~
+ ~ (2m + l)B.mlm ~ ai ~ (2j + 1) m=0
+ Z B.m(2rn + 11I"
o:3
i=0
j =0
871"
i=0
× BqljC3(m, j) + 8B.111] + ~ - alKl ,
c:(3
[11]
× ai ~ (2j + 1)BJ~C,(m, j)] [9] j=0
where Journal o f Colloid and Interface Science, Vol. 93, No. 1, May 1983
98
G L E N D I N N I N G A N D RUSSEL
next section to reconsider the limit of thin double layers.
TABLE I N u m b e r of Terms for Four-Place Convergence
Constant potential
C O R R E C T I O N TO THE D E R J A G U I N APPROXIMATION
Constantcharge
ar (r-2)
l
2
lO
20
1
2
10
20
0.1-0.5 0.5-1.0 1.0-2.0 2.0-3.0 3.0-4.0 4.0-5.0 >5.0
11 8 6 6 5 5 4
13 9 9 7 7 6 4
21 14 14 14 13 12 9
25 22 17 17 16 15 15
16 8 7 6 5 3 3
22 10 10 7 6 5 4
42 27 17 15 12 12 11
44 38 27 21 18 16 11
C3(n, m)
=
2(m + 1)(m + 2)m (2m+3)(2m+ 1)' n=m+
1
2 m ( m + 1 ) ( r n - 1) (2m + l)(2m 1) ' n
1
I
O,
The Derjaguin approximation for the electrostatic force between charged spheres with thin double layers has been widely exploited but only recently (2, 8) tested to determine the range of aK for which it is numerically accurate. In this section we derive the next term in a regular perturbation expansion for ar >> 1 to gauge the error and to provide a check for our numerical results. This extends the work of Ohshima et al. (8) to include the case of constant surface charge. Expressing [1] in cylindrical coordinates centered at the midpoint between the spheres with the rescaling Z ~
m-
n4=m++_l.
The numerical evaluation of the forces involved first truncating the linear systems [4] and [5] and using Gaussian elimination to determine the coefficients a~ and then summing [9] or [11]. Double precision was used routinely with a few computations in quadruple precision to establish that round-off error was insignificant. The m a x i m u m number representable on the computer ( 10 75) limited the computations to aK < 35. Recursion formulas determine the Kn accurately and rapidly from K0 and K1 but the rapid growth of round-off error dictated that the exact series representation be used for In. For aK ~ O(1) 10 terms in the expansion suffice for all but the smallest separations. For dimensionless separations less than 0.1 up to 44 terms were required to converge the force to within 0.1% (Table I). The results are most meaningful when compared with the limiting cases noted in the introduction. Hence we digress in the Journal of Colloid and Interface Science, Vol. 93, No. 1, May 1983
z/aK
r ~ r/(aK) 1/2 produces 024
Oz ~
1 1 0
~b
aK r Or
04
[12]
F--
Or
with 04
0
Oz
at
z
0
and either
4=1 or
a~ nr a~ nz-- + -
Oz
1
aK Or
on the sphere surface defined by
+
r7.
where 2ho is the m i n i m u m separation and the scaled normals are
nz = - ( 1 + rZ]-l/ZaKj r4 [ = r-1 nr t (LTl
1,2 "
POISSON-BOLTZMANN EQUATION For ag >> 1 we expand as
0 [1~] -~K
02@2 Oz 2
@ = @1 + @2+ . . . aK
lr2
h(r) = ho + ~
0@1
r Or r O~
0@2 -0 Oz
_lr 4
10
@2--
at
z=0
0@1 @2 = - & - Oz
+ 8 aK + " ' "
[16]
or
=hi+h2+ aK -
99
-
,
.
,
0%
n~=-I
r2 +2aK + "'"
nr=r+
..-
(hi - ho){~z' + 2 0@1 1 -- ~ (hi
[13]
and collapse the second b o u n d a r y conditions onto z = hi through the Taylor series
-
ho)
02@1~ ~ - ~ at
z =
hi.
The solutions are: (a) constant potential cosh z cosh hi
@1 -- - -
sinh z cosh hi
@(h) = @1(hl)
@2 = f z - -
1{ 0@1 } + - - @2(hl)+h2 (hi) + - • • aK &-z
[
1
1
@1 tanh hi f h l + ~ (hi - ho) 2 ,
-
[17]
with
0@ 0@1 1 10@2 O~ (h) = O z (hi) + --aK ( Oz (hi)
f = tanh h i + ( h i - ho)(1 - 2 tanh2h0; (b) constant charge
+ (hi -~ h°)2 02@10z 2 (hi)} + • • ".
[141
Then the potential is determined by O(1)
02@1 Oz 2
cosh z @1 sinh hi sinh z @2 = z sinh h~ {ctanh hi + (hi - ho)
@1 = 0
× (1 - 2 ctanh 2 hi)}
0@1 - 0
at
Oz
z=0
× ctanh hi{ 1 + h i ctanh hi 2 hi(hi - ho) sinh 2 hi
@1 = 1
or
0@1 Oz
-
¢2( 0 , O)
1
sinh z sinh hi
at
.....
z=h~,
[15]
1 (h~ - ho2)} 2
[18]
In both cases the correction for finite aK reduces the potential in the gap between the particles, e.g.,
(ho2 tanh 2 ho,
const potential
-~l(O'O)~,ctanh" ho(1 + ho tanh ho),
const charge [19]
<0. Journal of Colloid and Interface Science,
Vol.93, No. 1, May 1983
100
GLENDINNING
A N D RUSSEL
3°L 2£ I 2.2
uY 1.4 1.0 0
12
.6
1,8
3D
2.4
oK ( r - 2 )
FIG. 2. Ratio of linear superposition force to exact force for constant potential. (-- • --) aK = 0.1, (. . . .
) aK = 0 . 5 ,
(
) aK =
1.0, ( ....
) aK =
The force acting on either particle [6] can be converted through use of the divergence theorem to an integral over the symmetry plane of the form
10.0.
with the upper sign pertaining to constant potential and the lower to constant charge. Here e-2h0
Fo = 2~-aK 1 _ e -2h°
F=Trfo~I(O~b~2+aKlp2} is the Derjaguin result. Note that the limiting forms = aK
+0
1~2
aK
0:(F
+ (h,- ho,t
l j; dh,. [2o] h o ~
The correction term consists of a negative effect of the reduced potential within the gap, which lowers the local osmotic pressure, and a positive contribution from the Maxwell stress arising from the O(aK) -1/2 radial gradients in ~,. Hence the correction can be of either sign. The final result for the force is
F=Fa 1
(aK)_i 1 + e -2h° 1 e_4h °
× -t-
T- 1 + h o
l__+_e-2h°
6e2ho(1 ± e-2ho)2 ln(1 + e 2ho)l}
Journal of Colloid and Interface Science,
h0~:
I 0.102, const potential
)
F-d- 1 a K ~
, -1/ho, const charge
-I
1
aK~+--h0 -2
[221
indicate that the Derjaguin approximation errs seriously in the case of constant charge at small separations. And in that limit the entire expansion breaks down when ho < (aK)-l since the correction becomes comparable to the leading order term. The results for constant potential conform with those of Ohshima et al. (8). DISCUSSION OF RESULTS
[21]
Vol. 93, No. 1, May 1983
The linear superposition approximation for the force has the dimensionless form (3)
PO1SSON-BOLTZMANN
EQUATION
10t
6,ol o.of'!/ ', \
L.L~ Cu~O
20 \\')', .6
o
t2
1.8
2.4
30
ex(r-2) FIG. 3. R a t i o o f linear superposition force to exact force for c o n s t a n t charge. ( - - . - - ) ) aK 1.0, ( - - - - - ) aK = 2.0, ( . . . . . ) aK = 5.0, ( . . . . ) aK = 10.0.
F-
(1 + aKR) R2
x exp{--aK(R
I
2)}
const potential 1 (1 + aK) 2, ~ const charge J
[231 with the difference between the two cases arising solely from the scaling of the force. The curves in Figs. 2 and 3 illustrate the deviation of this far-field limit from the exact
aK = 0.1,
result at small separations. The approximation becomes exact except at contact for aK ~ 0 but the error increases with increasing ar at all separations. The force ratio is inverted between the two figures since linear superposition overestimates the force for interactions at constant potential but underestimates the effect at fixed surface charge density. Figure 4 compares the exact theory with the Derjaguin approximation for the constant potential case. For aK > 1 the agreement
3.0
2.4
1.8
u_ ~
t,2
r
0
1.0
I
I
I
2.0
3.0
4.0
5.0
Gx(r-2)
FIG. 4. R a t i o o f D e r j a g u i n force to exact force for c o n s t a n t potential. ( = 2.0, ( . . . . ) aK = 10.0, ( - - - ) aK = 20.0.
) a~ = 1.0, ( - - - - - )
aK
Journal of Colloid and Interface Science, Vol. 93, No. 1, May 1983
102
GLENDINNING
AND
RUSSEL
o
-I
-2
I I
0
I 2
I 3
I 4
az (r-2) FIG. 5. C o r r e c t i o n t o D e r j a g u i n f o r c e f o r c o n s t a n t p o t e n t i a l . (.... ) a~ = 10.0, ( . . . . . ) aK = 2 0 . 0 , a n d ( ) Eq. [36].
is quite good. The disparity increases with separation but only becomes significant for aK ~ O(1) since the forces decay exponentially. The plot of the correction to the Derjaguin result scaled by aK in Fig. 5 conforms to the trends predicted by the O(1/aK) perturbation theory in the previous section. For aK > 20 the first-order perturbation theory becomes exact.
) aK =
~J
.J
.J
~
O/ 0
~
I
----~
2
1.0, ( - - - - - )
aK = 2 . 0 ,
The situation for approach at constant charge shown in Fig. 6 is quite different. The Derjaguin approximation becomes accurate for O(1) separations only for aK > 30 and always errs substantially at small separations. Figure 7 confirms this error to be in qualitative accord with the O(1/aK) correction derived in the previous section, with quantitative agreement indicated for aK > 20. The
.J
II•Zz•-z_--
(
.J
----
I
I
I
I
I
I
3
4
5
6
7
8
aK(r-2)
FIG. 6. R a t i o o f D e r j a g u i n f o r c e t o e x a c t f o r c e f o r c o n s t a n t c h a r g e . ( - - ) (.... ) aK = 10.0, ( - - - ) aK = 2 0 . 0 .
Journal of Colloid and Interface Science, Vol. 93, No. 1, May 1983
aK = 1.0, ( - - - - - )
aK = 2 . 0 ,
.
~.'i .......
-_ . . . . . . . .
- .....
7 ..........
i
I
0
I
I
I
I
I
2
3
4
oK{r-2) ) aK
FIG. 7. C o r r e c t i o n to D e r j a g u i n force for c o n s t a n t charge. (. . . . ) aK = 10.0, ( - - - ) au = 20,0, a n d ( ) Eq. [36].
1.0, (-----) aK = 2.0,
/!
10 iI // // iI ii iI
6
/
/
i/ oK
LSA
0
~
~
2
4
ox(r -2) iO B
8
6
E
4
LSA 2
0 o
4
', aK (r-2)
FIG. 8. V a l i d i t y o f results for (A) c o n s t a n t p o t e n t i a l a n d (B) c o n s t a n t charge. W i t h i n each region, the label d e n o t e s which, if any, a p p r o x i m a t e f o r m for the force is w i t h i n 10% o f the e x a c t value. 103 Journal of Colloid and Interface Science, Vol. 93, No. I, May 1983
104
GLENDINNING AND RUSSEL
problem at small separations arises from the increasing magnitude of radial gradients in potential, i.e., the term on the fight-hand side of [12]. As that term becomes large the Derjaguin approximation, and the regular perturbation expansion of the previous section, breaks down. At the same time, however~ the potential becomes sufficiently large to render invalid the linearization on which the entire analysis is based. Figures 8A and B summarize the regions of applicability for the various theories for constant potential and constant charge, respectively. The boundaries correspond to approximately a 10% error for the approximate theories. At constant potential the approximations cover much of the domain, but not at constant charge since the Derjaguin result only becomes accurate for aK > 30. Available numerical solutions (2) define in part the effect of nonlinearity on the conclusions presented here. For constant poten-
Journal of Colloidand InterfaceScience, VoL93, No. 1, May 1983
tial conditions the error in the force is ~ 10% at e~ko/kTS4 for aK = 10, but the criteria for validity of the aK >> 1 approximation appear to be unaffected. The situation for the constant charge case, where the error could be much greater, remains to be defined. REFERENCES 1. McCartney, L, N., and Levine, S., J. Colloid Interface Sci. 30, 345 (1969). 2. Chan, B. K. C., and Chan, D. Y. C., preprint. 3. Bell, G. M., Levine, S., and McCartney, L. N., J. Colloid Interface Sci. 33, 335 (1970). 4. Hogg, R., Healy, T. W., and Fuerstenau, D. W., Trans. Faraday Soc. 62, 1638 (1966). 5. Bell, G. M., and Peterson, G. C., J. Colloidlnterface Sci. 41, 542 (1972). 6. Langbein, D., "Theory of van der Waals Attraction." Springer, Berlin, 1974. 7. Marcelja, S., Mitchell, D. J., Ninham, B. W., and Sculey, M. J. J. Chem. Soc. Faraday Trans. H 73, 630 (1977). 8. Ohshima, H., Chan, D. Y. C., Healy, T. W., and White, L. R., J. Colloid Interface Sci. 90, 17 (1982).