Exact solutions for hydrodynamic interactions of two squirming spheres

We provide exact solutions of the Stokes equations for a squirming sphere close to a no-slip surface, both planar and spherical, and for the interactions between two squirmers, in three dimensions. These allow the hydrodynamic interactions of swimming microscopic organisms with confining boundaries, or each other, to be determined for arbitrary separation and, in particular, in the close proximity regime where approximate methods based on point singularity descriptions cease to be valid. We give a detailed description of the circular motion of an arbitrary squirmer moving parallel to a no-slip spherical boundary or flat free surface at close separation, finding that the circling generically has opposite sense at free surfaces and at solid boundaries. While the asymptotic interaction is symmetric under head-tail reversal of the swimmer, in the near field microscopic structure can result in significant asymmetry. We also find the translational velocity towards the surface for a simple model with only the lowest two squirming modes. By comparing these to asymptotic approximations of the interaction we find that the transition from near- to far-field behaviour occurs at a separation of about two swimmer diameters. These solutions are for the rotational velocity about the wall normal, or common diameter of two spheres, and the translational speed along that same direction, and are obtained using the Lorentz reciprocal theorem for Stokes flows in conjunction with known solutions for the conjugate Stokes drag problems, the derivations of which are demonstrated here for completeness. The analogous motions in the perpendicular directions, i.e. parallel to the wall, currently cannot be calculated exactly since the relevant Stokes drag solutions needed for the reciprocal theorem are not available.


I. INTRODUCTION
Swimming microorganisms do not live in an infinite, unbounded fluid domain, but instead inhabit complex geometries confined by fluid interfaces and solid boundaries, and populated by other organisms and passive particles. Much, if not most, of the rich variety of behaviour that is seen [50,62] cannot be explained outside of the context of confinement. The most basic interaction is of a single swimmer with a boundary or object; even in such cases we see striking behaviour, such as the 'waltzing' of a pair of Volvox colonies [25]. Flagellated bacteria such as Escherichia coli and Vibrio alginolyticus are known to trace out circles [1,21,58] near solid boundaries and, remarkably, if the boundary is replaced by a free surface the direction of rotation changes [20,49]. An important effect is the attraction of swimmers to boundaries, noted by Rothschild [83] for bull spermatozoa but subsequently also seen in bacteria [2,29], which is thought to contribute to the navigation of sperm cells in the female reproductive system [18] and plays a fundamental role in biofilm formation at surfaces [72]. Swimmers adhere to surfaces due to strong lubrication forces, causing catalytic selfpropelled rods to orbit colloidal spheres [91] and the microalga Chlamydomonas reinhardtii to orbit cylindrical posts until flagellar beating detaches it [11]. As the number of interacting components increases so does the complexity of the behaviour, resulting in the phase behaviours seen in dense suspensions of bacteria, including long-range orientational order [10], and the formation of largescale turbulent structures [22,27] and stable spiral vortices [97]. The latter can only occur in confinement. An understanding of these phenomena is important for the design of microfluidic systems, such as devices to direct swimmers [18] and harness them for mass transport [47,95], and to extract mechanical work from their activity [19]; and to construct biomimetic artificial swimmers [26,33,38,75] which may fulfil a number of nanotechnological and medical roles.
The interaction of swimming microorganisms with each other and their environment is a complex combination of several factors, including biology, such as the taxis which allows them to move in search of food and tolerable living conditions, and to aggregate and form patterns through chemical signalling and quorum sensing [9,77]; and hydrodynamics and physical contact. It has been shown that the scattering of swimmers from planes [44] and posts [11] arises from physical contact of the flagellae with the surface, so that hydrodynamics is not the dominant contributor to the phenomena in these cases. Nonetheless, the fact that there is physical contact with the surface emphasises that any hydrodynamic effects have to be considered in this contact regime. Other cases are less clear-cut; for instance, the typical density profile of a suspension of swimmers close to a wall has been reproduced both by considering hydrodynamics [2] and Brownian motion combined with collisions [53].
The importance of hydrodynamics to the interactions of swimmers means exact solutions to the Stokes equations are desirable. In the case of a single swimmer in an unbounded domain several such solutions exist, notably for the motion of a single axisymmetric squirmer [4,54], later generalised to non-axisymmetric slip velocities [73], and for the motion of a 'treadmilling' spheroidal [52] or toroidal swimmer [51,80,93], as well as a two-dimensional analogue for a squirming cylinder [3] or waving sheet [92]. The squirmer solutions have been used to find the advection of tracer particles due to a squirmer [81]. Dropping instead to two dimensions a number of additional solutions in confinement become available using conformal mapping techniques, such as the motion of an active cylinder near a planar or concave boundary [13,14,74], or under a free surface [15]. However, in three dimensions, hydrodynamic interactions have only been calculated by approximate methods. For instance, Ishikawa et al. [39] find the far-field interactions between two squirmers by considering multipole expansions of stresslets, and the near-contact interactions using lubrication theory. In many cases point singularity methods are valuable. The image systems calculated by Blake & Chwang [5] allow the interaction with walls to be found [87,99], recovering the experimentally observed attraction to walls [2] and the swimming in circles close to surfaces [49,74]. Other phenomena relevant to low Reynolds number swimming that have been successfully treated using point singularities include flagellar beating for feeding [37] and resulting in synchronisation [8]. A comparison of real flow fields with point-singularity approximations shows them to be in good agreement even close to the organism in several cases [23,24].
We find exact solutions for the axisymmetric translation and rotation of a squirmer in the presence of a spherical or planar boundary. These are valid at any separation, both in the far-field where point singularity solutions are accurate and in the contact limit of vanishing separation, where such approximate solutions are not accurate. They also account for any type of squirming motion and not simply the lowest order modes considered by Ishikawa et al. [39] and that point singularity descriptions are restricted to. These solutions are obtained using the Lorentz reciprocal theorem for Stokes flows [35], first applied to calculate the motion of individual swimmers by Stone & Samuel [90] and recently extended to a many-body setting [74]. In this form of the reciprocal theorem, the stress tensor associated with the Stokes drag on an object of the same shape as the swimmer serves as the integration kernel to extract the speed and angular frequency of the swimming from the slip velocity. This can be viewed as a specialisation of the boundary-element method [79], simplified by the requirement that a swimmer be free of net forces and torques. The simplicity of the Stokes drag solution on a sphere means this calculation is straightforward for a single spherical microorganism, although since the full hydrodynamic solution of a squirming sphere exists [4,54,73] the only advantage of the reciprocal theorem is computational convenience. Indeed, the swimming speed found for an active sphere self-propelling by means of a metachronal wave on its surface by Stone & Samuel [90] had been derived by other means not long previously by Ehlers et al. [28]. Nevertheless the simplicity of the calculation means it has become a standard tool in the active matter literature [34,88]. More recently there have been interesting extensions of the reciprocal theorem and other related integral theorems to cases such as propulsion by the Marangoni effect [63] and self-propulsion through viscoelastic and non-Newtonian fluids [48] where direct solutions are not so readily available.
The reciprocal theorem is immediately applicable to swimmers in confined fluid domains. This has been exploited in two dimensions to find the motion of squirming [13] and self-diffusiophoretic cylinders close to walls [14] using as a conjugate solution the Stokes drag on a cylinder in the half-space [43], and it has been noted that the reciprocal theorem may be used to find the motion of any number of swimmers [74]. More generally, the reciprocal theorem may be used to find the full hydrodynamics for a given active problem, provided the existence of an appropriate conjugate solution: If the Green's function for the Stokes equations in a particular confined geometry is known, the reciprocal theorem extracts the flow due to activity on the boundaries. Following such an approach, Michelin & Lauga [65] found the fluid flux through an active pipe by relation to the flow solution for a no-slip channel. Furthermore, an approximate integration kernel for a swimmer in a given geometry may be constructed using one of the many existing flows for point forces; relevant examples include the solution for point singularities near walls [5], outside spheres [37], between two plates [55] or in a cylindrical pipe [56].
By following this approach, we are able to find exact solutions for swimmer interactions by using exact solutions for the conjugate Stokes drag problem. Such solutions are available for the Stokes drag on a pair of spheres, or of a single sphere close to a planar wall or fluid interface. The symmetries of the geometry mean there are two independent directions, namely the common diameter of the two spheres and any axis perpendicular to this, and the solution consists of translation and rotation in each of these directions, so that the general motion separates into four components that can be treated individually. The axisymmetric rotation was solved exactly by Jeffery [41], and has since been supplemented by the closely related solution for rotation of a sphere beneath a planar fluid interface [7,70]. The solution for axisymmetric translation was given by Stimson & Jeffery [89] and found to be in remarkably good agreement with experiment [36]. The special case of sedimentation of a sphere towards a solid plane was subsequently given a more detailed analysis both in the limit of large separation [6] and of contact [12], with the latter giving a comparison to results obtained from lubrication theory. Several attempts to find the non-axisymmetric motions and rotations have been unable to give the solution in a closed form; the problem is reduced to a system of difference equations, of which an analytic solution has not been found. Nevertheless it is possible to compute the flow to any degree of accuracy [17, 30-32, 67, 68, 71].
All of these results rely on the use of bispherical coordinates, in which any configuration of two convex or concave spherical boundaries, as well as the intermediate limit of a plane, is an isosurface. This coordinate system greatly simplifies the imposition of boundary conditions; furthermore, since it is conformally equivalent to spherical coordinates, Laplace's equation is separable [40], allowing a general solution to the Stokes equations to be written down [42]. Another notable application of Stimson & Jeffery's solution to swimmer problems is to study the hydrodynamics of catalytic dimers. Catalytic dimers are artificial self-propelled particles composed of a pair of chemically active colloidal beads powered by self-diffusiophoresis [84]; their simplicity facilitates manufacture and allows experiments involving many interacting units [94]. Using the reciprocal theorem together with bispherical coordinates, Popescu et al. [78] and Michelin & Lauga [64] calculated exact expressions for the propulsion speed of catalytic dimers, and were able to discuss optimisation of their swimming speed through changes to the relative sizes and separation of the two beads. Bispherical coordinates also provide a way to calculate hydrodynamic interactions of two spherical objects, such as a sphere sedimenting against a plane [12]; here we investigate their use in calculating interactions driven by force-free swimming, for which the reciprocal theorem is ideally suited [74]. In this way Mozaffari et al. [66] and Sharifi-Mood et al. [85] found the interaction of spherical self-diffusiophoretic particles with each other and with planar boundaries, finding that the chemical interaction is dominant over the hydrodynamic and usually results in repulsion, except where coverage of the chemically active site over the swimmers is large. Their consideration of non-axisymmetric components of motion necessitated numerical solution.
We outline the equations of viscous flow and discuss the use of the reciprocal theorem to obtain exact solutions for interactions in §II. In §III we review the Stokes drag solutions of Jeffery [41] and Stimson & Jeffery [89] that we use with the reciprocal theorem. This section is included for reference and may be skipped if desired. The main results of this paper are contained in §IV, where we find the motion of a squirmer interacting with a passive spherical boundary. The contribution to the motion from the azimuthal squirming coefficients is found explicitly for all orders and is shown for a model organism driven by a rotating cap, while a simple extension to case of interaction with a planar free surface discussed in §IV D. The meridional and radial squirming coefficients do not at present have a general form for the interaction valid at all orders; in §IV E we calculate the interaction due to the lowest two orders of these modes. Finally we discuss the results obtained and possible extensions to the work presented here in §V.

II. STOKES FLOWS AND THE RECIPROCAL THEOREM
The motions of a collection of swimmers or active particles can be determined from the Lorentz reciprocal theorem for Stokes flows. Specifically, for a collection of N force and torque-free swimmers generating motion through active surface slip velocities u s i on their boundaries ∂D i , i = 1, . . . , N , their translational speedsŨ i and rotationsΩ i are given by [74,90] Here,n is the unit outward normal to the fluid domain and σ ↔ is the stress tensor of a conjugate Stokes flow solution for the same set of particles acted upon by forces F i and torques T i , and with no slip boundary conditions. Thus the fundamental object in application of the reciprocal theorem to swimmer problems is the normal stress of the Stokes drag problem, where p is the pressure and µ the viscosity, whose integral against the slip velocities yields the swimmer motion. For instance, in the case of a single spherical swimmer the classic flow solution for the Stokes drag on a single sphere may be used to calculate the swimming speed of a squirmer [90]. Of course, for a spherical squirmer the full flow field, in addition to the swimmer motion, may be calculated directly [4,54,73], without significant additional effort, and gives more information. In the case of hydrodynamic interactions between two objects, even as simple as two squirming spheres, a direct solution for the full flow is not available.
However, full solutions are available for the Stokes drag of two spheres under conditions of axisymmetry, both for rotational and translational motion. These may be used in the reciprocal theorem to deduce the corresponding axisymmetric interactions of an arbitrary pair of squirming spheres, or of a single squirming sphere with a spherical, or planar, boundary. The solution is founded upon having an expression for the stress tensor for a conjugate Stokes drag problem. The work of Jeffery [41] and Stimson & Jeffery [89], as well as subsequent extensions and generalisations [6, 7, 12, 17, 45, 60, 67-69, 71, 76], gives the flow for these problems, from which the stress can be computed directly. However, to keep our work self-contained and in a consistent notation, b. c. It should be noted that the reciprocal theorem may also be used to calculate approximate expressions for the motion of a squirmer near a wall [16,74]: since the leading-order flow about a sedimenting sphere is a Stokeslet, the well-known solution for a Stokeslet and for a rotlet near a wall due to [5] may be used to construct stress tensors which can obtain the translational and rotational motion respectively, to third and fourth order in the swimmer's size; the approximation is improved by also including the stress tensor derived from Blake's solution for the source-dipole near a wall. The results given are identical to those found by matched-asymptotics [16], but by approximating the integration kernel rather than the swimming stroke the slip velocity may be kept completely general [74]. This approach is obviously extensible to other geometries: for instance, [37] has given the solution for a Stokeslet outside a sphere, which would allow the interactions between two swimmers to be determined approximately.

III. BISPHERICAL COORDINATES AND CONJUGATE SOLUTIONS
We record in this section the solutions to Stokes drag problems involving two spheres that we will use with the reciprocal theorem to obtain exact swimmer hydrodynamics. The reader who is primarily interested in these applications may safely skip to §IV and only refer back as necessary.

A. Bispherical coordinates
Problems involving two spheres, such as we consider here, are naturally treated by employing a bispherical coordinate system. If z + iρ is a complex coordinate on a Cartesian grid, the bipolar coordinate grid ξ + iη is defined by 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 = ±R. Finally, a rotation about the z axis gives an azimuthal coordinate φ, which coincides for bispherical and ordinary cylindrical coordinates. Surfaces of constant ξ are non-intersecting spheres centred on −R coth ξ, with radius r = R |cosech ξ|. We consider the fluid to be the region ξ 2 < ξ < ξ 1 , where ξ 1 is taken to be positive. The choice of ξ 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 two spheres, while the intermediate case ξ 2 = 0 represents the half-space, as shown in figure 1(c). In what follows we will use both cylindrical coordinates (z, ρ, φ) and bispherical coordinates (ξ, η, φ), and denote by W the conformal factor R/(cosh ξ − cos η) that appears frequently.

B. Coaxial rotation
The general solution for axisymmetric azimuthal flows was given by Jeffery [41], together with a number of specific examples, including the coaxial rotation of two spheres, and was subsequently expanded upon by Kanwal [45] to consider objects which do not intersect the axis of symmetry. For an axisymmetric, purely azimuthal Stokes flow, the fluid velocity takes the form u = u φ (z, ρ)e φ . 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 (∇ 2 −ρ −2 )u φ = 0, of which the general solution in bispherical coordinates is The constants c l and d l are determined from the boundary conditions of solid-body rotation of the sphere ξ = ξ 1,2 with angular velocity ω 1,2 about their common diameter, The force per area on the spheres' surfaces is purely azimuthal, and integrating it gives the torques These expressions are uniformly convergent since they are bounded independently of ξ 1,2 : taking the first term in (7) as an example, using that sinh nξ > n sinh ξ for n > 0 and ξ > 0 we have that for ξ 2 < 0, ξ 1 > 0. The uniform convergence of the remaining terms is demonstrated similarly, although note that (8) diverges as ξ 2 → ξ 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 ξ 1 and ξ 2 both tend to zero in a way that preserves the ratio r 1 /r 2 , Then, since (7) and (8) are uniformly convergent we may interchange the limit-taking with the summation, finding that the torques converge to the finite values where ζ(s, q) ≡ ∞ n=0 (n + q) −s is the Hurwitz zeta function and ζ(s, 1) ≡ ζ(s) the Riemann zeta function. When one sphere encloses the other (ξ 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 (1) reducing to T 1 (Ω 1 −Ω 2 ). It is natural to take the concave boundary to set the frame of reference. In the intermediate limit of a plane (ξ 2 = 0) the two solutions (11) and (12) coincide and as with the enclosed system the torque on the wall is equal and opposite to the torque on the finite-sized sphere.

C. Translation along the common axis
The translational motion of two spheres along their common diameter was first studied by [89] and subsequently adapted for the special case of a sphere sedimenting towards a planar surface [6,12]. 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 [17,60,[67][68][69]71], 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 (u − 1 2µ px), 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 [89]. 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 ξ = ξ 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 ψ that satisfies the biharmonic equation, ∇ 4 ψ = 0. For an axisymmetric flow, ψ may be chosen to be purely azimuthal, ψ = ψe φ , and then the general solution written in the form The boundary conditions that determine the real constants A l , B l , C l , D l are that the axial flow should equal the constant translation speeds V 1,2 e z of the two spheres ξ = ξ 1,2 , and that the radial flow u ρ should vanish on their surfaces. These conditions may be combined into the statement ∇(ρψ − 1 2 ρ 2 V j ) ξ=ξ j = 0 for j = 1, 2 [89], giving four equations to determine the four unknown coefficients, that reduce to the 2 × 2 block-diagonal form and Inversion is straightforward, and the explicit forms of the coefficients A l , B l , C l , D l are shown in appendix A.
In bispherical coordinates the flow is given by which allows the pressure at a point x to be calculated as where p ∞ 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 [6,89], 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 1 l , rather than in Gegenbauer polynomials C −1/2 n+1 . Furthermore the final scheme we arrive at for determining the coefficients A l , B l , C l , D l presents the 4 × 4 problem in block diagonal form, which we have not seen in the previous literature.
Finally, using (16) and (17) a stress tensor may be constructed and integrated by parts over the spheres to give the hydrodynamic drag force When the fluid has a finite volume or in the limiting case of the half-space (ξ 2 ≥ 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 [5], which is obtained from our solution in the limit Since the coefficients A l , B l , C l , D l have exponential decay in ξ 1 and ξ 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 [12] and we adapt their method here.

A. 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 [1,49], and can also shed light on the hydrodynamics of a daughter colony of Volvox inside its parent [25], or the contact interaction of swimmers with passive particles [98]. Many of these have been studied asymptotically using leading-order point-singularity descriptions [2,87], 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 [39], whose motion is given bỹ where the stress tensor is the corresponding to the Stokes drag problems given in § III. 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 α and its squirming motion is described in terms of a local orthonormal basis {s r , s θ , s φ } and polar coordinate system (θ s , φ s ) relative to this direction, as shown in figure 2. Its slip velocity may be decomposed into squirming modes [54] as where V n (x) ≡ −2P 1 n (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 [4,54]. The addition of the azimuthal modes C n [73] allows for axial rotation of the swimmer, as seen in several real microorganisms, with rotation speed Ω free = −C 1 about the axisymmetry axis.
To perform the integral in (19) it is convenient to express the swimmer's slip velocity in bispherical coordinates. Defining the angle β as the rotation angle about e φ between the cylindrical (e z , e ρ ) and bispherical (e ξ , e η ) basis vectors, the swimmer polar angle is found using the spherical cosine law to be, cos θ s = cos α cos β + sin α sin β cos φ, The slip velocity involves Legendre polynomials in cos θ s , which are expanded using the addition theorem for Legendre functions [61,86], P n (cos α cos β+ sin α sin β cos φ) = P n (cos α)P n (cos β) As the stress tensor in (19) is axisymmetric the integral over φ only affects the slip velocity components. Hence it is convenient to perform the φ-integral first and define a vector of the resulting azimuthally averaged slip velocity components, which is contracted against the stress tensor and integrated over η 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 α 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 α), as found using point-singularity models of swimmer interactions with walls [16,74,87].
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

B. 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 [1,49]. 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 [16,49,57,74], 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 (19) using the slip velocity (25) and the stress corresponding to axisymmetric rotation, where c i , d i are as given in (5) and (6). The factor of V l (cos β) may be written as a polynomial of order l − 1 in W , V l (cos β) = 2 sinh ξ sin η l(l + 1) where the coefficients w n (ξ) are determined using any of the various series representations of Legendre polynomials [96], so that (27) becomes (W/R) 1 2 is the generating function for Legendre polynomials [96]. Successive differentiations of the generating function give the identity which reduces (27) 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 while for interaction with a tracer it is where the torques are given by (7) and (8). To findΩ 1 , ω 2 must be chosen so that T 2 = 0. Conversely, choosing ω 2 such that T 1 = 0 allows the tracer motionΩ 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 §IV C 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 §III B. For a shell this gives and for a tracerΩ It can be readily verified that these coincide for r 1 /r 2 → 0 with the valuẽ 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 θ 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 as depicted schematically in figure 4(a). The slip velocity within the cap region is Ω c , which is balanced by a counter-rotation of the body, Ω 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 and the counter-rotation Ω b required to cancel out the free rotation is When θ 0 = π/2 we have that Ω b = Ω c , as expected on symmetry grounds. E. coli has a body counter-rotation measured to be on the order of one-tenth the rotation of its flagellar bundle, with large variation between specimens [59], and inversion of (38) gives an appropriate value of approximately 0.28π for θ 0 , which we idealise as π/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 α), 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 higherorder 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 α = π/4 the swimmer's rotation changes sense as it approaches the wall. When cos α = 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.

C. Asymptotics at large separation
The tendency of flagellated bacteria such as Escherichia coli and Vibrio alginolyticus to follow distinctive circular trajectories near boundaries [21] has resulted in considerable theoretical work on the rotation of a swimmer normal to a planar surface [20,49,57,74]. 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 [87]), results in rotation by a mechanism analogous to the turning of a tank. The availability of the exact solutions, (31) and (32), for this component of motion allow extension to further geometries than those appearing in the literature.
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 [73]. Focusing on the l = 1 contribution we recognise the sum in (31) and (32) as the torque, (8) and (7) respectively. Hence the contribution to the rotation from this squirming mode is −C 1 cos α for any separation and in any configuration. This is precisely the same as the rotation found asymptotically [16,73], 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, [57] 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 leadingorder order behaviour of a swimmer with arbitrary azimuthal slip velocity near a wall, with Ω 1 = C 2 (r 1 /2d) 4 P 2 (cos α)/5 [74]. Figure 3(a) shows that this behaviour agrees with our exact solution (31) up to a separation of about a squirmer diameter. An explicit form for the rotation near a wall is obtained from eq. (31) by setting ξ 2 = 0,Ω 2 = 0 and ξ 1 = log(d/r 1 + d 2 /r 2 1 − 1), a.
c.  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 θ 0 = π/4 and α = π/4 (black), cos α = 3 −1/2 (blue) and cos α = 0.783 (red). (c) Orientation-dependence of rotation near a wall as a function of distance for θ 0 = π/4. Rotation is normalised by the α-average, Ω α = l (l + 1 2 ) −1 π 0 dα sin α P l (cos α)Ω. (d) Orientation-dependence of normalised rotation near a free surface as a function of distance for θ 0 = π/4. 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.
The asymptotic rotation of a swimmer in the presence of a tracer may be calculated using the leading-order forms ξ 1 ∼ log(r 1 /d) and ξ 2 ∼ log(r 2 /d), giving a decay of This may be understood in terms of multipole reflections [46]. 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 crossover 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 crossover 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.

D. 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 [49] and subsequently observed experimentally [20] 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μ [70], the flow may be found explicitly by supposing an ansatz of the form (4) in each phase and matching flow and stress across the boundary. Then the torque on the sphere is where Λ = (µ−μ)/(µ+μ). When Λ = −1 the empty phase is infinitely viscous and corresponds a noslip wall; instead, when Λ = +1 the boundary is a free surface. Since the torque (7) corresponding to the two-sphere solution when r 1 = r 2 = r is it can be seen that the result for a free surface is recovered when ω 2 = ω 1 [7], while the result for a no-slip wall is given by ω 2 = −ω 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, (31). 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 (32) by substituting ξ 2 = −ξ 1 and ω 2 = ω 1 , givingΩ 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Ω 1 ∝ 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 (35) and for a free surfacẽ 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 [54]).

E. Translation
Phenomena such as the tendency of swimming microorganisms to aggregate at surfaces [2] 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 [82], causes it to align parallel or normal to the wall, respectively [87], while the inclusion of a source dipole ensures selfpropulsion [24]. By adopting the reciprocal theorem the behaviour due to any slip velocity may, in principle, be found [74]; 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 a.
b. c. translational motion is analogous to the calculation for rotation shown in §IV B, 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.
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 [39] 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 lowestorder image singularity will be a stresslet. In the far-field these point-singularity descriptions are sufficient to fully characterise the generic behaviour [87]. For the special case of interaction with a wall an explicit asymptotic estimate [16,74] of the swimming speed is available as 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. 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 crossover to far-field behaviour that is well-described by point-singularity models [16,39,74,87] 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 [37] 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 ∼ 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 [46]. The crossover between the two types of behaviour is shown in figure 6(c).
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 (44). 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 around 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 [13].
The hydrodynamic force in the conjugate problem is divergent in the near-field. This divergence has a leading part going as F i ∼ ξ −2 i , but the representation of the force as an infinite sum contains a harmonically divergent subleading term. By approximating the sum as an integral Cox a. b.
Collision trajectories from exact solutions compared to approximate results of [74]. (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 pointsingularity description is shown as a grey dotted line.
U c , may be calculated analytically and depends on the relative sizes as 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 [74], and where the maximum speed of a self-propulsive swimmer approaches its free-space value in proportion to the excluded area.

V. DISCUSSION
We have found exact expressions for the axisymmetric translation and rotation of a spherical squirmer close to convex, planar or concave no-slip boundary, as well as the axisymmetric rotation beneath a free surface, by making use of the Lorentz reciprocal theorem and the known Stokes drag solutions in these geometries. This covers the hydrodynamics at all separations, including at contact in the case of rotation, and for arbitrary squirming motion. The near-field regime where separations are comparable to the swimmer size or smaller is the regime of greatest relevance to many experimental settings and our exact solution provides rigorous, generic insight. In particular, while the radial and meridional squirming modes show the same asymptotic behaviour, in the near field, at separations smaller than a couple of swimmer diameters, they are markedly different, with the former giving a divergent interaction strength and the latter a hard repulsion. Azimuthal squirming results in the circling behaviour near boundaries seen in flagellated bacteria and our results describe this situation in some detail. The experimentally reported reversal of orbit direction at a free surface is found to be a generic effect. At large separations, the exact solution reproduces results found previously from asymptotic calculations using point singularity approximations of swimmers and also generalises these to interactions of squirmers with spherical boundaries and tracer particles.
Our solution is founded upon the reciprocal theorem for swimmer problems [90] and appears to be the first application of this method to deduce exact solutions that are not currently available by any other method. Given the widespread significance of hydrodynamic interactions, with confining surfaces and with other organisms, to swimmer motion, there are obvious merits to developing applications of this technique in other settings. For instance, we have only been able to provide a partial solution to the interaction of two swimmers, as the non-axisymmetric components of the motion have not been determined. This is because the solution is founded upon the reciprocal theorem and requires the corresponding Stokes drag problem to be solved. For the non-axisymmetric Stokes drag of two spheres, there is, at present, no exact closed-form solution, although there is a scheme in terms of a set of difference equations that could be solved numerically to any desired degree of accuracy. Such an approach would allow the full hydrodynamic interaction of an arbitrary pair of squirmers to be computed, although not in closed form. Furthermore the large range of validity of the approximate far-field solutions here compared with our exact results indicates that asymptotic estimates are valuable and there is merit to pursuing an approximate approach to find the non-axisymmetric behaviour. This may be done, for instance, by constructing an approximate stress tensor using the solution for a Stokeslet outside a sphere [37].

ACKNOWLEDGMENTS
We are grateful to Tom Machon, Marco Polin and especially George Rowlands for fruitful discussions. This work was partially supported by the UK EPSRC through Grant No. A.MACX.0002. The real coefficients A l , B l , C l , D l appearing in (13) are found by inverting the solution schemes (14)- (15).
Using these expressions the explicit form of the force on two spheres may be found using (18).