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 Reference Hinch1991). 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 Reference Fraenkel1969), 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 Reference Van Dyke1964). 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 Reference Pozrikidis1992)

(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 Reference Saffman and Delbrück1975). 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 Reference Kropinski, Ward and Keller1995.) 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 Reference Moon and Spencer1988):

(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 Reference Shintani, Umemura and Takano1983)

(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. Reference Kropinski, Ward and Keller1995). 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 Reference Bussell, Koch and Hammer1992; Bussell, Hammer & Koch Reference Bussell, Hammer and Koch1994; Dodd et al. Reference Dodd, Hammer, Sangani and Koch1995; Singh et al. Reference Singh, Sangani, Daniel and Koch2019). 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 (Reference Saffman1976) – 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 Reference Wakiya1975). 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 Reference Power1993). 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.