Long-range two-dimensional hydrodynamic interaction between a pair of mutually repellent disks

Abstract

While the problem governing Stokes flow about a single particle that is subject to an external force is ill posed in two dimensions (the ‘Stokes paradox’), the related problem of two mutually repellent particles is well posed. Motivated by self-assembly phenomena in thin viscous membranes, we consider this problem in the limit of remote particles. Such limits are typically handled in the literature using reflection techniques, which provide successive approximations to the mutual hydrodynamic interactions. Since their starting point is a single particle in an unbounded fluid domain, these techniques are futile in the present two-dimensional problem. We show how this apparent contradiction is resolved via use of singular perturbations. We obtain a two-term approximation for the velocity acquired by circular disks, considering both rigid and free particle surfaces. We also illustrate our perturbation scheme for elliptic disks, deriving a renormalised single-particle velocity. The utility of our asymptotic scheme is illustrated in the general problem of hydrodynamic interaction between a cluster of remote disks.

1. Introduction

In analysing the motion of molecular probes in liquid membranes (Saffman & Delbrück 1975), it is convenient (Goodrich 1969; Evans & Sackmann 1988) to treat the membrane as a zero-thickness surface with a Boussinesq–Scriven rheology (Scriven 1960). The hydrodynamic problem governing the translation of a single probe was addressed by Saffman (1976), who modelled the probe as a membrane-trapped circular disk. Saffman focused upon the limit of small probes (e.g. proteins), where viscous forces in the liquid substrate bounded by the membrane are presumably negligible (large Boussinesq numbers). While this may appear to result in a convenient two-dimensional (2-D) set-up, the hydrodynamic problem does not admit a solution: all possible flow fields that are compatible with the no-slip condition diverge logarithmically at large distances from the particle. This is the well-known Stokes paradox.

Saffman (1976) showed that in the case of a membrane-bound probe, the Stokes paradox is resolved by the incorporation of the substrate stresses, which enter the leading-order balance at large distances. (This regularisation mechanism typically ‘overrides’ the familiar resolution via fluid inertia which takes place at even larger distances.) The resulting two-scale problem was solved using matched asymptotic expansions, eventually leading to Saffman's celebrated mobility formula.

There are several related problems that are well posed even when the dynamic effect of the liquid substrate is completely neglected. One obvious such case, considered briefly by Saffman (1976), is that of a bounded membrane; another case involves curved membranes (e.g. spherical and cylindrical) that close on themselves (Henle & Levine 2010). Other related configurations, motivated by heat and mass transport problems, involve periodic (Hasimoto 1959; Sangani & Acrivos 1982) and random (Sangani & Yao 1988) arrays of cylindrical obstacles.

As the aforementioned examples illustrate genuine 2-D Stokes flow problems under force fields, they are of fundamental interest beyond their practical importance. We here address another such fundamental scenario, namely that of two identical membrane-bound disks that repel each other via equal and opposite central forces. The quantity of interest is the velocity attained by the disks. As in Saffman (1976), we focus at the limit of large Boussinesq numbers, Since there is no external force, the resulting 2-D hydrodynamic problem is well posed.

This problem of mutual interaction is pertinent to the understanding of the manner by which electrostatic repulsion between interface-bounded colloidal particles results in the formation of ordered phases (Bresme & Oettel 2007). The literature abounds with both experimental measurements (Aveyard et al. 2002; Wirth, Furst & Vermant 2014) and theoretical models (Hurd 1985; Frydel, Dietrich & Oettel 2007) of the repulsive forces, but it appears that a proper description of the hydrodynamic response for a given interaction is still lacking. With the Stokes equations depending only upon the instantaneous configuration of the particle pair, the hydrodynamic description is actually indifferent to the particular form of the repulsive force. This allows us to address a generic interaction that satisfies Newton's third law. For convenience, we assume here the standard case of a central force.

Since colloid repulsion on interfaces is long-ranged, we seek an approximation for well-separated particles. Beyond practical interest, this limit entails an apparent contradiction. In the standard transport literature (Leal 2007), remote hydrodynamic interactions are typically handled using reflection methods where, starting with a solution for a single particle, successive approximations are constructed for the impact of one particle on its neighbour. These methods are inapplicable in the present problem: first, because the single-particle solution is not unique, and second, because the logarithmic divergence results in a scheme where successive terms only become larger in magnitude.

We here solve this singular problem using the method of matched asymptotic expansions (Hinch 1991). Thus we separately analyse the flow in an ‘inner’ region, on the particle scale, and an ‘outer’ region, on the separation scale. Since the condition of velocity attenuation does not apply at the inner scale, no Stokes paradox emerges. The velocity field, however, is not unique, the freedom being associated with the unknown particle speed. On the outer scale, the two particles appear as singular point forces in opposite directions, ensuring the trivial satisfaction of the attenuation condition. Asymptotic matching between the inner and outer solutions provides the requisite particle speed.

Our robust scheme can handle various interfacial conditions. This is important due to uncertainties regarding the proper boundary condition that is applicable in interface-trapped colloids (Saffman & Delbrück 1975). Considering the simplest geometry of circular particles, the scheme is applied to both rigid boundaries, where the velocity satisfies a no-slip condition, and free boundaries, where it satisfies a shear-free condition. In both cases, we go beyond leading-order calculation and obtain the asymptotic correction to the particle speed. We further illustrate the applicability of the asymptotic method to particles of non-circular cross-section, analysing the interaction between disks of elliptic cross-section. Finally, we generalise our asymptotic scheme to illuminate the hydrodynamic interaction between remote particles.

The paper is arranged as follows. In the next section we formulate the problem governing the hydrodynamic interaction between two circular disks. In § 3, we describe the asymptotic structure in the singular limit of remote particles. The leading-order analysis is carried out in § 4, with asymptotic corrections being derived subsequently in § 5. The corresponding analysis for shear-free disks is given in § 6. The handling of non-circular shapes is illustrated in § 7 for the case of elliptic disks. Generalisations to a cluster of non-identical disks are presented in § 8, which may also serve as a useful recapitulation. We conclude in § 9.

2. Problem formulation

Two circular disks (radius $\ell$) are suspended in a viscous film (surface viscosity $\mu$). The instantaneous distance between their centres is $2\ell /\lambda$ ($\lambda <1$). The disks experience a mutual repulsion, represented by a central force of magnitude $F$ (with $F<0$ corresponding to mutual attraction). Our interest is in the velocity acquired by the disks.

We employ a dimensionless formulation, where length, force and velocity variables are normalised by $\ell$, $F$ and $F/\mu$, respectively. We use the Cartesian coordinates $(x',y')$ with the respective unit vectors denoted by $(\hat {\boldsymbol {\imath }},\hat {\boldsymbol {\jmath }})$. The $x'$-axis passes through the particles’ line of centres, with the origin midway between them (see figure 1). For concreteness, we denote the disks centred at $x'=\pm 1/\lambda$ as the ‘$\pm$’ disks. The forces on the $\pm$ disks are $\pm \hat {\boldsymbol {\imath }}$. Due to the problem symmetry, the velocities acquired by the ‘$\pm$’ disks are in the $x'$-direction; they are denoted by $\pm {\mathcal {U}}\hat {\boldsymbol {\imath }}$, respectively. The velocity ${\mathcal {U}}$, a function of $\lambda$, is the quantity of interest.

Figure 1. Geometry and coordinate systems.

The velocity field $\boldsymbol {u}'$ is governed by the continuity and Stokes equations, the no-slip condition

(2.1)\begin{equation} \text{$\boldsymbol{u}'={\pm} {\mathcal{U}}\hat{\boldsymbol{\imath}}\quad $ at the boundaries of the `${\pm}$' disks},\end{equation}

and the decay condition

(2.2)\begin{equation} \lim_{{x'}^2+{y'}^2\to\infty}\boldsymbol{u}'=\boldsymbol{0}. \end{equation}

The velocity ${\mathcal {U}}$, appearing in (2.1), is determined by the requirement that

(2.3)\begin{equation} \text{the hydrodynamic forces on the `${\pm}$' disks are ${\mp}\hat{\boldsymbol{\imath}}$, respectively.}\end{equation}

It is evident that $\hat {\boldsymbol {\imath }}\boldsymbol {\cdot }\boldsymbol {u}'$ is an odd function of $x'$ and an even function of $y'$, while $\hat {\boldsymbol {\jmath }}\boldsymbol {\cdot }\boldsymbol {u}'$ is an even function of $x'$ and an odd function of $y'$. In particular,

(2.4)\begin{equation} \left. \begin{array}{@{}l@{}} \hat{\boldsymbol{\imath}}\boldsymbol{\cdot}\boldsymbol{u}' = 0 \\[2pt] \hat{\boldsymbol{\jmath}}\boldsymbol{\cdot}\displaystyle\dfrac{\partial{\boldsymbol{u}'}}{\partial{x'}}=0 \end{array} \right\} \quad \text{at } x'=0.\end{equation}

In what follows, we directly exploit the aforementioned symmetry about the $y'$-axis. Thus we solve the problem only for $x'>0$, using (2.4) as a symmetry condition and abandoning the no-slip condition on the ‘$-$’ disk. We note that the symmetry about the $x'$-axis is trivially satisfied in the problem formulation, and need not be imposed.

3. Matched asymptotic expansions

Our interest is in the remote limit,

(3.1)\begin{equation} \lambda\ll1. \end{equation}

As already explained, this limit cannot be resolved by perturbing about the problem of a single disk, as the latter is ill posed. Instead, we employ the method of matched asymptotic expansions (Hinch 1991). Thus we utilise separate asymptotic expansions at the inner region, on the scale of the ‘+’ disk, and the outer region, on the scale of the disks’ separation.

The inner Cartesian coordinates are defined by

(3.2a,b)\begin{equation} x = x'-1/\lambda, \quad y = y'. \end{equation}

The inner position vector $\hat {\boldsymbol {\imath }} x+\hat {\boldsymbol {\jmath }} y$ is denoted by $\boldsymbol {x}$. We also employ polar coordinates $(\rho,\phi )$, defined by (see figure 1)

(3.3a,b)\begin{equation} x=\rho\cos\phi,\quad y=\rho\sin\phi. \end{equation}

The associated unit vectors are denoted by $(\hat {\boldsymbol {e}}_\rho,\hat {\boldsymbol {e}}_\phi )$, respectively. Note that $\rho = |\boldsymbol {x}|$. The velocity in the inner region, $\boldsymbol {u}=\hat {\boldsymbol {e}}_\rho u+\hat {\boldsymbol {e}}_\phi v$, is governed by the continuity and Stokes equations, the no-slip condition (cf. (2.1))

(3.4)\begin{equation} \left. \begin{array}{@{}l@{}} u = {\mathcal{U}}\cos\phi \\[2pt] v ={-}{\mathcal{U}}\sin\phi \end{array} \right\} \quad \text{at } \rho=1,\end{equation}

and the force requirement (cf. (2.3))

(3.5)\begin{equation} \text{hydrodynamic force on the `+' disk}={-}\hat{\boldsymbol{\imath}}.\end{equation}

We note that neither the decay condition (2.2) nor the symmetry constraint (2.4) applies at the inner region. Rather, the inner solution is required to match the outer one. On the other hand, the symmetry about the $x$-axis, inherited from the exact problem, must now be enforced explicitly. In terms of the polar components, it reads

(3.6a,b)\begin{equation} u(\rho,-\phi) = u(\rho,\phi), \quad v(\rho,-\phi) ={-}v(\rho,\phi).\end{equation}

For convenience we employ the stream function $\psi$, defined by

(3.7a,b)\begin{equation} u = \frac{1}{\rho}\,\frac{\partial{\psi}}{\partial{\phi}}, \quad v ={-} \frac{\partial{\psi}}{\partial{\rho}}, \end{equation}

whereby the continuity equation is trivially satisfied. The stream function is biharmonic. It must satisfy the symmetry condition (cf. (3.6))

(3.8)\begin{equation} \psi(\rho,-\phi) ={-}\psi(\rho,\phi), \end{equation}

where, with no loss of generality, $\psi$ is taken as zero at $y=0$ ($\phi =0,{\rm \pi}$). The no-slip condition (3.4) then gives

(3.9)\begin{equation} \left. \begin{array}{@{}l@{}} \psi = {\mathcal{U}} \sin\phi \\[2pt] \displaystyle\dfrac{\partial{\psi}}{\partial{\rho}} ={\mathcal{U}}\sin\phi \end{array} \right\} \quad \text{at } \rho=1.\end{equation}

The Cartesian coordinates $(X,Y)$ in the outer region are defined by

(3.10a,b)\begin{equation} X = \lambda x', \quad Y = \lambda y', \end{equation}

corresponding to length normalisation by $\ell /\lambda$. The position vector is $\boldsymbol {X}=\hat {\boldsymbol {\imath }} X+\hat {\boldsymbol {\jmath }} Y$. Note that $\boldsymbol {X}=\lambda \boldsymbol {x}'$ and $\boldsymbol {X}-\hat {\boldsymbol {\imath }}=\lambda \boldsymbol {x}$. The velocity field in the outer region is denoted by $\boldsymbol {U}$. It is governed by the continuity and Stokes equations; in addition, it satisfies the symmetry condition (cf. (2.4))

(3.11)\begin{equation} \left.\begin{array}{@{}l@{}} \hat{\boldsymbol{\imath}}\boldsymbol{\cdot}\boldsymbol{U} = 0 \\[2pt] \hat{\boldsymbol{\jmath}}\boldsymbol{\cdot}\displaystyle\dfrac{\partial{\boldsymbol{U}}}{\partial{X}}=0 \end{array} \right\}\quad \text{at } X=0,\end{equation}

and the decay condition (cf. (2.2))

(3.12)\begin{equation} \lim_{X^2+Y^2\to\infty}\boldsymbol{U}=\boldsymbol{0}. \end{equation}

Neither the no-slip condition (recall (2.1)) nor the force condition (3.5) applies directly in the outer region.

We employ the following asymptotic expansion for the inner velocity field,

(3.13)\begin{equation} \boldsymbol{u}(\boldsymbol{x};\lambda) \sim \boldsymbol{u}_0(\boldsymbol{x};\lambda) + \lambda\,\boldsymbol{u}_1(\boldsymbol{x};\lambda)+\lambda^2\,\boldsymbol{u}_2(\boldsymbol{x};\lambda) + \cdots; \end{equation}

it induces a comparable expansion of $\psi$. Following Fraenkel's warning (Fraenkel 1969), we do not separate asymptotic orders by logarithms of the expansion parameter. Thus we allow the respective coefficients to depend ‘weakly’ on $\lambda$ through its logarithm. This enables the use of the Van Dyke matching rule (Van Dyke 1964). Given condition (3.4), the disk velocity must possess the expansion

(3.14)\begin{equation} {\mathcal{U}}(\lambda) \sim \mathcal {\mathcal{U}}_0(\lambda) + \lambda\,{\mathcal{U}}_1(\lambda) + \lambda^2\,{\mathcal{U}}_2(\lambda) + \cdots, \end{equation}

where, again, we allow for a weak dependence upon $\lambda$ (see indeed (4.7)).

The counterpart of (3.13) in the outer region is

(3.15)\begin{equation} \boldsymbol{U}(\boldsymbol{X};\lambda) \sim \boldsymbol{U}_0(\boldsymbol{X};\lambda) + \lambda\,\boldsymbol{U}_1(\boldsymbol{X};\lambda)+\lambda^2\,\boldsymbol{U}_2(\boldsymbol{X};\lambda) + \cdots. \end{equation}

4. Leading-order solutions

The force condition (3.5) implies a fluid-velocity term that diverges as $\ln \rho$ at large $\rho$ (a 2-D Stokeslet); this is the well-known Stokes paradox. The leading-order solution in the inner region is obtained by requiring that the divergence rate at infinity is not worse. (That is not the case in higher asymptotic orders, of course.) The symmetry constraint (3.6a,b) then implies a superposition of a Stokeslet, a uniform stream, and an irrotational doublet (Pozrikidis 1992)

(4.1)\begin{equation} \boldsymbol{u}_0 = \frac{1}{4{\rm \pi}}\,(-{\boldsymbol{\mathsf{I}}}\ln\rho+\hat{\boldsymbol{e}}_\rho\hat{\boldsymbol{e}}_\rho)\boldsymbol{\cdot}\hat{\boldsymbol{\imath}} + c_0 \hat{\boldsymbol{\imath}} + \frac{d_0}{2{\rm \pi}\rho^2}\,(-{\boldsymbol{\mathsf{I}}} + 2 \hat{\boldsymbol{e}}_\rho\hat{\boldsymbol{e}}_\rho) \boldsymbol{\cdot}\hat{\boldsymbol{\imath}},\end{equation}

where the magnitude of the Stokeslet has been set by requirement (3.5). The no-slip condition (2.1) gives $d_0=-1/4$ and $c_0={\mathcal {U}}_0-1/8{\rm \pi}$, whereby

(4.2)\begin{equation} \boldsymbol{u}_0 = {\mathcal{U}}_0 \hat{\boldsymbol{\imath}} + \frac{1}{4{\rm \pi}}\left(-{\boldsymbol{\mathsf{I}}}\ln\rho+\hat{\boldsymbol{e}}_\rho\hat{\boldsymbol{e}}_\rho - \frac{1}{2}\,{\boldsymbol{\mathsf{I}}}\right)\boldsymbol{\cdot}\hat{\boldsymbol{\imath}} - \frac{1}{8{\rm \pi}\rho^2}\,(-{\boldsymbol{\mathsf{I}}} + 2 \hat{\boldsymbol{e}}_\rho\hat{\boldsymbol{e}}_\rho) \boldsymbol{\cdot}\hat{\boldsymbol{\imath}}. \end{equation}

The associated stream function is therefore given by

(4.3)\begin{equation} \psi_0 =\frac{1}{4{\rm \pi}} \left[ \left(4{\rm \pi}{\mathcal{U}}_0-\frac{1}{2}\right)\rho - \rho(\ln\rho-1) - \frac{1}{2\rho}\right] \sin\phi. \end{equation}

The velocity ${\mathcal {U}}_0$ cannot be determined from the inner analysis. To find it, we note that

(4.4)\begin{equation} \boldsymbol{u}_0 \sim{-} \frac{1}{4{\rm \pi}}\,\hat{\boldsymbol{\imath}} \ln\rho +{\mathcal{U}}_0\hat{\boldsymbol{\imath}} +\frac{1}{4{\rm \pi}}\left( \hat{\boldsymbol{e}}_\rho\hat{\boldsymbol{e}}_\rho- \frac{1}{2}\,{\boldsymbol{\mathsf{I}}}\right)\boldsymbol{\cdot}\hat{\boldsymbol{\imath}} +O(\rho^{{-}2}) \quad \text{as } \rho\to\infty,\end{equation}

and follow by considering the outer region.

The leading-order outer flow is set uniquely by constraints (3.11) and (3.12) together with the need to match (4.4). It therefore consists of a Stokeslet at $\boldsymbol {X}=\hat {\boldsymbol {\imath }}$, with a mirror image at $\boldsymbol {X}=-\hat {\boldsymbol {\imath }}$:

(4.5)\begin{align} \boldsymbol{U}_0 &= \frac{1}{4{\rm \pi}}\left[-{\boldsymbol{\mathsf{I}}}\ln |\boldsymbol{X}-\hat{\boldsymbol{\imath}}| +\frac{(\boldsymbol{X}-\hat{\boldsymbol{\imath}})(\boldsymbol{X}-\hat{\boldsymbol{\imath}})}{|\boldsymbol{X}-\hat{\boldsymbol{\imath}}|^2}\right]\boldsymbol{\cdot}\hat{\boldsymbol{\imath}}\nonumber\\ &\quad -\frac{1}{4{\rm \pi}}\left[-{\boldsymbol{\mathsf{I}}}\ln |\boldsymbol{X}+\hat{\boldsymbol{\imath}}| +\frac{(\boldsymbol{X}+\hat{\boldsymbol{\imath}})(\boldsymbol{X}+\hat{\boldsymbol{\imath}})}{|\boldsymbol{X}+\hat{\boldsymbol{\imath}}|^2}\right]\boldsymbol{\cdot}\hat{\boldsymbol{\imath}}. \end{align}

Rewriting in terms of the inner coordinates and expanding for small $\lambda$ gives

(4.6)\begin{align} \boldsymbol{U}_0 &\sim \frac{1}{4{\rm \pi}}\,[-\hat{\boldsymbol{\imath}}\ln(\lambda\rho) +\hat{\boldsymbol{e}}_\rho\cos\phi + \hat{\boldsymbol{\imath}}(\ln2-1)]\nonumber\\ &\quad - \frac{\lambda\rho}{8{\rm \pi}}\,(\hat{\boldsymbol{e}}_\rho-2\hat{\boldsymbol{\imath}}\cos\phi) -\frac{\lambda^2\rho^2}{32{\rm \pi}}\, (3\hat{\boldsymbol{\imath}}\cos2\phi-2\hat{\boldsymbol{e}}_\rho\cos\phi) + \cdots. \end{align}

Imposing $\operatorname {ord}(1)$–$\operatorname {ord}(1)$ Van Dyke matching using (4.1) and (4.6) readily gives

(4.7)\begin{equation} {\mathcal{U}}_0 = \frac{1}{4{\rm \pi}}\left(\ln\frac{2}{\lambda}-\frac{1}{2}\right).\end{equation}

Given (3.1), ${\mathcal {U}}_0$ is positive, as expected.

5. Asymptotic corrections

The only term of $\boldsymbol {u}_0$ that is not accounted for in (4.5) is the doublet, which decays as $\rho ^{-2}$. Moreover, the $\operatorname {ord}(\lambda )$ term in expansion (4.6) does not incorporate a uniform flow. These observations suggest that $\boldsymbol {U}_1\equiv \boldsymbol {0}$ and then ${\mathcal {U}}_1=0$. To determine $\boldsymbol {u}_1$, we note that $\operatorname {ord}(\lambda )$–$\operatorname {ord}(\lambda )$ Van Dyke matching using (4.6) implies

(5.1)\begin{equation} \boldsymbol{u}_1 \sim \frac{\rho}{8{\rm \pi}}\,(\hat{\boldsymbol{e}}_\rho-2\hat{\boldsymbol{\imath}}\cos\phi) \quad \text{as } \rho\to\infty, \end{equation}

or equivalently,

(5.2)\begin{equation} \psi_1 \sim \frac{\rho^2}{16{\rm \pi}}\sin2\phi \quad \text{as } \rho\to\infty.\end{equation}

The correction $\psi _1$ is governed by that condition, together with the symmetry constraint (3.8), the force-free requirement (recall (3.5)) and the condition (recall (3.9))

(5.3)\begin{equation} \psi_1 =\frac{\partial{\psi_1}}{\partial{\rho}} =0 \quad \text{at } \rho=1.\end{equation}

The unique biharmonic function that satisfies these conditions is

(5.4)\begin{equation} \psi_1 = \frac{\rho^2-2+\rho^{{-}2}}{16{\rm \pi}}\sin2\phi.\end{equation}

Performing $\operatorname {ord}(\lambda )$–$\operatorname {ord}(\lambda ^2)$ Van Dyke matching implies that as $|\boldsymbol {X}-\hat {\boldsymbol {\imath }}|\to 0$,

(5.5)\begin{align} \boldsymbol{U}_2 &\sim{-}\frac{1}{8{\rm \pi}}\left[ -\frac{{\boldsymbol{\mathsf{I}}}}{|\boldsymbol{X}-\hat{\boldsymbol{\imath}}|^2} + 2\,\frac{(\boldsymbol{X}-\hat{\boldsymbol{\imath}})(\boldsymbol{X}-\hat{\boldsymbol{\imath}})}{|\boldsymbol{X}-\hat{\boldsymbol{\imath}}|^4}\right] \boldsymbol{\cdot}\hat{\boldsymbol{\imath}} \nonumber\\ &\quad - \frac{1}{4{\rm \pi}}\, (\hat{\boldsymbol{\imath}}\hat{\boldsymbol{\imath}}-\hat{\boldsymbol{\jmath}}\hat{\boldsymbol{\jmath}}) \boldsymbol{:}\frac{(\boldsymbol{X}-\hat{\boldsymbol{\imath}})(\boldsymbol{X}-\hat{\boldsymbol{\imath}})(\boldsymbol{X}-\hat{\boldsymbol{\imath}})}{|\boldsymbol{X}-\hat{\boldsymbol{\imath}}|^4}, \end{align}

where the first term is associated with the doublet in (4.1), and the second term (wherein $\boldsymbol {:}$ denotes tensor contraction) is associated with the stresslet in (5.4). The Stokes flow $\boldsymbol {U}_2$ that satisfies (5.5) together with (3.11)–(3.12) is simply a superposition of (5.5) with its mirror image:

(5.6)\begin{align} \boldsymbol{U}_2 &= \frac{1}{8{\rm \pi}}\left[ \frac{{\boldsymbol{\mathsf{I}}}}{|\boldsymbol{X}-\hat{\boldsymbol{\imath}}|^2} - 2\,\frac{(\boldsymbol{X}-\hat{\boldsymbol{\imath}})(\boldsymbol{X}-\hat{\boldsymbol{\imath}})}{|\boldsymbol{X}-\hat{\boldsymbol{\imath}}|^4} -\frac{{\boldsymbol{\mathsf{I}}}}{|\boldsymbol{X}+\hat{\boldsymbol{\imath}}|^2} + 2\,\frac{(\boldsymbol{X}+\hat{\boldsymbol{\imath}})(\boldsymbol{X}+\hat{\boldsymbol{\imath}})}{|\boldsymbol{X}+\hat{\boldsymbol{\imath}}|^4}\right]\boldsymbol{\cdot}\hat{\boldsymbol{\imath}} \nonumber\\ &\quad - \frac{1}{4{\rm \pi}}\, (\hat{\boldsymbol{\imath}}\hat{\boldsymbol{\imath}}-\hat{\boldsymbol{\jmath}}\hat{\boldsymbol{\jmath}}) \boldsymbol{:} \left[ \frac{(\boldsymbol{X}-\hat{\boldsymbol{\imath}})(\boldsymbol{X}-\hat{\boldsymbol{\imath}})(\boldsymbol{X}-\hat{\boldsymbol{\imath}})}{|\boldsymbol{X}-\hat{\boldsymbol{\imath}}|^4} +\frac{(\boldsymbol{X}+\hat{\boldsymbol{\imath}})(\boldsymbol{X}+\hat{\boldsymbol{\imath}})(\boldsymbol{X}+\hat{\boldsymbol{\imath}})}{|\boldsymbol{X}+\hat{\boldsymbol{\imath}}|^4}\right]. \end{align}

Note that expansion of (5.6) about $\boldsymbol {X}=\hat {\boldsymbol {\imath }}$ produces, inter alia, an $\operatorname {ord}(1)$ uniform flow term. This suggests that ${\mathcal {U}}_2$ is non-zero.

We can now calculate $\boldsymbol {u}_2$ and ${\mathcal {U}}_2$. Performing $\operatorname {ord}(\lambda ^2)$–$\operatorname {ord}(\lambda ^2)$ Van Dyke matching, we obtain

(5.7)\begin{equation} \boldsymbol{u}_2 \sim{-}\frac{\rho^2}{32{\rm \pi}}\, (3\hat{\boldsymbol{\imath}}\cos2\phi-2\hat{\boldsymbol{e}}_\rho\cos\phi) - \frac{3}{32{\rm \pi}}\,\hat{\boldsymbol{\imath}} \quad \text{as } \rho\to\infty, \end{equation}

where the first term follows from the last term in (4.6), and the second term is the aforementioned uniform flow. In terms of $\psi _2$, (5.7) reads

(5.8)\begin{equation} \psi_2 \sim \frac{\rho^3}{64{\rm \pi}}\,(\sin\phi-\sin3\phi) - \frac{3\rho}{32{\rm \pi}}\sin\phi \quad \text{as } \rho\to\infty.\end{equation}

The most general biharmonic function that satisfies (5.8) together with the symmetry constraint (3.8) and the force-free requirement (which excludes terms proportional to $\rho \ln \rho$) is

(5.9)\begin{equation} \psi_2 = \left(\frac{\rho^3}{64{\rm \pi}} - \frac{3\rho}{32{\rm \pi}} + \frac{d_2}{2{\rm \pi}\rho}\right)\sin\phi - \frac{1}{64{\rm \pi}}\left(\rho^3 + \frac{g_2}{\rho} + \frac{h_2}{\rho^3}\right)\sin3\phi. \end{equation}

Applying the conditions (cf. (3.9))

(5.10)\begin{equation} \left. \begin{array}{@{}l@{}} \psi_2 = {\mathcal{U}}_2 \sin\phi \\[2pt] \displaystyle \dfrac{\partial{\psi_2}}{\partial{\rho}} ={\mathcal{U}}_2\sin\phi \end{array} \right\} \quad \text{at } \rho=1\end{equation}

gives ${\mathcal {U}}_2 = -1/16{\rm \pi}$ (as well as $d_2=1/32$, $g_2=-3$ and $h_2=2$).

We conclude that

(5.11)\begin{equation} {\mathcal{U}} \sim \frac{1}{4{\rm \pi}}\left(\ln\frac{2}{\lambda}-\frac{1}{2}\right)- \frac{\lambda^2}{16{\rm \pi}} + O(\lambda^3) \quad \text{for } \lambda\ll1. \end{equation}

6. Shear-free interface

As another illustration of the present methodology, we consider now the case where the particle boundaries are shear-free surfaces – an alternative model of membrane-trapped colloids (Saffman & Delbrück 1975). The no-slip condition (3.4) is then replaced by the shear-free condition

(6.1)\begin{equation} \left. \begin{array}{@{}l@{}} u = {\mathcal{U}}\cos\phi \\[2pt] \displaystyle\dfrac{\partial{ v}}{\partial{\rho}}-v ={\mathcal{U}}\sin\phi \end{array} \right\} \quad \text{at } \rho=1,\end{equation}

or equivalently (cf. (3.9)),

(6.2)\begin{equation} \left. \begin{array}{@{}l@{}} \psi = {\mathcal{U}} \sin\phi \\[2pt] \displaystyle\dfrac{\partial{^2\psi}}{\partial{\rho^2}} - \dfrac{\partial{\psi}}{\partial{\rho}} ={-} {\mathcal{U}} \sin\phi \end{array} \right\} \quad \text{at } \rho=1.\end{equation}

At leading order, we find that (cf. (4.3))

(6.3)\begin{equation} \psi_0 =\frac{\rho}{4{\rm \pi}}\,(4{\rm \pi}{\mathcal{U}}_0-\ln\rho) \sin\phi.\end{equation}

The associated velocity field is (cf. (4.2))

(6.4)\begin{equation} \boldsymbol{u}_0 = {\mathcal{U}}_0 \hat{\boldsymbol{\imath}} + \frac{1}{4{\rm \pi}}\,(-{\boldsymbol{\mathsf{I}}}\ln\rho+\hat{\boldsymbol{e}}_\rho\hat{\boldsymbol{e}}_\rho - {\boldsymbol{\mathsf{I}}})\boldsymbol{\cdot}\hat{\boldsymbol{\imath}} . \end{equation}

Clearly, the leading-order outer solution (4.5) is unaltered; in particular, (4.6) remains valid. Here, $\operatorname {ord}(1)$–$\operatorname {ord}(1)$ Van Dyke matching gives

(6.5)\begin{equation} {\mathcal{U}}_0 = \frac{1}{4{\rm \pi}}\ln\frac{2}{\lambda}.\end{equation}

Comparing with (4.7), a bubble drifts faster than a rigid particle. This is hardly surprising, as it would experience a smaller drag for a given velocity.

In calculating $\boldsymbol {u}_1$, we see that the matching condition (5.2) still holds. The correction $\psi _1$ is governed by that condition, together with the symmetry constraint (3.8), the force-free requirement (recall (3.5)) and the slip condition (recall (6.2))

(6.6)\begin{equation} \left. \begin{array}{@{}l@{}} \psi_1 = 0 \\[2pt] \displaystyle\dfrac{\partial{^2\psi_1}}{\partial{\rho^2}} - \dfrac{\partial{\psi_1}}{\partial{\rho}}= 0 \end{array} \right\} \quad \text{at } \rho=1.\end{equation}

The unique biharmonic function that satisfies these conditions is (cf. (5.4))

(6.7)\begin{equation} \psi_1 = \frac{\rho^2-1}{16{\rm \pi}}\sin2\phi.\end{equation}

In the absence of a doublet in (6.3), $\operatorname {ord}(\lambda )$–$\operatorname {ord}(\lambda ^2)$ Van Dyke matching implies that (cf. (5.5))

(6.8)\begin{equation} \boldsymbol{U}_2 \sim{-} \frac{1}{8{\rm \pi}}\, (\hat{\boldsymbol{\imath}}\hat{\boldsymbol{\imath}}-\hat{\boldsymbol{\jmath}}\hat{\boldsymbol{\jmath}}) \boldsymbol{:}\frac{(\boldsymbol{X}-\hat{\boldsymbol{\imath}})(\boldsymbol{X}-\hat{\boldsymbol{\imath}})(\boldsymbol{X}-\hat{\boldsymbol{\imath}})}{|\boldsymbol{X}-\hat{\boldsymbol{\imath}}|^4} \quad \text{as } |\boldsymbol{X}-\hat{\boldsymbol{\imath}}|\to 0. \end{equation}

Thus $\boldsymbol {U}_2$ consists only of stresslets (cf. (5.6)):

(6.9)\begin{equation} \boldsymbol{U}_2 ={-} \frac{1}{8{\rm \pi}}\, (\hat{\boldsymbol{\imath}}\hat{\boldsymbol{\imath}}-\hat{\boldsymbol{\jmath}}\hat{\boldsymbol{\jmath}}) \boldsymbol{:} \left[ \frac{(\boldsymbol{X}-\hat{\boldsymbol{\imath}})(\boldsymbol{X}-\hat{\boldsymbol{\imath}})(\boldsymbol{X}-\hat{\boldsymbol{\imath}})}{|\boldsymbol{X}-\hat{\boldsymbol{\imath}}|^4} +\frac{(\boldsymbol{X}+\hat{\boldsymbol{\imath}})(\boldsymbol{X}+\hat{\boldsymbol{\imath}})(\boldsymbol{X}+\hat{\boldsymbol{\imath}})}{|\boldsymbol{X}+\hat{\boldsymbol{\imath}}|^4}\right].\end{equation}

We can now calculate $\boldsymbol {u}_2$ and ${\mathcal {U}}_2$. Performing $\operatorname {ord}(\lambda ^2)$–$\operatorname {ord}(\lambda ^2)$ Van Dyke matching (cf. (5.8)) gives here

(6.10)\begin{equation} \psi_2 \sim \frac{\rho^3}{64{\rm \pi}}\,(\sin\phi-\sin3\phi) - \frac{\rho}{16{\rm \pi}}\sin\phi \quad \text{as } \rho\to\infty, \end{equation}

where the first term follows from (4.6), and the second term follows from expanding (6.9). The most general biharmonic function that satisfies (6.10) as well as the symmetry constraint (3.8) and force-free requirement is (cf. (5.9))

(6.11)\begin{equation} \psi_2 = \left(\frac{\rho^3}{64{\rm \pi}} - \frac{\rho}{16{\rm \pi}} + \frac{d_2}{2{\rm \pi}\rho}\right)\sin\phi - \frac{1}{64{\rm \pi}}\left(\rho^3 + \frac{g_2}{\rho} + \frac{h_2}{\rho^3}\right)\sin3\phi. \end{equation}

Applying the conditions (cf. (6.2))

(6.12)\begin{equation} \left. \begin{array}{@{}l@{}} \psi_2 = {\mathcal{U}}_2 \sin\phi \\[2pt] \displaystyle\dfrac{\partial{^2\psi_2}}{\partial{\rho^2}} - \dfrac{\partial{\psi_2}}{\partial{\rho}} ={-} {\mathcal{U}}_2 \sin\phi \end{array} \right\} \quad \text{at } \rho=1\end{equation}

yields ${\mathcal {U}}_2 = -1/16{\rm \pi}$ (incidentally, just as in the case of a rigid boundary), $d_2=1/16$, $g_2=-3$ and $h_2=2$. We conclude that

(6.13)\begin{equation} {\mathcal{U}} \sim \frac{1}{4{\rm \pi}}\ln\frac{2}{\lambda}- \frac{\lambda^2}{16{\rm \pi}} + O(\lambda^3) \quad \text{for } \lambda\ll1. \end{equation}

7. Elliptic disks

In principle, the present asymptotic scheme may be applied to particles of non-circular cross-section. The leading-order outer problem is unaffected, whereby (4.5) remains valid. The leading-order inner problem, on the other hand, depends upon the specific disk shape. (In that sense, the present scenario is somewhat reminiscent of a small-Reynolds-number analysis of a single elliptic particle, see Kropinski, Ward & Keller 1995.) To simplify the calculation, we now impose the force constraint (2.3) in that inner problem via the far-field specification (cf. (4.4))

(7.1)\begin{equation} \psi_0 \sim{-}\frac{\rho\ln\rho}{4{\rm \pi}}\sin\phi \quad \text{as } \rho\to\infty.\end{equation}

We here illustrate the generalisation to non-circular shapes, considering disks of elliptic cross-section. We first consider the symmetric case where one of the ellipse axes is aligned with the line of centres. The associated (dimensional) semi-axis is denoted by $a$; the other semi-axis is denoted by $b$. We choose the length scale $\ell$ as $(a+b)/2$. For simplicity, we consider here only rigid disks, where the no-slip condition applies.

In allowing for all possible ratios $a/b$, we need to address separately the cases $a>b$ and $a< b$.

7.1. Case $a>b$

For $a>b$, the dimensional distance from the origin to the focal points is $c=\sqrt {a^2-b^2}$. The associated dimensionless distance is

(7.2)\begin{equation} \tilde c=2\,\sqrt{\frac{a-b}{a+b}}. \end{equation}

It is natural to employ the elliptic coordinates $(\xi,\eta )$ defined by (Moon & Spencer 1988):

(7.3a,b)\begin{equation} x = \tilde c \cosh\xi\cos\eta, \quad y = \tilde c \sinh\xi\sin\eta. \end{equation}

The curves $\xi =\textrm {const.}$ are ellipses. In particular, the disk boundary is $\xi =\xi ^*$, where

(7.4)\begin{equation} \tilde c = 2\,\mathrm{e}^{-\xi^*} . \end{equation}

Note that the eccentricity $\mathcal {E}=c/a$ is given by

(7.5)\begin{equation} \mathcal{E} = \operatorname{sech}\xi^*. \end{equation}

We also note that at large distances, where $\xi$ is large,

(7.6a,b)\begin{equation} \rho \sim \mathrm{e}^{\xi-\xi^*}, \quad \phi \sim \eta.\end{equation}

In prescribing the boundary conditions governing the flow problem, we temporarily consider a co-moving reference frame, where the disk is stationary. The boundary conditions on the stream function in that frame, say $\tilde \psi$, are

(7.7)\begin{equation} \tilde\psi=\frac{\partial{\tilde\psi}}{\partial{\xi}}=0 \quad \text{at } \xi = \xi^*. \end{equation}

The most general biharmonic function that satisfies (7.7) together with the symmetry (3.8) and does not diverge more rapidly than (7.1) is (Shintani, Umemura & Takano 1983)

(7.8)\begin{equation} \tilde\psi = \tilde D \{ (\xi-\xi^*)\sinh\xi - \sinh\xi^*\cosh\xi^*\sinh\xi +\sinh^2\xi^* \cosh\xi\}\sin\eta. \end{equation}

The appropriate stream function in the laboratory reference frame is

(7.9)\begin{equation} \psi_0=\tilde\psi +{\mathcal{U}}_0 y. \end{equation}

Imposing the force constraint (7.1) using (7.4)–(7.6), we find that

(7.10)\begin{equation} \tilde D ={-}\frac{\mathrm{e}^{-\xi^*}}{4{\rm \pi}}, \end{equation}

whereby

(7.11)\begin{equation} \psi_0 \sim{-} \frac{\rho}{4{\rm \pi}}\,( \ln\rho - \mathrm{e}^{-\xi^*}\sinh\xi^* -4{\rm \pi}{\mathcal{U}}_0)\sin\phi \quad \text{as } \rho\to\infty. \end{equation}

This refinement of (7.1) constitutes the large-$\rho$ leading-order behaviour of $\psi _0$, in the aforementioned convention that abides by Fraenkel's warning.

7.2. Case $a< b$

For $a< b$, we define $c=\sqrt {b^2-a^2}$. The associated dimensionless distance is

(7.12)\begin{equation} \tilde c=2\,\sqrt{\frac{b-a}{b+a}}. \end{equation}

The elliptic coordinates are now defined by (cf. (7.3))

(7.13a,b)\begin{equation} y = \tilde c \cosh\xi\cos\eta, \quad -x = \tilde c \sinh\xi\sin\eta.\end{equation}

Relation (7.4) remains valid. The eccentricity, now defined as $\mathcal {E}=c/b$, is again given by (7.5). Here, at large distances (cf. (7.6)),

(7.14a,b)\begin{equation} \rho \sim \mathrm{e}^{\xi-\xi^*}, \quad \phi \sim {\rm \pi}/2 +\eta.\end{equation}

We note that conditions (7.7) still hold. The most general biharmonic function that satisfies (7.7) together with symmetry (3.8) and does not diverge more rapidly than a Stokeslet is (cf. (7.8))

(7.15)\begin{equation} \tilde\psi = \tilde D \{ (\xi-\xi^*)\cosh\xi + \sinh\xi^*\cosh\xi^*\cosh\xi -\cosh^2\xi^* \sinh\xi\}\cos\eta. \end{equation}

The laboratory-frame stream function $\psi _0$ is still given by (7.9). Considering large distances using (7.14), we find that (7.10) remains valid, whereby (cf. (7.11))

(7.16)\begin{equation} \psi_0 \sim{-} \frac{\rho}{4{\rm \pi}}\,( \ln\rho - \mathrm{e}^{-\xi^*}\cosh\xi^*-4{\rm \pi}{\mathcal{U}}_0)\sin\phi \quad \text{as } \rho\to\infty. \end{equation}

7.3. Particle velocity

With the leading-order inner flow available (up to the particle velocity), we can now employ asymptotic matching to obtain that velocity. Thus, making use of (7.11) and (7.16) as well as expansion (4.6) of the outer flow, we obtain, using Van Dyke matching,

(7.17)\begin{equation} {4{\rm \pi}}{\mathcal{U}}_0 =\ln\frac{1}{\lambda} + \mathcal{W}, \end{equation}

where

(7.18)\begin{equation} {\mathcal{W}} = \ln2 - \begin{cases} \mathrm{e}^{-\xi^*}\sinh\xi^*, & a>b, \\ \mathrm{e}^{-\xi^*}\cosh\xi^*, & a< b. \end{cases} \end{equation}

The case of a circle ($a=b$), where $\mathcal {E}=0$, is represented by the limit $\xi ^*\to \infty$; see (7.5). In that limit, we find ${\mathcal {W}} =\ln 2-1/2$, in agreement with (4.7). In the limit of a line segment ($\mathcal {E}=1$), where $\xi ^*=0$, we obtain

(7.19)\begin{equation} {\mathcal{W}} = \begin{cases} \ln 2, & b=0, \\ \ln2-1, & a=0. \end{cases} \end{equation}

Unsurprisingly, the velocity of a segment that is aligned with the line of centres is larger than that of a segment perpendicular to it.

Result (7.17) represents the particle velocity as a superposition of an ‘interaction velocity’, which depends only upon the normalised separation $1/\lambda$, and a ‘renormalised’ velocity, which depends only upon the eccentricity $\mathcal {E}$ via (7.5). The singularity of the interaction velocity in the limit $\lambda \to 0$ is a manifestation of the Stokes paradox: the single-particle problem is ill posed. In figure 2, we present the renormalised velocity as a function of eccentricity for both $a>b$ and $a< b$ (cf. Kropinski et al. 1995). Recall that $\mathcal {E}=\sqrt {|a^2-b^2|}/{\max (a,b)}$.

Figure 2. Renormalised velocity as a function of eccentricity, obtained from (7.18) using (7.5), for both $a>b$ and $a< b$.

7.4. Uniform mobility expression and extension to non-aligned ellipses

Making use of (7.2) and (7.4), we find that for $a>b$, $\mathrm {e}^{-\xi ^*}\sinh \xi ^*$ is simply given by $b/(a+b)$. Similarly, using (7.12) and (7.4), we find that for $a< b$, $\mathrm {e}^{-\xi ^*}\cosh \xi ^*$ is also given by $b/(a+b)$. We conclude that (7.17)–(7.18) may be combined to the simple formula

(7.20)\begin{equation} {4{\rm \pi}}{\mathcal{U}}_0 = \ln\frac{2}{\lambda} -\frac{b}{a+b}, \end{equation}

valid for all $a$ and $b$. For future reference, we note that the large-$\rho$ approximations (7.11) and (7.16) may be combined to

(7.21)\begin{equation} \psi_0 \sim{-} \frac{\rho}{4{\rm \pi}} \left( \ln\rho - \frac{b}{a+b} -4{\rm \pi}{\mathcal{U}}_0\right)\sin\phi \quad \text{as } \rho\to\infty. \end{equation}

The associated velocity field is (cf. (4.4))

(7.22)\begin{equation} \boldsymbol{u}_0 \sim {\mathcal{U}}_0\hat{\boldsymbol{\imath}} +\frac{1}{4{\rm \pi}}\left[ (-{\boldsymbol{\mathsf{I}}}\ln\rho + \hat{\boldsymbol{e}}_\rho\hat{\boldsymbol{e}}_\rho)\boldsymbol{\cdot}\hat{\boldsymbol{\imath}} - \frac{a}{a+b}\,\hat{\boldsymbol{\imath}}\right] \quad \text{as } \rho\to\infty.\end{equation}

It is now easy to generalise for the case of a non-aligned ellipse, for which the velocity $\boldsymbol {\mathcal {U}}_0$ of the ‘+’ ellipse is not necessarily directed along the $x'$-axis. Denoting by $\hat {\boldsymbol {e}}_a$ a unit vector parallel to the $2a$-axis, and by $\hat {\boldsymbol {e}}_b$ a unit vector parallel to the $2b$-axis, it is evident from symmetry arguments that the generalisation of (7.22) is

(7.23)\begin{equation} \boldsymbol{u}_0 \sim \boldsymbol{\mathcal{U}}_0 +\frac{1}{4{\rm \pi}}\,(-{\boldsymbol{\mathsf{I}}}\ln\rho + \hat{\boldsymbol{e}}_\rho\hat{\boldsymbol{e}}_\rho - {\boldsymbol{\mathsf{M}}})\boldsymbol{\cdot}\hat{\boldsymbol{\imath}} \quad \text{as } \rho\to\infty,\end{equation}

wherein

(7.24)\begin{equation} {\boldsymbol{\mathsf{M}}} = \frac{a}{a+b}\,\hat{\boldsymbol{e}}_a\hat{\boldsymbol{e}}_a + \frac{b}{a+b}\,\hat{\boldsymbol{e}}_b\hat{\boldsymbol{e}}_b. \end{equation}

Matching with (4.6) therefore gives

(7.25)\begin{equation} {4{\rm \pi}}\boldsymbol{\mathcal{U}}_0= \left(\ln\frac{2}{\lambda} -1\right)\hat{\boldsymbol{\imath}} +{\boldsymbol{\mathsf{M}}}\boldsymbol{\cdot}\hat{\boldsymbol{\imath}}. \end{equation}

Writing $\boldsymbol {\mathcal {U}}_0={\mathcal {U}}_0\hat {\boldsymbol {\imath }} + {\mathcal {V}}_0\hat {\boldsymbol {\jmath }}$, we obtain

(7.26a,b)\begin{align} {4{\rm \pi}}{\mathcal{U}}_0=\ln\frac{2}{\lambda} -1 + \frac{a}{a+b} \cos^2\alpha + \frac{b}{a+b}\sin^2\alpha, \quad {4{\rm \pi}} {\mathcal{V}}_0=\frac{a-b}{a+b} \cos\alpha\sin\alpha,\end{align}

wherein $\alpha$ is the angle from the $x'$-axis to the $2a$-axis, reckoned positive in the anticlockwise direction.

8. Generalisation to a collection of interacting particles

At this stage, our work left unanswered some natural questions to ask, such as how the result would change if the disks are not identical, if the forces are not centrally symmetric, or indeed if there are multiple disks. It turns out that these questions can be addressed, at least to leading order, by considering a generalisation to $N$ well-separated particles, which could be of arbitrary shape, size, orientation and surface properties (i.e. no-slip or shear-free interface).

In this general scenario, a dimensionless notation is non-beneficial, so we resort to a dimensional description using the position vector $\boldsymbol {r}$. The outer velocity field, driven by Stokeslets, is simply

(8.1)\begin{equation} \sum_{n=1}^N {\boldsymbol{\mathsf{J}}} (\boldsymbol{r}-\boldsymbol{r}^{(n)})\boldsymbol{\cdot} \boldsymbol{F}^{(n)}, \end{equation}

wherein $\boldsymbol {r}^{(n)}$ is the position of the centroid of the $n$th particle,

(8.2)\begin{equation} {\boldsymbol{\mathsf{J}}} (\boldsymbol{r}) = \frac{1}{4{\rm \pi}\mu}\left(-{\boldsymbol{\mathsf{I}}} \ln\frac{|\boldsymbol{r}|}{\mathcal{R}} + \frac{\boldsymbol{r}\boldsymbol{r}}{|\boldsymbol{r}|^2}\right) \end{equation}

is the Oseen–Burgers tensor, where a length scale $\mathcal {R}$ (e.g. a characteristic separation distance) is introduced, and $\boldsymbol {F}^{(n)}$ is the force on the $n$th particle due to its (non-hydrodynamic) interaction with its neighbours. Since Newton's third law necessitates

(8.3)\begin{equation} \sum_{n=1}^N\boldsymbol{F}^{(n)}=\boldsymbol{0}, \end{equation}

approximation (8.1) is actually independent of the arbitrary scale $\mathcal {R}$. Moreover, it follows from (8.3) that (8.1) approaches zero as $|\boldsymbol {r}|\to \infty$, as required (cf. (2.2)).

To approximate (8.1) near one of the particles, say particle ‘1’, we write $\boldsymbol {r}=\boldsymbol {r}^{(1)}+\boldsymbol {\rho }$ and expand for small $|\boldsymbol {\rho }|$ to obtain

(8.4)\begin{equation} \frac{1}{4{\rm \pi}\mu}\left(-{\boldsymbol{\mathsf{I}}} \ln\frac{|\boldsymbol{\rho}|}{\mathcal{R}} + \hat{\boldsymbol{e}}_\rho\hat{\boldsymbol{e}}_\rho\right)\boldsymbol{\cdot}\boldsymbol{F}^{(1)} + \sum_{n=2}^N {\boldsymbol{\mathsf{J}}} (\boldsymbol{r}^{(1)}-\boldsymbol{r}^{(n)})\boldsymbol{\cdot} \boldsymbol{F}^{(n)}, \end{equation}

wherein $\hat {\boldsymbol {e}}_\rho = \boldsymbol {\rho }/|\boldsymbol {\rho }|$, consistently with our earlier notation. To obtain the rectilinear velocity $\boldsymbol {V}^{(1)}$ of that ‘test particle’, we need a far-field approximation of the leading-order particle-scale (‘inner’) flow field. This generally necessitates the solution of the inner flow about a translating particle. (In general, the test particle also acquires an angular velocity to satisfy the torque-free condition; the flow associated with it, however, decays algebraically with distance from the particle, and does not affect the requisite approximation.)

For a no-slip circular disk, say of radius $\ell$, this approximation is simply (cf. (4.4))

(8.5)\begin{equation} \boldsymbol{V}^{(1)} + \frac{1}{4{\rm \pi}\mu}\left(-{\boldsymbol{\mathsf{I}}} \ln\frac{|\boldsymbol{\rho}|}{\ell} + \hat{\boldsymbol{e}}_\rho\hat{\boldsymbol{e}}_\rho - \frac{1}{2}\,{\boldsymbol{\mathsf{I}}}\right) \boldsymbol{\cdot} \boldsymbol{F}^{(1)}.\end{equation}

Asymptotic matching with (8.4) yields

(8.6)\begin{equation} \boldsymbol{V}^{(1)} = \frac{1}{4{\rm \pi}\mu} \left(\ln\frac{\mathcal{R}}{\ell}+ \frac{1}{2}\right) \boldsymbol{F}^{(1)} + \sum_{n=2}^N {\boldsymbol{\mathsf{J}}} (\boldsymbol{r}^{(1)}-\boldsymbol{r}^{(n)})\boldsymbol{\cdot} \boldsymbol{F}^{(n)}.\end{equation}

In the case of a shear-free disk, (8.5) is replaced by (cf. (6.4))

(8.7)\begin{equation} \boldsymbol{V}^{(1)} + \frac{1}{4{\rm \pi}\mu}\left(-{\boldsymbol{\mathsf{I}}} \ln\frac{|\boldsymbol{\rho}|}{\ell} + \hat{\boldsymbol{e}}_\rho\hat{\boldsymbol{e}}_\rho - {\boldsymbol{\mathsf{I}}}\right) \boldsymbol{\cdot} \boldsymbol{F}^{(1)},\end{equation}

whereby matching with (8.4) gives

(8.8)\begin{equation} \boldsymbol{V}^{(1)} = \frac{1}{4{\rm \pi}\mu} \left(\ln\frac{\mathcal{R}}{\ell}+ 1\right) \boldsymbol{F}^{(1)} + \sum_{n=2}^N {\boldsymbol{\mathsf{J}}} (\boldsymbol{r}^{(1)}-\boldsymbol{r}^{(n)})\boldsymbol{\cdot} \boldsymbol{F}^{(n)}. \end{equation}

It may appear as though (8.6) and (8.8) represent a superposition of a single-particle solution and advection by the other particles. This interpretation, however, is somewhat misleading given the dependence upon the arbitrary scale $\mathcal {R}$, which cancels out only when all terms are added together. Nonetheless, it is evident from the above matching procedure that the consideration of a more complicated test particle would retain the structure of (8.6) and (8.8), where the ‘single-particle’ velocity is modified but the ‘advection’ term remains unaltered. It follows that the leading-order velocity of a particle depends (through the single-particle term) on its own shape, size, orientation and surface properties (no-slip or shear-free), but not on the details of the other particles beyond their position and the interaction forces that they experience.

For the case of two particles with separation $2\boldsymbol {s}=\boldsymbol {r}^{(1)}-\boldsymbol {r}^{(2)}$, Newton's third law (8.3) gives $\boldsymbol {F}\stackrel {\mathrm {def}}=\boldsymbol {F}^{(1)}=-\boldsymbol {F}^{(2)}$. For a no-slip surface, we obtain from (8.6) the rectilinear velocity $\boldsymbol {V}\stackrel {\mathrm {def}}=\boldsymbol {V}^{(1)}$ of the test particle as

(8.9)\begin{equation} \boldsymbol{V}=\frac{1}{4{\rm \pi}\mu} \left[\left(\ln\frac{2|\boldsymbol{s}|}{\ell} +\frac{1}{2}\right)\boldsymbol{F} - \frac{\boldsymbol{s}\boldsymbol{s}}{|\boldsymbol{s}|^2} \boldsymbol{\cdot}\boldsymbol{F}\right], \end{equation}

which for $\boldsymbol {F}\|\boldsymbol {s}$ and $\boldsymbol {F}\bot \boldsymbol {s}$ reduces to

(8.10a,b)\begin{equation} \boldsymbol{V}=\frac{\boldsymbol{F}}{4{\rm \pi}\mu} \left(\ln\frac{2|\boldsymbol{s}|}{\ell} -\frac{1}{2}\right), \quad \boldsymbol{V}=\frac{\boldsymbol{F}}{4{\rm \pi}\mu} \left(\ln\frac{2|\boldsymbol{s}|}{\ell} +\frac{1}{2}\right),\end{equation}

respectively, with (8.10a) in agreement with (4.7). For shear-free surfaces, we obtain from (8.8) that

(8.11)\begin{equation} \boldsymbol{V}=\frac{1}{4{\rm \pi}\mu} \left[\left(\ln\frac{2|\boldsymbol{s}|}{\ell} +1\right)\boldsymbol{F} - \frac{\boldsymbol{s}\boldsymbol{s}}{|\boldsymbol{s}|^2} \boldsymbol{\cdot}\boldsymbol{F}\right], \end{equation}

which for $\boldsymbol {F}\|\boldsymbol {s}$ and $\boldsymbol {F}\bot \boldsymbol {s}$ reduces to

(8.12a,b)\begin{equation} \boldsymbol{V}=\frac{\boldsymbol{F}}{4{\rm \pi}\mu} \ln\frac{2|\boldsymbol{s}|}{\ell}, \quad \boldsymbol{V}=\frac{\boldsymbol{F}}{4{\rm \pi}\mu} \left(\ln\frac{2|\boldsymbol{s}|}{\ell} +1\right),\end{equation}

respectively, with (8.12a) in agreement with (6.5). Note that both (8.9) and (8.11) are indeed independent of $\mathcal {R}$, that the velocity of a shear-free disk is larger than that of a no-slip disk, and that the velocity in the transverse case $\boldsymbol {F}\bot \boldsymbol {s}$ is larger than that in the longitudinal case $\boldsymbol {F}\|\boldsymbol {s}$.

We may also consider the case where the test particle is an elliptic disk of major axis $a$ and minor axis $b$, with corresponding unit vectors $\hat {\boldsymbol {e}}_a$ and $\hat {\boldsymbol {e}}_b$. Using the analysis of § 7, the inner flow at large distances from the disk is (cf. (7.23))

(8.13)\begin{equation} \boldsymbol{V}^{(1)} +\frac{1}{4{\rm \pi}\mu} \left(-{\boldsymbol{\mathsf{I}}}\ln\frac{2|\boldsymbol{\rho}|}{a+b} + \hat{\boldsymbol{e}}_\rho\hat{\boldsymbol{e}}_\rho - {\boldsymbol{\mathsf{M}}}\right)\boldsymbol{\cdot}\boldsymbol{F}^{(1)} \quad \text{as } \rho\to\infty .\end{equation}

Asymptotic matching with (8.4) yields

(8.14)\begin{equation} \boldsymbol{V}^{(1)} = \frac{1}{4{\rm \pi}\mu} \left({\boldsymbol{\mathsf{I}}}\ln\frac{2\mathcal{R}}{a+b}+ {\boldsymbol{\mathsf{M}}}\right) \boldsymbol{\cdot}\boldsymbol{F}^{(1)} + \sum_{n=2}^N {\boldsymbol{\mathsf{J}}} (\boldsymbol{r}^{(1)}-\boldsymbol{r}^{(n)})\boldsymbol{\cdot} \boldsymbol{F}^{(n)}.\end{equation}

For the case of two disks with separation $2\boldsymbol {s}=\boldsymbol {r}^{(1)}-\boldsymbol {r}^{(2)}$, we obtain (cf. (8.9))

(8.15)\begin{equation} \boldsymbol{V}= \frac{1}{4{\rm \pi}\mu} \left({\boldsymbol{\mathsf{I}}}\ln\frac{4|\boldsymbol{s}|}{a+b} +{\boldsymbol{\mathsf{M}}} - \frac{\boldsymbol{s}\boldsymbol{s}}{|\boldsymbol{s}|^2}\right) \boldsymbol{\cdot}\boldsymbol{F}.\end{equation}

9. Concluding remarks

The present calculation has been motivated by the apparent conflict between the Stokes paradox and the observation that the mutual interaction between two particles in a viscous film is well-posed in 2-D Stokes flow, however large the separation between them. It is based upon the observation that standard reflection methods are not directly applicable to the problem, which requires instead the systematic use of matched asymptotic expansions. The asymptotic scheme developed herein can handle various interfacial conditions at the particle boundaries. Since the Stokes equations involve only the instantaneous configuration, the present scheme is applicable when the repulsive force depends upon the mutual separation.

In generalising the analysis of circular particles to non-circular shapes, we observe that the leading-order outer flow is unaffected by the shape. We illustrate the effect of shape by considering elliptic disks. The particle velocity is a superposition of an interaction term, which depends only upon the pair separation, and a renormalised term that depends only upon the ellipse eccentricity. By allowing for arbitrary angle between the ellipse axes and the line of centres, central forces may result in particle motion perpendicular to that line; see indeed (7.26b). This net drift in the absence of any net external force resembles self-propulsion. The indifference of the leading-order outer flow to the particle geometry and surface properties has been further exploited in extending the asymptotic analysis to a collection of non-identical disks.

It is important to emphasise the difference between the present analysis and earlier investigations of particle interaction in membranes (Bussell, Koch & Hammer 1992; Bussell, Hammer & Koch 1994; Dodd et al. 1995; Singh et al. 2019). These investigations considered particle motion that is driven by an external force, similarly to Saffman's analysis of a single particle. Since the Stokes paradox persists at these problems, they also require – following Saffman (1976) – the incorporation of small substrate viscosity.

In the case of disks of circular cross-section, it may be possible to obtain an exact solution of the general problem of arbitrary $\lambda$ using appropriate Bessel–Fourier expansions in bipolar coordinates (Wakiya 1975). Such an exact solution, however, does not necessarily provide qualitative insight regarding the asymptotic behaviour at the singular limit of remote disks (or, for that matter, in the other extreme of near contact). This requisite insight is provided by the present asymptotic scheme. For non-circular disks, the need for an asymptotic approach is even more paramount. Indeed, given the unavailability of orthogonal coordinate systems that can handle two non-circular disks, the exact analysis of the interaction between such disks requires numerical simulations (Power 1993). These, of course, are ill-suited for handling the scale disparity involved in the well-separated configuration. It is exactly at that limit where matched asymptotic expansions become indispensable.