Exact solutions for hydrodynamic interactions of two squirming spheres

3.1 Bispherical coordinates

Problems involving two spheres, such as we consider here, are naturally treated by employing a bispherical coordinate system. If $z+\text{i}\unicode[STIX]{x1D70C}$ is a complex coordinate on a Cartesian grid, the bipolar coordinate grid $\unicode[STIX]{x1D709}+\text{i}\unicode[STIX]{x1D702}$ is defined by

(3.1a,b ) $$\begin{eqnarray}\unicode[STIX]{x1D709}+\text{i}\unicode[STIX]{x1D702}=\ln \left[{\displaystyle \frac{z+\text{i}\unicode[STIX]{x1D70C}-R}{z+\text{i}\unicode[STIX]{x1D70C}+R}}\right],\quad z+\text{i}\unicode[STIX]{x1D70C}=-R{\displaystyle \frac{(\sinh \unicode[STIX]{x1D709}-\text{i}\sin \unicode[STIX]{x1D702})}{(\cosh \unicode[STIX]{x1D709}-\cos \unicode[STIX]{x1D702})}},\end{eqnarray}$$

where $R$ is a positive real number. This can be thought of as a stereographic projection of the lines of latitude and longitude on a sphere about a point on the equator, as demonstrated in figure 1(a), with the poles mapping to two symmetric points, $z=\pm R$ . Finally, a rotation about the $z$ axis gives an azimuthal coordinate $\unicode[STIX]{x1D719}$ , which coincides for bispherical and ordinary cylindrical coordinates. Surfaces of constant $\unicode[STIX]{x1D709}$ are non-intersecting spheres centred on $-R\,\text{coth}\,\unicode[STIX]{x1D709}$ , with radius $r=R\,|\text{cosech}\,\unicode[STIX]{x1D709}|$ . We consider the fluid to be the region $\unicode[STIX]{x1D709}_{2}<\unicode[STIX]{x1D709}<\unicode[STIX]{x1D709}_{1}$ , where $\unicode[STIX]{x1D709}_{1}$ is taken to be positive. The choice of $\unicode[STIX]{x1D709}_{2}$ then defines the geometry: if it is positive the fluid is the region between two nested spherical boundaries and if it is negative the fluid is external to the two spheres, while the intermediate case $\unicode[STIX]{x1D709}_{2}=0$ represents the half-space, as shown in figure 1(c). In what follows we will use both cylindrical coordinates $(z,\unicode[STIX]{x1D70C},\unicode[STIX]{x1D719})$ and bispherical coordinates $(\unicode[STIX]{x1D709},\unicode[STIX]{x1D702},\unicode[STIX]{x1D719})$ , and denote by $W$ the conformal factor $R/(\cosh \unicode[STIX]{x1D709}-\cos \unicode[STIX]{x1D702})$ that appears frequently.

Figure 1. (a) Stereographic projection of gridlines on a globe about a pole gives a polar grid (below), while a projection about an equatorial point gives a bipolar grid. (b) Conventions of the bispherical coordinate system $(\unicode[STIX]{x1D709},\unicode[STIX]{x1D702},\unicode[STIX]{x1D719})$ used in this work, related to the cylindrical basis $(z,\unicode[STIX]{x1D70C},\unicode[STIX]{x1D719})$ . The $\unicode[STIX]{x1D719}$ coordinate coincides for the two coordinate systems.

3.2 Coaxial rotation

The general solution for axisymmetric azimuthal flows was given by Jeffery (Reference Jeffery1915), together with a number of specific examples, including the coaxial rotation of two spheres, and was subsequently expanded upon by Kanwal (Reference Kanwal1961) to consider objects which do not intersect the axis of symmetry. For an axisymmetric, purely azimuthal Stokes flow, the fluid velocity takes the form $\boldsymbol{u}=u_{\unicode[STIX]{x1D719}}(z,\unicode[STIX]{x1D70C})\boldsymbol{e}_{\unicode[STIX]{x1D719}}$ . The pressure is constant everywhere and may be taken to be equal to zero without loss of generality, so the flow satisfies the scalar equation $(\unicode[STIX]{x1D6FB}^{2}-\unicode[STIX]{x1D70C}^{-2})u_{\unicode[STIX]{x1D719}}=0$ , of which the general solution in bispherical coordinates is

(3.2) $$\begin{eqnarray}u_{\unicode[STIX]{x1D719}}=W^{-1/2}\mathop{\sum }_{l=0}^{\infty }\left[c_{l}\,\text{e}^{(l+1/2)(\unicode[STIX]{x1D709}-\unicode[STIX]{x1D709}_{2})}+d_{l}\,\text{e}^{-(l+1/2)(\unicode[STIX]{x1D709}-\unicode[STIX]{x1D709}_{1})}\right]\,P_{l}^{1}(\cos \unicode[STIX]{x1D702}).\end{eqnarray}$$

The constants $c_{l}$ and $d_{l}$ are determined from the boundary conditions of solid-body rotation of the sphere $\unicode[STIX]{x1D709}=\unicode[STIX]{x1D709}_{1,2}$ with angular velocity $\unicode[STIX]{x1D714}_{1,2}$ about their common diameter

(3.3) $$\begin{eqnarray}\displaystyle & \displaystyle c_{l}=(2R)^{3/2}\mathop{\sum }_{n=0}^{\infty }\left(\unicode[STIX]{x1D714}_{2}\text{e}^{-(l+1/2)\left(2(n+1)(\unicode[STIX]{x1D709}_{1}-\unicode[STIX]{x1D709}_{2})+|\unicode[STIX]{x1D709}_{2}|\right)}-\unicode[STIX]{x1D714}_{1}\text{e}^{-(l+1/2)\left((2n+1)(\unicode[STIX]{x1D709}_{1}-\unicode[STIX]{x1D709}_{2})+|\unicode[STIX]{x1D709}_{1}|\right)}\right), & \displaystyle\end{eqnarray}$$

(3.4) $$\begin{eqnarray}\displaystyle & \displaystyle d_{l}=(2R)^{3/2}\mathop{\sum }_{n=0}^{\infty }\left(\unicode[STIX]{x1D714}_{1}\text{e}^{-(l+1/2)\left(2(n+1)(\unicode[STIX]{x1D709}_{1}-\unicode[STIX]{x1D709}_{2})+|\unicode[STIX]{x1D709}_{1}|\right)}-\unicode[STIX]{x1D714}_{2}\text{e}^{-(l+1/2)\left((2n+1)(\unicode[STIX]{x1D709}_{1}-\unicode[STIX]{x1D709}_{2})+|\unicode[STIX]{x1D709}_{2}|\right)}\right). & \displaystyle\end{eqnarray}$$

The force per area on the spheres’ surfaces is purely azimuthal, and integrating it gives the torques

(3.5) $$\begin{eqnarray}\displaystyle T_{1,2} & = & \displaystyle 8\unicode[STIX]{x03C0}\unicode[STIX]{x1D707}r_{1,2}^{3}\mathop{\sum }_{n=0}^{\infty }\left\{\unicode[STIX]{x1D714}_{1,2}{\displaystyle \frac{\sinh ^{3}\unicode[STIX]{x1D709}_{1,2}}{\sinh ^{3}\left(n(\unicode[STIX]{x1D709}_{1,2}-\unicode[STIX]{x1D709}_{2,1})+\unicode[STIX]{x1D709}_{1,2}\right)}}\right.\nonumber\\ \displaystyle & & \displaystyle \left.-\,\unicode[STIX]{x1D714}_{2,1}{\displaystyle \frac{\sinh ^{3}\unicode[STIX]{x1D709}_{1,2}}{\sinh ^{3}(n+1)(\unicode[STIX]{x1D709}_{1,2}-\unicode[STIX]{x1D709}_{2,1})}}\right\},\quad \unicode[STIX]{x1D709}_{2}<0<\unicode[STIX]{x1D709}_{1},\end{eqnarray}$$

(3.6) $$\begin{eqnarray}\displaystyle T_{1,2} & = & \displaystyle 8\unicode[STIX]{x03C0}\unicode[STIX]{x1D707}r_{1}^{3}(\unicode[STIX]{x1D714}_{1,2}-\unicode[STIX]{x1D714}_{2,1})\mathop{\sum }_{n=0}^{\infty }{\displaystyle \frac{\sinh ^{3}\unicode[STIX]{x1D709}_{1}}{\sinh ^{3}(n(\unicode[STIX]{x1D709}_{1}-\unicode[STIX]{x1D709}_{2})+\unicode[STIX]{x1D709}_{1})}},\quad 0\leqslant \unicode[STIX]{x1D709}_{2}<\unicode[STIX]{x1D709}_{1}.\end{eqnarray}$$

These expressions are uniformly convergent since they are bounded independently of $\unicode[STIX]{x1D709}_{1,2}$ : taking the first term in (3.5) as an example, using that $\sinh n\unicode[STIX]{x1D709}>n\sinh \unicode[STIX]{x1D709}$ for $n>0$ and $\unicode[STIX]{x1D709}>0$ we have that

(3.7) $$\begin{eqnarray}{\displaystyle \frac{\sinh ^{3}\unicode[STIX]{x1D709}_{1}}{\sinh ^{3}\left(n(\unicode[STIX]{x1D709}_{1}-\unicode[STIX]{x1D709}_{2})+\unicode[STIX]{x1D709}_{1}\right)}}<{\displaystyle \frac{\sinh ^{3}\unicode[STIX]{x1D709}_{1}}{\sinh ^{3}(n+1)\unicode[STIX]{x1D709}_{1}}}\leqslant {\displaystyle \frac{1}{(n+1)^{3}}},\end{eqnarray}$$

for $\unicode[STIX]{x1D709}_{2}<0,\,\unicode[STIX]{x1D709}_{1}>0$ . The uniform convergence of the remaining terms is demonstrated similarly, although note that (3.6) diverges as $\unicode[STIX]{x1D709}_{2}\rightarrow \unicode[STIX]{x1D709}_{1}$ .

A limit that is of particular interest for us in considering the near-field hydrodynamics of swimmers is that of vanishing separation where the spheres touch. In this limit $\unicode[STIX]{x1D709}_{1}$ and $\unicode[STIX]{x1D709}_{2}$ both tend to zero in a way that preserves the ratio $r_{1}/r_{2}$ ,

(3.8) $$\begin{eqnarray}\unicode[STIX]{x1D709}_{2}\sim \text{sgn}(\unicode[STIX]{x1D709}_{2}){\displaystyle \frac{r_{1}}{r_{2}}}\unicode[STIX]{x1D709}_{1}.\end{eqnarray}$$

Then, since (3.5) and (3.6) are uniformly convergent we may interchange the limit taking with the summation, finding that the torques converge to the finite values

(3.9) $$\begin{eqnarray}\displaystyle & \displaystyle T_{1,2}\rightarrow {\displaystyle \frac{8\unicode[STIX]{x03C0}\unicode[STIX]{x1D707}r_{1}^{3}r_{2}^{3}}{(r_{1}+r_{2})^{3}}}\left\{\unicode[STIX]{x1D714}_{1,2}\,\unicode[STIX]{x1D701}\left(3,\left(1+\frac{r_{1,2}}{r_{2,1}}\right)^{-1}\right)-\unicode[STIX]{x1D714}_{2,1}\,\unicode[STIX]{x1D701}(3)\right\},\quad \unicode[STIX]{x1D709}_{2}<0<\unicode[STIX]{x1D709}_{1}, & \displaystyle\end{eqnarray}$$

(3.10) $$\begin{eqnarray}\displaystyle & \displaystyle T_{1,2}\rightarrow {\displaystyle \frac{8\unicode[STIX]{x03C0}\unicode[STIX]{x1D707}r_{1}^{3}(\unicode[STIX]{x1D714}_{1,2}-\unicode[STIX]{x1D714}_{2,1})}{\left(1-{\displaystyle \frac{r_{1}}{r_{2}}}\right)^{3}}}\unicode[STIX]{x1D701}\left(3,\left(1-\frac{r_{1}}{r_{2}}\right)^{-1}\right),\quad 0\leqslant \unicode[STIX]{x1D709}_{2}<\unicode[STIX]{x1D709}_{1}, & \displaystyle\end{eqnarray}$$

where $\unicode[STIX]{x1D701}(s,q)\equiv \sum _{n=0}^{\infty }(n+q)^{-s}$ is the Hurwitz zeta function and $\unicode[STIX]{x1D701}(s,1)\equiv \unicode[STIX]{x1D701}(s)$ the Riemann zeta function. When one sphere encloses the other ( $\unicode[STIX]{x1D709}_{2}>0$ ) the torques on the two spheres are equal and opposite, meaning that only relative motion can be deduced from the reciprocal theorem, the left-hand side of (2.1) reducing to $T_{1}(\tilde{\unicode[STIX]{x1D6FA}}_{1}-\tilde{\unicode[STIX]{x1D6FA}}_{2})$ . It is natural to take the concave boundary to set the frame of reference. In the intermediate limit of a plane ( $\unicode[STIX]{x1D709}_{2}=0$ ) the two solutions (3.9) and (3.10) coincide and as with the enclosed system the torque on the wall is equal and opposite to the torque on the finite-sized sphere.

3.3 Translation along the common axis

The translational motion of two spheres along their common diameter was first studied by Stimson & Jeffery (Reference Stimson and Jeffery1926) and subsequently adapted for the special case of a sphere sedimenting towards a planar surface (Brenner Reference Brenner1961; Cox & Brenner Reference Cox and Brenner1967). The approach adopted involves the introduction of a streamfunction to solve the continuity equation by construction and then the Stokes equations reduce to a fourth-order operator acting on the streamfunction. Later studies (Dean & O’Neill Reference Dean and O’Neill1963; O’Neill Reference O’Neill1964; O’Neill & Stewartson Reference O’Neill and Stewartson1967; O’Neill & Majumdar Reference O’Neill and Majumdar1970a ,Reference O’Neill and Majumdar b ; Majumdar & O’Neill Reference Majumdar and O’Neill1977), motivated in part by a desire to extend to translations perpendicular to the common axis, follow the opposite approach: the Stokes equation is first solved by constructing a harmonic vector from an appropriate combination of the flow and the pressure ( $\boldsymbol{u}-(1/2\unicode[STIX]{x1D707})p\boldsymbol{x}$ ), the coefficients of which are then to be determined from boundary conditions and the imposition of incompressibility, resulting in a set of second-order difference equations that unfortunately proves analytically intractable. We quote the solution here in a form that is a minor variation of that given by Stimson & Jeffery (Reference Stimson and Jeffery1926). We also generalise to arbitrary translational speeds $V_{1}$ and $V_{2}$ and to the full range of geometries allowed by the bispherical coordinate system.

The domain has non-trivial cohomology in degree 2, associated with expansions or contractions of each of the two spherical boundary surfaces $\unicode[STIX]{x1D709}=\unicode[STIX]{x1D709}_{1,2}$ , such as might occur for two small gas bubbles. Neglecting these motions, the fluid velocity can be written as the curl of a vector $\unicode[STIX]{x1D74D}$ that satisfies the biharmonic equation, $\unicode[STIX]{x1D6FB}^{4}\unicode[STIX]{x1D74D}=\pmb{{\it 0}}$ . For an axisymmetric flow, $\unicode[STIX]{x1D74D}$ may be chosen to be purely azimuthal, $\unicode[STIX]{x1D74D}=\unicode[STIX]{x1D713}\boldsymbol{e}_{\unicode[STIX]{x1D719}}$ and then the general solution may be written in the form

(3.11) $$\begin{eqnarray}W^{-1/2}\unicode[STIX]{x1D713}=\mathop{\sum }_{l=1}^{\infty }\left({\mathcal{A}}_{l}\text{e}^{(l+3/2)\unicode[STIX]{x1D709}}+{\mathcal{D}}_{l}\text{e}^{-(l+3/2)\unicode[STIX]{x1D709}}+{\mathcal{B}}_{l}\text{e}^{(l-1/2)\unicode[STIX]{x1D709}}+{\mathcal{C}}_{l}\text{e}^{-(l-1/2)\unicode[STIX]{x1D709}}\right)P_{l}^{1}(\cos \unicode[STIX]{x1D702}).\end{eqnarray}$$

The boundary conditions that determine the real constants ${\mathcal{A}}_{l},{\mathcal{B}}_{l},{\mathcal{C}}_{l},{\mathcal{D}}_{l}$ are that the axial flow should equal the constant translation speeds $V_{1,2}\,\boldsymbol{e}_{z}$ of the two spheres $\unicode[STIX]{x1D709}=\unicode[STIX]{x1D709}_{1,2}$ and that the radial flow $u_{\unicode[STIX]{x1D70C}}$ should vanish on their surfaces. These conditions may be combined into the statement $\unicode[STIX]{x1D735}(\unicode[STIX]{x1D70C}\unicode[STIX]{x1D713}-(1/2)\unicode[STIX]{x1D70C}^{2}V_{j})_{\unicode[STIX]{x1D709}=\unicode[STIX]{x1D709}_{j}}=\pmb{{\it 0}}$ for $j=1,2$  (Stimson & Jeffery Reference Stimson and Jeffery1926), giving four equations to determine the four unknown coefficients, that reduce to the $2\times 2$ block-diagonal form

(3.12) $$\begin{eqnarray}\displaystyle & & \displaystyle \hspace{-10.00002pt}\left[\begin{array}{@{}cc@{}}(\text{e}^{(l+3/2)\unicode[STIX]{x1D709}_{1}}+\text{e}^{(l+3/2)\unicode[STIX]{x1D709}_{2}}) & (\text{e}^{(l-1/2)\unicode[STIX]{x1D709}_{1}}+\text{e}^{(l-1/2)\unicode[STIX]{x1D709}_{2}})\\ (2l+3)(\text{e}^{(l+3/2)\unicode[STIX]{x1D709}_{1}}-\text{e}^{(l+3/2)\unicode[STIX]{x1D709}_{2}}) & (2l-1)(\text{e}^{(l-1/2)\unicode[STIX]{x1D709}_{1}}-\text{e}^{(l-1/2)\unicode[STIX]{x1D709}_{2}})\end{array}\right]\!\left[\begin{array}{@{}c@{}}{\mathcal{A}}_{l}+{\mathcal{D}}_{l}\text{e}^{-(l+3/2)(\unicode[STIX]{x1D709}_{1}+\unicode[STIX]{x1D709}_{2})}\\ {\mathcal{B}}_{l}+{\mathcal{C}}_{l}\text{e}^{-(l-1/2)(\unicode[STIX]{x1D709}_{1}+\unicode[STIX]{x1D709}_{2})}\end{array}\right]\nonumber\\ \displaystyle & & \displaystyle \hspace{-10.00002pt}\quad ={\displaystyle \frac{R}{\sqrt{2}}}\left[\begin{array}{@{}c@{}}\left({\displaystyle \frac{\text{e}^{-(l+3/2)|\unicode[STIX]{x1D709}_{1}|}}{(2l+3)}}-\text{e}^{-(l-1/2)|\unicode[STIX]{x1D709}_{1}|}/(2l-1)\right)V_{1}+\left({\displaystyle \frac{\text{e}^{-(l+3/2)|\unicode[STIX]{x1D709}_{2}|}}{(2l+3)}}-{\displaystyle \frac{\text{e}^{-(l-1/2)|\unicode[STIX]{x1D709}_{2}|}}{(2l-1)}}\right)V_{2}\\ 2V_{1}\text{e}^{-(l+1/2)|\unicode[STIX]{x1D709}_{1}|}\sinh \unicode[STIX]{x1D709}_{1}-2V_{2}\text{e}^{-(l+1/2)|\unicode[STIX]{x1D709}_{2}|}\sinh \unicode[STIX]{x1D709}_{2}\end{array}\right]\nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$

and

(3.13) $$\begin{eqnarray}\displaystyle & & \displaystyle \!\!\left[\begin{array}{@{}cc@{}}(\text{e}^{(l+3/2)\unicode[STIX]{x1D709}_{1}}-\text{e}^{(l+3/2)\unicode[STIX]{x1D709}_{2}}) & (\text{e}^{(l-1/2)\unicode[STIX]{x1D709}_{1}}-\text{e}^{(l-1/2)\unicode[STIX]{x1D709}_{2}})\\ (2l+3)(\text{e}^{(l+3/2)\unicode[STIX]{x1D709}_{1}}+\text{e}^{(l+3/2)\unicode[STIX]{x1D709}_{2}}) & (2l-1)(\text{e}^{(l-1/2)\unicode[STIX]{x1D709}_{1}}+\text{e}^{(l-1/2)\unicode[STIX]{x1D709}_{2}})\end{array}\right]\left[\begin{array}{@{}c@{}}{\mathcal{A}}_{l}-{\mathcal{D}}_{l}\text{e}^{-(l+3/2)(\unicode[STIX]{x1D709}_{1}+\unicode[STIX]{x1D709}_{2})}\\ {\mathcal{B}}_{l}-{\mathcal{C}}_{l}\text{e}^{-(l-1/2)(\unicode[STIX]{x1D709}_{1}+\unicode[STIX]{x1D709}_{2})}\end{array}\right]\nonumber\\ \displaystyle & & \displaystyle \!\!\quad ={\displaystyle \frac{R}{\sqrt{2}}}\left[\begin{array}{@{}c@{}}\left({\displaystyle \frac{\text{e}^{-(l+3/2)|\unicode[STIX]{x1D709}_{1}|}}{(2l+3)}}-{\displaystyle \frac{\text{e}^{-(l-1/2)|\unicode[STIX]{x1D709}_{1}|}}{(2l-1)}}\right)V_{1}-\left({\displaystyle \frac{\text{e}^{-(l+3/2)|\unicode[STIX]{x1D709}_{2}|}}{(2l+3)}}-{\displaystyle \frac{\text{e}^{-(l-1/2)|\unicode[STIX]{x1D709}_{2}|}}{(2l-1)}}\right)V_{2}\\ 2V_{1}\text{e}^{-(l+1/2)|\unicode[STIX]{x1D709}_{1}|}\sinh \unicode[STIX]{x1D709}_{1}+2V_{2}\text{e}^{-(l+1/2)|\unicode[STIX]{x1D709}_{2}|}\sinh \unicode[STIX]{x1D709}_{2}\end{array}\right]\end{eqnarray}$$

Inversion is straightforward, and the explicit forms of the coefficients ${\mathcal{A}}_{l},{\mathcal{B}}_{l},{\mathcal{C}}_{l},{\mathcal{D}}_{l}$ are shown in appendix A.

In bispherical coordinates, the flow is given by

(3.14) $$\begin{eqnarray}\displaystyle \boldsymbol{u}=\unicode[STIX]{x1D735}\times \unicode[STIX]{x1D74D} & = & \displaystyle \left[W^{-1/2}{\displaystyle \frac{1}{\sin \unicode[STIX]{x1D702}}}\unicode[STIX]{x2202}_{\unicode[STIX]{x1D702}}\left(\sin \unicode[STIX]{x1D702}\,W^{-1/2}\unicode[STIX]{x1D713}\right)-{\displaystyle \frac{3\sin \unicode[STIX]{x1D702}}{2R}}W^{1/2}W^{-1/2}\unicode[STIX]{x1D713}\right]\boldsymbol{e}_{\unicode[STIX]{x1D709}}\nonumber\\ \displaystyle & & \displaystyle -\,\left[W^{-1/2}\unicode[STIX]{x2202}_{\unicode[STIX]{x1D709}}\left(W^{-1/2}\unicode[STIX]{x1D713}\right)-{\displaystyle \frac{3\sinh \unicode[STIX]{x1D709}}{2R}}W^{1/2}W^{-1/2}\unicode[STIX]{x1D713}\right]\boldsymbol{e}_{\unicode[STIX]{x1D702}},\end{eqnarray}$$

which allows the pressure at a point $\boldsymbol{x}$ to be calculated as

(3.15) $$\begin{eqnarray}p(\boldsymbol{x})=p_{\infty }+\int _{\infty }^{\boldsymbol{x}}\text{d}\boldsymbol{l}\boldsymbol{\cdot }\unicode[STIX]{x1D735}p=p_{\infty }+\unicode[STIX]{x1D707}\int _{\infty }^{\boldsymbol{x}}\text{d}\boldsymbol{l}\boldsymbol{\cdot }\unicode[STIX]{x1D6FB}^{2}\boldsymbol{u},\end{eqnarray}$$

where $p_{\infty }$ is its asymptotic value. There is no dependence on the path of integration as the domain is simply connected.

The solution we have presented is equivalent to those given previously (Stimson & Jeffery Reference Stimson and Jeffery1926; Brenner Reference Brenner1961), although it is not identical because the manner in which we have solved the continuity equation differs slightly. A consequence is the expansion of the vector potential in associated Legendre polynomials $P_{l}^{1}$ , rather than in Gegenbauer polynomials $C_{n+1}^{-1/2}$ . Furthermore the final scheme we arrive at for determining the coefficients ${\mathcal{A}}_{l},{\mathcal{B}}_{l},{\mathcal{C}}_{l},{\mathcal{D}}_{l}$ presents the $4\times 4$ problem in block diagonal form, which we have not seen in the previous literature.

Finally, using (3.14) and (3.15) a stress tensor may be constructed and integrated by parts over the spheres to give the hydrodynamic drag force

(3.16) $$\begin{eqnarray}F_{1,2}=\pm 4\unicode[STIX]{x03C0}\unicode[STIX]{x1D707}\sqrt{2R}\mathop{\sum }_{l=1}^{\infty }l(l+1)\left\{\begin{array}{@{}cl@{}}-({\mathcal{A}}_{l}+{\mathcal{B}}_{l}), & \quad \unicode[STIX]{x1D709}_{1,2}\geqslant 0,\\ ({\mathcal{C}}_{l}+{\mathcal{D}}_{l}), & \quad \unicode[STIX]{x1D709}_{1,2}<0.\end{array}\right.\end{eqnarray}$$

When the fluid has a finite volume or in the limiting case of the half-space ( $\unicode[STIX]{x1D709}_{2}\geqslant 0$ ), the net force on the fluid is zero since the contributions from the two boundaries are equal and opposite. This is also seen for a point force in the half-space (Blake & Chwang Reference Blake and Chwang1974), which is obtained from our solution in the limit $\unicode[STIX]{x1D709}_{1}\rightarrow \infty$ , $\unicode[STIX]{x1D709}_{2}\rightarrow 0$ , with $R$ held constant.

Since the coefficients ${\mathcal{A}}_{l},{\mathcal{B}}_{l},{\mathcal{C}}_{l},{\mathcal{D}}_{l}$ have exponential decay in $\unicode[STIX]{x1D709}_{1}$ and $\unicode[STIX]{x1D709}_{2}$ the sums converge rapidly for large separation; however, as the separation between the spheres vanishes the forces diverge. This limit has been treated in detail by Cox & Brenner (Reference Cox and Brenner1967) and we adapt their method here.

4.1 Squirming

The results of the previous section for Stokes drag of two spheres allow a variety of axisymmetric swimmer motions to be determined via the reciprocal theorem. Despite the absence of expressions for non-axisymmetric motions this is enough to, for instance, give an exact description of the circular motion of microorganisms such as E. coli close to planar boundaries (Berg Reference Berg2000; Lauga et al. Reference Lauga, Diluzio, Whitesides and Stone2006), and can also shed light on the hydrodynamics of a daughter colony of Volvox inside its parent (Drescher et al. Reference Drescher, Leptos, Tuval, Ishikawa, Pedley and Goldstein2009), or the contact interaction of swimmers with passive particles (Wu & Libchaber Reference Wu and Libchaber2000). Many of these have been studied asymptotically using leading-order point-singularity descriptions (Berke et al. Reference Berke, Turner, Berg and Lauga2008; Spagnolie & Lauga Reference Spagnolie and Lauga2012), but using the exact solutions for Stokes drag we are able to describe the behaviour for arbitrarily small separation and arbitrary squirming motions. Since the reciprocal theorem and the Stokes equations are linear, it suffices to calculate the interaction of a swimmer and passive sphere (Ishikawa et al. Reference Ishikawa, Simmonds and Pedley2006), whose motion is given by

(4.1) $$\begin{eqnarray}\left.\tilde{U} _{1}F_{1}+\tilde{U} _{2}F_{2}+\tilde{\unicode[STIX]{x1D6FA}}_{1}T_{1}+\tilde{\unicode[STIX]{x1D6FA}}_{2}T_{2}=-\int _{0}^{2\unicode[STIX]{x03C0}}\text{d}\unicode[STIX]{x1D719}\int _{0}^{\unicode[STIX]{x03C0}}W^{2}\sin \unicode[STIX]{x1D702}\,\text{d}\unicode[STIX]{x1D702}\,(\boldsymbol{u}^{s}\boldsymbol{\cdot }\boldsymbol{\unicode[STIX]{x1D70E}}\boldsymbol{\cdot }\hat{\boldsymbol{n}})\right|_{\unicode[STIX]{x1D709}=\unicode[STIX]{x1D709}_{1}},\end{eqnarray}$$

where the stress tensor is the corresponding to the Stokes drag problems given in § 3.

We consider a swimmer of radius $r_{1}$ centred a perpendicular distance $d$ away from the surface of a passive sphere of radius $r_{2}$ , which may be convex, flat or concave and we respectively term the ‘tracer’, ‘wall’ or ‘shell’. The swimmer approaches the surface at an angle $\unicode[STIX]{x1D6FC}$ and its squirming motion is described in terms of a local orthonormal basis $\{\boldsymbol{s}_{r},\boldsymbol{s}_{\unicode[STIX]{x1D703}},\boldsymbol{s}_{\unicode[STIX]{x1D719}}\}$ and polar coordinate system $(\unicode[STIX]{x1D703}_{s},\unicode[STIX]{x1D719}_{s})$ relative to this direction, as shown in figure 2. Its slip velocity may be decomposed into squirming modes (Lighthill Reference Lighthill1952) as

(4.2) $$\begin{eqnarray}\boldsymbol{u}^{s}=\mathop{\sum }_{n\geqslant 1}\left[A_{n}P_{n}(\cos \unicode[STIX]{x1D703}_{s})\,\boldsymbol{s}_{r}+B_{n}V_{n}(\cos \unicode[STIX]{x1D703}_{s})\,\boldsymbol{s}_{\unicode[STIX]{x1D703}}+r_{1}C_{n}V_{n}(\cos \unicode[STIX]{x1D703}_{s})\,\boldsymbol{s}_{\unicode[STIX]{x1D719}}\right],\end{eqnarray}$$

where $V_{n}(x)\equiv -2P_{n}^{1}(x)/n(n+1)$ and $A_{n}$ , $B_{n}$ and $C_{n}$ are real coefficients with the units of velocity. The free swimming speed, asymptotically far from the surface, is $U_{free}=(2B_{1}-A_{1})/3$  (Lighthill Reference Lighthill1952; Blake Reference Blake1971b ). The addition of the azimuthal modes $C_{n}$  (Pak & Lauga Reference Pak and Lauga2014) allows for axial rotation of the swimmer, as seen in several real microorganisms, with rotation speed $\unicode[STIX]{x1D6FA}_{free}=-C_{1}$ about the axisymmetry axis.

Figure 2. A spherical swimmer of radius $r_{1}$ , located a perpendicular distance $d$ away from the surface of a shell, wall or tracer of radius $r_{2}$ . The swimmer surface is parametrised by the coordinate $(\unicode[STIX]{x1D703}_{s},\unicode[STIX]{x1D719}_{s})$ and the swimmer approaches the passive sphere at an angle $\unicode[STIX]{x1D6FC}$ to the common diameter. The $\boldsymbol{e}_{z}$ axis points into the wall. Inset: the cylindrical and bispherical bases on the swimmer’s surface are related by a rotation of angle $\unicode[STIX]{x1D6FD}$ about $\boldsymbol{e}_{\unicode[STIX]{x1D719}}$ . The head–tail axis defining the swimmer spherical basis is denoted $\boldsymbol{s}$ .

To perform the integral in (4.1) it is convenient to express the swimmer’s slip velocity in bispherical coordinates. Defining the angle $\unicode[STIX]{x1D6FD}$ as the rotation angle about $\boldsymbol{e}_{\unicode[STIX]{x1D719}}$ between the cylindrical ( $\boldsymbol{e}_{z},\boldsymbol{e}_{\unicode[STIX]{x1D70C}}$ ) and bispherical ( $\boldsymbol{e}_{\unicode[STIX]{x1D709}},\boldsymbol{e}_{\unicode[STIX]{x1D702}}$ ) basis vectors,

(4.3) $$\begin{eqnarray}\cos \unicode[STIX]{x1D6FD}\equiv -\boldsymbol{e}_{z}\boldsymbol{\cdot }\boldsymbol{e}_{\unicode[STIX]{x1D709}}={\displaystyle \frac{1-\cosh \unicode[STIX]{x1D709}\cos \unicode[STIX]{x1D702}}{\cosh \unicode[STIX]{x1D709}-\cos \unicode[STIX]{x1D702}}}=\cosh \unicode[STIX]{x1D709}-{\displaystyle \frac{\sinh ^{2}\unicode[STIX]{x1D709}}{\cosh \unicode[STIX]{x1D709}-\cos \unicode[STIX]{x1D702}}},\end{eqnarray}$$

the swimmer polar angle is found using the spherical cosine law to be,

(4.4) $$\begin{eqnarray}\cos \unicode[STIX]{x1D703}_{s}=\cos \unicode[STIX]{x1D6FC}\cos \unicode[STIX]{x1D6FD}+\sin \unicode[STIX]{x1D6FC}\sin \unicode[STIX]{x1D6FD}\cos \unicode[STIX]{x1D719},\end{eqnarray}$$

as illustrated in figure 2, and the surface unit vectors are given by the transformation

(4.5) $$\begin{eqnarray}\left[\begin{array}{@{}c@{}}\boldsymbol{s}_{r}\\ \boldsymbol{s}_{\unicode[STIX]{x1D703}}\\ \boldsymbol{ s}_{\unicode[STIX]{x1D719}}\end{array}\right]=\left[\begin{array}{@{}ccc@{}}-1 & 0 & 0\\ 0 & {\displaystyle \frac{\unicode[STIX]{x2202}_{\unicode[STIX]{x1D6FD}}\cos \unicode[STIX]{x1D703}_{s}}{\sin \unicode[STIX]{x1D703}_{s}}} & {\displaystyle \frac{\sin \unicode[STIX]{x1D6FC}\sin \unicode[STIX]{x1D719}}{\sin \unicode[STIX]{x1D703}_{s}}}\\ 0 & {\displaystyle \frac{\sin \unicode[STIX]{x1D6FC}\sin \unicode[STIX]{x1D719}}{\sin \unicode[STIX]{x1D703}_{s}}} & -{\displaystyle \frac{\unicode[STIX]{x2202}_{\unicode[STIX]{x1D6FD}}\cos \unicode[STIX]{x1D703}_{s}}{\sin \unicode[STIX]{x1D703}_{s}}}\end{array}\right]\left[\begin{array}{@{}c@{}}\boldsymbol{e}_{\unicode[STIX]{x1D709}}\\ \boldsymbol{ e}_{\unicode[STIX]{x1D702}}\\ \boldsymbol{e}_{\unicode[STIX]{x1D719}}\end{array}\right].\end{eqnarray}$$

The slip velocity involves Legendre polynomials in $\cos \unicode[STIX]{x1D703}_{s}$ , which are expanded using the addition theorem for Legendre functions (Sneddon Reference Sneddon1956; Maleček & Nádeník Reference Maleček and Nádeník2001),

(4.6) $$\begin{eqnarray}\displaystyle & & \displaystyle P_{n}(\cos \unicode[STIX]{x1D6FC}\cos \unicode[STIX]{x1D6FD}+\sin \unicode[STIX]{x1D6FC}\sin \unicode[STIX]{x1D6FD}\cos \unicode[STIX]{x1D719})=P_{n}(\cos \unicode[STIX]{x1D6FC})P_{n}(\cos \unicode[STIX]{x1D6FD})\nonumber\\ \displaystyle & & \displaystyle \quad +\,2\mathop{\sum }_{m=1}^{n}{\displaystyle \frac{(n-m)!}{(n+m)!}}P_{n}^{m}(\cos \unicode[STIX]{x1D6FC})P_{n}^{m}(\cos \unicode[STIX]{x1D6FD})\cos m\unicode[STIX]{x1D719}.\end{eqnarray}$$

As the stress tensor in (4.1) is axisymmetric the integral over $\unicode[STIX]{x1D719}$ only affects the slip velocity components. Hence it is convenient to perform the $\unicode[STIX]{x1D719}$ -integral first and define a vector of the resulting azimuthally averaged slip velocity components,

(4.7) $$\begin{eqnarray}\displaystyle \hspace{-10.00002pt}\langle \boldsymbol{u}^{s}\rangle _{\unicode[STIX]{x1D719}} & \equiv & \displaystyle \langle \boldsymbol{u}^{s}\boldsymbol{\cdot }\boldsymbol{e}_{\unicode[STIX]{x1D709}}\rangle _{\unicode[STIX]{x1D719}}\boldsymbol{e}_{\unicode[STIX]{x1D709}}+\langle \boldsymbol{u}^{s}\boldsymbol{\cdot }\boldsymbol{e}_{\unicode[STIX]{x1D702}}\rangle _{\unicode[STIX]{x1D719}}\boldsymbol{e}_{\unicode[STIX]{x1D702}}+\langle \boldsymbol{u}^{s}\boldsymbol{\cdot }\boldsymbol{e}_{\unicode[STIX]{x1D719}}\rangle _{\unicode[STIX]{x1D719}}\boldsymbol{e}_{\unicode[STIX]{x1D719}}\nonumber\\ \displaystyle \hspace{-10.00002pt} & = & \displaystyle -\mathop{\sum }_{n\geqslant 1}P_{n}(\cos \unicode[STIX]{x1D6FC})\left[A_{n}P_{n}(\cos \unicode[STIX]{x1D6FD})\boldsymbol{e}_{\unicode[STIX]{x1D709}}+B_{n}V_{n}(\cos \unicode[STIX]{x1D6FD})\boldsymbol{e}_{\unicode[STIX]{x1D702}}-r_{1}C_{n}V_{n}(\cos \unicode[STIX]{x1D6FD})\boldsymbol{e}_{\unicode[STIX]{x1D719}}\right],\end{eqnarray}$$

which is contracted against the stress tensor and integrated over $\unicode[STIX]{x1D702}$ to give the motion. The radial and meridional modes $A_{n}$ and $B_{n}$ cannot drive axisymmetric rotation, since the normal stress corresponding to axisymmetric rotation is purely azimuthal; similarly, the azimuthal modes $C_{n}$ cannot drive axisymmetric translation. The dependence of the motion on the swimmer’s orientation $\unicode[STIX]{x1D6FC}$ is a purely geometric factor for each order of squirming mode, and at large separations where higher-order modes may be neglected the orientation dependence is simply $P_{2}(\cos \unicode[STIX]{x1D6FC})$ , as found using point-singularity models of swimmer interactions with walls (Spagnolie & Lauga Reference Spagnolie and Lauga2012; Davis & Crowdy Reference Davis and Crowdy2015; Papavassiliou & Alexander Reference Papavassiliou and Alexander2015).

The contributions from the tangential modes, $B_{n},C_{n}$ , are evaluated straightforwardly (albeit tediously) using orthogonality of Legendre polynomials. The radial modes, $A_{n}$ , pick up a contribution from the pressure, which may be rewritten in terms of the flow by integrating by parts using the identity

(4.8) $$\begin{eqnarray}W^{2}\sin (\unicode[STIX]{x1D702})P_{n}(\cos \unicode[STIX]{x1D6FD})\equiv -{\displaystyle \frac{R^{2}}{n(n+1)}}\unicode[STIX]{x2202}_{\unicode[STIX]{x1D702}}\left[{\displaystyle \frac{\sin \unicode[STIX]{x1D702}}{\sinh ^{2}\unicode[STIX]{x1D709}}}\unicode[STIX]{x2202}_{\unicode[STIX]{x1D702}}P_{n}(\cos \unicode[STIX]{x1D6FD})\right].\end{eqnarray}$$

4.2 Rotation

In this section we calculate explicitly the rotational motion of a squirmer close to a surface. Combined with self-propulsion parallel to the surface this rotation results in circling behaviour, which has been observed experimentally for flagellated bacteria such as E. coli and Vibrio alginolyticus in close proximity to a planar boundary (Berg Reference Berg2000; Lauga et al. Reference Lauga, Diluzio, Whitesides and Stone2006). The effect is highly local, with the gap between the bacterium and the wall typically much smaller than the size of the bacterium itself. While point-singularity methods predict such behaviour just as a result of the $C_{2}$ mode and indeed agree that it should be strongly localised close to the surface, with an inverse-fourth dependence on the separation (Lauga et al. Reference Lauga, Diluzio, Whitesides and Stone2006; Lopez & Lauga Reference Lopez and Lauga2014; Davis & Crowdy Reference Davis and Crowdy2015; Papavassiliou & Alexander Reference Papavassiliou and Alexander2015), higher-order modes can be expected to play an important role at such small gap widths.

The induced rotation is calculated by performing the integral (4.1) using the slip velocity (4.7) and the stress corresponding to axisymmetric rotation,

(4.9) $$\begin{eqnarray}\displaystyle & & \displaystyle \tilde{\unicode[STIX]{x1D6FA}}_{1}T_{1}+\tilde{\unicode[STIX]{x1D6FA}}_{2}T_{2}=-2\unicode[STIX]{x03C0}\unicode[STIX]{x1D707}\mathop{\sum }_{l\geqslant 1}\mathop{\sum }_{i=1}^{\infty }r_{1}C_{l}P_{l}(\cos \unicode[STIX]{x1D6FC})\int _{0}^{\unicode[STIX]{x03C0}}\sin \unicode[STIX]{x1D702}\,\text{d}\unicode[STIX]{x1D702}\,P_{i}^{1}(\cos \unicode[STIX]{x1D702})V_{l}(\cos \unicode[STIX]{x1D6FD})\nonumber\\ \displaystyle & & \displaystyle \quad \times \left(W^{1/2}\left(i+\frac{1}{2}\right)\left(c_{i}\,\text{e}^{(i+1/2)(\unicode[STIX]{x1D709}_{1}-\unicode[STIX]{x1D709}_{2})}-d_{i}\right)+W^{3/2}{\displaystyle \frac{3\sinh \unicode[STIX]{x1D709}}{2R}}\left(c_{i}\,\text{e}^{(i+1/2)(\unicode[STIX]{x1D709}_{1}-\unicode[STIX]{x1D709}_{2})}+d_{i}\right)\right)\bigg|_{\unicode[STIX]{x1D709}=\unicode[STIX]{x1D709}_{1}},\nonumber\\ \displaystyle & & \displaystyle\end{eqnarray}$$

where $c_{i},d_{i}$ are as given in (3.3) and (3.4). The factor of $V_{l}(\cos \unicode[STIX]{x1D6FD})$ may be written as a polynomial of order $l-1$ in $W$ ,

(4.10) $$\begin{eqnarray}\displaystyle V_{l}(\cos \unicode[STIX]{x1D6FD}) & = & \displaystyle {\displaystyle \frac{2\sinh \unicode[STIX]{x1D709}\sin \unicode[STIX]{x1D702}}{l(l+1)}}{\displaystyle \frac{W}{R}}P_{l}^{\prime }\left(\cosh \unicode[STIX]{x1D709}-\sinh ^{2}\unicode[STIX]{x1D709}{\displaystyle \frac{W}{R}}\right)\nonumber\\ \displaystyle & \equiv & \displaystyle {\displaystyle \frac{2\sin \unicode[STIX]{x1D702}}{l(l+1)}}{\displaystyle \frac{W}{r_{1}}}\mathop{\sum }_{n=0}^{l-1}w_{n}(\unicode[STIX]{x1D709})\left({\displaystyle \frac{W}{R}}\right)^{n},\end{eqnarray}$$

where the coefficients $w_{n}(\unicode[STIX]{x1D709})$ are determined using any of the various series representations of Legendre polynomials (Whittaker & Watson Reference Whittaker and Watson1996), so that (4.9) becomes

(4.11) $$\begin{eqnarray}\displaystyle & & \displaystyle \tilde{\unicode[STIX]{x1D6FA}}_{1}T_{1}+\tilde{\unicode[STIX]{x1D6FA}}_{2}T_{2}=-2\unicode[STIX]{x03C0}\unicode[STIX]{x1D707}\mathop{\sum }_{l\geqslant 1}\mathop{\sum }_{n=0}^{l-1}{\displaystyle \frac{2R^{3/2}}{l(l+1)}}w_{n}(\unicode[STIX]{x1D709})\mathop{\sum }_{i=1}^{\infty }C_{l}P_{l}(\cos \unicode[STIX]{x1D6FC})\nonumber\\ \displaystyle & & \displaystyle \quad \times \left(\left(i+\frac{1}{2}\right)\left(c_{i}\,\text{e}^{(i+1/2)(\unicode[STIX]{x1D709}_{1}-\unicode[STIX]{x1D709}_{2})}-d_{i}\right)\int _{0}^{\unicode[STIX]{x03C0}}\sin \unicode[STIX]{x1D702}\,\text{d}\unicode[STIX]{x1D702}\,P_{i}^{1}(\cos \unicode[STIX]{x1D702})\left({\displaystyle \frac{W}{R}}\right)^{n+3/2}\sin \unicode[STIX]{x1D702}\right.\nonumber\\ \displaystyle & & \displaystyle \quad \left.\left.+\,{\displaystyle \frac{3\sinh \unicode[STIX]{x1D709}}{2}}\left(c_{i}\,\text{e}^{(i+1/2)(\unicode[STIX]{x1D709}_{1}-\unicode[STIX]{x1D709}_{2})}+d_{i}\right)\!\int _{0}^{\unicode[STIX]{x03C0}}\!\sin \unicode[STIX]{x1D702}\,\text{d}\unicode[STIX]{x1D702}\,P_{i}^{1}(\cos \unicode[STIX]{x1D702})\!\left({\displaystyle \frac{W}{R}}\right)^{n+5/2}\sin \unicode[STIX]{x1D702}\right)\right|_{\unicode[STIX]{x1D709}=\unicode[STIX]{x1D709}_{1}}\!\!.\end{eqnarray}$$

$(W/R)^{1/2}$ is the generating function for Legendre polynomials (Whittaker & Watson Reference Whittaker and Watson1996). Successive differentiations of the generating function give the identity

(4.12) $$\begin{eqnarray}\left({\displaystyle \frac{W}{R}}\right)^{n+3/2}\sin \unicode[STIX]{x1D702}=\sqrt{2}{\displaystyle \frac{(-2)^{n+1}}{(2n+1)!!}}\mathop{\sum }_{m=1}^{\infty }P_{m}^{1}(\cos \unicode[STIX]{x1D702})\left[{\displaystyle \frac{1}{\sinh \unicode[STIX]{x1D709}}}\unicode[STIX]{x2202}_{\unicode[STIX]{x1D709}}\right]^{n}\text{e}^{-(m+1/2)|\unicode[STIX]{x1D709}|},\end{eqnarray}$$

which reduces (4.9) to a pair of integrals over orthogonal associated Legendre polynomials. Performing these integrals and resumming the results we find that near a concave shell or wall the motion is

(4.13) $$\begin{eqnarray}(\tilde{\unicode[STIX]{x1D6FA}}_{1}-\tilde{\unicode[STIX]{x1D6FA}}_{2})T_{1}=-8\unicode[STIX]{x03C0}\unicode[STIX]{x1D707}\unicode[STIX]{x1D714}_{1}r_{1}^{3}\mathop{\sum }_{l\geqslant 1}C_{l}P_{l}(\cos \unicode[STIX]{x1D6FC})\mathop{\sum }_{n=0}^{\infty }{\displaystyle \frac{\sinh ^{3}\unicode[STIX]{x1D709}_{1}\sinh ^{l-1}n(\unicode[STIX]{x1D709}_{1}-\unicode[STIX]{x1D709}_{2})}{\sinh ^{l+2}\left(n(\unicode[STIX]{x1D709}_{1}-\unicode[STIX]{x1D709}_{2})+\unicode[STIX]{x1D709}_{1}\right)}},\end{eqnarray}$$

while for interaction with a tracer it is

(4.14) $$\begin{eqnarray}\displaystyle & & \displaystyle \hspace{-30.00005pt}\tilde{\unicode[STIX]{x1D6FA}}_{1}T_{1}+\tilde{\unicode[STIX]{x1D6FA}}_{2}T_{2}=8\unicode[STIX]{x03C0}\unicode[STIX]{x1D707}r_{1}^{3}\mathop{\sum }_{l\geqslant 1}C_{l}P_{l}(\cos \unicode[STIX]{x1D6FC})\nonumber\\ \displaystyle & & \displaystyle \hspace{-30.00005pt}\quad \times \,\mathop{\sum }_{n=0}^{\infty }\left[\unicode[STIX]{x1D714}_{2}{\displaystyle \frac{\sinh ^{3}\unicode[STIX]{x1D709}_{1}\sinh ^{l-1}\left(n(\unicode[STIX]{x1D709}_{1}-\unicode[STIX]{x1D709}_{2})-\unicode[STIX]{x1D709}_{2}\right)}{\sinh ^{l+2}(n+1)(\unicode[STIX]{x1D709}_{1}-\unicode[STIX]{x1D709}_{2})}}-\unicode[STIX]{x1D714}_{1}{\displaystyle \frac{\sinh ^{3}\unicode[STIX]{x1D709}_{1}\sinh ^{l-1}n(\unicode[STIX]{x1D709}_{1}-\unicode[STIX]{x1D709}_{2})}{\sinh ^{l+2}\left(n(\unicode[STIX]{x1D709}_{1}-\unicode[STIX]{x1D709}_{2})+\unicode[STIX]{x1D709}_{1}\right)}}\right],\end{eqnarray}$$

where the torques are given by (3.5) and (3.6). To find $\tilde{\unicode[STIX]{x1D6FA}}_{1}$ , $\unicode[STIX]{x1D714}_{2}$ must be chosen so that $T_{2}=0$ . Conversely, choosing $\unicode[STIX]{x1D714}_{2}$ such that $T_{1}=0$ allows the tracer motion $\tilde{\unicode[STIX]{x1D6FA}}_{2}$ to be found. These expressions are exact, for any separation, any axisymmetric slip velocity and any of the geometries covered by bispherical coordinates. We will show in § 4.3 that this reproduces the rotational hydrodynamic interactions that have been determined previously using asymptotic methods, such as minimal reflections of point singularities. First, however, we examine the limit of small separation, where the swimmers approach contact and the hydrodynamic interactions are strongest.

The limit of vanishing separation is treated using dominated convergence as described in § 3.2. For a shell this gives

(4.15) $$\begin{eqnarray}\tilde{\unicode[STIX]{x1D6FA}}_{1}-\tilde{\unicode[STIX]{x1D6FA}}_{2}\rightarrow -\mathop{\sum }_{l\geqslant 1}C_{l}P_{l}(\cos \unicode[STIX]{x1D6FC})\mathop{\sum }_{k=0}^{l-1}(-1)^{k}\binom{l-1}{k}\left(1-{\displaystyle \frac{r_{1}}{r_{2}}}\right)^{-k}{\displaystyle \frac{\unicode[STIX]{x1D701}\left(3+k,\left(1-{\displaystyle \frac{r_{1}}{r_{2}}}\right)^{-1}\right)}{\unicode[STIX]{x1D701}\left(3,\left(1-{\displaystyle \frac{r_{1}}{r_{2}}}\right)^{-1}\right)}}\end{eqnarray}$$

and for a tracer

(4.16) $$\begin{eqnarray}\displaystyle & & \displaystyle \tilde{\unicode[STIX]{x1D6FA}}_{1}\rightarrow -\mathop{\sum }_{l\geqslant 1}C_{l}P_{l}(\cos \unicode[STIX]{x1D6FC})\mathop{\sum }_{k=0}^{l-1}(-1)^{k}\binom{l-1}{k}\left(1+{\displaystyle \frac{r_{1}}{r_{2}}}\right)^{-k}\nonumber\\ \displaystyle & & \displaystyle \quad \times \left[{\displaystyle \frac{\unicode[STIX]{x1D701}\left(3,\left(1+{\displaystyle \frac{r_{2}}{r_{1}}}\right)^{-1}\right)\unicode[STIX]{x1D701}\left(3+k,\left(1+{\displaystyle \frac{r_{1}}{r_{2}}}\right)^{-1}\right)-\unicode[STIX]{x1D701}(3)\unicode[STIX]{x1D701}(3+k)}{\unicode[STIX]{x1D701}\left(3,\left(1+{\displaystyle \frac{r_{2}}{r_{1}}}\right)^{-1}\right)\unicode[STIX]{x1D701}\left(3,\left(1+{\displaystyle \frac{r_{1}}{r_{2}}}\right)^{-1}\right)-\unicode[STIX]{x1D701}(3)^{2}}}\right].\end{eqnarray}$$

It can be readily verified that these coincide for $r_{1}/r_{2}\rightarrow 0$ with the value

(4.17) $$\begin{eqnarray}\tilde{\unicode[STIX]{x1D6FA}}_{1}\rightarrow -\mathop{\sum }_{l\geqslant 1}C_{l}P_{l}(\cos \unicode[STIX]{x1D6FC})\mathop{\sum }_{k=0}^{l-1}(-1)^{k}\binom{l-1}{k}{\displaystyle \frac{\unicode[STIX]{x1D701}(3+k)}{\unicode[STIX]{x1D701}(3)}},\end{eqnarray}$$

representing the rotation of a squirmer touching a no-slip wall.

To illustrate the near-field behaviour that can be found exactly using the reciprocal theorem, we give a specific example of a swimmer whose slip velocity is an azimuthal circulation within a polar-cap region of opening angle $\unicode[STIX]{x1D703}_{0}$ . Although crude, this provides a squirmer representation of a rotating flagellar bundle, and counter-rotating cell body. Explicitly, we take the slip velocity to be

(4.18) $$\begin{eqnarray}\boldsymbol{u}^{s}=\left\{\begin{array}{@{}cc@{}}\unicode[STIX]{x1D6FA}_{c}r_{1}\sin \unicode[STIX]{x1D703}_{s}\boldsymbol{s}_{\unicode[STIX]{x1D719}},\quad & 0<\unicode[STIX]{x1D703}_{s}<\unicode[STIX]{x1D703}_{0},\\ -\unicode[STIX]{x1D6FA}_{b}r_{1}\sin \unicode[STIX]{x1D703}_{s}\boldsymbol{s}_{\unicode[STIX]{x1D719}},\quad & \unicode[STIX]{x1D703}_{0}<\unicode[STIX]{x1D703}_{s}<\unicode[STIX]{x03C0},\end{array}\right.\end{eqnarray}$$

as depicted schematically in figure 4(a). The slip velocity within the cap region is $\unicode[STIX]{x1D6FA}_{c}$ , which is balanced by a counter-rotation of the body, $\unicode[STIX]{x1D6FA}_{b}$ , chosen so that the coefficient $C_{1}=0$ to remove any free rotation and focus on the effects of interactions. The squirming coefficients are given by

(4.19) $$\begin{eqnarray}C_{l}=-{\displaystyle \frac{(2l+1)}{4}}\left[\unicode[STIX]{x1D6FA}_{c}\int _{0}^{\unicode[STIX]{x1D703}_{0}}\text{d}\unicode[STIX]{x1D703}\sin ^{2}\unicode[STIX]{x1D703}P_{l}^{1}(\cos \unicode[STIX]{x1D703})-\unicode[STIX]{x1D6FA}_{b}\int _{\unicode[STIX]{x1D703}_{0}}^{\unicode[STIX]{x03C0}}\text{d}\unicode[STIX]{x1D703}\sin ^{2}\unicode[STIX]{x1D703}P_{l}^{1}(\cos \unicode[STIX]{x1D703})\right],\end{eqnarray}$$

and the counter-rotation $\unicode[STIX]{x1D6FA}_{b}$ required to cancel out the free rotation is

(4.20) $$\begin{eqnarray}\unicode[STIX]{x1D6FA}_{b}=\unicode[STIX]{x1D6FA}_{c}{\displaystyle \frac{(2+\cos \unicode[STIX]{x1D703}_{0})}{(2-\cos \unicode[STIX]{x1D703}_{0})}}\tan ^{4}\left({\displaystyle \frac{\unicode[STIX]{x1D703}_{0}}{2}}\right).\end{eqnarray}$$

When $\unicode[STIX]{x1D703}_{0}=\unicode[STIX]{x03C0}/2$ we have that $\unicode[STIX]{x1D6FA}_{b}=\unicode[STIX]{x1D6FA}_{c}$ , as expected on symmetry grounds. E. coli has a body counter-rotation measured to be of the order of one-tenth the rotation of its flagellar bundle, with large variation between specimens (Magariyama, Sugiyama & Kudo Reference Magariyama, Sugiyama and Kudo2001), and inversion of (4.20) gives an appropriate value of approximately $0.28\unicode[STIX]{x03C0}$ for $\unicode[STIX]{x1D703}_{0}$ , which we idealise as $\unicode[STIX]{x03C0}/4$ .

The dependence of the swimmer’s rotation on its orientation at large distances is given by the slowest-decaying squirming mode, $C_{2}$ , and hence by $P_{2}(\cos \unicode[STIX]{x1D6FC})$ , which is head–tail symmetric. However in the near field there may be significant asymmetry in the orientation dependence which could persist for relatively large separations. Figure 4(c) shows how the orientation dependence changes for distances up to 100 times the swimmer radius, for a model E. coli interacting with a no-slip wall. A comparison between the interaction with a no-slip wall of a spherical-cap swimmer calculated using all modes up to $C_{100}$ , and an equivalent squirming sphere with only the dominant far-field $C_{2}$ mode, is shown in figure 4(b) and further illustrates the importance of including higher-order modes in calculating near-field interactions: at separations of the order of the swimmer’s size the effect of including the higher-order modes can be dramatic. In the case of $\unicode[STIX]{x1D6FC}=\unicode[STIX]{x03C0}/4$ the swimmer’s rotation changes sense as it approaches the wall. When $\cos \unicode[STIX]{x1D6FC}=3^{-1/2}$ the contribution of the $C_{2}$ mode is identically zero since $P_{2}(3^{-1/2})=0$ ; however there is still motion driven by higher-order modes of non-negligible magnitude.

4.3 Asymptotics at large separation

The tendency of flagellated bacteria such as E. coli and Vibrio alginolyticus to follow distinctive circular trajectories near boundaries (DiLuzio et al. Reference DiLuzio, Turner, Mayer, Garstecki, Weibel, Berg and Whitesides2005) has resulted in considerable theoretical work on the rotation of a swimmer normal to a planar surface (Lauga et al. Reference Lauga, Diluzio, Whitesides and Stone2006; Di Leonardo et al. Reference Di Leonardo, Dell’Arciprete, Angelani and Iebba2011; Lopez & Lauga Reference Lopez and Lauga2014; Papavassiliou & Alexander Reference Papavassiliou and Alexander2015). This phenomenon has been explained by noting that the rotation of the body and counter-rotation of the tail gives a flow field resembling a rotlet dipole and, when aligned parallel to a wall (as extensile swimmers generically do (Spagnolie & Lauga Reference Spagnolie and Lauga2012)), results in rotation by a mechanism analogous to the turning of a tank. The availability of the exact solutions, (4.13) and (4.14), for this component of motion allow extension to further geometries than those appearing in the literature.

Figure 3. The rotational speed, $\tilde{\unicode[STIX]{x1D6FA}}$ , due to the $C_{2}$ squirming mode, in units of $C_{2}/P_{2}(\cos \unicode[STIX]{x1D6FC})$ , as a function of $d$ . (a) Near a no-slip (black) and free (red) planar boundary, compared to the $d^{-4}$ decay predicted by approximate models (grey dashed). The rotation near a free surface has the opposite sense to that near a solid boundary. (b) Inside a shell of radius 1.2 (blue), 1.5 (grey), 2 (black) and 4 (red), and the wall limit (dashed). Inset: behaviour at small separation. (c) Near a tracer of radius 0.5 (orange), 1 (grey), 2 (black) and 10 (blue), and the wall limit (black dashed).

The flow field generated by the $C_{l}$ contribution to the slip velocity has an asymptotic decay of $d^{-(l+2)}$ in an unbounded domain (Pak & Lauga Reference Pak and Lauga2014). Focusing on the $l=1$ contribution we recognise the sum in (4.13) and (4.14) as the torque, (3.6) and (3.5) respectively. Hence the contribution to the rotation from this squirming mode is $-C_{1}\cos \unicode[STIX]{x1D6FC}$ for any separation and in any configuration. This is precisely the same as the rotation found asymptotically (Pak & Lauga Reference Pak and Lauga2014; Davis & Crowdy Reference Davis and Crowdy2015), and demonstrates that this mode corresponds only to self-rotation and does not result in interaction. Therefore, the slowest-decaying contribution to the normal rotation of a squirmer due to interactions is from $C_{2}$ , which represents a rotlet dipole. Here we discuss the interaction of this squirming mode with a passive sphere as a leading-order behaviour which is generic for all swimmers.

Using dimensional analysis, Lopez & Lauga (Reference Lopez and Lauga2014) argued that a squirmer circling parallel to a wall has an asymptotic angular frequency decaying no slower than $d^{-4}$ ; this was confirmed by considering the flow induced by a rotlet dipole near a wall, and indeed has been found to be the leading-order behaviour of a swimmer with arbitrary azimuthal slip velocity near a wall, with $\tilde{\unicode[STIX]{x1D6FA}}_{1}=C_{2}(r_{1}/2d)^{4}P_{2}(\cos \unicode[STIX]{x1D6FC})/5$  (Papavassiliou & Alexander Reference Papavassiliou and Alexander2015). Figure 3(a) shows that this behaviour agrees with our exact solution (4.13) up to a separation of approximately a squirmer diameter. An explicit form for the rotation near a wall is obtained from equation (4.13) by setting $\unicode[STIX]{x1D709}_{2}=0$ , $\tilde{\unicode[STIX]{x1D6FA}}_{2}=0$ and $\unicode[STIX]{x1D709}_{1}=\log (d/r_{1}+\sqrt{d^{2}/r_{1}^{2}-1})$ , and it is found that the rotation induced by interaction with the wall goes as $d^{-(l+2)}$ for $C_{l}$ . The behaviour in a shell, shown in figure 3(b), also has this form since the separation is always smaller than the radius of curvature of the shell. Note that the rotation of the swimmer when it is precisely in the centre is zero by symmetry and changes sense as the swimmer crosses between hemispheres.

Figure 4. The behaviour of a ‘spherical-cap’ type swimmer near a wall, calculated using squirming modes up to $C_{100}$ . (a) Schematic of the swimmer. (b) Near-field discrepancy between exact solution (solid) and asymptotic $C_{2}$ mode behaviour (dashed) for swimmer with $\unicode[STIX]{x1D703}_{0}=\unicode[STIX]{x03C0}/4$ and $\unicode[STIX]{x1D6FC}=\unicode[STIX]{x03C0}/4$ (black), $\cos \unicode[STIX]{x1D6FC}=3^{-1/2}$ (blue) and $\cos \unicode[STIX]{x1D6FC}=0.783$ (red). (c) Orientation dependence of rotation near a wall as a function of distance for $\unicode[STIX]{x1D703}_{0}=\unicode[STIX]{x03C0}/4$ . Rotation is normalised by the $\unicode[STIX]{x1D6FC}$ -average,  $\langle \tilde{\unicode[STIX]{x1D6FA}}\rangle _{\unicode[STIX]{x1D6FC}}=\sum _{l}(l+1/2)^{-1}\int _{0}^{\unicode[STIX]{x03C0}}\text{d}\unicode[STIX]{x1D6FC}\,\sin \unicode[STIX]{x1D6FC}\,P_{l}(\cos \unicode[STIX]{x1D6FC})\tilde{\unicode[STIX]{x1D6FA}}$ . (d) Orientation dependence of normalised rotation near a free surface as a function of distance for $\unicode[STIX]{x1D703}_{0}=\unicode[STIX]{x03C0}/4$ .

The asymptotic rotation of a swimmer in the presence of a tracer may be calculated using the leading-order forms $\unicode[STIX]{x1D709}_{1}\sim \log (r_{1}/d)$ and $\unicode[STIX]{x1D709}_{2}\sim \log (r_{2}/d)$ , giving a decay of

(4.21) $$\begin{eqnarray}\tilde{\unicode[STIX]{x1D6FA}}_{1}\sim C_{2}P_{2}(\cos \unicode[STIX]{x1D6FC}){\displaystyle \frac{r_{1}^{4}r_{2}^{5}}{d^{9}}}.\end{eqnarray}$$

This may be understood in terms of multipole reflections (Kim & Karrila Reference Kim and Karrila2013). At large separation the swimmer’s motion is driven by the flow reflected in the tracer, which has the leading behaviour of a stresslet since the tracer must remain force free. Dimensional analysis suggests that the reflected flow at the swimmer should have a strength going as $d^{-6}$ , and therefore a vorticity of $d^{-7}$ , but for this case of an axisymmetric, azimuthal flow the leading reflected flow is identically zero, so the rotation is driven by a vorticity of $d^{-9}$ . Figure 3(c) shows a cross-over to this behaviour when the separation exceeds the radius of the tracer. In the near field the passive sphere resembles as a wall and we see a $d^{-4}$ dependence of the swimmer’s rotational speed. For the passive sphere (not shown) the cross-over is not seen, and the asymptotic interaction is $d^{-4}$ , equal to the asymptotic vorticity generated by a rotlet dipole; thus, the dominant effect of the $C_{2}$ squirming mode is the motion of the tracer, and by superposition two squirmers with this slip velocity would tend to move each other more than themselves.

4.4 Rotation near a free surface

An interesting application of the exact solution presented above is to find the rotation of a squirmer close to a free surface. It has been hypothesised (Lauga et al. Reference Lauga, Diluzio, Whitesides and Stone2006) and subsequently observed experimentally (Di Leonardo et al. Reference Di Leonardo, Dell’Arciprete, Angelani and Iebba2011) that the circular trajectories of E. coli near a free surface have the opposite sense to those near a no-slip wall. Using the known hydrodynamic solution for the rotation of a sphere beneath the interface between two fluid phases we find that both cases of rotation near a wall and a free surface may be described as image systems, using the two-sphere solution presented previously. This allows the swimming close to such boundaries to be found and compared without further calculation, and we find that the change of direction depending on the type of boundary is generic and explained by these image systems.

If a sphere rotates beneath the flat interface between the fluid that contains it, and another fluid of viscosity $\tilde{\unicode[STIX]{x1D707}}$  (O’Neill & Ranger Reference O’Neill and Ranger1979), the flow may be found explicitly by supposing an ansatz of the form (3.2) in each phase and matching flow and stress across the boundary. Then the torque on the sphere is

(4.22) $$\begin{eqnarray}T_{1}=8\unicode[STIX]{x03C0}\unicode[STIX]{x1D707}\unicode[STIX]{x1D714}_{1}r_{1}^{3}\mathop{\sum }_{n=0}^{\infty }(-\unicode[STIX]{x1D6EC})^{n}{\displaystyle \frac{\sinh ^{3}\unicode[STIX]{x1D709}_{1}}{\sinh ^{3}(n+1)\unicode[STIX]{x1D709}_{1}}},\end{eqnarray}$$

where $\unicode[STIX]{x1D6EC}=(\unicode[STIX]{x1D707}-\tilde{\unicode[STIX]{x1D707}})/(\unicode[STIX]{x1D707}+\tilde{\unicode[STIX]{x1D707}})$ . When $\unicode[STIX]{x1D6EC}=-1$ the empty phase is infinitely viscous and corresponds a no-slip wall; instead, when $\unicode[STIX]{x1D6EC}=+1$ the boundary is a free surface. Since the torque (3.5) corresponding to the two-sphere solution when $r_{1}=r_{2}=r$ is

(4.23) $$\begin{eqnarray}T_{1}=8\unicode[STIX]{x03C0}\unicode[STIX]{x1D707}r^{3}\mathop{\sum }_{n=0}^{\infty }\left[\unicode[STIX]{x1D714}_{1}{\displaystyle \frac{\sinh ^{3}\unicode[STIX]{x1D709}_{1}}{\sinh ^{3}(n+1)\unicode[STIX]{x1D709}_{1}}}-(\unicode[STIX]{x1D714}_{1}+\unicode[STIX]{x1D714}_{2}){\displaystyle \frac{\sinh ^{3}\unicode[STIX]{x1D709}_{1}}{\sinh ^{3}(2n+2)\unicode[STIX]{x1D709}_{1}}}\right],\end{eqnarray}$$

it can be seen that the result for a free surface is recovered when $\unicode[STIX]{x1D714}_{2}=\unicode[STIX]{x1D714}_{1}$  (Brenner Reference Brenner1964), while the result for a no-slip wall is given by $\unicode[STIX]{x1D714}_{2}=-\unicode[STIX]{x1D714}_{1}$ . Hence a rotating sphere near a free surface has as its image system a corotating sphere which decreases the torque compared to the free-space value, while near a wall the image system is an antirotating sphere which increases the torque.

The rotation near a free surface may then be calculated exactly using the reciprocal theorem and compared to our expressions for squirming near a wall, (4.13). Although the activity of the squirmer generates tangential flows on the interface, since the stress in the conjugate problem is zero there is no contribution to the reciprocal theorem from an integral over the free surface and an expression for the rotation is obtained immediately from (4.14) by substituting $\unicode[STIX]{x1D709}_{2}=-\unicode[STIX]{x1D709}_{1}$ and $\unicode[STIX]{x1D714}_{2}=\unicode[STIX]{x1D714}_{1}$ , giving

(4.24) $$\begin{eqnarray}\tilde{\unicode[STIX]{x1D6FA}}_{1}T_{1}=-8\unicode[STIX]{x03C0}\unicode[STIX]{x1D707}r_{1}^{3}\unicode[STIX]{x1D714}_{1}\mathop{\sum }_{l=2}^{\infty }C_{l}P_{l}(\cos \unicode[STIX]{x1D6FC})\mathop{\sum }_{n=0}^{\infty }(-1)^{n}{\displaystyle \frac{\sinh ^{3}\unicode[STIX]{x1D709}_{1}\sinh ^{l-1}n\unicode[STIX]{x1D709}_{1}}{\sinh ^{l+2}(n+1)\unicode[STIX]{x1D709}_{1}}}.\end{eqnarray}$$

In both cases of a wall and a free surface the leading far-field contribution from the $l$ -th squirming mode is equal and opposite, with a strength $\tilde{\unicode[STIX]{x1D6FA}}_{1}\propto d^{-(l+2)}$ . Hence the nature of the boundary determines the sense of rotation generically in the asymptotic limit. In the contact limit, which for a wall is given by (4.17) and for a free surface

(4.25) $$\begin{eqnarray}\tilde{\unicode[STIX]{x1D6FA}}_{1}\rightarrow -{\displaystyle \frac{4}{3}}\mathop{\sum }_{l=2}^{\infty }C_{l}P_{l}(\cos \unicode[STIX]{x1D6FC})\mathop{\sum }_{k=0}^{l-1}(-1)^{k}\binom{l-1}{k}\left(1-2^{-(k+2)}\right){\displaystyle \frac{\unicode[STIX]{x1D701}(3+k)}{\unicode[STIX]{x1D701}(3)}},\end{eqnarray}$$

we find that the contribution to the rotation from each squirming mode is smaller at a free surface than at a wall, and the ratio of these contributions has a faster-than-exponential decay with increasing $l$ , indicating that higher-order effects due to microscopic details of a swimmer are less important at a free surface than at a wall. This can be seen in figure 4(d), which gives the orientation dependence of the same spherical-cap swimmer as considered before near a free surface; compared to the analogous trace for rotation near a wall, figure 4(c), the one for a free surface has the opposite sign and is smoother.

The reciprocal theorem relies on the swimmer problem and the conjugate Stokes drag being defined in the same region. Here we have assumed that the free surface does not deform in either solution, so that this region is the half-space with an embedded sphere in both cases; however, there is no reason not to expect deformation, particularly in the close proximity regime and furthermore, it cannot be assumed that this deformation of the surface will be the same for a swimmer and a dragged sphere. Nevertheless, if the deformation of the surface is sufficiently small it may be treated in a linear fashion by projecting onto the plane and expressing as a slip velocity in both solutions (much as Lighthill’s deforming squirmer has its activity projected onto the surface of a sphere for determination of the swimming speed (Lighthill Reference Lighthill1952)).

4.5 Translation

Phenomena such as the tendency of swimming microorganisms to aggregate at surfaces (Berke et al. Reference Berke, Turner, Berg and Lauga2008) have been neatly explained by modelling a swimmer as a collection of point singularities. The sign of a swimmer’s stresslet, which characterises it as extensile or contractile (Ramaswamy Reference Ramaswamy2010), causes it to align parallel or normal to the wall, respectively (Spagnolie & Lauga Reference Spagnolie and Lauga2012), while the inclusion of a source dipole ensures self-propulsion (Drescher et al. Reference Drescher, Goldstein, Michel, Polin and Tuval2010). By adopting the reciprocal theorem the behaviour due to any slip velocity may, in principle, be found (Papavassiliou & Alexander Reference Papavassiliou and Alexander2015); while here the conjugate solution used restricts us to axisymmetric motions, we have the freedom to generalise to curved surfaces. The actual calculation for the translational motion is analogous to the calculation for rotation shown in § 4.2, but is rather more involved and will not be shown explicitly; explicit expressions are given in appendix B. The general result for an arbitrary squirming mode has not yet been found and each contribution must be calculated separately. Instead, we will attempt to describe the behaviour using a few illustrative examples.

Figure 5. The interactions of the $A_{1}$ (dashed) and $B_{1}$ (solid) modes. (a) The speed $\tilde{U}$ of a squirmer, in units of the free swimming speed $U_{free}$ , as a function of $d$ inside a shell of radius 1.2 (blue), 1.5 (grey), 2 (black), 4 (red). The wall limits are shown as dotted lines. (b) The speed difference $\unicode[STIX]{x0394}\tilde{U} =|U_{c}-U_{free}|$ at the centre of a shell as a function of shell radius $r_{2}$ , showing an excluded volume dependence. (c) $\tilde{U}$ as a function of $d$ near a tracer of radius 0.5 (orange), 1 (grey), 2 (black) and 10 (blue).

We consider the first few translational squirming modes, $A_{1}$ , $B_{1}$ , $A_{2}$ and $B_{2}$ . The first two of these set the self-propulsive speed in free space and asymptotically resemble source dipoles. $A_{2}$ and $B_{2}$ give the asymptotic stresslet of the swimmer (Ishikawa et al. Reference Ishikawa, Simmonds and Pedley2006) and while they generate no motion in an unbounded domain they are of fundamental importance in the interactions of the swimmer with boundaries, since both the swimmer and the boundaries must remain force free and the lowest-order image singularity will be a stresslet. In the far field these point-singularity descriptions are sufficient to fully characterise the generic behaviour (Spagnolie & Lauga Reference Spagnolie and Lauga2012). For the special case of interaction with a wall an explicit asymptotic estimate (Davis & Crowdy Reference Davis and Crowdy2015; Papavassiliou & Alexander Reference Papavassiliou and Alexander2015) of the swimming speed is available as

(4.26) $$\begin{eqnarray}{\displaystyle \frac{\text{d}d}{\text{d}t}}={\displaystyle \frac{1}{3}}(2B_{1}-A_{1})-{\displaystyle \frac{1}{5}}(B_{2}-A_{2})\left({\displaystyle \frac{r_{1}}{2d}}\right)^{2}P_{2}(\cos \unicode[STIX]{x1D6FC});\end{eqnarray}$$

thus, asymptotically, $A_{1}$ and $B_{1}$ contribute behaviour that differs only in a numerical factor, while behaviour due to $A_{2}$ and $B_{2}$ is distinguishable only by a sign change.

Figure 6. The speed $\tilde{U}$ of a squirmer due to modes $A_{2}$ and $B_{2}$ , normalised so that $A_{2}P_{l}(\cos \,\unicode[STIX]{x1D6FC})=B_{2}P_{l}(\cos \,\unicode[STIX]{x1D6FC})=1$ , as a function of $d$ . Dashed line is $A_{2}$ , solid line is $B_{2}$ . (a) Interaction with a no-slip wall. Dotted grey line is point-singularity approximation. (b) Interaction with a shell of radius 1.2 (blue), 1.5 (grey), 2 (black) and 4 (red). $\tilde{U} =0$ in the centre of the shell and the motion is equal and opposite in the other hemisphere. (c) Interaction with a tracer of radius 0.5 (orange), 1 (grey), 2 (black) and 10 (blue), with the red dotted line demonstrating the wall limit.

Figures 5 and 6 show, respectively, the normal translation speed of a swimmer induced by the $A_{1}$ and $B_{1}$ , and the $A_{2}$ and $B_{2}$ modes respectively, in interaction with a passive concave or convex sphere. It can be seen that the far-field equivalence of $A_{1}$ and $B_{1}$ , and $A_{2}$ and $B_{2}$ , also holds for the concave and convex geometries. These figures indicate that the cross-over to far-field behaviour that is well described by point-singularity models (Ishikawa et al. Reference Ishikawa, Simmonds and Pedley2006; Spagnolie & Lauga Reference Spagnolie and Lauga2012; Davis & Crowdy Reference Davis and Crowdy2015; Papavassiliou & Alexander Reference Papavassiliou and Alexander2015) occurs at very small separations, of the order of a few swimmer diameters.

The $A_{2}$ and $B_{2}$ squirming modes generate an asymptotic flow field of $d^{-2}$ , while the propulsive modes $A_{1}$ and $B_{1}$ give a flow field decaying as $d^{-3}$ . Hence mixing is dominated by the swimmers’ dipoles, and the speed of a passive tracer has a dependence of $d^{-2}$ , by Fáxen’s law, until the separation becomes small and higher-order effects become important. Since $A_{2}$ and $B_{2}$ do not drive self-propulsion, the motion of the swimmer resulting from these modes is due to reflected flow in the boundary of the passive sphere. At separations smaller than the tracer’s radius of curvature the leading order of the reflected flow is equal to that of the flow itself, and gives rise to a local $d^{-2}$ behaviour. As the separation increases the finite size of the tracer becomes important. Higdon (Reference Higdon1979) gives the image system for a force dipole in a fixed, finite-sized sphere as the sum of a Stokeslet with leading-order strength proportional to $d^{-2}$ and a dipole with strength ${\sim}d^{-3}$ . Thus, if the passive sphere were fixed the leading-order reflection would go as $d^{-3}$ , but as it is free to move in such a way as to cancel any force acting on it, we see $d^{-5}$ . This dependence may also be calculated using a second-order multipole expansion, in which case the leading-order motion of the swimmer is driven by the reflected stresslet inside the tracer (Kim & Karrila Reference Kim and Karrila2013). The cross-over between the two types of behaviour is shown in figure 6(c).

Figure 7. Collision trajectories from exact solutions compared to approximate results of Papavassiliou & Alexander (Reference Papavassiliou and Alexander2015). (a) Collision trajectories for slip velocity given by first-order modes $A_{1}$ (red) and $B_{1}$ (black). (b) Collision trajectories for slip velocity given by second-order modes $A_{2}$ (red) and $B_{2}$ (black). The trajectory predicted by point-singularity description is shown as a grey dotted line.

The availability of exact solutions for the motion due to these squirming modes means it is possible to calculate explicit trajectories in time, albeit only for motion along the common diameter of two spheres. Hence we consider a head-on collision of the swimmer with a wall, which allows comparison to analogous trajectories calculated by integrating the approximate results (4.26). This is shown in figure 7, and it can be seen that the trajectories differ very little between radial and tangential modes (red and black lines respectively), and the corresponding asymptotic approximation (grey dashed line). These trajectories become distinct only at very small separation, again of the order of a swimmer diameter.

In contrast, the near-field behaviour due to radial and tangential slip is rather different. It can be seen from figures 5 and 6 that the tangential modes  $B_{1}$ and  $B_{2}$ result in a swimming speed which goes to zero as contact with the surface is approached; this results in collision taking a longer time than predicted by a point-singularity. The radial modes $A_{1}$ and $A_{2}$ result in acceleration to a finite speed as contact is approached; this results from the incompatibility of the boundary conditions of no-slip and radial flow on two touching surfaces. The consequences of this can be seen in figure 7, where the radial modes $A_{1}$ and $A_{2}$ result in collisions in finite time while the tangential modes $B_{1}$ and $B_{2}$ cause deceleration close to contact. Although we have been unable to explicitly integrate the expressions of the swimming speed to determine whether physical contact occurs within finite time or not, we note that the exact solution for a swimming disc with tangential slip collides with a wall in infinite time, as may be verified using the equations of motion calculated by Crowdy (Reference Crowdy2011).

The hydrodynamic force in the conjugate problem is divergent in the near field. This divergence has a leading part going as $F_{i}\sim \unicode[STIX]{x1D709}_{i}^{-2}$ , but the representation of the force as an infinite sum contains a harmonically divergent subleading term. By approximating the sum as an integral Cox & Brenner (Reference Cox and Brenner1967) found that this subleading divergence may be expressed as ${\sim}\log \unicode[STIX]{x1D709}_{i}$ as $\unicode[STIX]{x1D709}_{i}\rightarrow 0$ . When calculating swimmer motions using the reciprocal theorem the force appears as a denominator, with the numerator given by the integral of the slip velocity against the conjugate stress tensor. This numerator also diverges, although no faster than the force; specifically, the divergence is the same as for the force for the integrals involving $A_{1}$ and $A_{2}$ , and one power of $\unicode[STIX]{x1D709}$ slower for those with $B_{1}$ and $B_{2}$ . This is enough to explain the behaviour shown in figures 5–7 and a detailed analysis of the subleading terms in the vein of Cox & Brenner (Reference Cox and Brenner1967) is unnecessary.

An interesting consequence of the near-field distinction between radial and tangential slip is the behaviour of a swimmer inside a small shell, with a radius smaller than the threshold for cross-over to asymptotic behaviour. When a $B_{1}$ swimmer is inside such a shell the swimming speed is attenuated by the presence of the boundaries, and by symmetry attains a maximum value, $U_{c}$ , in the centre of the shell. Conversely, a swimmer with $A_{1}$ activity has an increased speed due to the interactions, as a result of the divergent interaction of the radial modes near boundaries. A comparison is shown for a variety of shell sizes in figure 5(a). The speed at the centre of the shell, $U_{c}$ , may be calculated analytically and depends on the relative sizes as

(4.27) $$\begin{eqnarray}U_{c}=U_{free}-{\displaystyle \frac{5}{3}}(A_{1}+B_{1})\left({\displaystyle \frac{r_{1}}{r_{2}}}\right)^{3}\left({\displaystyle \frac{1-\left({\displaystyle \frac{r_{1}}{r_{2}}}\right)^{2}}{1-\left({\displaystyle \frac{r_{1}}{r_{2}}}\right)^{5}}}\right).\end{eqnarray}$$

Thus as the shell becomes large, $U_{c}$ approaches the free swimming speed in proportion to the volume of fluid displaced by the swimmer, see figure 5(b). The same occurs in two dimensions where an exact result is available for the swimming of an active disc inside a circular boundary (Papavassiliou & Alexander Reference Papavassiliou and Alexander2015) and where the maximum speed of a self-propulsive swimmer approaches its free-space value in proportion to the excluded area.