Electrostatic Force between Two Dielectric Particles in Electrorheological Fluids: beyond Spherical Particles

In an effort to increase the shear yield stress of dielectric electrorheological fluids, we focus on the electrostatic force of different forms of particles in a dielectric polarization model. By solving Laplace’s equation and applying the multiple image method and the finite element method, the analytical and numerical solutions of the electrostatic force of a two-sphere structure have been studied. The results suggest that when the dielectric mismatch factor is large and when the positions of the two spheres are nearly in contact with each other, most of the analytical solutions either overor underestimate the force. Additionally, the structure of particles beyond the spherical form is considered. Three example cases are studied to shed light on how different geometries of particles may affect the electrostatic force, thereby influencing the shear yield stress of the fluid.


Introduction
Electrorheological (ER) fluids are regarded as smart fluids, whose resistance to flow (from a liquid to a solid gel, and back) can be quickly (with response times on the order of milliseconds [1]) and dramatically altered by an applied electric field. Several applications have been proposed, including for clutches, brakes, shock absorbers, hydraulic valves, polishing, and artificial muscles [2][3][4][5][6][7][8][9]. In this work, we consider dielectric electrorheological (DER) fluids, which are made of micrometer-sized dielectric particles in an electrically insulating liquid. Compared to that of magnetorheological (MR) fluids, the yield stress of DER fluids is much lower. The shear yield stress of a general DER material is typically of several kPa when subjected to an electric field of 1 MATEC Web of Conferences [10], while the typical shear yield stress of MR fluid could vary from 2-3 kPa with no magnetic field to 50-100 kPa with an applied magnetic field of 1 240 A m   [11]. The DER fluid has therefore found limited use in the mechanical industry. Equation Chapter 1 Section 1 Our aim is to increase its yield stress, which, at the micro level, means to increase the interparticle force (mainly the electrostatic force, and other forces are discussed in Refs. [12][13][14]) according to the dielectric polarization model (DPM), where the electrical response is governed by linear electrostatics. Recently, dielectric effects have gained increased attention in other fields. For instance, E. Luijten et al. investigated polarization effects in charge-stabilized dielectric colloidal suspensions. They introduced an efficient technique to resolve polarization charges in dynamical dielectric geometries and demonstrated that dielectric effects qualitatively alter the predicted self-assembled structures [15,16].
In this work, we ignore temperature effects and entropic effects, i.e., the first and the third term in the expression of the free energy W [17,18] where N is the total number of dielectric particles in the DER fluid; B k is the Boltzmann constant; T is the thermodynamic temperature, which is usually the room temperature; and are, respectively, the local electric field and electric displacement with and without the dielectric particles. S is the configuration entropy, which can be well-approximated by where  is the volume available to the many-particle system and the thermal de Broglie wavelength is  , which is very small compared to the size of the particle. Generally, the particle diameters range from 0.1 to 100 μm [13]. In fact, the electrostatic energy usually dominates the thermal and entropic effects [19], i.e., the electrostatic energy of a particle is much In the dipole DPM, the interparticle force depends on the electric dipole moment. Thus, it has inspired us to focus on the hierarchical structure of the particles (stick, lamella, needle-like, etc.) assembled at different scales. Charges polarized in an electric field tend to concentrate at certain positions of the particles surface so that the charges can be used efficiently to increase the shear yield stress. By optimizing the hierarchical structure combination, a high performance electrorheological fluid can be efficiently fabricated. To the best of our knowledge, work of such nature has not yet been performed. Various investigations, such as of particle patterns, arrangement structures, surface/core dielectric properties and interface structures, could be conducted to characterize the fluids. With a better characterization of electrorheological fluids, they may find wider applications in various domains, such as aviation, noise control, and automobile engineering [20].
In previous papers, the separation of variables and multiple image methods are used to determine the analytical formula of the electrostatic force between two spherical particles in the DPM [21][22][23]. In this paper, by proposing to investigate the problem via the finite element method (FEM) approach, we are no longer confined by the regular shape of the particles and the other parameters that we study.
In the first part, we considered the rudimentary case of the model, i.e., the point-dipole (PD) model. Then, in the following section, we consider models beyond the PD case. Next, we proceed to consider three special cases beyond the basic form as a demonstration, which provides several ideas to increase the shear yield stress in the experiment. Eq uatio n Se c tion (Ne xt) 2 Solution in the point-dipole model Based on Sect. 4.4 of Ref. [24], we derive the formula of the electrostatic force in the point-dipole model in a completely rigorous way. As shown in Figure 1, we consider two uncharged dielectric spheres and , of radius a , having identical dielectric constants 0 p , fixed in place in a liquid, with a dielectric constant 0 f , far from the electrodes so that it appears to be in a uniform electric field 0Ê  E z . We assume and as two By separation of variables and azimuthal symmetry, the general solution, called the harmonic function, can be written in the form of Legendre polynomials, and the solution of this problem can be found using the method of undetermined coefficients. We note the dielectric mismatch factor / cos , where the Clausius-Mossotti factor is 1 2 Then, we determine the relation between the potential and the dipole moment Eventually, we obtain the expression of the force F that exerts on and the dimensionless force E is the scalar quantity of the external uniform electric field, and 2 b  is commonly observed for low to moderate field strengths. At large field strengths, b can decrease below 2 over a range in / p f   (3-500) [13], which can be explained by the nonlinear conduction [25].

Solution beyond the point-dipole model
Beyond the dipole approximation, to solve the problem of two dielectric spheres completely, the exact analytical solution can be calculated in bispherical coordinates, which is a coordinate system in which the separation of variables is possible [21].E q uatio n Se c tio n (Ne xt) A much easier method is to apply the multiple image method to attain the multipole-induced-dipole (MID) model. The nomenclature is similar to that of the pointdipole model. In this case, however, we consider and with different radii a and b and different Clausius-Mossotti factors a  and b  , respectively. Note d , the distance between the centers of two spheres. Additionally, we no longer assume that the two particles are independent of each other. In Ref. [ Therefore, the expression of the dipole moment can be written in a unified form p are scalar quantities of the dipole moments in the point-dipole model 3 3 The hyperbolic angle parameter  satisfies The force between the spheres is given by The direction is along the line joining the centers of the spheres. The positive value represents repulsion, and the negative value represents attraction.

MATEC Web of Conferences
When n=1, it is the point-dipole case , which is the same as Eq. (2.10) with 0 Considering the n=1 and n=2 terms, we obtain the expression in the dipole-induced-dipole (DID) model.
The ratio  of the forces between the DID and the PD models can be written as More terms can be considered, but the calculation becomes much more complex. The multipole-induceddipole model includes all the terms. After extensive calculations, we obtain the expression of the electrostatic force in the MID model Inspired by the PD model, the feature parameter  could be extended as and MID F can be generalized to the projection of the electrostatic force between two spheres in the direction of  r . Several limitations of Eq. (3.8) warrant comment. As it is an approximation in the MID model, more details will be discussed in the following part.

Preliminary validation works for the solution beyond the spherical form
So far, by confining ourselves to the hypothesis of an ideal spherical particle, there is a difference between the theoretical and experimental results: what about other forms of particle ? By varying the boundary Ω  , the case is transferred to a Dirichlet problem for Laplace's equation, whose analytical solution is difficult to calculate by adding correction terms as in the case of a two-sphere problem. An alternative approach is to use a numerical solution; the finite element method (FEM) is investigated for its feasibility in approximating the solution for problems of such nature. Many researchers have considered cases beyond a typical two-sphere model. For instance, Wen et al. considered the coated sphere structure, and shown that the use of an effective dielectric constant concept based on the first-principles method can yield quantitative predictions [19]. Tao et al. adopted an infinite chain of dielectric particles [22], while Ahn et al. investigated nsphere clusters [26]. J. G. Cao et al. considered the microstructure evolution from a random structure that consists of spherical particles to chains and then to stable lamellar patterns [27]. However, to the best of the authors' knowledge, there are few studies that consider models beyond that of spherical particles. Hence, before studying the non-spherical case, we start from the abovementioned spherical model, the simplest case, to verify its numerical solution accuracy when compared with the analytical solution.Equation :

Model building
To ensure the accuracy of the numerical solution, axisymmetric structures are used in the model. As shown in Figure 3, to create an environment of a uniform

Method
The center of is fixed along the longitudinal axis, and is fixed at the origin. According to the principle of virtual work and Newton's third law, the scalar quantity of the force F that exerts on can be written as 3 3 where i n denotes the normal to the boundary i  , whose direction is from i  to 0  , and the dielectric mismatch Therefore, to determine the force, we simply calculate the distribution of the potential, where we apply a volume-discretization method, namely, the finite element method (FEM), due to its adaptability to complex geometries and its ease of handling discontinuous gradients of a variable. The basic steps of the FEM are shown in Figure 4. Based on the three-dimensional FEM in Sect. 5.3 of Ref. [28], we calculate the distribution of the potential, and therefore the electrostatic energy and the force.
To solve the problem beyond the regular form of particles and pay more attention to the post-processing, we use an unstructured grid for the FEM mesh. Tetrahedrons are employed as the cell shapes. The Bowyer-Watson algorithm [29,30] is used to compute the Delaunay triangulations. The parameters of the tetrahedral mesh generated using the Delaunay tessellation method are described below. The maximum element growth rate is 1.3. The curvature factor is 0.2. The average element quality is 0.7704. The subdivision parameter on the x,y,z-direction scale is mesh  for 1,2  and 1 for 0  . Figure 5 shows the FEM unstructured mesh used for the numerical calculation.
Solving the sparse matrix equation, the conjugate gradient method and the smoothed aggregation algebraic multigrid (SA-AMG) method are used [31], where the convergence criterion (the relative error after each iteration compared to that of the previous one) is set to be 3 1 10   to ensure the result converges well.

Electrostatic energy
The electrostatic energy of a single spherical particle in a uniform electric field can be calculated according to   In the simple-sphere case,  is defined to be 2. As shown in Figure 6, for a series of test cases, the maximum relative error is 4 5.54 10   with 9.5

 
, which suggests the validity of the simulation.
As for the case of two spherical particles, the calculation becomes complex. In Figure 7, the  Figure 8 presents the interaction energy of one particle W for the double-sphere case as a function of the spacing d at  =1, 2, 3, 5, and 10. Using Eq. (3.3), we obtain the curves shown in Figure 9.

Electrostatic force for the double-sphere case
We compare the numerical simulation result with the

Example 1: Earth-and-moon structure
In the first example, we consider two spheres with different radii. The sum of the volumes of the spheres is fixed and equals the volume of two spheres with the same radii a . Note  is the ratio of the two radii. The radii can then be written as F decreases sharply, which reveals that to obtain the biggest force, the particles should be as identical as possible. The difference in size will result in great changes in the electrostatic force. Now, we consider the sub-structures that are uniformly distributed across a sphere. To create the model, we consider N number of points distributed uniformly on a sphere, of which several approximate solutions have been proposed [32][33][34]. Consider 20 N  . Since the regular dodecahedron has 20 vertices, its circumscribed sphere, the one that touches the regular dodecahedron at all vertices, can be considered as the base sphere. As shown in Figure 11, the Cartesian coordinates of the 20 vertices of a regular dodecahedron where the radius of its circumscribed sphere is 3 can be defined as follows: In the following examples, we let

Example 2: Litchi structure
Consider two dielectric basic spheres and of identical radii a , with a distance d between them, on which there are 20 small spheres of radii cp a  whose centers are situated on each of the 20 vertices of a regular dodecahedron as defined above for each basic sphere. We perform the rotation transformation for one of the basic sphere; see Figure 12.    needs to be studied. Figure 13 presents the electrostatic energy 2  and approximately 2.08 for the two-sphere and litchi structures, respectively. It is noted that in both cases, when the two particles are increasingly near to each other, there is a significant increase in the force F , which tends to a finite value. The force at the limit position of the litchi structure is 129% that of the two-sphere case, which indicates that the litchi structure can help increase the electrostatic force.    In this example, by replacing small spheres with spikes, we create a sea urchin structure, which is illustrated in Figure 15. , respectively. As seen in Figure 16, unlike the litchi structure, the sea urchin structure contributes little to the electrostatic force and reduces the force: the force at the limit position of the sea urchin structure is 69% that of the two-sphere case. This is mainly because although this structure seems to be more closely integrated, in order to meet the assumption that no deformation occurs when the two particles are in contact with each other, the relative distance / d a   is quite limited (here  should be larger than 3.35), which makes the force smaller than that in the two-sphere case where 2   .

Conclusions
In this paper, we derive formulae of the electrostatic force between two spherical dielectric particles in the PD, DID, and MID models whose results are Eqs. (2.10), (3.7), and (3.8), respectively. The study of the convergence suggests that in the case of two identical spheres in MID models, the simplified version of Eq. (3.8), i.e., Eq. (4.2), is inapplicable when either of the following occurs (i) 50 . The numerical simulation based on the finite element method (FEM) is then necessary, by which we conclude the following: • In the n-sphere case, the electrostatic energy per sphere at the limit position is a function of  . When n is fixed, with the increase of  , there exists a maximum energy.
• Compared to the case of having different particle sizes, the force between identical particles has the greatest force, which suggests making the particles as identical as possible in the experiment to increase the shear yield stress.
• Sub-structures uniformly grown on spherical particles can increase the force, but the effect depends on the balance of  and the form of the sub-structures. A litchi structure is preferred.
Compared to the analytical methods, the FEM is a more efficient method to calculate the electrostatic force between particles, regardless of its form and relative permittivity. Thanks to this method, we could conceive a