Huada
Lian
^{a} and
Jian
Qin
*^{b}
^{a}Department of Materials Science and Engineering, Stanford University, Stanford, California 94305, USA
^{b}Department of Chemical Engineering, Stanford University, Stanford, California 94305, USA. E-mail: jianq@stanford.edu

Received
18th September 2017
, Accepted 27th October 2017

First published on 27th October 2017

Evaluating the interaction energy between charged dielectric spheres in close contact is difficult because interfacial charge polarization becomes increasingly stronger as separation decreases. The limit of high dielectric permittivity is amenable to theoretical treatment because the particles become equipotential upon close contact. For particles with low permittivity, the constitutive equation for the induced surface charges demanded by the boundary conditions is derived in this work, based on which the contact energy of symmetric cases is evaluated. The known result for conducting spheres is recovered as a limiting case.

## Design, System, ApplicationSuccessful design of the structural stability of molecular systems needs to be based on a thorough understanding of the inter-component interactions. This work focuses on the strong, nearly singular polarization interaction between charged dielectric spheres in close contact. The polarization interaction is caused by the surface charges accumulated at the interfaces between the particles and medium, which exhibit distinct dielectric permittivities. The interaction was known to be strong and singular for conducting spheres. We consider dielectric spheres and develop an analytical formalism that reveals a similar type of singularity and allows for calculation of polarization energy at small separations. For two spheres of equal size and charges, we evaluate the cohesive energy and show that the result includes the conducting case as a special example. |

(1) |

Despite extensive efforts for dielectric particles,^{5–7} no analogous results on contact energy have been reported. When dielectric spheres have large separation, the polarization energy can be calculated by approximating the induced polarization charges as dipoles or, as separation decreases, as a sum of multipoles. A systematic approach based on multiple-scattering formalism has been developed to calculate the polarization energy of many polarizable spheres, which carry monopolar, dipolar or multipolar charges.^{8–10} This approach expresses the total electrostatic energy of dielectric spheres carrying charges of ideal multipoles as a sum of two-body, three-body and other interactions. The two-body terms are normal Coulomb interactions. The three-body terms are interactions between two spheres mediated by one polarizable surface. The four-body terms involve two polarizable surfaces, and the higher order terms are constructed analogously. Each of these terms can be efficiently evaluated by using image lines.^{8,10} At a scaling level, the two-body terms are proportional to the inverse of average separation R = ρ^{−1/3} between spheres, where ρ is the number density of spheres. The three-body terms are smaller by a factor (a/R)^{−3}. Each higher order term is smaller by an additional factor (a/R)^{−6}.

Therefore, it is expected that the three-body terms are sufficient for dilute systems or for spheres at large separation. As average separation increases, however, higher order terms become increasingly more relevant. Upon close contact, essentially all the higher order polarization interactions are needed, which is the origin of the singular behavior identified for contact energy between conducting spheres. For dielectric spheres, since the Derjaguin approximation^{4} cannot be applied to obtain the contact charge—the dielectrics are not equipotential—the contact singularity cannot be readily identified. Preliminary numerical calculation based on multiple-scattering formalism has revealed the difficulty of resolving this contact energy: in one case, the number of higher order polarization terms increases with reduced separation;^{9} in another, greater spatial resolution is needed.^{5}

This work aims at calculating the contact energy between two charged dielectric spheres. Instead of relying on the perturbative multiple-scattering approach for electrostatic potential, the boundary conditions demanded by the induced surface charges are derived using a bispherical coordinate system. Then, the potential and energy are calculated from the induced charges for varying inter-particle separations. The expressions for elastance coefficients are derived. The contact energy for spheres with arbitrary dielectric permittivity ε_{in} is evaluated. It is demonstrated that, in the conducting limit, i.e., ε_{in} → ∞, the contact energy reduces to that of Maxwell's result (eqn (1)). The dependence on permittivity allows one to estimate the cohesive energy for dielectric dimers carrying symmetric charges.

The particle configuration is sketched in Fig. 1. The dielectric permittivity inside and outside the two spheres are ε_{in} and ε_{out}, respectively. We are interested in the variation of polarization energy with gap distance d. At a large distance, the energy is the total of self-polarization energy of two isolated spheres, , which is caused by the interaction between the charges Q_{1} and Q_{2} and the induced surface charges on their own surfaces. At a reduced distance, the normal Coulomb energy, E_{∞} = Q_{1}Q_{2}/(4πε_{out}R_{12}), is needed. In the limit of close contact, d → 0, the mutual polarization interaction becomes relevant.^{8} Evaluating the energy and electrostatic potential requires solving the boundary value problem for Poisson's equation explicitly.^{11}

The boundary condition is conveniently formulated using the bispherical coordinate system,^{12} parameterized by (σ, τ, ϕ). Here, 0 ≤ σ ≤ π is the angle spanned between two foci located at (0, 0, −b) and (0, 0, +b). The term −∞ < τ ≡ ln(d_{−}/d_{+}) < ∞ is defined based on the ratio of distance d_{−} to foci (0, 0, −b) and distance d_{+} to foci (0, 0, +b), so that τ > 0 in half space z > 0, τ < 0 in half space z < 0, and τ = 0 in the x–y plane. The azimuthal angle is defined by 0 ≤ ϕ < 2π.

The surfaces of the two spheres in Fig. 1 have a constant τ value: τ = ζ > 0 and τ = −η < 0. The coordinates of points on the surface satisfy x^{2} + y^{2} + (z − b coth τ)^{2} = (b csch τ)^{2}. The two spheres approach each other when ζ and η both approach zero. The center of sphere τ is (0, 0, b coth τ) in Cartesian and (0, 2τ, 0) in bispherical coordinate systems. The scale factors for coordinates σ, τ and ϕ are

(2) |

[2h_{σ}h_{σ1}(cosh(τ − τ_{1}) − cosγ)]^{1/2}, | (3) |

By introducing the surface charge densities ρ_{ζ} and ρ_{η} into the two spheres, the electrostatic potential ϕ(r) at position r can be written as a sum:

(4) |

Here, ϕ_{0}(r) is the potential of source charges,

(5) |

(6) |

Here, ε_{r} ≡ ε_{in}/ε_{out} is the ratio of dielectric permittivities inside and outside the particle. It is clear that ρ_{η} can be interpreted as (1 − ε_{r})n·∇ϕ(r(S_{η})) and ρ_{ζ} as (1 − ε_{r})n·∇ϕ(r(S_{ζ})).

Eqn (4) is a formal solution for Poisson's equation valid at all separations. We are particularly interested in the limit when separation d = R_{12} − a_{1} − a_{2} is smaller compared to both radii, a_{1} and a_{2}. In the bispherical coordinate system, this can be achieved by setting b → 0. At the same time, to keep a_{1} = b/sinhζ and a_{2} = b/sinhη constant, b needs to be of the same order as sinhζ ∼ ζ or sinhη ∼ η. In the limit of close contact, we have ζ = b/a_{1} and η = b/a_{2}. The distance between the two surfaces can be written as , in which . Using a as the unit of length, we find that b/a is asymptotically small and that b/a ≃ (d/a)^{1/2}.

ε_{in}E_{ς}(−) = ε_{out}E_{ς}(+), ε_{in}E_{η}(−) = ε_{out}E_{η}(+). | (7) |

The electrical field at the boundary of sphere ζ includes contributions E^{(0)}_{ςς} from the source charge in the ζ particle, E^{(0)}_{ςη} from the source charge in the η particle, E^{(s)}_{ςς} from the induced charge in the ζ particle, and E^{(s)}_{ςη} from the induced charge in the η particle. The electrical field at the boundary of sphere η can be decomposed similarly. Thus, we have E_{ς} = E^{(0)}_{ςς} + E^{(0)}_{ςη} + E^{(s)}_{ςς} + E^{(s)}_{ςη} and E_{η} = E^{(0)}_{ης} + E^{(0)}_{ηη} + E^{(s)}_{ης} + E^{(s)}_{ςς}. The self-polarization term E^{(0)}_{ςς} equals inside sphere ζ (at a_{1} − 0) and outside sphere ζ (at a_{1} + 0). The self-polarization term E^{(0)}_{ηη} equals inside sphere η and outside sphere η. They both satisfy the required boundary conditions in eqn (7) alone. On the other hand, E^{(0)}_{ςη} and E^{(s)}_{ςη} are continuous across the ζ surface, and E^{(s)}_{ςς} is discontinuous; similar arguments hold for the η surface. So the boundary conditions can also be written as:

(8) |

The subscripts “out” and “in” remind us from which side the field is evaluated.

Eqn (8) serves as the basis for solving the surface charges. The source terms E^{(0)}_{ςη} and E^{(0)}_{ης} can be trivially evaluated. The remaining induced terms can all be expressed as integrals over the surface charges. For instance, the potential at an arbitrary point r = (σ, τ, ϕ) generated by the surface charge ρ_{ζ} can be written as:

(9) |

Here, an auxiliary charge density, _{ζ} ≡ ρ_{ζ}/(coshζ − cosσ′)^{3/2}, has been introduced in order to simplify the notation. The factor cosγ ≡ cosσcosσ′ + sinσsinσ′cos(ϕ − ϕ′) is the angle cosine between unit vectors oriented along (σ, ϕ) and (σ′, ϕ′).

The straightforward differentiation applied to the potential yields the electric field along the direction −

(10) |

The term in the bracket is continuous across the normal direction of the surface since the logarithmic divergence arising from the weak singularity at the denominator does not depend on whether surface ζ is approached from the inside or the outside. So we denote it as E

It is convenient to introduce the relative dielectric discontinuity, , and divide both sides by (coshζ − cosσ)

(11) |

Apart from the dependence on free charges, the equation above is identical to eqn (31) of ref. 12.

The analogous result for particle η is obtained by exchanging η and ζ:

(12) |

Eqn (11) and (12) are the boundary conditions that will be solved to obtain the induced charge densities.

The expansion coefficients are computed by projections, e.g.,

By substituting the above mode expansion coefficients to eqn (11) and (12), and using the expansions and , we obtain the following algebraic relation for the coefficients from the boundary conditions on sphere ζ,

(13) |

Here, the coefficients are defined by g_{} ≡ e^{−(+1/2)(ζ+η)}, , and An analogous expression can be found from the boundary conditions on sphere η. For the uniaxial case of interest here, only the m = 0 terms survive. So we drop the explicit reference to index m and adopt the convention that all the coefficients are evaluated at m = 0. Then, the constitutive equation can be written as:

(14) |

The coefficients C_{′} and I_{′} are reduced to and . The corresponding constitutive equation for particle η is obtained by exchanging A and B and by replacing U with V:

(15) |

In the conducting limit, ε → 1, the coefficient relation can be further simplified to:

In the far field regime, ζ → ∞ and η → ∞, whereas at the same time ζ/sinhζ → a_{1} and η/sinhη → a_{2}. The cross terms all vanish because the coupling g_{} decays to zero.

In the near field regime, ζ → 0, η → 0, b/ζ → a_{1} and b/η → a_{2}; the small parameters ζ, η, and b are of the same order of magnitude. In this limit, the coefficient I_{′} is of the order . To solve the unknown A_{} and B_{}, the mode coefficients U_{} and V_{} from the fields produced by free charges are needed. The numerical results for the conducting limit obtained from such an approach are illustrated in Fig. 2 and compared to the prediction of ref. 3.

Fig. 2 Comparison of surface charge densities on two spheres with equal radius Q_{ζ} = 2Q_{η} and separation distance R/a = 3 with known results for conducting spheres.^{3} The thick curves are the total densities and the dashed lines denote the free charge density, both of which have units of Q_{ζ}/4πa^{2}. The bispherical coordinate σ is related to θ by . The upper sign corresponds to the right half space and the lower sign corresponds to the left half space. |

(16) |

(17) |

The coefficients V

Setting ζ = η and A

(18) |

Once the coefficients A_{} and B_{} have been solved, the polarization energy can be calculated by using:

(19) |

By writing coefficients A_{} and B_{} explicitly in terms of charges Q_{ζ} and Q_{η}, the elastance coefficients can be obtained from the summation of entries in array T.

Then, the contact energy for dielectric spheres with symmetric charges can be calculated as follows. First, the interaction energy for a given dielectric permittivity, which includes the pairwise Coulomb energy and the mutual polarization energy, with the self-polarization energy omitted, is calculated at several separations, close to the contact value 2, which is the last term given by eqn (18) and (19). Then, the separation dependence is fitted and extrapolated to the close-contact value. This procedure is demonstrated in Fig. 3, for two spheres with identical radii and charges Q_{1} = Q_{2} = 1. Two families of data points are calculated from the cases with strong polarization, i.e., ε → ±1. In both cases, the fitting is based on quadratic polynomials, and a smooth extrapolation towards R/a = 2 can be made.

The extrapolated values for the contact energy are plotted against the values of ε in Fig. 4. Three limits are worth noting. (1) At ε = 0, the dielectric interface and polarization contribution vanishes. So the total energy is dominated by the normal Coulomb interaction, which is 0.5 on this normalized scale. (2) At ε = 1, the behavior approaches that of conducting spheres, as demonstrated in Fig. 2. The value for the contact energy obtained from the Maxwell expression (eqn (1)), which is 1/ln(2) − 1 ≃ 0.44, is recovered. (3) At ε = −1, the behavior of weak dielectrics, with ε_{in} ≪ ε_{out}, is recovered. This is the typical case for an aqueous solution of colloidal particles. In this limit, the polarization energy becomes independent of the internal dielectric permittivity, and is determined by the medium permittivity. From the known result^{13} for the electrostatic potential ϕ_{w} produced by a point charge outside dielectric particles with a low value of permittivity, the polarization energy can be estimated to be (1/2)(Q_{1}ϕ_{1} + Q_{2}ϕ_{2}), in which ϕ_{1} and ϕ_{2} are the electrostatic potential ϕ_{w} produced by charges Q_{1} and Q_{2}, and evaluated at the location of Q_{1} and Q_{2}, which is a correction to the self-polarization energy. Taking the limit R/a = 2, it is straightforward to show that Q_{1}ϕ_{1} = Q_{2}ϕ_{2} = 1/3 + ln(3/4 = 0.05) on the normalized scale. Adding this value to the Coulomb energy gives an estimate of 0.55, a value very close to that obtained from the numerical approach.

The formalism developed, although only applied to spherical geometries, can be readily modified to treat the cases involving flat interfaces, by setting ζ = 0 or η = 0, or to charged cylinders, by modifying the metric h_{ϕ} related to the azimuthal degree of freedom. To adapt the formalism to arbitrary source charge distribution, the same type of multipole expansion can be applied, and eqn (13) is the most general form of constitutive equation. This formalism, however, only works for continuum media and breaks down when the inter-particle separation becomes comparable to the average separation between discrete charges on the surface. The dielectric permittivity of particles need not be identical, as is clearly shown by the expression for the boundary conditions. Further, the surface polarization can be induced by externally applied and oscillating fields, as demonstrated in ref. 12 in the treatment of van der Waals interactions.

The numerical data was only provided for the symmetric case. The contact energy for the asymmetric case such as Q_{1} = 1 and Q_{2} = −1 cannot be obtained by applying the extrapolation in Fig. 3, because the contact charge accumulates rapidly at small separations, which leads to a nearly singular behavior for the separation dependence of polarization energy, as demonstrated in ref. 4. For the conducting case, this leads to a dependence. For the dielectric cases, this dependence is weakened; no analogous expression based on the Derjaguin approximation can be readily obtained, since the dielectrics are not equipotential. Yet, in the absence of an analytical expression for the strength of this singularity, any numerical extrapolation will fail at sufficiently small separation. We thus leave the clarification of this issue to future work.

- J. Maxwell, A Treatise on Electricity and Magnetism, Clarendon Press, 1873 Search PubMed .
- A. Russell, Proc. Phys. Soc., London, 1922, 35, 10 CrossRef .
- J. Lekner, Proc. R. Soc. A, 2012, 468, 2829 CrossRef .
- J. Qin, N. W. Krapf and T. A. Witten, Phys. Rev. E, 2016, 93, 022603 CrossRef PubMed .
- K. Barros, D. Sinkovits and E. Luijten, J. Chem. Phys., 2014, 140, 64903 CrossRef PubMed .
- K. Barros and E. Luijten, Phys. Rev. Lett., 2014, 113, 1–5 CrossRef PubMed .
- Z. Gan, S. Jiang, E. Luijten and Z. Xu, SIAM J. Sci. Comput., 2016, 38, 375–395 CrossRef .
- J. Qin, J. J. de Pablo and K. F. Freed, J. Chem. Phys., 2016, 145, 124903 CrossRef PubMed .
- J. Qin, J. Li, V. Lee, H. Jaeger, J. J. de Pablo and K. F. Freed, J. Colloid Interface Sci., 2016, 469, 237–241 CrossRef CAS PubMed .
- K. Gustafson, G. Xu, K. F. Freed and J. Qin, J. Chem. Phys., 2017, 147, 064908 CrossRef PubMed .
- W. R. Smythe, Static and Dynamic Electricity, McGraw-Hill Book Company, Inc., 1950 Search PubMed .
- D. Mitchell and B. Ninham, J. Chem. Phys., 1972, 56, 1117–1126 CrossRef CAS .
- A. P. dos Santos, A. Bakhshandeh and Y. Levin, J. Chem. Phys., 2011, 135, 044124 CrossRef PubMed .

This journal is © The Royal Society of Chemistry 2018 |