Manipulating Interactions between Dielectric Particles with Electric Fields: A General Electrostatic Many-Body Framework

We derive a rigorous analytical formalism and propose a numerical method for the quantitative evaluation of the electrostatic interactions between dielectric particles in an external electric field. This formalism also allows for inhomogeneous charge distributions, and, in particular, for the presence of pointlike charges on the particle surface. The theory is based on a boundary integral equation framework and yields analytical expressions for the interaction energy and net forces that can be computed in linear scaling cost, with respect to the number of interacting particles. We include numerical results that validate the proposed method and show the limitations of the fixed dipole approximation at small separation between interacting particles. The proposed method is also applied to study the stability and melting of ionic colloidal crystals in an external electric field.

To begin with, we assume that the external harmonic potential we consider satifies Φ ext ∈ H 1 loc (R 3 ) with the associated external electric field E ext := −∇Φ ext ∈ L 2 loc (R 3 ), where L 2 loc (R 3 ) and H 1 loc (R 3 ) denote the spaces of locally square integrable functions and locally square integrable functions with locally square integrable first derivatives, respectively. Next, we emphasise that, as is common in the mathematical literature, the solution to the PDE (4), i.e., the perturbed electrostatic potential Φ, is typically understood as an element of the space H 1 (Ω − ) ∪ H 1 (Ω + ) and is therefore not, in general, continuous. Strictly speaking therefore, the transmission conditions in Equation (4) must be understood in the sense of so-called Dirichlet and Neumann traces in the Sobolev spaces H 1 2 (∂Ω) and H − 1 2 (∂Ω) respectively. A detailed description of trace operators and fractional Sobolev spaces is beyond the scope of this article and can, for instance, be found in. 1 Concerning the mapping properties of the single layer potential and boundary operators, it can be shown that for any s ∈ R, the mapping S extends as a bounded linear map from the Sobolev space H s (∂Ω) to H s+ 3 /2 loc (R 3 ) and the operator V extends as an invertible, bounded linear map from H s (∂Ω) to H s+1 (∂Ω) (see, e.g., 1 for a concise exposition on Sobolev spaces and for precise definitions and properties of the single layer potential). "Local" versions of the single layer potential and boundary operators which we have used frequently in this article are formally defined as follows: For each i ∈ {1, . . . , N }, we have In addition, we have used extensively in this article, the so-called Dirichlet-to-Neumann map, denoted DtN. Mathematically, the map DtN : H s (∂Ω) → H s−1 (∂Ω), s ∈ R is defined as follows: Given some boundary function λ ∈ H s (∂Ω), let u λ denote the harmonic extension of λ in Ω − .
Then DtNλ ∈ H s−1 (∂Ω) is the normal derivative (more precisely, the Neumann trace) of u λ on the boundary ∂Ω. We emphasise that in contrast to the single layer potential and boundary operator, the DtN map is a purely local operator, i.e., for any λ ∈ H s (∂Ω), DtNλ| ∂Ω i depends only on λ| ∂Ω i .
Concerning the regularity of solutions to the BIE (6), we recall from Equation (2)   The set of spherical harmonics is dense in L 2 (S 2 ) and is therefore well-suited for the choice of basis functions in the Galerkin discretisation of BIE (9).

Definition (Approximation Spaces)
Let max ∈ N be a discretisation parameter. First, on each sphere ∂Ω i , i = 1, . . . , N we define a local approximation space W max (∂Ω i ) as where we introduced for notational convenience the basis functions Y i m : Next, we define the global approximation space W max as

Mathematical Proofs of Theorems 2.1 and 2.2
In this section we provide proofs of Theorems 2.1 and 2.2 from Section 2.4. For technical reasons, it is useful to begin with the proof of Theorem 2.2. This result shows that the definition of the interaction energy that we have provided in this article using quantities of interest from the integral equation (6) is consistent with the electric field-based definition of the interaction energy as derived directly from the PDEs (3) and (4). Throughout this section, we will use the notation and setting introduced in Sections 2.1, 2.2 and 2.4.

Proof of Theorem 2.2
Let j ∈ {1, . . . , N } and let B r be an open ball large enough so that Ω − ⊂ B r . We begin by defining precisely E j j , i.e., the electric field produced only due to the sphere ∂Ω j in the absence of both the external field E ext as well as the other spheres. Maxwell's equations imply that E j j := −∇Φ j j where the self-potential Φ j j satisfies the PDE (c.f., Equation (4)) where we remind the reader that σ s, j := σ s | ∂Ω j and σ p, j := σ p | ∂Ω j .
Next, to aid the subsequent exposition, we define the auxiliary quantity We may now use simple algebra and the fact that Φ tot = Φ+Φ ext (see Section 2.1) to deduce Next, we recall from the PDEs (3) and (4) that Φ is harmonic on Ω − ∪ Ω + , Φ ext is harmonic on R 3 , and Φ j j is harmonic on Ω j ∪ (R 3 \ Ω j ). Therefore we can appeal to Green's first identity to simplify the above integrals as Recalling the interface conditions from the PDEs (4) and (1), we can further simplify several of these integral as where we remind the reader that σ ext = −(κ − κ 0 )∂ n Φ ext . Using the fact that λ, λ ext and λ j j are the restrictions on the spheres of the potentials Φ, Φ ext , and Φ j j respectively, we can deduce Comparing this final expression with Equation (2) allows us to deduce the required result (28).
Next, we will prove Theorem 2.1 which shows that Definition (24) of the approximate electrostatic forces is consistent with the usual notion in the chemistry literature of the forces as the negative sphere-centered gradients of the electrostatic interaction energy. In order to present a concise and well-structured proof, we will first prove two lemmas.
Proof: Let i ∈ {1, . . . , N } be fixed. A simple application of the product rule yields that Using the fact that both σ s and σ p are independent of changes in the locations {x i } N i =1 of the sphere centres locations, 2 we further obtain that Finally, it is straightforward to see that in fact so that we obtain the expression Consequently, it remains to compute the sphere-centred gradient of λ max . This is a slightly technical task so to aid the subsequent exposition, we first introduce some additional notation.
Equipped with the notation introduced above, we now take the gradient on both sides of Equation (4). Using the product rule together with the fact that the Dirichlet-to-Neumann map is independent of changes in the locations {x i } N i =1 of the sphere centers, we obtain that or equivalently, after collecting terms Next, recalling that ν max satisfies Equation (23), it is easy to deduce that with indices i ∈ {1, . . . , N }, ∈ {0, . . . , max } and |m| ≤ . We therefore conclude from Equation (5) that Recalling now the last term on the right-hand side of Equation (3), we deduce that Next, a direct calculation and comparison with the Galerkin discretisation (18) reveals that Using the definition of ν max as given by Equation (22), we obtain that Finally, a direct but tedious computation can be used to show that 2,3 Combining therefore the developments (6) and (7) Proof: Recall the notation λ ext := Φ ext | ∂Ω ∈ H Since Φ ext is harmonic in R 3 and therefore in particular on Ω i , it follows that we can write

Proof of Theorem 2.1
We are now ready to state the proof of Theorem 2.1. Before proceeding to the proof, let us simply remark that the relation (26) in Theorem 2.1 remains true if exact quantities are considered, i.e., if the force defined by (24) is built upon the exact induced charge ν being solution to the BIE (7) and where the energy corresponds to E int as defined by (27).
Let i ∈ {1, . . . , N } be fixed. By the definition of the approximate electrostatic interaction energy, we have .
We now simplify each of the terms (I), (II), (III), and (IV). First, we observe that the self energy term (IV) is defined entirely through functions that are independent of changes in the location of the center x i of the sphere ∂Ω i , even in the case j = i . This can be seen by noticing that σ s,i + σ p,i , λ ii max L 2 (∂Ω i ) remains constant as one displaces x i by any translation. Consequently, we obtain that (IV) ≡ 0.
The term (I) can be simplified using Lemmas 2.1 and 2.2 as .
(Using Lemma 2.2) (9) Next, we simplify the term (IB). Indeed, a direct calculation shows that where the second line follows from a similar calculation as done to obtain Equation (7).
In order to simplify the term (II), we again recall that the free charges σ s , σ p are independent of changes in the location of the center x i of the sphere ∂Ω i . Consequently, we obtain Therefore, using a calculation similar to the one used to obtain Equation (10), we deduce that where σ max s and σ max p are the best approximations in W max of σ s and σ p respectively.
Next, we attempt to simplify the term (III). A simple application of the product rule together with Lemma 2.2 yields that Once again, a direct calculation of the form used to obtain Equation (10) allows us to conclude that (III) = − σ ext , ∇ x i λ max ext L 2 (∂Ω i ) = σ max ext , E ext L 2 (∂Ω i ) .