2. Modelling approach

We consider a rigid spherical particle of radius $a$ attached to the interface ${\it\Sigma}_{12}$ between two immiscible fluids 1 and 2, as depicted in figure 1. The circular TCL is pinned to the particle surface ${\it\Sigma}_{1p}\cup {\it\Sigma}_{2p}$ . Since the particle radius is the characteristic length scale of the flow field, we may neglect undulations of the TCL, which are typically tens of nanometres in size (Stamou et al. Reference Stamou, Duschl and Johannsmann2000; Park & Furst Reference Park and Furst2011), when dealing with microparticles. For simplicity, we assume the Bond number $\mathit{Bo}:={\rm\Delta}{\it\rho}ga^{2}/{\it\sigma}_{12}$ , containing the maximum density difference ${\rm\Delta}{\it\rho}$ between the three phases, the gravitational acceleration $g$ as well as the fluid–fluid interfacial tension ${\it\sigma}_{12}$ , to be zero. The assumption $\mathit{Bo}=0$ leads to the equilibrium configuration shown in figure 1(a), characterized by a planar interface and an equilibrium contact angle ${\it\Theta}$ . Note that there is no conceptual difficulty in including a non-zero Bond number within the asymptotic model described below, since Lee et al. (Reference Lee, Chadwick and Leal1979) have demonstrated that the perturbation quantity may be a combination of the Bond number and the capillary number $\mathit{Ca}:={\it\mu}_{1}U/{\it\sigma}_{12}$ used in the present work. In the definition of $\mathit{Ca}$ , ${\it\mu}_{1}$ denotes the dynamic viscosity of fluid 1, while $U$ denotes the absolute value of the undisturbed flow velocity far from the particle when the situation is analysed in the particle’s rest frame (see figure 1 b). The force $\boldsymbol{F}$ driving the particle is supposed to act at the particle’s centre of mass and parallel to the plane given by the fluid interface in equilibrium. As discussed above, the plane containing the circular TCL has rotated by the trim angle ${\it\alpha}$ once the particle motion has reached a steady state (see figure 1 b). In the regime of particle Reynolds numbers $\mathit{Re}:=Ua/{\it\nu}_{1}$ much smaller than unity, where ${\it\nu}_{1}$ denotes the kinematic viscosity of fluid 1, incompressible Newtonian fluid flow is governed by the Stokes equation

and the continuity equation

where $\boldsymbol{u}$ is the velocity and $p$ is the pressure. In (2.1) and (2.2), indices 1 and 2 denoting the respective fluid phases have been omitted. On the particle surface, the flow field obeys the no-slip condition

while approaching the homogeneous velocity field $\boldsymbol{U}$ far from the particle. Across the fluid–fluid interface, which is assumed pure and therefore non-viscous (zero Boussinesq number), implying absence of Marangoni stresses, the velocity is continuous,

and we have the interfacial momentum balance equation

where $\unicode[STIX]{x1D644}$ and  $\unicode[STIX]{x1D64E}$ are the unit tensor and the rate-of-strain tensor, respectively, $\boldsymbol{n}_{12}$ is the fluid interface normal, and $\boldsymbol{{\rm\nabla}}^{s}:=(\unicode[STIX]{x1D644}-\boldsymbol{n}_{12}\otimes \boldsymbol{n}_{12})\boldsymbol{\cdot }\boldsymbol{{\rm\nabla}}$ . In (2.4) and (2.5), we have used the jump bracket notation

Note that the last part of (2.6) implies that the orientation of the normal vector with respect to the fluid phases has been specified in such a way that it points from fluid 2 into fluid 1. Concerning the particle having the outer normal vector $\boldsymbol{n}_{p}$ , steady state is reached when the balances of linear momentum

and angular momentum with respect to the particle’s centre of mass

are fulfilled. Note that the inertial contributions to the balances (2.7) and (2.8) have been neglected since the Reynolds number is much smaller than unity. In (2.7) and (2.8), $\unicode[STIX]{x1D64F}:=-p\unicode[STIX]{x1D644}+2{\it\mu}\unicode[STIX]{x1D64E}$ is the stress tensor, $\boldsymbol{r}_{\mathit{TCL}}$ is the position vector of a point on the TCL, parametrized by $s$ such that the vector $\partial \boldsymbol{r}_{\mathit{TCL}}/\partial s\times \boldsymbol{n}_{12}$ points away from the particle, and $\boldsymbol{r}_{p}$ is the position vector of a point on the particle surface. Recall that by definition the driving force $\boldsymbol{F}$ itself does not produce a torque on the particle as it acts on its centre of mass. We concentrate on the case of a low viscosity ratio between the two fluids, implying ${\it\mu}_{2}/{\it\mu}_{1}\rightarrow 0$ , because in that case the effects studied in this work are expected to be most pronounced, and also because explicit analytical results can be found. Upon introducing the scaling

and incorporating the definition of the jump bracket (2.6) and of the capillary number $\mathit{Ca}$ , we may write the interfacial momentum balance (2.5) as

where the stress contribution due to fluid 2 is omitted according to the assumption of low viscosity ratio, ${\it\mu}_{2}/{\it\mu}_{1}\rightarrow 0$ . The presence of the factor $\mathit{Ca}^{-1}$ on the right-hand side of (2.10) allows for a perturbation approach in which zeroth-order stresses cause a first-order deformation of the fluid interface, and correspondingly for higher orders. In the present work, we focus on the first-order deformation. An extension to higher orders is formally straightforward, but is not expected to yield much additional insight. For any quantity $f$ , the perturbation series for small capillary numbers is given by

Accordingly, the interfacial stress balance (2.10) becomes

where the planarity of the zeroth-order fluid interface, $\tilde{\boldsymbol{{\rm\nabla}}}^{s}\boldsymbol{\cdot }\boldsymbol{n}_{12}^{(0)}=0$ , has been used.

3. Flow field and stress distribution

As we are interested in the first-order curvature of the fluid interface represented by $\tilde{\boldsymbol{{\rm\nabla}}}^{s}\boldsymbol{\cdot }\boldsymbol{n}_{12}^{(1)}$ , it is obvious from (2.12) that the flow field of order $\mathit{Ca}^{(0)}$ is to be described. Physically, the case $\mathit{Ca}=0$ implies that the interfacial tension is large enough to counterbalance any jump in the normal stresses across the fluid interface without noticeable interfacial deformation. Moreover, because of the low viscosity ratio assumption, the stress vector at the fluid–fluid interface is perpendicular to the latter, similar to a symmetry boundary. The flow scenario within the half-space occupied by fluid 1 at $\mathit{Ca}=0$ and ${\it\mu}_{2}/{\it\mu}_{1}=0$ can thus be imitated by notionally removing the fluid interface and replacing fluid 2 by fluid 1. At the same time, the part of the rigid particle originally wetted by fluid 1 must be mirrored with respect to the planar fluid interface, making the latter a symmetry plane for the flow problem. Figure 2 provides a sketch of the two setups which are equivalent with respect to the flow of fluid 1. The procedure just described allows for a convenient treatment of the flow problem based on Stokes flow around a body moving parallel to its symmetry plane (usually termed asymmetric flow as opposed to axisymmetric flow). In the case of a spherical particle with an equilibrium contact angle ${\it\Theta}$ (see figure 1 a) of $90^{\circ }$ , the flow field of order $\mathit{Ca}^{0}$ is given by the well-known Stokes flow problem past a sphere (Happel & Brenner Reference Happel and Brenner1983).

Figure 2.

Equivalent problems for describing the flow field within fluid 1 at zero capillary number and zero viscosity ratio. (a) Original scenario including the planar fluid–fluid interface and the spherical particle. (b) Transformed scenario with symmetry plane (dash-dotted) and two fused equal spheres; the control volume (dashed) utilized in the calculation of the torque is also shown.

Figure 3 defines the spherical coordinate system used in the subsequent asymptotic analysis of the flow field. Herein, the $x$ -axis points downwards into the region originally occupied by fluid 1. The depicted configuration corresponds to ${\it\varepsilon}>0$ and therefore ${\it\Theta}<90^{\circ }$ , as shown in figure 1(a). Clearly, we have the relation

The undisturbed fluid infinitely remote from the particle flows into the negative $z$ direction when $U$ is positive, that is $\boldsymbol{U}=-U\boldsymbol{e}_{z}$ . A formal solution to the Stokes flow problem around two fused spheres with arbitrary ${\it\varepsilon}$ was developed by Zabarankin (Reference Zabarankin2007). However, that solution is not suited to the present work since we intend to derive explicit analytical expressions capturing the main effects, while evaluating the results by Zabarankin (Reference Zabarankin2007) requires the solution of a Fredholm integral equation. Therefore, we choose an asymptotic approach developed by Brenner (Reference Brenner1964) (see also Happel & Brenner Reference Happel and Brenner1983). The method is based on a representation of the particle shape (in our case, the fused spheres) as a power series in a geometric perturbation quantity (in our case, ${\it\varepsilon}$ ) by means of surface spherical harmonics. In this way, the no-slip boundary condition (2.3) on the generally non-spherical particle surface may be approximated by a more complicated velocity boundary condition on a spherical particle surface. Upon fitting the general solution of the Stokes equation (2.1), as given by Lamb (Reference Lamb1932), in terms of solid spherical harmonics, to the boundary condition projected onto the spherical surface written in terms of surface spherical harmonics, the specific solution to the flow problem can be determined. Such an asymptotic approach implies that we focus on a situation in which the contact angle ${\it\Theta}$ is not too far from $90^{\circ }$ .

Figure 3.

Spherical coordinates and the direction of the undisturbed flow at infinity.

Consulting figures 3 and 1(a), we see that the coordinates of a point $\boldsymbol{r}_{p}=r_{p}\boldsymbol{e}_{r}=x_{p}\boldsymbol{e}_{x}+y_{p}\boldsymbol{e}_{y}+z_{p}\boldsymbol{e}_{z}$ on the particle–fluid interface obey the relation

in the Cartesian system, and consequently

in spherical coordinates and non-dimensional quantities. Within the scope of the present work, it is sufficient to restrict the analysis to first order in ${\it\varepsilon}$ . Additionally, in order to simplify notation, we usually omit higher-order terms in ${\it\varepsilon}$ , keeping in mind that the subsequent expressions containing ${\it\varepsilon}$ are valid only up to first order. The particle shape is therefore approximately given by

According to the method by Brenner (Reference Brenner1964), we need to expand the last term in (3.4) into a series of surface spherical harmonics,

Table 1.

Expansion coefficients $c_{km}$ for two fused spheres approximated to first order in ${\it\varepsilon}$ as calculated from (3.9); coefficients with $m<0$ and/or $m$ odd vanish.

where $f_{k}$ is a surface spherical harmonic of order $k$ defined by

In more detail, using the associated Legendre polynomials $P_{k}^{m}(\cos {\it\theta})$ , we write (Byerly Reference Byerly1893; Dassios & Vafeas Reference Dassios and Vafeas2001)

with the normalization coefficients

The expansion coefficients $c_{km}$ are given by

We have omitted the terms corresponding to $m<0$ , associated with $\sin (m{\it\varphi})$ , in (3.7) and (3.8) since the respective $c_{km}$ can be shown to vanish in the expansion (3.5). In table 1, the first few expansion coefficients are listed. Now that we know the spherical harmonics expansion of the $O({\it\varepsilon}^{})$ -approximation of two fused spheres, we turn to the calculation of the solution

and analogously for the pressure $\tilde{p}^{(0)}$ . Recalling that ${\it\varepsilon}=0$ describes the case of a particle with a contact angle of $90^{\circ }$ (or, equivalently, two fused spheres collapsed into a single sphere), it is obvious that

which is the solution to the classical problem of Stokes flow past a sphere (Happel & Brenner Reference Happel and Brenner1983). Regarding the $O({\it\varepsilon}^{})$ -term in (3.10), Brenner (Reference Brenner1964) has constructed the boundary condition

projected onto a single sphere. The formalism leading to the velocity and pressure fields $\tilde{\boldsymbol{u}}^{(0,1)}$ and $\tilde{p}^{(0,0)}$ can be found elsewhere (Brenner Reference Brenner1964; Happel & Brenner Reference Happel and Brenner1983). However, although that method can in principle be directly applied to the problem of fused spheres, corrections to certain formulae given in the above references are needed. If the formalism is applied consistently with Brenner (Reference Brenner1964) and Happel & Brenner (Reference Happel and Brenner1983), the results for the velocity and pressure fields do not solve the Stokes equation. For the special case $\boldsymbol{U}=-U\boldsymbol{e}_{z}$ discussed here, we will thus provide the corrected expressions using the same notation as Brenner (Reference Brenner1964). To be precise, equations (3.22) and (3.23) given by Brenner (Reference Brenner1964) need to be reformulated. The goal is to express the right-hand side of

through a sum of surface spherical harmonics to find the $\text{}_{k}Y_{n}^{(1)}$ . Since $\tilde{\boldsymbol{U}}=-\boldsymbol{e}_{z}$ in the present work, after some calculation we arrive at

replacing equation (3.22) by Brenner (Reference Brenner1964) and allowing us to read off

In conjunction with equations (2.12a–c) by Brenner (Reference Brenner1964) and Lamb’s solution (Brenner Reference Brenner1964, equations (2.13a,b)), (3.16) yields the expressions for $\tilde{\boldsymbol{u}}^{(0,1)}$ and $\tilde{p}^{(0,1)}$ . As those expressions are very lengthy, we waive writing them down at this point and go over to deriving the force and torque acting on the particle as well as the stress at the fluid–fluid interface. Due to the symmetry, the drag force $\boldsymbol{F}_{D}^{(0)}$ on the particle in figure 2(a) is simply half the drag force on the mirrored particle in figure 2(b). The latter can be calculated based on equation (2.20) by Brenner (Reference Brenner1964), which is well-suited to the perturbation approach used above, and is found to be

As expected, if ${\it\varepsilon}=0$ , the drag is half the Stokes drag of a sphere (Radoev et al. Reference Radoev, Nedjalkov and Djakovich1992; Danov et al. Reference Danov, Aust, Durst and Lange1995; Ally & Amirfazli Reference Ally and Amirfazli2010). Equally intuitive is the increase in drag with increasing ${\it\varepsilon}$ , that is, deeper immersion of the particle in fluid 1. We now turn to the torque acting on the particle, occurring in the balance of angular momentum (2.8) and therefore required to determine the trim angle ${\it\alpha}$ . Instead of directly calculating the torque on the particle, we write down the $y$ -component of the balance of angular momentum for the control volume enclosed by the interfaces ${\it\Sigma}_{12},{\it\Sigma}_{1p}$ and ${\it\Sigma}_{1\infty }$ shown in figure 2(b), which reads

where the vectors $\boldsymbol{n}^{(0)}$ are outer normals with respect to the control volume. Again, inertial terms have been neglected according to the assumption of small Reynolds numbers. The second integral on the right-hand side of (3.18) is equal to the negative of the viscous torque acting on the particle (from (2.8)) in the case ${\it\mu}_{2}/{\it\mu}_{1}=0$ . The third integral can be shown to vanish as the radius of the hemispherical interface ${\it\Sigma}_{1\infty }$ goes to infinity. We therefore only need to calculate the integral over the fluid–fluid interface ${\it\Sigma}_{12}$ whose integrand contains the stress distribution on ${\it\Sigma}_{12}$ (which is equal to the stress difference across the interface when ${\it\mu}_{2}/{\it\mu}_{1}=0$ ). Conveniently, the same stress distribution is required for the calculation of the $O(\mathit{Ca}^{})$ interfacial deformation further below. We have already noted that the fluid–fluid interface of order $\mathit{Ca}^{0}$ is a symmetry boundary to which the stress vector is perpendicular. Since for $\mathit{Ca}=0$ we see from figure 3 that the outer normal at ${\it\Sigma}_{12}$ is $\boldsymbol{e}_{{\it\varphi}}$ at ${\it\varphi}={\rm\pi}/2$ and $-\boldsymbol{e}_{{\it\varphi}}$ at ${\it\varphi}=3{\rm\pi}/2$ , only the component $\unicode[STIX]{x1D61B}_{{\it\varphi}{\it\varphi}}^{(0)}$ is non-zero. The latter is given by

at both ${\it\varphi}={\rm\pi}/2$ and ${\it\varphi}=3{\rm\pi}/2$ . Remarkably, the stress distribution has the form of an infinite series in negative powers of $\tilde{r}$ . As a consequence, we need to check the convergence of any quantity dependent on $\tilde{\unicode[STIX]{x1D61B}}_{{\it\varphi}{\it\varphi}}^{(0)}$ . Now we are in a position to evaluate

which will be used further below for finding the trim angle ${\it\alpha}$ .

4. Interfacial deformation and particle trim angle

In order to parametrize the shape of the fluid–fluid interface ${\it\Sigma}_{12}$ of order $\mathit{Ca}$ – resulting, within the perturbation approach, from the stress distribution of order $\mathit{Ca}^{0}$ calculated above – we introduce circular cylindrical coordinates $({\it\varrho},{\it\vartheta},{\it\chi})$ ,

allowing us to utilize the expression (3.19) for the stress distribution by setting $\tilde{r}\rightarrow \tilde{{\it\varrho}}$ and ${\it\theta}\rightarrow 2{\rm\pi}-{\it\vartheta}$ . Note that in the cylindrical coordinate system $0\leqslant {\it\vartheta}<2{\rm\pi}$ , in contrast to $0\leqslant {\it\theta}<{\rm\pi}$ for spherical coordinates, which does not introduce a sign error but helps us to avoid the distinction between the regions of positive and negative $y$ (that is, between ${\it\varphi}={\rm\pi}/2$ and ${\it\varphi}=3{\rm\pi}/2$ at $x=0$ ) necessary in the context of spherical coordinates (cf. (4.6) below). Within the cylindrical coordinate system, we measure the interfacial deformation through its undulation in the direction of the ${\it\chi}$ -axis,

By choosing the orientation of the normal vector $\boldsymbol{n}_{12}$ consistently with the jump condition (2.5), we obtain

and insert the expansion $h=h^{(1)}\mathit{Ca}+O(\mathit{Ca}^{2})$ . We are led to

as well as to the curvature

Figure 4.

Geometry of the TCL when both ${\it\alpha}$ and ${\it\varepsilon}$ are positive.

In order to evaluate the interfacial momentum balance (2.12), we need to specify the sign of the stress deforming the interface. Since the left-hand side of (2.12) is equal to the term $\tilde{\unicode[STIX]{x1D64F}}^{(0)}\boldsymbol{\cdot }\boldsymbol{n}_{12}^{(0)}$ and we are interested in the normal component of the stress, we may write (using $\boldsymbol{n}_{12}^{(0)}=\boldsymbol{e}_{{\it\chi}}=\boldsymbol{e}_{x}$ according to (4.4))

Insertion of (4.5), (4.6), and (3.19) – now in cylindrical coordinates – into the jump condition (2.12) and scalar multiplication with $\boldsymbol{n}_{12}^{(0)}$ yields

In addition, boundary conditions for $h$ are to be specified. At ${\it\varrho}\rightarrow \infty$ , obviously we have $h\rightarrow 0$ because the origin of the cylindrical coordinate system coincides with the centre of the disk enclosed by the TCL (see figure 4). A closer look is needed concerning the boundary condition on the TCL which has an elliptic shape with principal axes $2a\sqrt{1-{\it\varepsilon}^{2}}$ and $2a\sqrt{1-{\it\varepsilon}^{2}}\cos {\it\alpha}$ when projected onto the $yz$ -plane. The coordinates of a point $\boldsymbol{r}_{\mathit{TCL}}$ on the TCL can be found by first describing a circle of radius $a\sqrt{1-{\it\varepsilon}^{2}}$ embedded in the $yz$ -plane and, secondly, by rotating the coordinate system by an angle of ${\it\alpha}$ around the $y$ -axis. One obtains

Since the trim angle ${\it\alpha}$ depends on the particle velocity, it can be written as ${\it\alpha}={\it\alpha}^{(1)}\mathit{Ca}+O(\mathit{Ca}^{2})$ . In conjunction with $\sqrt{1-{\it\varepsilon}^{2}}=1+O({\it\varepsilon}^{2})$ and by omission of the terms $O({\it\varepsilon}^{2})$ , (4.8) is simplified to

Considering the parametrization (4.2) and expanding

we read off the boundary condition

It is now straightforward to find the solution to (4.7) consistent with the boundary conditions (4.11) and $\tilde{h}^{(1)}(\tilde{{\it\varrho}}\rightarrow \infty ,{\it\vartheta})=0$ , which is given by

Note that the trim angle ${\it\alpha}^{(1)}$ is still unknown and is to be determined from the angular momentum balance (2.8). Owing to the identical stress term $\partial \boldsymbol{r}_{\mathit{TCL}}/\partial s\times \boldsymbol{n}_{12}$ in (2.7) and (2.8), we can derive the expression for $\boldsymbol{F}$ in parallel to the one for ${\it\alpha}^{(1)}$ . The parametrization is chosen as $\text{d}\tilde{s}=\text{d}{\it\vartheta}$ (recall that the stress vector will be oriented outwards with respect to the particle). Since $\partial \boldsymbol{e}_{{\it\varrho}}/\partial {\it\vartheta}=\boldsymbol{e}_{{\it\vartheta}}$ , we may deduce from (4.9)

The fluid–fluid interfacial normal vector $\boldsymbol{n}_{12}$ , or – as a consequence of (4.4) – the derivatives of $\tilde{h}^{(1)}$ with respect to $\tilde{{\it\varrho}}$ and ${\it\vartheta}$ , need to be evaluated at the TCL. Analogously to (4.10), this means that the derivatives of $\tilde{h}^{(1)}$ in (4.4) are to be evaluated at $(1,{\it\vartheta})$ as long as we limit the analysis to first order in $\mathit{Ca}$ . As a result, we find

Obviously, the integral along the TCL occurring in (2.7) vanishes, which implies that interfacial tension does not contribute to the drag. The driving force $\boldsymbol{F}^{(0)}$ is therefore equal to the negative of the drag force $\boldsymbol{F}_{D}^{(0)}$ . In the case of the $y$ -component of (2.8), normalized by ${\it\sigma}_{12}a^{2}$ , a similar calculation using the result (4.14) yields

As has been discussed above in the context of (3.18), the viscous torque acting on the particle is equal to the negative of the viscous torque (3.24) acting on the fluid volume enclosed by the interfaces ${\it\Sigma}_{12},{\it\Sigma}_{1p}$ and ${\it\Sigma}_{1\infty }$ (recall figure 2 b):

The $y$ -component of (2.8) can now be solved by insertion of the expressions (4.15) and (4.16). We are led to the result for the particle trim angle to first order in $\mathit{Ca}$ ,

where the numerical value of 1.146 corresponds to the sum of the first nine terms. The convergence of expression (4.17) is demonstrated in figure 5. By means of the result for ${\it\alpha}^{(1)}$ we can write down the final formula describing the interfacial shape,

which is visualized in figure 6 for a capillary number of $\mathit{Ca}=0.3$ and an equilibrium contact angle of ${\it\Theta}=75^{\circ }$ . The lower half-space is occupied by fluid 1 and the particle is moving to the right of the picture. A white line piercing through the north and south poles of the sphere has been added to indicate the change in orientation, that is, the trim angle. From figure 6 it can be seen that the contact angle changes around the particle’s circumference. While the contact angle retains its equilibrium value at ${\it\vartheta}={\rm\pi}/2$ and $3{\rm\pi}/2$ , it deviates most strongly from that value at ${\it\vartheta}=0$ and ${\rm\pi}$ . From geometric considerations (see also figure 4), the deviation at ${\it\vartheta}=0$ is found to be

where again nine terms have been included in the summation. At ${\it\vartheta}={\rm\pi}$ , the contact angle difference ${\it\delta}{\it\Theta}$ has the same magnitude but the opposite sign. Figure 7 shows cross-sections of the interfacial shape and the TCL in the plane ${\tilde{y}}=0$ for three different contact angles ${\it\Theta}$ . The particle trim angle, indicated by the inclination of the TCL plane between $\tilde{z}=-1$ and $\tilde{z}=1$ , obviously increases when the particle plunges deeper into fluid 1 (larger ${\it\varepsilon}$ ) with its relatively high viscosity, which is an intuitive result.

Figure 5.

Convergence plot for the $O({\it\varepsilon}^{})$ -term in (4.17).

Figure 6.

Visualization of the $O(\mathit{Ca}^{1})$ -configuration for an equilibrium contact angle ${\it\Theta}$ of $75^{\circ }$ and a capillary number of $\mathit{Ca}=0.3$ . The trim angle can be estimated from the inclination of the white line. The gap between the particle surface and the fluid–fluid interface arises because the boundary condition is evaluated on a cylindrical surface according to (4.10).

Figure 7.

Cross-section of the interfacial shape $\tilde{h}=\tilde{h}^{(1)}\mathit{Ca}$ according to (4.18) and of the corresponding TCL for three different equilibrium contact angles ${\it\Theta}$ and a capillary number of $\mathit{Ca}=0.3$ : (a)  ${\it\Theta}=75^{\circ }$ , (b)  ${\it\Theta}=90^{\circ }$ , (c)  ${\it\Theta}=105^{\circ }$ . The inclination of the TCL is given by the particle trim angle ${\it\alpha}={\it\alpha}^{(1)}\mathit{Ca}$ from (4.17); the particle is moving in the $\tilde{z}$ -direction.

Figure 8.

Coordinate systems for two interacting particles $p1$ and $p2$ of separation $d$ viewed from the negative $x/{\it\chi}$ -direction; integration path $c=c_{I}\cup c_{II}\cup c_{III}\cup c_{IV}$ used in the evaluation of the capillary force.

5. Pair interaction via linear superposition approximation

Within the LSA, the interfacial deformation (4.18) around a single particle is used to construct the two-particle deformation necessary for studying capillary interactions via linear superposition. Since the method can only be expected to be accurate at large inter-particle separation $d\ll a$ , we truncate the infinite series in (4.18) after the first, asymptotically dominant, $\tilde{{\it\varrho}}^{-2}$ -term. Returning to dimensional quantities, but keeping the capillary number as an abbreviation, we may express the deformation $h_{pi}$ around a single particle $pi$ in the form

In (5.1), the cylindrical coordinate system $({\it\varrho}_{i},{\it\vartheta}_{i},{\it\chi})$ is centred at particle $i$ . Since we will later use only the coordinate system of particle 1, we set $({\it\varrho}_{1},{\it\vartheta}_{1},{\it\chi})=({\it\varrho},{\it\vartheta},{\it\chi})$ for simplicity (see figure 8). The notation $\mathit{Ca}_{pi}$ instead of $\mathit{Ca}$ accounts for a potentially different velocity between the particles (in magnitude and sign). This may apply to the case when the particles differ in the properties governing their reaction to the external field driving them. Figure 8 shows the two particles with their centres separated by a distance $d\gg a$ and with the line connecting their centres coinciding with the coordinate line ${\it\vartheta}={\it\beta}$ . In the course of the LSA, it is required to re-express the single-particle deformation $h_{p2}$ , which is given by (5.1) within the $({\it\varrho}_{2},{\it\vartheta}_{2},{\it\chi})$ -system, in terms of the coordinates $({\it\varrho},{\it\vartheta},{\it\chi})$ . From figure 8, we can geometrically derive the relations

so that we have

Regarding the interfacial deformation $h_{\mathit{LSA}}$ in the presence of two particles, the LSA is based on the assumption that

Physically, the dynamical deformation of the interface by the particles gives rise to an interaction force. In the present case of two moving particles, the interaction force is expected to be dependent on both the relative positions of the particles (described by ${\it\beta}$ ), their separation, and the respective capillary numbers $\mathit{Ca}_{pi}$ . To be sure, the boundary conditions on the particle surfaces are generally not fulfilled by the sum of the deformations according to (5.4). However, in the region halfway between the particles, the approximation is supposed to be accurate when the separation $d$ is large. Therefore, a reasonable path of integration ( $c=c_{I}\cup c_{II}\cup c_{III}\cup c_{IV}$ , see figure 8) for the interfacial tension force runs through that region and is closed at infinite distance from the particles. The interaction force $\boldsymbol{F}_{\mathit{int}}$ may be written as

where the parametrization of a point $\boldsymbol{r}_{c}$ on the integration path $c$ based on $s$ is chosen such that the integrand becomes a vector pointing away from particle 1. Prior to its evaluation, the integral in (5.5) can be simplified considerably with the aid of figure 8. First, we understand that the contributions along the paths $c_{II}$ and $c_{IV}$ mutually cancel because the interfacial deformation vanishes at infinite distance from the particles. We are left with the integrals along $c_{I}$ and $c_{III}$ , which diverge when evaluated separately since the interfacial tension force (more specifically, its line density) becomes constant at infinity. Therefore, both integrals need to be evaluated simultaneously. In effect, the contribution from $c_{III}$ can be taken into account by subtraction of a vector of constant orientation under the integral since, again, the interfacial force is constant at infinity. Choosing $\text{d}s=\text{d}{\it\vartheta}$ , the vector $\boldsymbol{r}_{c}$ describing a point on the integration path is given by

with $h_{\mathit{LSA}}$ according to (5.4). The normal vector, which is to be evaluated at $\boldsymbol{r}_{c}$ using $h_{\mathit{LSA}}$ , is calculated from (4.4). As discussed above, we need to subtract a vector

representing the contribution from $c_{III}$ . Note that with the parametrization $\text{d}s=\text{d}{\it\vartheta}$ the vector (5.7) is not constant as it would have been if $\text{d}s$ had been set equal to $\text{d}y$ , say. The foregoing expressions are valid for any choice of $\mathit{Ca}_{p1}$ and $\mathit{Ca}_{p2}$ . Within the following discussion, however, we focus on the most relevant case of two particles having identical velocity vectors $U\boldsymbol{e}_{z}$ and thus set $\mathit{Ca}=\mathit{Ca}_{p1}=\mathit{Ca}_{p2}$ . As the evaluation of (5.5) and of the analogous relation for the torque $\boldsymbol{M}_{\mathit{int}}$ with respect to the centre of particle 1 are straightforward from the above considerations, we may directly state the results

Figure 9.

Interaction map qualitatively describing the components of the interaction force (5.8) acting on particle 1 dependent on the relative particle positions: (a)  $y$ - and $z$ -components, (b)  $x$ -component. The correspondence between length of the vectors and the force magnitude differs between the subfigures.

in $(x,y,z)$ -coordinates, where $H:=(1/2+11{\it\varepsilon}/8)\mathit{Ca}$ . Recall that the particle motion due to the external force is aligned with the $z$ direction. The first term ${\sim}d^{-2}$ in (5.8) only results in vertical ( $x$ -)adjustment of the particle positions. On the other hand, the remaining two terms ${\sim}d^{-5}$ act in the plane of the undisturbed fluid interface. In figure 9 the force components are shown qualitatively in the form of a vector pinned at a position determined by the respective value of ${\it\beta}$ . The interaction map may be read as follows. For a given ${\it\beta}$ , the connecting line between the particles $p1$ and $p2$ intersects the dashed circle in figure 9. The vector plotted at the intersection point represents the force acting on particle 1 within the actual configuration. Figure 9(a) displays the $y$ - and $z$ -components of the interaction force while figure 9(b) shows its $x$ -component. The scaling of the forces differs between figures 9(a) and 9(b) for reasons of optimal visualization, owing to the different dependences of the force components on the separation $d$ (see (5.8)). In addition to the direction and the relative strength of the forces depending on the particle configuration, figure 9 also contains information on the torque resulting from the interaction. It can be shown that the forces sketched in the figure as pinned to the dashed circle do indeed provide a full representation of the force distribution, meaning that the vector between the centre of particle $p1$ and the respective intersection point on the dashed circle, namely

(see figure 8), equals the lever arm for the calculation of the torque on the particle by means of the interaction force. In other words, we have the relation

which can be easily verified via (5.8) and (5.9). As an example, consider the black dot in figure 9 symbolizing an intersection point. The force on particle 1 can be read off as pointing to the northeast of the picture while the resulting torque tends to rotate the particle in the anticlockwise direction because the line of action of the force does not run through the particle centre but passes south of it.

From (5.8) and, accordingly, figure 9 it can be seen that the interaction force vanishes when the particles move along a common line. When the second particle’s position is on the $y$ -axis, the force is attractive and acts along the line of centres. Between those extreme cases the interaction force changes continuously. The force component $F_{\Vert }$ along the line of centres

is generally attractive. Following the commonly accepted method of analysing interfacial deformation by means of a multipole expansion (Danov et al. (Reference Danov, Kralchevsky, Naydenov and Brenn2005, equation 2.16); Oettel & Dietrich Reference Oettel and Dietrich2008), we identify the interfacial shape (4.18) as dipolar. As a consequence, the interaction between particles driven by an external force, given by (5.8), may also be termed dipolar. In fact, the interaction force between two fixed electric dipoles with equal dipole moments lying within a common plane has the form (Israelachvili Reference Israelachvili2011, table 2.2, setting ${\it\theta}_{1}={\it\theta}_{2}={\rm\pi}/2-{\it\beta}$ and ${\it\phi}=0$ )

showing the same variation with ${\it\beta}$ as in the case of driven particles.

Since the particles move along the fluid–fluid interface and their motion causes fluid flow, the discussion of interactions between the particles is only complete if we comment on hydrodynamic interactions. For two equal spheres at low Reynolds number in an unbounded fluid acted upon by the same force, it is well known that the particle velocity vectors are always parallel (Russel et al. Reference Russel, Saville and Schowalter1989, p. 52; Happel & Brenner Reference Happel and Brenner1983, pp. 242–243). In that context, the particles’ direction of motion need not be parallel to the driving force vector and the inter-particle separation remains constant during the motion. In the case of particles moving along a fluid–fluid interface considered in the present work, we have already seen that the hydrodynamic problem at zero capillary number can be treated by means of an equivalent problem of a particle in an unbounded fluid (cf. figure 2). Therefore, the aspect of hydrodynamic interactions for the $O(\mathit{Ca}^{0})$ -problem is also equivalent to the interaction of particles in an unbounded fluid. The same reasoning has been followed by Vassileva et al. (Reference Vassileva, van den Ende, Mugele and Mellema2005). As a consequence, the hydrodynamic interaction between particles moving along a fluid interface by action of the same external driving force does not influence their separation, that is, the hydrodynamic force vectors acting on each particle are equal. We conclude that capillary effects govern the interaction between driven interfacial particles at least up to order $\mathit{Ca}$ .

In order to estimate the strength of the dipolar interaction force, varying as $d^{-5}$ according to (5.8), a comparison with another capillary interaction force is instructive, namely the quadrupolar force arising from an undulated TCL (Stamou et al. Reference Stamou, Duschl and Johannsmann2000; Danov et al. Reference Danov, Kralchevsky, Naydenov and Brenn2005; Domínguez et al. Reference Domínguez, Oettel and Dietrich2008). The latter also decays as $d^{-5}$ . Using the values given by Stamou et al. (Reference Stamou, Duschl and Johannsmann2000) in their estimation, we arrive at the quadrupolar interaction force $F_{\unicode[STIX]{x2723}}=48{\rm\pi}{\it\sigma}_{12}a^{-1}{\it\delta}x\,\tilde{d}^{-5}$ . Thereby, ${\it\delta}x$ is a characteristic undulation of the TCL. Park & Furst (Reference Park and Furst2011) report a value of ${\it\delta}x/a=0.03$ for the particles investigated, leading to $F_{\unicode[STIX]{x2723}}\approx 0.14\,{\it\sigma}_{12}a\tilde{d}^{-5}$ . On the other hand, the parallel component of the dipolar interaction force in the case ${\it\varepsilon}=0$ and ${\it\beta}={\rm\pi}/2$ is $F_{\Vert }\approx 21.21{\it\sigma}_{12}a\mathit{Ca}^{2}\tilde{d}^{-5}$ . Therefore, the two forces are equal at a capillary number of $\mathit{Ca}\approx 0.08$ corresponding to a particle velocity of $U\approx 1.7~\text{mm}~\text{s}^{-1}$ for the 1000 cSt silicone oil used by Ally & Amirfazli (Reference Ally and Amirfazli2010). If, of course, the TCL is less irregular, the dipolar force may dominate at an even smaller velocity. The same obviously holds for a liquid of higher viscosity.