9832
Langmuir 2005, 21, 9832-9842
Determination of the Dynamic Electrophoretic Mobility of a Spherical Colloidal Particle through a Novel Numerical Solution of the Electrokinetic Equations† Matthew A. Preston,‡ Ralph Kornbrekke,§ and Lee R. White*,‡ Center for Complex Fluids Engineering, Department of Chemical Engineering, Carnegie Mellon University, Pittsburgh, Pennsylvania 15213, and Lubrizol Corporation, 29400 Lakeland Boulevard, Wickliffe, Ohio 44092 Received March 1, 2005. In Final Form: July 10, 2005 The standard equations developed to describe the electrophoretic motion of a charged particle immersed in an electrolyte subjected to an oscillating electric field are solved numerically with a new technique suitable for stiff systems. The focus of this work is to use this solution to determine the dynamic particle mobility, one of several quantities that can be extracted from these equations. This solution is valid from low frequencies to indefinitely high frequencies and has no restriction on zeta potential, double-layer thickness, or electrolyte composition. The solution has been used to calculate the dynamic electrophoretic mobility of a particle for a wide range of double-layer thicknesses and zeta potentials. The solution agrees with analytic approximations obtained previously by other authors under the conditions of a thin double layer and low zeta potential. The results are also consistent with calculations valid at frequencies where the ion diffusion length extends a significant distance beyond the double layer as obtained by another numerical technique.
1. Introduction The measurement of high-frequency electrokinetic properties of colloidal suspensions has been advanced by the development of electroacoustic techniques.1,2 When an alternating electric field is applied to a colloidal system, the electrophoretic motion of the particles gives rise to an acoustic wave, of which the electrokinetic sonic amplitude (ESA) can be measured and used to compute the dynamic electrophoretic mobility µ,3 defined implicitly (for a spherical particle) by
U ) µ(ω)E
(1.1)
where U is the particle velocity and E is the applied electric field of radial frequency ω. The measurement of electrophoretic mobility is a noninvasive probe for particle characterization. Static electrophoresis experiments, with the use of an appropriate theory,4,5 for example, can be used in the determination of the zeta potential ζ and the surface conductance parameter. High-frequency electrophoresis experiments can yield more information than static measurements because the dynamic mobility is related not only to the electrical properties of the surface but also to the particle size a and the size distribution.3 The ability to carry out such characterization rests on the availability of a stable and accurate technique to solve the standard set of †
Part of the Bob Rowell Festschrift special issue. * Corresponding author. E-mail:
[email protected]. ‡ Carnegie Mellon University. § Lubrizol Corporation. (1) Oja, T.; Petersen, G. L.; Cannon, D. W. U.S. Patent 4,497,208, 1985. (2) O’Brien, R. W.; Cannon, D. W.; Rowlands, W. N. J. Colloid Interface Sci. 1995, 173, 406. (3) O’Brien, R. W. J. Fluid Mech. 1988, 190, 71. (4) O’Brien, R. W.; White, L. R. J. Chem. Soc., Faraday Trans. 2 1978, 74, 1607. (5) Mangelsdorf, C. S.; White, L. R. J. Chem. Soc., Faraday Trans. 1990, 86, 2859.
electrokinetic equations.6 O’Brien3 developed an analytic approximation for the dynamic mobility of a single particle in the case of a thin double layer, and Mangelsdorf and White7 and Ohshima8 both found a general expression under the conditions of arbitrary double-layer thicknesses and low zeta potentials. A full numerical treatment of the equations by Mangelsdorf and White9 improved upon earlier calculations,10 but this technique is limited by numerical difficulties at high frequencies. Ennis and White11 developed a technique for numerical extension into the high-frequency regime. The numerical technique of Mangelsdorf and White9 contains the implicit assumption that the functions describing the ion dynamics and the perturbed electrical potential decay on a length scale that is greater than the double-layer thickness. As a result, this technique is not applicable to the study of nonaqueous colloidal suspensions, which typically have thick double layers due to the very weak ionization of the charge-carrying micelles.12 In addition, aqueous colloidal suspensions with very low electrolyte concentration or any suspension subjected to high-frequency electric fields requires an alternate solution method in order to determine the dynamic mobility. More recently, other authors have published solutions to the electrokinetic equations but have focused on the determination of quantities other than the dynamic mobility. Hill et al.13 used a finite difference method to study steady electrophoresis of polymer-coated particles. (6) Hunter, R. J. Foundations of Colloid Science; Clarendon Press: Oxford, U.K., 1989; Vol. 2. (7) Mangelsdorf, C. S.; White, L. R. J. Colloid Interface Sci. 1993, 160, 275. (8) Ohshima, H. J. Colloid Interface Sci. 1996, 179, 431. (9) Mangelsdorf, C. S.; White, L. R. J. Chem. Soc., Faraday Trans. 2 1992, 88, 3567. (10) DeLacey, E. H. B.; White, L. R. J. Chem. Soc., Faraday Trans. 2 1981, 77, 2007. (11) Ennis, J.; White, L. R. J. Colloid Interface Sci. 1996, 178, 460. (12) Morrison, I. D. Colloids Surf., A 1993, 71, 1. (13) Hill, R. J.; Saville, D. A.; Russel, R. B. J. Colloid Interface Sci. 2003, 258, 56.
10.1021/la050555d CCC: $30.25 © 2005 American Chemical Society Published on Web 09/20/2005
Electrophoretic Mobility of a Colloidal Particle
Langmuir, Vol. 21, No. 22, 2005 9833
The group used the same method to study the frequencydependent electrokinetic equations by determining the polarizability and complex conductivity of particles with polyelectrolyte14 and neutral15 coatings. Zhou et al. (using the standard method of Mangelsdorf and White) calculated the electric polarizability of a dielectric particle.16 Although the model and technique described in this work are capable of handling an arbitrary electrolyte, to demonstrate its applicability to nonaqueous dispersions the model is simplified to represent an electrolyte consisting of two monovalent ions of equal limiting conductance. This simplification is based on the work of authors who have proposed charging mechanisms including the adsorption of micelles to the particle surface followed by proton transfer17 and two-body interactions of neutral micelles. These mechanisms result in large charge carriers of essentially equal diameters and drag coefficients. The primary contribution of this work is the determination of the dynamic electrophoretic mobility through a numerical solution of the electrokinetic equations that does not rely on any assumptions regarding the relative sizes of the length scales. The solution contains no restriction on the frequency range, electrolyte concentration, or magnitude of the zeta potential. The equations are formulated as a boundary-value problem and solved using a collocation method.18
The other characteristic distances in the problem are the ion diffusion length and the hydrodynamic penetration depth, which are incorporated as ratios to the doublelayer thickness as
γj2 )
γ2 )
iωλj 2
)
κ kBT iωF0 2
)
κ η0
- layer decay length [double ] ion diffusion length
Lφj(x) + γj2[φj(x) - Y(x)] )
[
]
h(x) dy dφj - mj z dx j dx x j ) 1, 2,..., N (2.1)
N
LY(x) )
βj exp[-zjy(x)][Y(x) - φj(x)] ∑ j)1
L[Lh(x) + γ2h(x)] ) -
1 dy
(2.2)
N
∑βj exp[-zjy(x)]φj(x)
x dxj)1
(2.3)
where the second-order differential operator L is given by
L)
2 d 2 d2 - 2 + 2 x dx dx x
(2.4)
κ2 )
e2 r0kBT
2
(2.7) The parameters mj and βj, the ionic drag coefficients λj, and the mobility scaling factor µ0 are given by
mj )
2λjµ0 e zj2n∞j
βj )
N
∑
(2.8)
(2.9)
zk2n∞k
k)1
λj )
µ0 )
NAVe2 Λ∞j
2 r0kBT 3 eη0
(2.10)
(2.11)
The valence zj, the bulk concentration n∞j , and the limiting ionic conductance Λ∞j are given for a general electrolyte of N ionic species. The solvent is characterized by the relative dielectric constant r, temperature T, viscosity η0, and density F0. The fundamental electrical charge is written as e, kB is Boltzmann’s constant, NAV is Avogadro’s number, 0 is the dielectric permittivity in a vacuum, and i is x-1. The scaled zeta potential ζˆ ) eζ/kBT enters the problem through the scaled electrostatic potential y(x), which is determined by solving the PoissonBoltzmann equation:
d2y dx2
+
2 dy x dx
N
+
βj exp[-zjy(x)] ) 0 ∑ j)1 y(κa) ) ζˆ
80 y(x) 9 xf∞
and x is the radial distance, scaled by the Debye length κ-1
(2.6)
double - layer decay length [hydrodynamic penetration depth]
2. Method Governing Equations. The theory describing the frequency-dependent behavior of a spherical colloidal particle immersed in an electrolyte subjected to an electric field has been presented elsewhere (ref 8, for example). The system that describes the ion dynamics functions φj, perturbed potential function Y, and hydrodynamics function h can be written as
2
(2.12)
The boundary conditions on eqs 2.1-2.3 at the particle surface are9
N
zj2n∞j ∑ j)1
(2.5)
(14) Hill, R. J.; Saville, D. A.; Russel, W. B. J. Colloid Interface Sci. 2003, 263, 478. (15) Hill, R. J.; Saville, D. A.; Russel, W. B. J. Colloid Interface Sci. 2003, 268, 230. (16) Zhou, H.; Preston, M. A.; Tilton, R. D.; White, L. R. J. Colloid Interface Sci. 2005, 285, 845. (17) Pugh, R. J.; Matsunaga, T.; Fowkes, F. M. Colloids Surf. 1983, 7, 183. (18) Ascher, U.; Christiansen, J.; Russell, R. D. COLSYS; Netlib, http://www.netlib.org.
dφj | )0 dx x)κa
(2.13)
p dY |x)κa Y(κa) ) 0 dx rκa
(2.14)
h(κa) )
dh | )0 dx x)κa
(2.15)
The boundary conditions far from the particle surface are
9834
Langmuir, Vol. 21, No. 22, 2005
Preston et al.
determined by solving the homogeneous forms of eqs 2.12.3:
Lχ + A‚χ ) 0
(2.16)
L[Lh(x) + γ2h(x)] ) 0
(2.17)
The solution to these equations can be written as9
()
χk 98 x + C0 xf∞
2
κa
N
+
x
CjRkj f (xR )(x) ∑ j)1
The primary motivation of this work is to generate a numerical solution of this system under conditions that are applicable to nonaqueous systems. As discussed above, the charge carriers in these systems are of equal size. Thus, the drag coefficient is the same for each ion and can be written as
λI ) 6πη0RI
and by eq 2.10, the limiting conductance and the chargecarrier radius RI are then related by
j
k ) 1, 2,..., N + 1 (2.18)
[ ] [ ]
RI )
x 2B 8 µ(ω) + 2 + CN+1 f (γ)(x) h9 xf∞ 3 3x
(2.19)
(γ) dh 98 µˆ (ω) 1 - 4B + C df N+1 dx 3 3x3 dx
(2.20)
xf∞
µ(ω) µ0
2 1 NAVe Λ∞I 6πη0
(2.28)
The solution conductivity K is
where Cj (j ) 0, 1, ..., N+1) represents asymptotic coefficients to be determined along with the scaled electrophoretic mobility:
µˆ (ω) )
(2.27)
(2.21)
In these equations
r0kBT K ) (κa)2 λIa2 Under these conditions,
(
γI2 0 - γI2 A ) 0 γI2 - γI2 1 /2 1/2 - 1
(2.29)
)
(2.30)
and the eigenvalues and eigenvectors of A are 3
B)
(κa) (Fp - F0) 3F0
(2.22)
()
for a particle with density Fp and
φ1(x) φ2(x) ‚ χ(x) ) ‚ φN(x) Y(x)
(2.23)
The integration constants Cj can be determined (as described below) along with the scaled electrophoretic mobility after eqs 2.1-2.3 are solved subject to boundary conditions on the particle surface (eqs 2.13-2.15) and far from the particle (eqs 2.18-2.20). The function f (β)(x) satisfies
Lf (β)(x) + β2f (β)(x) ) 0 f (β)(κa) ) 1
(2.24)
The solution of system 2.24 can be written as
f (β)(x) )
( )( x0 x
2
)
1 - iβx exp[iβ(x - x0)] (2.25) 1 - iβx0
Every derivative of f (β) can be written as the product of a rational function and the same function, so the ratio of any function f (xRj) to its derivative is written as
rj(x) )
x [ 1
2 Rj
ix -
1
]
xRj
(2.26)
The matrix R consists of the N + 1 column eigenvectors of A. The N + 1 eigenvalues of A are Rj, where R0 ) 0.
R ) {0, γI2, γI2 - 1}
(2.31)
1 - 1 γI2 R) 1 1 γI2 1 0 1
(2.32)
(
)
Standard Solution Method. The standard numerical scheme employed by O’Brien and White4 for the static mobility problem was adapted by Mangelsdorf and White9 to solve the equations described above. The governing equations are solved as an initial value problem while systematically setting a single asymptotic coefficient equal to 1 while the others are set equal to zero. The function values are recorded at the surface of the particle, and once each homogeneous problem and the two inhomogeneous problems are solved, the true coefficients are determined by solving a matrix equation relating the computed surface values to the surface boundary conditions. Mangelsdorf and White noted that they were unable to compute a solution for small κa systems at high frequency. They attributed this difficulty to the fact that the ion diffusion length becomes the dominant length scale, giving rise to stiffness in the differential equations. Ennis and White11 pointed out that the method relies on the derivation of an asymptotic form for the perturbed ion densities that extend beyond the double layer. Figure 1 shows the ratio of the double-layer thickness to the ion diffusion length for two typical systems as a function of frequency. When this ratio is greater than 1, the ion diffusion length extends beyond the double-layer thickness, and calculations above this frequency will be increasingly time-consuming and ultimately inaccurate and unstable. Under these circumstances, the residual quantities in the asymptotic forms of each function cause the incorporation of an exponentially large error after a few steps using initial value methods. Hill et al.13 used a finite difference calculation to determine the electrophoretic mobility of a polymer-coated
Electrophoretic Mobility of a Colloidal Particle
Langmuir, Vol. 21, No. 22, 2005 9835
Conditions 2.13-2.15 and 2.18-2.20 (with the appropriate derivatives) are then rewritten as the following 15 equations (5 at the surface and 5 far from the surface with their derivatives) in terms of 5 unknown constants:
dΦ( | ) -1 dx x)κa p dΨ | Ψ(κa) ) dx x)κa rκa p mI d2y -1 µ ˆ (ω) r 3γI2 dx2
( )
Figure 1. Ratio of the double-layer thickness to the ion diffusion length for (A) system I and (B) system II with κa ) 1.
colloidal particle with disparate length scales. Here, we use a related technique, which is described in more detail below. The code is suitable for solving boundary-value problems but requires a reformulation of the system. It is necessary to write the boundary conditions, both at the surface and far from the particle, free of any undetermined parameters. Problem Reformulation. The N + 3 boundary conditions on the particle surface (eqs 2.13-2.15) and the N + 3 conditions far from the particle (eqs 2.18-2.20) make up the 2N + 6 conditions for the system (eqs 2.1-2.3) of the same order. Those far from the particle are written in terms of N + 3 undetermined parameters (µˆ (ω) and Cj, j ) 0,1,..., N + 1). The software used to solve this system uses a boundary-value collocation method that requires us to reformulate each of the conditions by eliminating the unknown constants in favor of linear relationships between the functions and their derivatives. Consequently, by taking linear combinations of the functions χk, h, and h′ and their derivatives, we end up with 2N + 6 conditions specified in terms of known values in addition to N + 3 equations for the unknown constants. The unknown constants can be determined a postiori from calculated function values. For the case of an electrolyte consisting of two monovalent ions of equal limiting conductance, we can simplify the boundary conditions by defining new functions. The transformation shown below in eq 2.33 removes the growing parts of the solution. This simplifies the boundary conditions and improves the computation speed.
x)κa
| ]
p dy rκa dx
x)κa
κa 3 dH 1 | ) - µˆ (ω) dx x)κa 3 κa 2 Φ((xmax) ) C0 - C1 f (xR1)(xmax) + xmax
(2.38) (2.39) (2.40)
( )
γI2C2 f (xR2)(xmax) (2.41) dΦ( (κa)2 df (xR1) |x)xmax ) -2C0 | - C1 + 3 dx dx x)xmax x max
( )
Ψ(xmax) ) C0
κa xmax
2
df (xR2) | (2.42) γI2C2 dx x)xmax + C2 f (xR2)(xmax) mI dy µˆ (ω) 2 |x)κa (2.43) 3γI dx
(κa)2 df (xR2) dΨ |x)xmax ) -2C0 | + C 2 dx dx x)xmax x 3 max
mI d2y µˆ (ω) 2 2|x)κa (2.44) 3γI dx B 2 H(xmax) ) µˆ (ω) 3 x
2 max
+ C3 f (γ)(xmax)
dH B 4 | ) - µˆ (ω) dx x)xmax 3 x
(2.45)
df (γ) | + C3 (2.46) dx x)xmax
B d2f (γ) d2H | ) 4µ ˆ (ω) + C |x)xmax (2.47) x)x 3 max dx2 xmax4 dx2
mI dy µˆ (ω) Y(x) ) Ψ(x) + x + 3γI2 dx x h(x) ) Η(x) + µˆ (ω) 3
-
H(κa) ) - µˆ (ω)
3 max
φ((x) ) Φ((x) + x
(2.33)
This transformation results in eqs 2.1-2.3 being represented as the following order 10 system:
LΦ((x) + γI2[Φ((x) - Ψ(x)] ) H(x) dy dΦ( dy - mI (2.34) ( ( dx dx x dx 1 LΨ(x) - Ψ(x)coshy(x) ) - [Φ+(x)e-y(x) + Φ-(x)e+y(x)] 2 (2.35)
[
[ |
(2.37)
]
L[LH(x) + γ2H(x)] ) 1 dy dy [Φ (x)e-y(x) + Φ-(x)e+y(x)] - coshy(x) (2.36) 2x dx + dx
d3H d3f (γ) B | ) -16µ ˆ (ω) + C |x)xmax x)xmax 3 dx3 xmax5 dx3
(2.48)
By introducing transformation 2.33, µˆ (ω) is now incorporated into the surface boundary conditions. To write these conditions free of any unknown constants, eqs 2.38 and 2.39 are combined to replace eq 2.38
p dΨ |x)κa Ψ(κa) ) dx rκa
( )
[ |
p mI d2y -1 + H(κa) r κaγI2 dx2
x)κa
-
| ]
p dy rκa dx
x)κa
(2.49)
9836
Langmuir, Vol. 21, No. 22, 2005
Preston et al.
and eqs 2.39 and 2.40 are combined to replace eq 2.39
H(xmax) )
{ [ H′′(x)
Η(κa) -
dΗ | κa ) 0 dx x)κa
(2.50)
There is now one equation (2.40) at the particle surface that contains µˆ (ω). This condition can be replaced by one far from the particle surface at the expense of increasing the order of the system by 1. This is done by introducing an auxiliary function v(x) that satisfies
dv )0 dx
v(κa) ) H(κa)
(2.52)
Because v is a constant, condition 2.40 at the surface is replaced by the following condition far from the particle:
κa )0 3
v(xmax) + µˆ (ω)
(2.53)
More algebra is required to eliminate all of the coefficients from the equations far from the particle. To do so, linear combinations are taken of eqs 2.41-2.48 and 2.54. These 11 equations can be written as 6 boundary conditions far from the surface (shown below) and 5 equations describing the unknown coefficients (shown later). The boundary conditions are written as relationships between each function and its derivatives:
Φ((xmax) )
[{[Φ ′(x) + Φ ′(x)][r (x)γ + 2x] 1 x - 2γ Ψ′(x)[r (x) + ]} 2 2(γ - 1) +
2
2
2
I
r1(x) ( [Φ+′(x) - Φ-′(x)] 2
{
[
- [4H′′(x) + xH′′′(x)] x4H′′(x)
}
[Φ+′(x) + Φ-′(x)] r2(x)γI2 +
]
(2.54)
[
{
x 2
]
}
]
{
[
2
[
(
C1 ) C2 )
{
1
[
]
{
C3 ) -
]
mI[2(γI2 - 1)y′(x) + (2r2 + γI2x)y′′(x)] 24B(γI2 - 1)γI2
(xR2)
2
]}
[
}
]]}
)
(2.59)
x)xmax
Φ-′(x) - Φ+′(x)
(2.60)
2f (xR1)′(x)
[
Φ+′(x) + Φ-′(x) - 2Ψ′(x) -
(
x4H′′(x)
6 6ix 3x2 ix3 - [4H′′(x) + xH′′′(x)] 4 - 3 - 2 + γ γ γ γ + x4H′′(x)
[
2(γI - 1)f ′(x) mIy′′(x) 6ix 3x2 ix3 6 (4H′′(x) + xH′′′(x)) 4 - 3 - 2 + + 2 γ 6BγI γ γ γ
x)xmax
1 x - 2γI Ψ′(x) r2(x) + 2 2(γ 2 - 1) I 2
{
6 κa 6ix 3x2 [4H′′(x) + xH′′′(x)] 4 - 3 - 2 + 12B γ γ γ ix3 + x4H′′(x) (2.58) γ x)xmax
Solution Method. Equations 2.34-2.36 and 2.51 subject to 2.37, 2.49-2.50, and 2.52 at the particle surface and 2.54-2.58 far from the surface form a system of coupled ordinary differential equations in complex variables of order 11. The system is solved using COLSYS,18 which is a general-purpose code for solving boundaryvalue ordinary differential equations by collocation. This code has been shown to be more effective than those based on shooting or initial-value techniques for some applications with complex behavior and has been used successfully in analyzing problems such as the vibration of curved beams,19 pump flows with interacting flow components,20 and material separation through simulated moving bed technology.21 The technique involves setting a sequence of meshes over the problem space until a specified tolerance is obtained. Between each interval on the mesh, a polynomial that approximates the solution is stored. As discussed above, expressions for the asymptotic coefficients are written as linear relationships of the derivatives of the functions, which can now be solved for numerically. These expressions are
]
mIy′′(x)(2r2 + x)
[
v(xmax) ) -
6 6ix 3x2 ix3 - 2 + + γ γ4 γ3 γ
24B(γI2 - 1)
x)xmax
1 {H′′(x)[4i + xγ] + ixH′′′(x)}x)xmax 3γ (2.57)
x4H′′(x)
Ψ(xmax) )
[{
H′(xmax) ) -
C0 )
I
2
I
]}
x3 Φ+′(x) + Φ-′(x) - 2γI2Ψ′(x) 2 4(κa) (γI - 1) mIy′′(x) 6 6ix 3x2 ix3 (4H′′(x) + xH′′′(x)) 4 - 3 - 2 + + 6B γ γ γ γ
2
-
[
(2.56)
(2.51)
subject to the condition
]
2 x x2 2ix ix2 + - 2 + H′′′(x) 6 3γ γ 6γ 2γ2
x)xmax
(2.55)
µˆ (ω) )
{
]]} }
1 [4H′′(x) + xH′′′(x)] xγ 2 4
[
)
(2.61)
x)xmax
(2.62) x)xmax
]
1 6ix 3x2 ix3 6 [4H′′(x) + xH′′′(x)] 4 - 3 - 2 + + 4B γ γ γ γ x4H′′(x)
}
x)xmax
(2.63)
Electrophoretic Mobility of a Colloidal Particle
Langmuir, Vol. 21, No. 22, 2005 9837
double-layer distortion into account was determined by O’Brien3:
µˆ (ω) ) ζˆ
1 - iγκa (1 + F(ω)) H-Γ
( (
) ( ) ) ( ) ) [ ] x (
3iω0r 3iω0p - 2λ 2K 2K F(ω) ) 3iω0r 3iω0p 2 1+ 2λ 2K 2K
(3.4)
1-
λ) Figure 2. Scaled static mobility as a function of scaled zeta potential for system I with κa ) (A) 1 and (B) 100. The solid lines are calculated using the theory of O’Brien and White (4) and the solid squares are the low-frequency points from the theory described in this paper. Table 1. System Properties Used to Generate Figures 1-14 property
system I
system II
particle radius (nm) temperature, solvent (K) dielectric constant, solvent (-) density, solvent (g mL-1) viscosity, solvent (cP) density, particle (g mL-1) dielectric constant, particle (-) limiting ionic conductance (cm2 Ω-1 (g equiv)-1) charge-carrier radius (nm)
100.0 298.15 78.54688 0.99704 0.8904 2.0 2.0 74.9
100.0 298.15 2.0 0.745 1.378 2.0 2.0 0.59
0.15
10.0
At a given frequency, the dynamic mobility can be determined by eq 2.63. As a consistency check, this quantity can also be calculated by
µˆ (ω) ) -
3 H(κa) κa
(2.64)
3. Results and Discussion The calculations have been performed on two systems. System I is a spherical particle immersed in an aqueous KCl solution, and system II is a carbon black sphere immersed in dodecane with a polyisobutylene succinimide dispersant. The properties of the suspensions are shown in Table 1. The static mobility (or the low-frequency dynamic mobility magnitude) has been calculated as a function of zeta potential by the standard O’Brien and White method4 and the method described above and is shown in Figure 2. As discussed in ref 4, the mobility increases monotonically at κa ) 1, but at κa ) 100, a maximum in mobility occurs at about ζˆ ) 6. Solution in Low/High κa Limits. The calculations based on the method described above are compared to several analytic approximations. The solution in the limit of low κa can be written as
1 2 µˆ (ω) ) ζˆ 3 H-Γ 1 H ) 1 - iγκa - (γκa)2 3 Γ)
(3.1) (3.2)
2
2γ B 3 κa
(3.3)
A solution in the limit of high κa that takes the effect of
3λIµ0 3 1 ζˆ 1+ exp 2 κa e 2
(3.5)
(3.6)
Figures 3 and 4 show a comparison of the numerical calculation of dynamic mobility at high (100) κa and low (0.1) κa to the analytic approximations in those limits for systems I and II. In the frequency range shown, there is strong agreement between the calculated and approximate solutions in these limits. One significant deviation in both systems is that the high-frequency phase angle calculated from O’Brien’s approximation3 is -45° (via the asymptote), while the phase angle calculated from the full theory is -90° (via the asumptote). (The resultant velocity amplitude is zero when the applied electric field is at its extreme values and vice versa.) A comparison of the two systems shows very little qualitative difference at low κa; however in the thin doublelayer case, a rapid drop in mobility magnitude is observed in the 1-10 MHz range of system II but not in system I. In this frequency range, the bulky charge carriers do not respond as quickly as the smaller aqueous ones. At even higher frequencies, the response of the double layer becomes electrolyte-independent. Aqueous Solvent (System I). The calculated dynamic mobility for system I is shown at increasing zeta potentials for three different values of κa in Figures 5-9. Figure 5 shows the dynamic mobility at κa ) 1. As shown in Figure 2, the low-frequency mobility magnitude increases monotonically as the zeta potential is increased and approaches a value just less than 4. The high-frequency mobility magnitude increases monotonically as well, but at moderate frequencies (from 0.1 to 10 MHz), a mobility maximum (in zeta potential) is apparent. There is strong agreement between the phase angle calculated by the low-zeta approximation, as shown in Figure 5b. The open squares in Figure 5 also show the greatest frequency at which the dynamic mobility can be calculated by ref 9. At this point, the assumption that the ion diffusion length extends beyond the double-layer thickness is no longer tenable. Under the conditions of Figure 5, the two length scales extend an equivalent distance at 0.06 MHz, and the double layer extends to nearly 4 times the ion diffusion length at 1 MHz. The increase in mobility magnitude with frequency has been explained by Gibb and Hunter22 by the corresponding phase lead, as seen in Figure 5b, which has a more pronounced effect with increasing zeta potential. The electric field that is felt by the moving particle is the result of the applied electric field and a back field caused by the double-layer distortion. When the double-layer is thick, (19) Tarnopolskaya, T.; Hoog, F. R. D.; Fletcher, N. H.; Sound, J. Vibration 1999, 228, 69. (20) O ¨ ztekin, A.; Seymour, B. R.; Varley, E. Stud. Appl. Math. 2001, 107, 1. (21) Lea˜o, C. P.; Rodrigues, A. E. Comput. Chem. Eng. 2004, 28, 1725. (22) Gibb, S. E.; Hunter, R. J. J. Colloid Interface Sci. 2000, 224, 99.
9838
Langmuir, Vol. 21, No. 22, 2005
Preston et al.
Figure 3. (a) Scaled dynamic mobility magnitude and (b) phase angle as functions of frequency for system I with ζˆ ) 2 and κa ) (A) 0.1 and (B) 100. The solid lines are calculated using the theory described in this article. The dotted lines are calculated using (A) the low-κa approximation (eq 3.1) and (B) the theory of O’Brien.3
Figure 4. (a) Scaled dynamic mobility magnitude and (b) phase angle as functions of frequency for system II with ζˆ ) 2 and κa ) (A) 0.1 and (B) 100. The solid lines are calculated using the theory described in this article. The dotted lines are calculated using (A) the low-κa approximation (eq 3.1) and (B) the theory of O’Brien.3
Figure 5. (a) Scaled dynamic mobility magnitude and (b) phase angle as functions of frequency for system I with κa ) 1. The arrows point in the direction of increasing zeta potential as ζˆ is increased from 1 to 10 (in increments of 1). The dotted line in b is calculated using the low-zeta analytic approximation.7,8 The open squares indicate the maximum frequency at which the dynamic mobility can be calculated using the method described in Mangelsdorf and White.9
at certain frequencies, the distortion lags behind the applied field, and the field felt by the particle results in a phase lead over the applied field. As the frequency is increased even further, the local field results in a phase lag, which ultimately continues to decrease to its limit and the mobility magnitude decreases to zero. The dynamic mobility at κa ) 10 is shown in Figures 6 and 7 over the same range of increasing zeta potential. In Figure 6, the low-frequency mobility magnitude increases monotonically (with zeta) to a maximum value, and then in Figure 7, the low-frequency mobility mag-
nitude decreases with an increase in zeta. The zeta potential at which the low-frequency mobility maximum appears (ζˆ ) 5) is consistent with ref 4. The mobility increases monotonically with zeta potential at all frequencies in Figure 6 (ζˆ e 5), but a maximum in zeta potential is apparent at higher frequencies in Figure 7 (ζˆ g 5). As in the case where κa ) 1, at significantly high frequencies (f g 100 MHz) the mobility magnitude increases monotonically with increasing zeta again. The dynamic mobility of the particle in system I under the conditions of a thin double layer (κa ) 100) is shown
Electrophoretic Mobility of a Colloidal Particle
Langmuir, Vol. 21, No. 22, 2005 9839
Figure 6. (a) Scaled dynamic mobility magnitude and (b) phase angle as functions of frequency for system I with κa ) 10. The arrows point in the direction of increasing zeta potential as ζˆ is increased from 1 to 5 (in increments of 1). The dotted line in b is calculated using the low-zeta analytic approximation.7,8
Figure 7. (a) Scaled dynamic mobility magnitude and (b) phase angle as functions of frequency for system I with κa ) 10. The arrows point in the direction of increasing zeta potential as ζˆ is increased from 5 to 10 (in increments of 1).
Figure 8. (a) Scaled dynamic mobility magnitude and (b) phase angle as functions of frequency for system I with κa ) 100. The arrows point in the direction of increasing zeta potential as ζˆ is increased from 1 to 6 (in increments of 1). The dotted line in b is calculated using the low-zeta analytic approximation.7,8
in Figures 8 and 9. The low-frequency mobility magnitude increases monotonically for ζˆ e 6 (Figure 8a) and decreases monotonically for ζˆ g 6 (Figure 9a). These results are predicted by ref 4 and shown in Figure 2. The most significant difference from those solutions at lower values of κa is the appearance of an additional feature above 10 MHz, which is most apparent in Figure 9b at the largest zeta potentials. A decrease in the magnitude of the phase lag causes an increase in the mobility magnitude. Just as the phase lead at around 1 MHz for the κa ) 10 system was responsible for the reversing of order in the zeta-ordered series, this feature
depicts the same order reversal. The effect is that ultimately, for ζˆ g 6, the mobility magnitude increases monotonically at a given frequency above 1 GHz. Hydrocarbon Solvent (System II). The calculations above have been reproduced under the conditions of system II, and these are shown in Figures 10-14. Many of the features seen in these graphs are similar to those seen in system I; however, the bulky size and sluggish nature of the charge carriers in system II tends to exaggerate many of the features and shed some light on their nature. Figure 10 depicts the dynamic mobility as the zeta potential is increased. Just as in system I, the mobility
9840
Langmuir, Vol. 21, No. 22, 2005
Preston et al.
Figure 9. (a) Scaled dynamic mobility magnitude and (b) phase angle as functions of frequency for system I with κa ) 100. The arrows point in the direction of increasing zeta potential as ζˆ is increased from 6 to 10 (in increments of 1).
Figure 10. (a) Scaled dynamic mobility magnitude and (b) phase angle as functions of frequency for system II with κa ) 1. The arrows point in the direction of increasing zeta potential as ζˆ is increased from 1 to 10 (in increments of 1). The dotted line in b is calculated using the low-zeta analytic approximation.7,8 The open squares indicate the maximum frequency at which the dynamic mobility can be calculated using the method described in Mangelsdorf and White.9
Figure 11. (a) Scaled dynamic mobility magnitude and (b) phase angle as functions of frequency for system II with κa ) 10. The arrows point in the direction of increasing zeta potential as ζˆ is increased from 1 to 5 (in increments of 1). The dotted line in b is calculated using the low-zeta analytic approximation.7,8
magnitude increases monotonically at all low and high frequencies. The mobility maximum (in frequency) appears as a plateau rather than a sharp peak, but at high frequencies (f > 0.1 GHz), the dynamic mobility appears to be independent of the drag coefficient of the charge carrier. The frequency at which the dynamic mobility can be calculated by ref 9 for system II is about 2 orders of magnitude lower than for system I. The ion diffusion length and the double-layer thickness extend an equivalent distance at about 0.5 kHz, and at 3.1 kHz, the double
layer extends nearly 2.5 times further than the doublelayer thickness. One major difference between system I and system II is apparent for the calculations at κa ) 10. The dynamic mobility for system II under this condition is shown in Figures 11 and 12. At low zeta potentials (ζˆ e 2), rather than seeing an increase in the mobility magnitude (or a constant magnitude) at moderate frequencies, a mobility decrease can be seen in Figure 11a. The phase angle in Figure 11b shows that in this range of zeta potential, where a phase lead was seen in system I, a phase lag appears.
Electrophoretic Mobility of a Colloidal Particle
Langmuir, Vol. 21, No. 22, 2005 9841
Figure 12. (a) Scaled dynamic mobility magnitude and (b) phase angle as functions of frequency for system II with κa ) 10. The arrows point in the direction of increasing zeta potential as ζˆ is increased from 5 to 10 (in increments of 1).
Figure 13. (a) Scaled dynamic mobility magnitude and (b) phase angle as functions of frequency for system II with κa ) 100. The arrows point in the direction of increasing zeta potential as ζˆ is increased from 1 to 6 (in increments of 1). The dotted line in b is calculated using the low-zeta analytic approximation.7,8
Figure 14. (a) Scaled dynamic mobility magnitude and (b) phase angle as functions of frequency for system II with κa ) 100. The arrows point in the direction of increasing zeta potential as ζˆ is increased from 6 to 10 (in increments of 1).
This phase lag decreases in magnitude to zero, and then the angle decreases again to -90° (via the asymptote) in the same way. After the decrease and increase in phase angle, the behavior is qualitatively the same as in system I as well as in system II at all zeta potentials. This effect is repeated in the case of κa ) 100 (Figures 13 and 14) at higher zeta potentials (ζˆ e 6). From Figure 13b, it can be noted that the increasing phase lag approaches a limit with decreasing zeta potential that is not predicted by the low zeta potential approximations. 4. Conclusions A numerical solution of the electrokinetic equations has been developed and used to calculate the dynamic elec-
trophoretic mobility of a particle immersed in an aqueous KCl solution and a nonaqueous solution with large charge carriers over a wide range of zeta potential and doublelayer thickness. The most common cases of those discussed above are the aqueous, thin double-layer system (system I, κa ) 100) and the nonaqueous, thick double-layer system (system II, κa ) 1). The method of solving the electrokinetic equations described in this article is stable over wide ranges of double-layer thickness and zeta potential, as shown in the calculations presented in this work. Because the asymptotic forms are written to include all of the terms that change on all of the relevant length scales, the solution
9842
Langmuir, Vol. 21, No. 22, 2005
is valid from the low-frequency regime to indefinitely high frequencies (shown here to 10 GHz). Double-layer distortion is responsible for the complex phase lead/phase lag behavior in the dynamic mobility. The distortion is strongly frequency-dependent, and its nature is dependent on particle characteristics as well as electrolytes. The model described in this article can be used with electroacoustic measurements to determine particle char-
Preston et al.
acteristics, such as zeta potential and particle size, of a colloidal system where previously published techniques are not adequate to describe the electrokinetic behavior. The equations used can also be adapted to study nonsymmetric electrolytes. In forthcoming work, this model will be used to analyze the behavior of such systems. A program performing this calculation is available from us. LA050555D
