1. Introduction
The problem of determining the fluid flow around slender bodies with approximately circular cross-sections in a zero-Reynolds-number environment has been central to many biological and physical contexts, such as flagellate motility (Gillies et al. Reference Gillies, Cannon, Green and Pacey2009; Smith etal. Reference Smith, Gaffney, Blake and Kirkman-Brown2009) and soft deformable sensors (Roper et al. Reference Roper, Dreyfus, Baudry, Fermigier, Bibette and Stone2006; Guglielmini et al. Reference Guglielmini, Kushwaha, Shaqfeh and Stone2012). Since the advent of their study, there have been significant advances in both the asymptotic and the numerical methods used to analyse such bodies, many of which exploit a small slenderness ratio 
$\epsilon$. This quantity is defined to be the ratio of the maximum cross-sectional radius of the body and the half-length of the centreline, invariably taken to be much less than unity. This assumption holds for many biological and physical applications, for example the slender flagella of mammalian spermatozoa (Cummins & Woodall Reference Cummins and Woodall1985), with slender-body theories thus having wide-ranging potential applications in the physical and biological sciences. An early theoretical study by Keller & Rubinow (Reference Keller and Rubinow1976) considered a locally cylindrical object with a curvilinear centreline, calculating solutions near to and far from the surface of the body. This led to an integral equation relating the flow field to the force per unit length on the body, enabling iterative calculations of the force per unit length on a moving slender body and, ultimately, the accompanying flow field. Analogous calculations in less general settings have also been performed by Cox (Reference Cox1970) and Lighthill (Reference Lighthill1976), each representing marked improvements over the resistive coefficients of Hancock (Reference Hancock1953) and Gray & Hancock (Reference Gray and Hancock1955). Each of these early methods, however, had sizeable inaccuracies at the body endpoints, where surface velocities are predicted to be infinite due to the use of purely cylindrical geometries. The theory of Johnson (Reference Johnson1980) addressed such issues, offering an improved slender-body theory with algebraic accuracy and reduced end errors.
In this work by Johnson (Reference Johnson1980), and indeed in what will follow in this study, an ansatz for the fluid velocity field is motivated by the solution of Chwang & Wu (Reference Chwang and Wu1975) for the Stokes flow around prolate spheroids, with the flow field represented by an integral of Stokeslet and potential dipole singularities. These integrals range over the centreline between two effective foci, rather than the entire body length, which, as long as the ends of the body are prolate spheroidal caps, circumvents pathological end errors. Johnson's analysis splits the slender body into a central section and two end regions, with the weighting of the potential dipole changing between them. In this theory, arclength-dependent cross-sectional radius functions, denoted here by 
$\eta (s)$ for arclength 
$s$, can be prescribed along the body, but must satisfy the condition that their derivative 
$\mathrm {d}\eta /\mathrm {d}s$ vanishes at, and only at, the centre of the slender body. This restriction of Johnson's theory prohibits the modelling, for example, of many biologically realistic shapes, including tapered flagella, asymmetric cells, and wave-like or ribbed structures. Indeed, such shapes remain absent even from the scope of recent works that have expanded slender-body theories to non-circular cross-sections (Koens & Lauga Reference Koens and Lauga2016; Borker & Koch Reference Borker and Koch2019). Overcoming this significant limitation on the scope of slender theories represents the primary objective of this study, aiming to broaden the range of biological and physical systems that may be justifiably modelled with slender-body theory.
The simulation of more-complex slender geometries is already possible using methods of greater computational complexity, for example the highly accurate but computationally expensive boundary element method (Pozrikidis Reference Pozrikidis2002). A commonplace approach, representing a balance between efficiency, accuracy and flexibility, is the method of regularised Stokeslets, which has seen widespread use in the study of both slender and non-slender bodies (Gillies et al. Reference Gillies, Cannon, Green and Pacey2009; Cortez & Nicholas Reference Cortez and Nicholas2012; Olson et al. Reference Olson, Lim and Cortez2013; Ishimoto & Gaffney Reference Ishimoto and Gaffney2018), though there remains a notable weakness associated with this regularised methodology. Whilst the freedom in the choice of mollifier used has been explored extensively by Zhao et al. (Reference Zhao, Lauga and Koens2019), absent from the literature is a rigorous justification for any particular choice of regularisation parameter, with studies typically opting to use a constant characteristic length scale. In this work, we will aim to address this gap in current understanding, identifying choices of regularisation parameter that will enforce a commonly applied boundary condition on the surface of a slender body to leading order, subject to a flow-field ansatz.
Hence, in this study, we will derive and present a slender-body theory that is applicable to a wider range of geometries than have previously been permissible in such theories, exploiting regularised singularities and giving an appropriate choice of regularisation parameter that encodes the body shape. Representing a stark contrast to previous uses of regularised singularities, here we will aim to satisfy a leading-order no-slip boundary condition on the body surface, with the regularisation parameter explicitly accounting for variations in cross-sectional radius with arclength. For a special but physically pertinent class of slender bodies that are rigid with a straight centreline, we additionally derive a surprising analytic solution, which in this case circumvents the need for numerical solution methods and enables the direct evaluation of regularised Stokeslet weightings. We will evidence the validity of our slender-body theory by explicit numerical evaluation of the surface velocity for a range of non-uniform slender bodies. In doing so we will exploit the smooth integrand in order to make use of non-specialised quadratures, representing a significant improvement in computational complexity over existing slender-body theories that make use of singular integral kernels.
2. The slender-body problem
2.1. Describing the slender body
Consider a three-dimensional slender body with centreline 
$\boldsymbol {\xi }(s)$, where 
$s \in [-L,L]$ is an arclength parameter along the centreline and 
$L$ is the half-length of the object. We assume local axisymmetry about the centreline, so that the cross-section at any point is circular and in a plane transverse to the centreline at arclength 
$s$. Hence the shape is wholly defined by the centreline and the non-negative cross-sectional radius 
$\eta (s)$, with this latter function having zeros at 
$s={\pm }L$. The assumption of slenderness is captured via the aspect ratio 
$\epsilon$, defined here as
\begin{equation} \epsilon = \frac{\max_{s\in[-L,L]}\{\eta(s)\}}{L}, \end{equation}
and taken to be significantly less than unity, i.e. 
$\epsilon \ll 1$.
In an inertial reference frame, which will remain consistent throughout this work, the instantaneous fluid velocity field at the point 
$\boldsymbol {x}$ is denoted 
$\boldsymbol {u}( \boldsymbol {x})$ and governed by the familiar dimensional Stokes equations
\begin{equation} \boldsymbol{\nabla} \boldsymbol{\cdot} \boldsymbol{u} = 0, \quad \boldsymbol{0} = - \boldsymbol{\nabla} p + \mu {\nabla}^{2} \boldsymbol{u}, \end{equation}
where 
$p( \boldsymbol {x})$ is the accompanying pressure field and 
$\mu$ is the dynamic viscosity of the fluid, here assumed constant. We will impose the no-slip boundary condition on the surface of the slender body, though we will first introduce further notation in order to write this succinctly.
Denoting unit tangent, normal and binormal vectors to the curved centreline of the slender body as 
$\boldsymbol {e}_t(s)$, 
$\boldsymbol {e}_n(s)$ and 
$\boldsymbol {e}_b(s)$, respectively, the orthonormal triad 
$\{ \boldsymbol {e}_t, \boldsymbol {e}_n, \boldsymbol {e}_b\}$ satisfies the standard Frenet–Serret relations
\begin{equation} \boldsymbol{e}_t(s) = \frac{\partial \boldsymbol{\xi}}{\partial s}, \quad \frac{\partial \boldsymbol{e}_t}{\partial s} = \kappa(s) \boldsymbol{e}_n(s), \quad \boldsymbol{e}_b(s) = \boldsymbol{e}_t(s) \mathbin{\times} \boldsymbol{e}_n(s), \end{equation}
where 
$\kappa (s)$ is the centreline curvature. We additionally define a radial unit vector 
$\boldsymbol {e}_r$, embedded in a transverse cross-section to the centreline, in terms of a cross-sectional angle 
$\phi$ as
\begin{equation} \boldsymbol{e}_r(s,\phi) = \boldsymbol{e}_n(s) \cos \phi + \boldsymbol{e}_b(s) \sin \phi. \end{equation}Throughout, we will assume that the slender body is inextensible and unshearable, so that material cross-sections are always perpendicular to the centreline tangent (Antman Reference Antman2005). Under this assumption, as exemplified in figure 1, we parameterise the surface of the body explicitly as
\begin{equation} \boldsymbol{X}(s,\phi) = \boldsymbol{\xi}(s) + \eta(s) \boldsymbol{e}_r(s,\phi), \end{equation}
and write 
$\boldsymbol {U}(s,\phi )$ for the surface velocity in the inertial frame at this point 
$\boldsymbol {X}(s,\phi )$. This velocity may be decomposed as
\begin{equation} \boldsymbol{U}(s,\phi) = \boldsymbol{V}(s) + \boldsymbol{\Omega}(s) \mathbin{\times} \eta(s) \boldsymbol{e}_r(s,\phi), \end{equation}
where 
$\boldsymbol {V}(s)$ and 
$\boldsymbol {\Omega }(s)$ are the translational and angular velocity components, respectively, with the latter measured about the centreline point 
$\boldsymbol {\xi }(s)$. Here 
$\boldsymbol {V}(s)$ represents the velocity of the centreline, which may vary with arclength due to non-uniform background flows or centreline motion such as that in a beating flagellum. Similarly, the angular velocity 
$\boldsymbol {\Omega }(s)$ also varies with arclength in general, though we will later see that these effects, and indeed the angular velocity in its entirety, will be subleading with respect to the body slenderness parameter.
Describing a slender body with circular cross-section of arclength-dependent radius 
$\eta (s)$, centreline 
$\boldsymbol {\xi }(s)$ and cross-sectional angle 
$\phi$. Points on the surface of the body are parameterised as 
$\boldsymbol {X}(s,\phi )$, where 
$s\in [-L,L]$ is an arclength parameter.
We may now succinctly formulate the no-slip boundary condition on the surface of the body. We impose that the fluid velocity on the surface of the slender body must exactly equal the surface velocity, a condition that may be stated simply as
\begin{equation} \boldsymbol{u}( \boldsymbol{X}(s,\phi)) = \boldsymbol{U}(s, \phi). \end{equation}2.2. Non-dimensionalisation
We non-dimensionalise lengths via the half-length, 
$L$, writing 
$s=L\hat {s}$ for dimensionless arclength variable 
$\hat {s}\in [-1,1]$ and explicitly decomposing the cross-sectional radius as
\begin{equation} \eta(s) = L \hat{\eta}(\hat{s}) = \epsilon L \tilde{\eta}(\hat{s}) \end{equation}
for dimensionless 
$\hat {\eta },\tilde {\eta }$. Recalling that the cross-sectional radius is at most 
$\epsilon L$ over the slender body, where we assume 
$\epsilon \ll 1$, we have 
$\tilde {\eta }\in [0,1]$ by construction with the upper bound tight.
All velocities are non-dimensionalised via a typical velocity scale, 
$U$, with the expression for the surface velocity (2.6) given in non-dimensional form as
\begin{equation} \skew2\hat{\boldsymbol{U}}(\hat{s},\phi) = \skew2\hat{ \boldsymbol{V}}(\hat{s}) + \epsilon \hat{ \boldsymbol{\Omega}}(\hat{s}) \mathbin{\times} \tilde{\eta}(\hat{s}) \boldsymbol{e}_r(\hat{s},\phi). \end{equation}
Throughout we adopt the standard viscous flow pressure field scaling 
$\mu U/L$, and henceforth will exclusively consider dimensionless quantities and thus drop the somewhat cumbersome 
$\hat {\cdot }$ and 
$\tilde {\cdot }$ notation for dimensionless variables. For completeness, we note that the dimensionless form of (2.5) is now
\begin{equation} \boldsymbol{X}(s,\phi) = \boldsymbol{\xi}(s) + \epsilon \eta(s) \boldsymbol{e}_r(s,\phi), \end{equation}from which the slenderness of the body is evident.
2.3. Leading-order problem
Assuming that the translational and angular velocities of the slender body are each 
$O(1)$ as 
$\epsilon \to 0$, the dimensionless surface velocity given in (2.9) reduces to
\begin{equation} \boldsymbol{U}(s,\phi) = \boldsymbol{V}(s) + O(\epsilon) \end{equation}so that the leading-order boundary condition on the body may be written as
\begin{equation} \boldsymbol{u}( \boldsymbol{X}(s,\phi)) = \boldsymbol{V}(s) + O(\epsilon) \end{equation}
for 
$\phi \in [0,2{\rm \pi} )$ and 
$s\in [-1,1]$. Thus, at leading order, the surface velocity of the slender body is equal to the velocity of the centreline for each 
$s$, and henceforth we will consider only this leading-order problem. We do, however, note that, by linearity of the Stokes equations, a solution of this reduced problem corresponds to a leading-order solution of the full motion of the slender body, and is exact for bodies undergoing pure translational motion.
3. A regularised theory
3.1. Johnson's approach
We now outline the methodology of Johnson (Reference Johnson1980), and will later seek to adapt the approach to admit a wider range of cross-sectional radius functions. Johnson begins by posing an ansatz for the leading-order velocity field based on the exact singularity representation for the flow field around a translating prolate spheroid given by Chwang & Wu (Reference Chwang and Wu1975). Fundamental to this expression are well known singularities of Stokes flow, the Stokeslet, here written 
$\boldsymbol{\mathsf{S}}^{0}$, and the potential dipole, denoted 
$\boldsymbol{\mathsf{D}}^{0}$, defined as
\begin{equation} \boldsymbol{\mathsf{S}}^{0}( \boldsymbol{x}, \boldsymbol{y}) = \frac{{{\boldsymbol{\mathsf{l}}}}}{\left\lvert { \boldsymbol{x}- \boldsymbol{y}}\right\rvert} + \frac{{{\boldsymbol{\mathsf{Q}}}}( \boldsymbol{x}, \boldsymbol{y})}{\left\lvert { \boldsymbol{x}- \boldsymbol{y}}\right\rvert^{3}}\,, \quad \boldsymbol{\mathsf{D}}^{0}( \boldsymbol{x}, \boldsymbol{y}) = -\frac{{{\boldsymbol{\mathsf{l}}}}}{\left\lvert { \boldsymbol{x}- \boldsymbol{y}}\right\rvert^{3}} + \frac{3{{\boldsymbol{\mathsf{Q}}}}( \boldsymbol{x}, \boldsymbol{y})}{\left\lvert { \boldsymbol{x}- \boldsymbol{y}}\right\rvert^{5}}, \end{equation}
respectively, where 
${{\boldsymbol{\mathsf{Q}}}}( \boldsymbol {x}, \boldsymbol {y}) = ( \boldsymbol {x}- \boldsymbol {y}) \otimes ( \boldsymbol {x}- \boldsymbol {y})$ for points 
$\boldsymbol {x}$ and 
$\boldsymbol {y}$ and 
${{\boldsymbol{\mathsf{l}}}}$ is the identity tensor. Johnson's analysis is limited initially to slender bodies with cross-sectional radius functions that approach prolate spheroids towards the ends of the body, i.e.
\begin{equation} \eta^{2}(s) \to (1-s^{2})(1+O(\epsilon^{2})) \quad \text{ for } s \to \pm 1, \end{equation}
though a further, more-significant restriction on 
$\eta$ later arises. By splitting the centreline into central and end regions, Johnson derives the uniformly valid, algebraically accurate flow field as
\begin{equation} \boldsymbol{u}( \boldsymbol{x}) \sim \int_{-e}^{e} \boldsymbol{\mathsf{S}}^{0}( \boldsymbol{x}, \boldsymbol{\xi}(s')) \boldsymbol{\alpha}(s') + W(s') \boldsymbol{\mathsf{D}}^{0}( \boldsymbol{x}, \boldsymbol{\xi}(s')) \boldsymbol{\beta}(s')~\mathrm{d}{s'}, \end{equation}
where 
$W$ is some scalar weight function and 
$\boldsymbol {\alpha }, \boldsymbol {\beta }$ are the unknown strengths of the distributions of singular Stokeslets and potential dipoles, respectively, each of which may vary along the centreline. The limits of integration in (3.3) are given by the quantity 
$e = \sqrt {1-\epsilon ^{2}}$, corresponding to the eccentricity of a prolate ellipsoid with minor axes 
$\epsilon$ and hence termed the effective eccentricity of the slender body. Significantly, these integration limits do not include the singular endpoints 
$s'={\pm }1$, ensuring that the velocity flow-field ansatz (3.3) is regular.
Upon enforcing the leading-order boundary condition from (2.12), Johnson notes that the resulting integral must be independent of the cross-sectional angle 
$\phi$ in order to give a well-defined slender-body theory, with no-slip giving 
$\boldsymbol {\alpha }(s)$ to leading order. This condition, after some manipulation, necessitates the relation
\begin{equation} \boldsymbol{\beta}(s) = \frac{(1-e^{2})}{2}\frac{s}{1-s^{2}} \frac{\eta(s)}{\eta'(s)} \, \boldsymbol{\alpha}(s), \end{equation}
here differing in sign to the original expression of Johnson (Reference Johnson1980) due to contrasting definitions of the potential dipole. This simple relation allows the resulting integral equation from the boundary condition to be cast in terms of a single unknown quantity, the Stokeslet strength, 
$\boldsymbol {\alpha }(s)$, which may then be solved for and the resulting flow field determined everywhere outside and on the slender body.
However, there is an additional inherent constraint on the prescribed cross-sectional radius function, in that the expression given in (3.4) must not attain singular values for 
$s\in [-e,e]$. As 
$\eta (s) > 0$ for 
$s \in (-1,1)$, noting in particular that this inequality is strict, 
$\eta '(s)$ may have a zero only at 
$s = 0$, necessarily accompanied by a suitable condition on the limit of 
$s/\eta '(s)$ for 
$s \in [-e,e]\setminus \{0\}$. This significantly restricts the range of permissible shapes allowed in the slender-body theory of Johnson (Reference Johnson1980), excluding for example the tapered flagella of mammalian sperm (Fawcett Reference Fawcett1970).
3.2. The prolate ellipsoid in uniform flow
Before moving to generalise Johnson's slender body theory, we state the analytic result of Chwang & Wu (Reference Chwang and Wu1975) that pertains to uniform flow around a straight translating prolate ellipsoid, which will be instructive in later analysis. With 
$\boldsymbol {V}_P$ denoting the constant background flow and 
$\boldsymbol {\xi }_P(s)$ denoting the straight centreline, the classical work of Chwang & Wu (Reference Chwang and Wu1975) shows that, for any 
$s\in [-1,1]$ and 
$\phi \in [0,2{\rm \pi} )$,
\begin{equation} \boldsymbol{V}_{P} = \int_{-e}^{e} \left[\boldsymbol{\mathsf{S}}^{0}( \boldsymbol{X}_{P}(s,\phi), \boldsymbol{\xi}_{P}(s')) - \frac{1-e^{2}}{2e^{2}}(e^{2}-s'^{2}) \boldsymbol{\mathsf{D}}^{0}( \boldsymbol{X}_P(s,\phi), \boldsymbol{\xi}_P(s'))\right] \boldsymbol{\alpha}_P~\mathrm{d}{s'} \end{equation}has exact solution given by
\begin{equation} \boldsymbol{\alpha}_P = \left[\dfrac{ \boldsymbol{e}_t\otimes \boldsymbol{e}_t}{-2e+(1+e^{2})\log{\dfrac{1+e}{1-e}}} + \dfrac{2({{\boldsymbol{\mathsf{l}}}} - \boldsymbol{e}_t\otimes \boldsymbol{e}_t)}{2e+(3e^{2}-1)\log{\dfrac{1+e}{1-e}}}\right] \boldsymbol{V}_{P}, \end{equation}
written with respect to the constant orthonormal triad 
$\{ \boldsymbol {e}_t, \boldsymbol {e}_n, \boldsymbol {e}_b\}$ and recalling 
$e=\sqrt {1-\epsilon ^{2}}$. Here, 
$\boldsymbol {X}_P(s,\phi ) = \boldsymbol {\xi }_P(s) + \epsilon \eta _P(s) \boldsymbol {e}_r(s,\phi )$ for 
$\eta _P(s)=\sqrt {1-s^{2}}$, corresponding to a prolate ellipsoid. Notably, as the uniform flow 
$\boldsymbol {V}_{P}$ is independent of the cross-sectional angle 
$\phi$, it is necessarily the case that the integral of (3.5) is also independent of 
$\phi$.
3.3. A regularised theory
Motivated by the exact solution for a prolate ellipsoid given above, we pose a new velocity-field ansatz based on regularised Stokeslets and potential dipoles, defined as
\begin{equation} \boldsymbol{u}( \boldsymbol{x}) = \int_{-e}^{e} \left[\boldsymbol{\mathsf{S}}^{\chi(s')}( \boldsymbol{x}, \boldsymbol{\xi}(s')) - \frac{1-e^{2}}{2e^{2}}(e^{2}-s'^{2}) \boldsymbol{\mathsf{D}}^{\chi(s')}( \boldsymbol{x}, \boldsymbol{\xi}(s'))\right] \boldsymbol{\alpha}(s')~\mathrm{d}{s'}, \end{equation}
where 
$\boldsymbol{\mathsf{S}}^{\chi (s')}$ and 
$\boldsymbol{\mathsf{D}}^{\chi (s')}$ are the regularised singularities of Cortez et al. (Reference Cortez, Fauci and Medovikov2005) and Ainley et al. (Reference Ainley, Durkin, Embid, Boindala and Cortez2008), respectively, given explicitly by
\begin{equation} \boldsymbol{\mathsf{S}}^{\chi(s')}( \boldsymbol{x}, \boldsymbol{y}) = \frac{(\left\lvert { \boldsymbol{x}- \boldsymbol{y}}\right\rvert^{2} + 2\chi(s')){{\boldsymbol{\mathsf{l}}}}}{(\left\lvert { \boldsymbol{x}- \boldsymbol{y}}\right\rvert^{2} + \chi(s'))^{3/2}} + \frac{{{\boldsymbol{\mathsf{Q}}}}( \boldsymbol{x}, \boldsymbol{y})}{(\left\lvert { \boldsymbol{x}- \boldsymbol{y}}\right\rvert^{2}+\chi(s'))^{3/2}}, \end{equation}
\begin{equation} \boldsymbol{\mathsf{D}}^{\chi(s')}( \boldsymbol{x}, \boldsymbol{y}) = -\frac{(\left\lvert { \boldsymbol{x}- \boldsymbol{y}}\right\rvert^{2} - 2\chi(s')){{\boldsymbol{\mathsf{l}}}}}{(\left\lvert { \boldsymbol{x}- \boldsymbol{y}}\right\rvert^{2} + \chi(s'))^{5/2}} + \frac{3{{\boldsymbol{\mathsf{Q}}}}( \boldsymbol{x}, \boldsymbol{y})}{(\left\lvert { \boldsymbol{x}- \boldsymbol{y}}\right\rvert^{2}+\chi(s'))^{5/2}}. \end{equation}
Defining 
$r^{2}( \boldsymbol {x}, \boldsymbol {y},s')=\left \lvert { \boldsymbol {x}- \boldsymbol {y}}\right \rvert ^{2} + \chi (s')$ for convenience, we write these regularised singularities as
\begin{equation} \boldsymbol{\mathsf{S}}^{\chi(s')}( \boldsymbol{x}, \boldsymbol{y}) = \frac{{{\boldsymbol{\mathsf{l}}}}}{r} + \frac{\chi(s'){{\boldsymbol{\mathsf{l}}}}}{r^{3}} + \frac{{{\boldsymbol{\mathsf{Q}}}}( \boldsymbol{x}, \boldsymbol{y})}{r^{3}}, \end{equation}
\begin{equation} \boldsymbol{\mathsf{D}}^{\chi(s')}( \boldsymbol{x}, \boldsymbol{y}) = -\frac{{{\boldsymbol{\mathsf{l}}}}}{r^{3}} + \frac{3\chi(s'){{\boldsymbol{\mathsf{l}}}}}{r^{5}} + \frac{3{{\boldsymbol{\mathsf{Q}}}}( \boldsymbol{x}, \boldsymbol{y})}{r^{5}}, \end{equation}
where we have suppressed the arguments of 
$r$ for brevity. The regularisation parameter 
$\chi (s')$ is a posteriori defined as
\begin{equation} \chi(s) = \epsilon^{2} [(1-s^{2}) - \eta^{2}(s)], \end{equation}
so that the degree of regularisation varies with the cross-sectional radius. It differs significantly from the regularisation parameters of previous studies (Cortez et al. Reference Cortez, Fauci and Medovikov2005; Gillies et al. Reference Gillies, Cannon, Green and Pacey2009; Cortez & Nicholas Reference Cortez and Nicholas2012; Olson et al. Reference Olson, Lim and Cortez2013; Ishimoto & Gaffney Reference Ishimoto and Gaffney2018) in that it depends on the arclength parameter 
$s$, and further may take negative values, though for consistency with previous works we retain its description as a regularisation parameter. Owing to the limits of integration in (3.7), the integral remains regular even if 
$\chi (s)=0$ for some 
$s\in [-e,e]$, whilst cases where 
$\chi (s)<0$ are more subtle and will be implicitly addressed in the subsequent analysis.
Herein fixing 
$s\in [-1,1]$, enforcing the leading-order boundary condition (2.12) at the surface of the slender body gives
\begin{align} \boldsymbol{V}(s) &= \int_{-e}^{e} \left[\boldsymbol{\mathsf{S}}^{\chi(s')}( \boldsymbol{X}(s,\phi), \boldsymbol{\xi}(s')) - \frac{1-e^{2}}{2e^{2}}(e^{2}-s'^{2}) \boldsymbol{\mathsf{D}}^{\chi(s')}( \boldsymbol{X}(s,\phi), \boldsymbol{\xi}(s'))\right] \boldsymbol{\alpha}(s')~\mathrm{d}{s'}\nonumber\\ &\quad + O(\epsilon), \end{align}
required to hold 
$\forall \phi \in [0,2{\rm \pi} )$. In the subsequent analysis we will denote this integrand by 
$\boldsymbol {f}{}(s,s',\phi )$, and in particular note its similarity to that of (3.5). In general, this integral equation must be solved numerically for the regularised Stokeslet density 
$\boldsymbol {\alpha }(s)$, notably with the necessary computation rendered trivial by the regular integrand. However, in order to show that it may indeed admit a leading-order solution, we will seek to show that the integral of (3.13) is independent of 
$\phi$ with errors algebraic in 
$\epsilon$, a necessary condition for a leading-order solution for 
$\boldsymbol {\alpha }(s)$ to be well defined.
Proposition The integral of (3.13)
$,$ written succinctly as 
$\int _{-e}^{e}\,\boldsymbol {f}{}(s,s',\phi )\,\mathrm{d}{s'},$ is independent of cross-sectional angle 
$\phi$ up to errors of 
$O(\epsilon ),$ subject to minimal smoothness conditions on the centreline and radius functions.
Proof. Following the approach of Keller & Rubinow (Reference Keller and Rubinow1976) and Johnson (Reference Johnson1980), we seek a uniformly valid expansion of 
$\boldsymbol {f}{}(s,s',\phi )$, retaining only leading-order terms as 
$\epsilon \to 0$ and yielding overall algebraic accuracy, and will similarly aim to show that its integral is in fact independent of 
$\phi$. The near-singular nature of the regularised Stokeslet and potential dipole suggest two regimes: an inner region where 
$s'-s=O(\epsilon )$, where the local regularisation dominates, and an outer region where 
$s'-s$ is 
$\operatorname {ord}(1)$, that is both 
$O(1)$ and not 
$o(1)$, and the cross-sectional radius function is subleading. When later seeking to match inner and outer expansions we should a priori consider the cases 
$s'<s$ and 
$s'>s$ separately. However, symmetries of the integrand render these solutions identical, hence, in what follows, we do not distinguish between these two cases.
3.3.1. Outer expansion
In the outer region, with 
$s'-s=\operatorname {ord}\!{(1)}$, we expand 
$r( \boldsymbol {X}(s,\phi ), \boldsymbol {\xi }(s'),s')$ to leading order in 
$\epsilon$ to give
\begin{equation} r( \boldsymbol{X}(s,\phi), \boldsymbol{\xi}(s'),s') = \underbrace{\left\lvert { \boldsymbol{\xi}(s) - \boldsymbol{\xi}(s')}\right\rvert}_{r_O(s,s',\phi)} + O(\epsilon), \end{equation}
denoting the leading-order outer expansion by 
$r_O$ and recalling that 
$\chi$ is 
$O(\epsilon ^{2})$. In doing so we are assuming that the centreline does not come close to self-intersection, so that 
$r_O$ as written is 
$\operatorname {ord}(1)$. We may now succinctly expand the integrand in the outer region of the domain of integration, first writing
\begin{equation} {{\boldsymbol{\mathsf{Q}}}}( \boldsymbol{X}(s,\phi), \boldsymbol{\xi}(s')) = \underbrace{[ \boldsymbol{\xi}(s)- \boldsymbol{\xi}(s')]\otimes[ \boldsymbol{\xi}(s)- \boldsymbol{\xi}(s')]}_{{{\boldsymbol{\mathsf{Q}}}}_O(s,s',\phi)} + O(\epsilon), \end{equation}giving the regularised Stokeslet and potential dipole as
\begin{align} \boldsymbol{\mathsf{S}}^{\chi(s')}(s,s',\phi) & =\underbrace{\frac{{{\boldsymbol{\mathsf{l}}}}}{r_O} + \frac{{{\boldsymbol{\mathsf{Q}}}}_O}{r_O^{3}}}_{\boldsymbol{\mathsf{S}}_O^{\chi(s')}} + O(\epsilon), \end{align}
\begin{align} \boldsymbol{\mathsf{D}}^{\chi(s')}(s,s',\phi) & =\underbrace{-\frac{{{\boldsymbol{\mathsf{l}}}}}{r_O^{3}} + \frac{3\kern-.5pt{{\boldsymbol{\mathsf{Q}}}}_O}{r_O^{5}}}_{\boldsymbol{\mathsf{D}}_O^{\chi(s')}}+O(\epsilon), \end{align}
respectively, again recalling that 
$\chi (s')=O(\epsilon ^{2})$ and omitting arguments for brevity. In particular, the potential dipole is 
$\operatorname {ord}{(1)}$ and premultiplied by a factor scaling with 
$1-e^{2}=\epsilon ^{2}$, and thus the potential dipole term is subleading in the outer expansion of the integrand. Furthermore, we note that 
${{\boldsymbol{\mathsf{Q}}}}_O$ and 
$r_O$ are independent of 
$\phi$, hence 
$\boldsymbol{\mathsf{S}}_O^{\chi (s')}$ is independent of 
$\phi$. Thus, the full outer expansion of the integrand, correct to errors of 
$O(\epsilon )$ and denoted by 
$\boldsymbol {f}{}_O(s,s',\phi )$, also does not depend on the cross-sectional angle 
$\phi$, with 
$\boldsymbol {f}{}_O$ being given by
\begin{equation} \boldsymbol{f}_O(s,s') = \left[\frac{{{\boldsymbol{\mathsf{l}}}}}{r_O} + \frac{{{\boldsymbol{\mathsf{Q}}}}_O}{r_O^{3}}\right] \boldsymbol{\alpha}(s'), \end{equation}
where 
$r_O$ and 
$Q_O$ each have their arguments of 
$s$ and 
$s'$ omitted for brevity.
3.3.2. Leading-order Prandtl matching term
We will later determine a uniformly valid composite solution at leading order, and thus will require the Prandtl matching term, which is given by the limit of the outer expansion 
$f_O(s,s')$ as 
$s'\to s$. Noting that
\begin{equation} r_O(s,s') \sim \left\lvert {s'-s}\right\rvert, \quad {{\boldsymbol{\mathsf{Q}}}}_O(s,s')\sim (s'-s)^{2} \boldsymbol{e}_t(s)\otimes \boldsymbol{e}_t(s), \end{equation}
as 
$s'\to s$, this limiting behaviour is readily found to be
\begin{equation} \boldsymbol{f}_M(s,s') = \left[\frac{{{\boldsymbol{\mathsf{l}}}}}{\left\lvert {s'-s}\right\rvert} + \frac{ \boldsymbol{e}_t(s)\otimes \boldsymbol{e}_t(s)}{\left\lvert {s'-s}\right\rvert}\right] \boldsymbol{\alpha}(s). \end{equation}
Here we have assumed that the derivative of 
$\boldsymbol {\alpha }(s')$ is 
$O(1)$ as 
$\epsilon \to 0$, which may be validated a posteriori. Inherited from 
$\boldsymbol {f}{}_O$, this matching term is also independent of 
$\phi$.
3.3.3. Inner expansion
In the inner region, with 
$s'-s=O(\epsilon )$, we introduce the inner variable 
$\sigma =(s'-s)/ \epsilon =O(1)$, and expand
\begin{align} r( \boldsymbol{X}(s,\phi), \boldsymbol{\xi}(s'),s') &= \left[\left\lvert { \boldsymbol{\xi}(s)- \boldsymbol{\xi}(s') + \epsilon\eta(s) \boldsymbol{e}_r(s,\phi)}\right\rvert^{2} + \chi(s')\right]^{1/2}\nonumber\\ &=\left[\left\lvert {-\epsilon\sigma \boldsymbol{e}_t(s) + \epsilon\eta(s) \boldsymbol{e}_r(s,\phi)}\right\rvert^{2} +\chi(s)+ O(\epsilon^{3}) \right]^{1/2}\nonumber\\ &=\epsilon\left[\sigma^{2} + \eta^{2}(s) + (1-s^{2}-\eta^{2}(s))\right]^{1/2}+ O(\epsilon^{2})\nonumber\\ &=\underbrace{\epsilon\left[\sigma^{2} + 1 - s^{2}\right ]^{1/2}}_{r_I(s,s',\phi)}+ \ O(\epsilon^{2})\,, \end{align}
having Taylor expanded 
$\boldsymbol {\xi }(s')$ and 
$\chi (s')$ about 
$s$ and noting that 
$\boldsymbol {e}_t(s)\cdot \boldsymbol {e}_r(s,\phi )=0$ by definition. In expanding we have assumed that the centreline curvature 
$\kappa$ is such that 
$\epsilon \kappa =o(1)$, in addition to 
$\mathrm {d}{\chi (s')} / \mathrm {d}{s'}=O(\epsilon ^{2})$, once more noting that 
$\chi (s')=O(\epsilon ^{2})$ by construction. With 
$\chi (s')$ dependent on the radius function via (3.12), this latter condition is satisfied by imposing the weak constraint 
$\eta \,\mathrm {d}{\eta } / \mathrm {d}{s'}=O(1)$, recalling that 
$\eta \in [0,1]$. Notably, we see from (3.21) that 
$r_I(s,s',\phi )$ is in fact independent of 
$\phi$, being purely a function of 
$s$ and 
$s'$.
As in the outer solution, we expand the outer product 
${{\boldsymbol{\mathsf{Q}}}}$ as
\begin{align} {{\boldsymbol{\mathsf{Q}}}}( \boldsymbol{X}(s,\phi), \boldsymbol{\xi}(s')) &= \left[ \boldsymbol{\xi}(s)- \boldsymbol{\xi}(s') + \epsilon\eta(s) \boldsymbol{e}_r(s,\phi)\right] \otimes \left[ \boldsymbol{\xi}(s)- \boldsymbol{\xi}(s') + \epsilon\eta(s) \boldsymbol{e}_r(s,\phi)\right] \nonumber\\ &=\underbrace{\epsilon^{2}\left[-\sigma \boldsymbol{e}_t(s) + \eta(s) \boldsymbol{e}_r(s,\phi)\right] \otimes \left[-\sigma \boldsymbol{e}_t(s) + \eta(s) \boldsymbol{e}_r(s,\phi)\right]}_{{{\boldsymbol{\mathsf{Q}}}}_I(s,s',\phi)} + \ O(\epsilon^{3})\,, \end{align}
noting that the only dependence on the variable of integration 
$s'$ is through the coefficients of 
$\boldsymbol {e}_t(s)$ in the multiplicands, here written in terms of the inner variable 
$\sigma$. The regularised Stokeslet and potential dipole are now given as
\begin{equation} \boldsymbol{\mathsf{S}}^{\chi(s')}(s,s',\phi)=\underbrace{\frac{{{\boldsymbol{\mathsf{l}}}}}{r_I} + \frac{\chi(s){{\boldsymbol{\mathsf{l}}}}}{r_I^{3}} + \frac{{{\boldsymbol{\mathsf{Q}}}}_I}{r_I^{3}}}_{\boldsymbol{\mathsf{S}}_I^{\chi(s')}} + \ O(1), \end{equation}
\begin{equation} \boldsymbol{\mathsf{D}}^{\chi(s')}(s,s',\phi)=\underbrace{-\frac{{{\boldsymbol{\mathsf{l}}}}}{r_I^{3}} + \frac{3\chi(s){{\boldsymbol{\mathsf{l}}}}}{r_I^{5}} + \frac{3\kern-.5pt{{\boldsymbol{\mathsf{Q}}}}_I}{r_I^{5}}}_{\boldsymbol{\mathsf{D}}_I^{\chi(s')}} + \ O(\epsilon^{-2}) \end{equation}
having suppressed the arguments on the right-hand side and again Taylor expanded 
$\chi (s')$ about 
$s$, recalling that these expansions are asymptotic only for 
$s'-s=O(\epsilon )$. The leading-order inner expansion of the integrand may then be written as
\begin{align} \boldsymbol{f}_I(s,s',\phi) = \left[\frac{{{\boldsymbol{\mathsf{l}}}}}{r_I} +\frac{\chi(s){{\boldsymbol{\mathsf{l}}}}}{r_I^{3}} + \frac{{{\boldsymbol{\mathsf{Q}}}}_I}{r_I^{3}} - \frac{1-e^{2}}{2e^{2}}(e^{2}-s'^{2})\left(-\frac{{{\boldsymbol{\mathsf{l}}}}}{r_I^{3}} + \frac{3\chi(s){{\boldsymbol{\mathsf{l}}}}}{r_I^{5}} + \frac{3\kern-.5pt{{\boldsymbol{\mathsf{Q}}}}_I}{r_I^{5}}\right)\right] \boldsymbol{\alpha}(s), \end{align}
expanding 
$\boldsymbol {\alpha }(s')\sim \boldsymbol {\alpha }(s)$ as in (3.20). With component-wise division implicit in the notation used here and below, we have 
$\boldsymbol {f}{}/ \boldsymbol {f}{}_I = 1+O(\epsilon )$ when 
$s'-s=O(\epsilon )$, and have chosen not to expand the coefficient of the regularised potential dipole for later algebraic convenience. Nonetheless, this coefficient is 
$O(\epsilon ^{2})$, so that all terms in 
$\boldsymbol {f}_I$ are of the same order in the inner region.
3.3.4. Uniformly valid expansion
As commonly considered in asymptotic analyses (Bender & Orszag Reference Bender and Orszag1999), we construct a uniformly valid composite expansion
\begin{equation} \boldsymbol{f}_C(s,s',\phi) = \boldsymbol{f}_O(s,s') + \boldsymbol{f}_I(s,s',\phi) - \boldsymbol{f}_M(s,s'), \end{equation}
with 
$\,\boldsymbol {f}{}/ \boldsymbol {f}{}_C = 1 + O(\epsilon )$ for all 
$s'\in [-e,e]$. We now proceed to show that 
$\int _{-e}^{e}\, \boldsymbol {f}_C(s,s',\phi )~\mathrm {d}{s'}$ is in fact independent of 
$\phi$, such that the leading-order boundary condition (3.13) may be satisfied. From (3.26) we see that we need only consider the integral of the inner part of the composite expansion,
\begin{equation} \int_{-e}^{e}\left[\frac{{{\boldsymbol{\mathsf{l}}}}}{r_I} +\frac{\chi(s){{\boldsymbol{\mathsf{l}}}}}{r_I^{3}} + \frac{{{\boldsymbol{\mathsf{Q}}}}_I}{r_I^{3}} - \frac{1-e^{2}}{2e^{2}}(e^{2}-s'^{2})\left(-\frac{{{\boldsymbol{\mathsf{l}}}}}{r_I^{3}} + \frac{3\chi(s){{\boldsymbol{\mathsf{l}}}}}{r_I^{5}} + \frac{3\kern-.5pt{{\boldsymbol{\mathsf{Q}}}}_I}{r_I^{5}}\right)\right] \boldsymbol{\alpha}(s)~\mathrm{d}{s'}. \end{equation}
However, the inner expansion has dependence on 
$\phi$ only through the outer product 
${{\boldsymbol{\mathsf{Q}}}}_I$, so we may consider only the reduced integral
\begin{align} & \int_{-e}^{e}\left[\frac{{{\boldsymbol{\mathsf{Q}}}}_I}{r_I^{3}} - \frac{1-e^{2}}{2e^{2}}(e^{2}-s'^{2})\frac{3\kern-.5pt{{\boldsymbol{\mathsf{Q}}}}_I}{r_I^{5}}\right] \boldsymbol{\alpha}(s)~\mathrm{d}{s'} \nonumber\\ &\quad = \int_{-(e+s)/\epsilon}^{(e-s)/\epsilon}\left[\frac{{{\boldsymbol{\mathsf{Q}}}}_I}{r_I^{3}} - \frac{\epsilon^{2}}{2e^{2}}(e^{2}-(s+\epsilon\sigma)^{2})\frac{3\kern-.5pt{{\boldsymbol{\mathsf{Q}}}}_I}{r_I^{5}}\right] \boldsymbol{\alpha}(s)\epsilon~\mathrm{d}{\sigma}, \end{align}
where we have rewritten the integral with respect to the inner variable 
$\sigma$. Careful consideration of the form of 
${{\boldsymbol{\mathsf{Q}}}}_I$ given in (3.22) reveals that its only dependence on the cross-sectional angle 
$\phi$ is through terms which are homogeneous degree zero or degree one in 
$\sigma$. Omitting factors that are constant with respect to the integration variable, we write the integrals of these terms explicitly as
\begin{equation} J_p = \int_{-(e+s)/\epsilon}^{(e-s)/\epsilon}\left[\frac{\sigma^{p}}{\left(\sigma^{2}+1-s^{2}\right)^{3/2}} - \frac{3}{2e^{2}}(e^{2} - (s+\epsilon\sigma)^{2})\frac{\sigma^{p}}{\left(\sigma^{2} + 1-s^{2}\right)^{5/2}}\right]~\mathrm{d}{\sigma} \end{equation}
for 
$p\in \{0,1\}$.
Whilst these integrals may be readily computed to give 
$J_0=J_1=0$, we can instead exploit the link between our theory and that of Chwang & Wu (Reference Chwang and Wu1975). Again noting the similarity between our leading-order boundary condition (3.13) and the exact solution for a prolate ellipsoid given in (3.5), we see that the analysis performed above can be applied to the prolate ellipsoid. Analogous to the above, conducting such an analysis on the integral of (3.5) highlights that any 
$\phi$ dependence at leading order can be present only in the quantity
\begin{equation} \int_{-(e+s)/\epsilon}^{(e-s)/\epsilon}\left[\frac{{{\boldsymbol{\mathsf{Q}}}}_P}{r_I^{3}} - \frac{\epsilon^{2}}{2e^{2}}(e^{2}-(s+\epsilon\sigma)^{2})\frac{3\kern-.5pt{{\boldsymbol{\mathsf{Q}}}}_P}{r_I^{5}}\right] \boldsymbol{\alpha}_P\epsilon~\mathrm{d}{\sigma}, \end{equation}where
\begin{equation} {{\boldsymbol{\mathsf{Q}}}}_P(s,s',\phi) = \epsilon^{2}\left[-\sigma \boldsymbol{e}_t(s) + \eta_P(s) \boldsymbol{e}_r(s,\phi)\right] \otimes \left[-\sigma \boldsymbol{e}_t(s) + \eta_P(s) \boldsymbol{e}_r(s,\phi)\right] \end{equation}
and 
$r_I$ is unmodified from (3.21) due to our choice of 
$\chi (s)$ collapsing to zero for the prolate spheroid. As above, the 
$\phi$ dependence in (3.30) is modulated by the integrals 
$J_0$ and 
$J_1$. Omitting detailed calculation, this 
$\phi$ dependence may be written as
\begin{equation} J_0 g(\phi) + J_1h(\phi), \end{equation}
where 
$g$ and 
$h$ are linearly independent functions of 
$\phi$. With the solution of Chwang & Wu (Reference Chwang and Wu1975) being necessarily independent of 
$\phi$, we therefore must have that 
$J_0=J_1=0$.
Hence, with 
$J_p$ not dependent on the radius function 
$\eta$ or on the centreline 
$\boldsymbol {\xi }$, we retain 
$J_0=J_1=0$ in the general case. Therefore, the integral of the composite expansion is independent of 
$\phi$ to 
$O(\epsilon )$ errors, and we may thus conclude that 
$\int _{-e}^{e}\ \boldsymbol {f}{}(s,s',\phi )~\mathrm {d}{s'}$ is also independent of 
$\phi$ with 
$O(\epsilon )$ error terms. □
3.4. Analytic solution for straight bodies in constant flow
Whilst in full generality the unknown Stokeslet density in (3.13) must be solved for numerically, for slender bodies with a constant surface velocity 
$\boldsymbol {V}(s)= \boldsymbol {V}$ and straight centreline, so that 
$\boldsymbol {\xi }(s')=s' \boldsymbol {e}_t$ for constant 
$\boldsymbol {e}_t$, a simple analytic solution for the leading-order density may be obtained.
Motivated by the solution for the prolate ellipsoid, we seek a solution with 
$\boldsymbol {\alpha }(s)= \boldsymbol {\alpha }$ constant, independent of arclength. Coupled with the assumption of a straight centreline, so that in particular 
$r_O(s,s')=\left \lvert {s'-s}\right \rvert +O(\epsilon )$ and is thus asymptotic to its limiting behaviour as 
$s'\to s$, given in (3.19), we see that the outer and matching terms of the composite expansion (3.26) cancel at leading order, leaving only
\begin{equation} \boldsymbol{f}_C(s,s',\phi) = \boldsymbol{f}_I(s,s',\phi) \end{equation}
as the uniformly valid composite expansion with 
$\boldsymbol {f}{}/ \boldsymbol {f}{}_C=1+O(\epsilon )$ for all 
$s'\in [-e,e]$. What remains as the surface velocity is simply
\begin{equation} \boldsymbol{V} = \int_{-e}^{e}\ \boldsymbol{f}_I(s,s',\phi)~\mathrm{d}{s'} + O(\epsilon), \end{equation}
with the integral being of the same form as that in (3.27). Motivated by the cancellation of regularised Stokeslet and potential dipole contributions that yielded 
$\phi$ independence in the previous analysis, we consider the terms explicitly involving 
$\chi (s)$ in the numerators of (3.27). With 
$s$ being constant with respect to the variable of integration, we may write these terms as
\begin{equation} \chi(s) \boldsymbol{\alpha}\int_{-e}^{e}\left[\frac{1}{r_I^{3}} - \frac{1-e^{2}}{2e^{2}}(e^{2}-s'^{2})\frac{3}{r_I^{5}}\right]~\mathrm{d}{s'}, \end{equation}
with the integrand being equal to that of 
$J_0$ when written in the inner variable 
$\sigma$. Hence, as 
$J_0=0$, the terms explicitly involving 
$\chi (s)$ contributions to the numerators vanish in the integral, leaving the surface-velocity ansatz as
\begin{equation} \boldsymbol{V} = \int_{-e}^{e}\left[\frac{{{\boldsymbol{\mathsf{l}}}}}{r_I} + \frac{{{\boldsymbol{\mathsf{Q}}}}_I}{r_I^{3}} - \frac{1-e^{2}}{2e^{2}}(e^{2}-s'^{2})\left(-\frac{{{\boldsymbol{\mathsf{l}}}}}{r_I^{3}} + \frac{3\kern-.5pt{{\boldsymbol{\mathsf{Q}}}}_I}{r_I^{5}}\right)\right] \boldsymbol{\alpha}~\mathrm{d}{s'} + O(\epsilon). \end{equation}
A simple direct comparison of the leading-order component of this integral equation with the formulation of Chwang & Wu (Reference Chwang and Wu1975) for a prolate ellipsoid, stated here in (3.5), highlights differences only within the outer product term 
${{\boldsymbol{\mathsf{Q}}}}_I$, noting that 
$r_I=r$ for the ellipsoid. For clarity, denoting terms relating to the prolate ellipsoid analogue as 
${{\boldsymbol{\mathsf{Q}}}}_P$ and 
$\eta _P(s)$ as above, we recapitulate
\begin{align} {{\boldsymbol{\mathsf{Q}}}}_I(s,s',\phi) & = \epsilon^{2}\left[-\sigma \boldsymbol{e}_t(s) + \eta(s) \boldsymbol{e}_r(s,\phi)\right] \otimes \left[-\sigma \boldsymbol{e}_t(s) + \eta(s) \boldsymbol{e}_r(s,\phi)\right], \end{align}
\begin{align} {{\boldsymbol{\mathsf{Q}}}}_P(s,s',\phi) & = \epsilon^{2}\left[-\sigma \boldsymbol{e}_t(s) + \eta_P(s) \boldsymbol{e}_r(s,\phi)\right] \otimes \left[-\sigma \boldsymbol{e}_t(s) + \eta_P(s) \boldsymbol{e}_r(s,\phi)\right]. \end{align}
From these expressions it is clear that the only difference between 
${{\boldsymbol{\mathsf{Q}}}}_I$ and 
${{\boldsymbol{\mathsf{Q}}}}_P$ is through the 
$\eta (s)$ and 
$\eta _P(s)$ terms, where we recall 
$\eta _P(s) = \sqrt {1-s^{2}}$ defines the shape of the prolate ellipsoid. However, employing the same approach as used above to establish 
$\phi$ independence of (3.28), recalling that 
$J_0=J_1=0$, all terms involving 
$\eta (s)$ and 
$\eta _P(s)$ in the numerators vanish upon integration, so that in fact there is no dependence of (3.36) on the shape 
$\eta (s)$.
Thus, having seen that all differing terms in the integrand provide no contribution to the integral at leading order, the relation between 
$\boldsymbol {\alpha }$ and 
$\boldsymbol {V}$ is precisely that between 
$\boldsymbol {\alpha }_P$ and 
$\boldsymbol {V}_P$ for the prolate ellipsoid, up to algebraic errors scaling with 
$\epsilon$ or less. Hence, for slender bodies with straight centrelines and a uniform surface velocity, taking 
$\boldsymbol {V}_P= \boldsymbol {V}$ we conclude that
\begin{equation} \boldsymbol{\alpha}(s) = \boldsymbol{\alpha}_P + O(\epsilon). \end{equation}4. Numerical verification
At first considering slender bodies with straight centrelines and a uniform surface velocity, 
$\boldsymbol {V}$, we begin by verifying that the boundary condition (3.13) is indeed satisfied with algebraic errors by 
$\boldsymbol {\alpha }(s)= \boldsymbol {\alpha }_P$, independent of cross-sectional angle 
$\phi$.
Taking 
$\boldsymbol {V}=[1,1,1]^{T}$ and rescaling so that 
$\left \lvert { \boldsymbol {V}}\right \rvert =1$ throughout, in figure 2 we report 
$\mathcal {E}=\max _{s,\phi }\left \lvert { \boldsymbol {V} - \boldsymbol {V}_{num}(s,\phi )}\right \rvert _{\infty }$ for a variety of slender bodies, where 
$\left \lvert {\,\cdot \,}\right \rvert _{\infty }$ denotes the infinity norm. Here, 
$\boldsymbol {V}_{num}(s,\phi )$ is the numerical value of the integral given in (3.13), and the maximum is approximated by computing 
$\left \lvert { \boldsymbol {V} - \boldsymbol {V}_{num}(s,\phi )}\right \rvert _{\infty }$ at 
$10^{4}$ sample points on the surface of the slender body. Integrals to evaluate surface velocities are computed numerically using the MATLAB® command quadv to a tolerance of 
$10^{-6}$, significantly less than the chosen 
$\epsilon$, which here is 
$\epsilon =10^{-2}$, though what follows holds in more generality. Figure 2 demonstrates the validity of our proposed slender-body theory, with maximal errors being 
$O(\epsilon )$ for a range of radius functions, and, reassuringly, less than numerical tolerance for the case of a prolate ellipsoid (see figure 2a). Remarkably, such low error is retained even when considering an example slender body with a cusp, as highlighted by figure 2(d,e), suggesting a robustness of this regularised slender-body theory to cases outside its analytical scope, though the largest errors in the boundary velocity occur at the cusps of these examples. We additionally validate the asymptotic scaling of the error in the boundary condition by repeating the above error calculations for a range of values of 
$\epsilon$, with the resulting 
$\epsilon$-dependent errors in the surface velocity shown in figure 3(a), verifying an approximate linear dependence on 
$\epsilon$. The precise functional forms of 
$\eta (s)$ corresponding to the shapes shown in figure 2 are given in appendix A, and the presented results are qualitatively unchanged when altering the prescribed uniform background velocity 
$\boldsymbol {V}$, with absolute errors scaling linearly with the magnitude of 
$\boldsymbol {V}$.
Maximum error 
$\mathcal {E}$ in the no-slip boundary condition for various slender bodies with a straight centreline and constant surface velocity, using 
$\boldsymbol {\alpha }(s)= \boldsymbol {\alpha }_P$. The characteristic shapes of the slender bodies, corresponding to the functional forms given in appendix A, are presented above, with the shapes shown stretched vertically for visual clarity. Here, with 
$\epsilon =10^{-2}$, we report the maximum error, as measured in the infinity norm, to a unit surface velocity, denoted 
$\mathcal {E}$, for each shape, showing errors that are 
$O(\epsilon )$. This level of error is retained even when considering bodies with a highly oscillatory radius function (c) or a cusp (d–e), with this latter observation suggesting a robustness of this methodology outside of its analytically derived scope. In panels (b) to (e), sections of body surface where 
$\left \lvert { \boldsymbol {V} - \boldsymbol {V}_{num}}\right \rvert _{\infty } \geq \mathcal {E}/10$ are highlighted with red circles, from which we note that the largest errors are present at the endpoints of the slender body and in regions of high curvature. No such points are highlighted in (a) as the boundary condition is satisfied to working precision everywhere on the prolate ellipsoid.
Error in the no-slip boundary condition as a function of slenderness 
$\epsilon$, and a comparison of flow fields with the boundary element method. (a) Maximum error 
$\mathcal {E}$ as a function of 
$\epsilon$ for the straight slender bodies shown in figure 2, excluding the prolate ellipsoid of figure 2(a) as it exhibits machine-precision errors independent of 
$\epsilon$. We see evidenced an approximate linear dependence of error on epsilon, with individual data points shown as a black dots. (b, c) Pointwise infinity-norm difference in flow fields computed via the proposed slender-body theory and the boundary element method, denoted 
$\left \lvert { \boldsymbol {u}- \boldsymbol {u}_B}\right \rvert _{\infty }$ for boundary element flow field 
$\boldsymbol {u}_B$. Having taken 
$\epsilon =10^{-2}$, we see in (b) a good agreement between the methodologies, serving as validation of the presented slender-body theory. Larger disparity in (c) is resultant of poor accuracy of the boundary element computation, limited by the overwhelming computational resources needed to accurately capture highly slender geometries. Taking straight centrelines, body shapes are as in figures 2(b) and 5, and vertical axes have been rescaled for visual clarity.
Further, we validate the flow fields as given by the proposed slender-body theory via direct comparison to those computed using the boundary element method summarised by Pozrikidis (Reference Pozrikidis2002), as implemented by Walker et al. (Reference Walker, Wheeler, Ishimoto and Gaffney2019). In these boundary element computations we represent the surface of the slender body by a discrete mesh containing 
$8\times 10^{4}$ flat triangles, linearly interpolating unknowns across mesh elements. Taking 
$\epsilon =10^{-2}$, figure 3(b) reports the pointwise infinity norm of the difference between the slender-body theory flow field 
$\boldsymbol {u}$ and that computed via the boundary element method around a sample slender body, with the latter flow being denoted 
$\boldsymbol {u}_B$. We observe agreement between the methodologies at approximately 
$O(\epsilon )$, serving as strong numerical validation of the analytically justified slender-body theory presented above. A further flow field comparison is presented in figure 3(c), though in the case of this particular slender body there is significant pointwise disparity between the flow fields computed by the two methodologies. However, with the significant differences localised around the most narrow region of the slender body, we find that the dominant contributions to the apparent error in fact stem from limitations of the boundary element computation. Indeed, owing to the very high surface curvatures of the slender body, the computational meshes needed to accurately capture the effects of the geometry on the flow around this particularly narrow segment of the slender body are more highly resolved than is computationally feasible in practice using this boundary element implementation, with further refinement of the result prohibited by the overwhelming computational resources required. This highlights a general weakness of gold-standard boundary element computations when applied to slender bodies, with refined meshes required in order to attain reasonable accuracy in the near-field of the most slender regions. Representing a significant advantage of slender theories over other computational methods, this limitation is not shared by slender-body theories such as that presented in this work, with accuracy not reliant on high-resolution computational meshes and body geometry analytically captured in the flow-field ansatz.
Considering slender bodies more generally, we now prescribe non-straight centrelines as well as non-trivial radius functions, which, in contrast to bodies with straight centrelines, necessitate finding a numerical solution of the integral equation (3.13) for the unknown 
$\boldsymbol {\alpha }(s)$. We discretise the unknown density 
$\boldsymbol {\alpha }(s)$ crudely into 100 piecewise-constant elements, resulting in a linear system of low dimension that may be readily inverted, with invertibility established numerically for each considered slender body, including those violating the assumptions of the asymptotic theory. We take 
$\phi =0$, though remark that the following results are qualitatively independent of both this choice and the number of elements used in the discretisation. Evaluating the error 
$\mathcal {E}$ in the boundary condition as defined above, here approximated at 
$10^{3}$ points along the body, we present in figure 4 the error for a number of sample centrelines, fixing the radius function to be that of figure 5 but with observations holding in more generality. We again see limited errors in the fluid velocity on the surface of the body in all examples satisfying the assumptions of our analysis, recovering component-wise differences that are 
$O(\epsilon )$ and validating our simple slender theory. However, candidate centrelines breaking the assumptions of body curvature do not give rise to such small errors, instead giving 
$\operatorname {ord}{(1)}$ variations in the boundary condition as shown in figure 4(d,e), demonstrating that curvature assumptions may not be relaxed. Whilst the examples shown here use plane curves as sample centrelines for visual clarity, the above remarks additionally hold for non-planar centrelines, as well as for a range of values of 
$\epsilon$. For the cases considered, we additionally note a slight increase in the computational requirements of the adaptive quadrature applied to non-straight bodies in comparison to straight slender bodies, though we retain the efficiency associated with slender-body theories.
Maximum error 
$\mathcal {E}$ in the no-slip boundary condition for various slender bodies with non-straight centrelines and constant surface velocity, having solved for 
$\boldsymbol {\alpha }(s)$ numerically. The centrelines of the slender bodies, corresponding to the functional forms given in appendix B, are presented above, each having the radius function corresponding to figure 2(b). Here with 
$\epsilon =10^{-2}$, for each shape we report the maximum error, as measured in the infinity norm, to a unit surface velocity, denoted 
$\mathcal {E}$, showing errors that are 
$O(\epsilon )$ in all cases except panels (d) and (e), which each violate assumptions on centreline curvature, suggesting that such assumptions are necessary for validity of our slender-body theory in practice. These plots have aspect ratio 1 : 1 and are independently scaled for visual clarity.
The characteristic shape generated by the radius function used in the evaluation of the presented slender-body theory on a non-straight body, with centrelines shown in figure 4, with the radius here shown stretched vertically and relative to a straight centreline for visual clarity. The body radius is proportional to 
$\sqrt {1-s^{2}}(1.1+\sin {{\rm \pi} s})$, where 
$s\in [-1,1]$.
5. Discussion
In this work, we have presented, justified and numerically verified a leading-order theory of non-uniform slender bodies in Stokes flow, accurate to errors algebraic in the slenderness parameter. Motivated by the methodology of Johnson (Reference Johnson1980) and the classical solutions of Chwang & Wu (Reference Chwang and Wu1975), we have significantly expanded the range of permissible radius functions considered by Johnson (Reference Johnson1980), enabling the future modelling of a variety of biologically and physically motivated slender bodies in flow, for example wave-like synthetic structures or the tapered flagella of mammalian spermatozoa. Having been afforded additional flexibility by this study and the complimentary work of Borker & Koch (Reference Borker and Koch2019), slender-body theory may also play a wider role in shape-optimisation studies (Keaveny & Shelley Reference Keaveny and Shelley2011), as it is computationally efficient in comparison with other standard methodologies, in particular those reliant on refined computational representations of body geometry, such as the boundary element method (PozrikidisReference Pozrikidis2002).
Previously, practical applications of many slender-body theories have been plagued with numerical complexity, requiring specialised quadratures or other approaches to overcome singularities in the integral kernels (Shelley & Ueda Reference Shelley and Ueda2000; Tornberg & Shelley Reference Tornberg and Shelley2004). Similar to the regularised slender-body theory of Cortez & Nicholas (Reference Cortez and Nicholas2012), our presented theory overcomes such issues by utilising the regularised singularities of Cortez (Reference Cortez2001), and additionally requires the evaluation of fluid velocities only on the surface of the slender body rather than on its centreline, enabling direct numerical evaluation of the integral kernel of (3.7). This regularised slender theory therefore gives rise to a computationally simple methodology that may be readily implemented in non-specialised frameworks to yield stable numerical solutions. Representing a simplification over more-complex mesh-based methodologies, body shape is encoded directly and explicitly in the regularisation parameter. This is of particular pertinence to bodies with straight centrelines, with this work having found that the leading-order solution for the regularised Stokeslet density in this case is given simply by that for a corresponding prolate ellipsoid. Hence, explorations of straight slender bodies with various cross-sectional radius functions do not necessitate repeated solution for the density 
$\boldsymbol {\alpha }(s)$, with only minor modifications to the integral kernel required in order to obtain the flow field via (3.7).
Moreover, having defined our choice of regularisation parameter in terms of the body shape, in the context of the flow-field ansatz the analysis in this work provides a rigorous relation between geometry and regularisation. Such a link has previously been lacking in studies making use of regularised singularities, with ad hoc links between the regularisation parameter and body radius typically being argued, with the former often taken to be spatially independent. Indeed, in the context of the flow-field ansatz posed in this study, such a choice of regularisation parameter does not correspond to a filament with constant radius, instead being characteristic of spheroidal slender bodies only. Thus, providing an explicit relationship between regularisation and geometry, this study represents a step towards understanding the effects of spatially varying regularisation parameters on the flow solutions obtained using regularised singularities.
Numerical explorations and validations of our approach suggested that, at least in the cases considered, some assumed smoothness properties of the cross-sectional radius function may be relaxed yet still give rise to a leading-order theory, though this phenomenon did not carry through to centrelines with high curvature or indeed cusps, which led to 
$\operatorname {ord}{(1)}$ errors in the boundary condition. This robustness to example cross-sectional radius functions that fall outside the proven scope of this methodology suggests potential applications to an even wider class of body geometries, though significant numerical and asymptotic analysis is required to establish such accuracy with confidence. Further analysis is also required to address questions of uniqueness inherent to regularised theories more generally, with the validity of our regularised slender-body theory evidenced here by direct comparison to a gold-standard boundary element methodology (Pozrikidis Reference Pozrikidis2002; Walker et al. Reference Walker, Wheeler, Ishimoto and Gaffney2019).
In summary, we have presented a simple regularised slender-body theory with algebraic accuracy, applicable to a wide range of slender filament geometries with circular cross-sections, which significantly expands upon existing methodologies. Furthermore, the presented theory overcomes numerical instabilities associated with many previous slender-body theories, making use of singularities with spatially dependent regularisation whilst retaining algebraic errors in the body slenderness ratio. The accuracy, simplicity and flexibility provided by this approach has the potential to facilitate extensive investigations of low-Reynolds-number motion, including numerous physical and biological applications. Finally, in a special but physically relevant case we have derived a leading-order solution to the corresponding flow problem, linking the flow about straight slender bodies to an elegant classical solution.
Acknowledgments
B.J.W. is supported by the UK Engineering and Physical Sciences Research Council (EPSRC), grant EP/N509711/1. K.I. is supported by JSPS-KAKENHI for Young Researchers (18K13456) and JST, PRESTO grant number JPMJPR1921. This publication is based, in part, on work supported by award no. KUK-C1-013-04, made by King Abdullah University of Science and Technology (KAUST). The research materials supporting this publication can be accessed at https://doi.org/10.5287/bodleian:xqzd9M45j.
Declaration of interests
The authors report no conflict of interest.
Appendix A. Radius functions of figure 2
The radius functions used in figure 2 have functional forms given by
(a)
$\sqrt{1-s^{2}} ,$(b)
$\sqrt{1-s^{2}}(1-0.1\cos{2{\rm \pi} s}) ,$(c)
$\sqrt{1-s^{2}}(1.1+\sin{9{\rm \pi} s}) ,$(d)
$(1-\left\lvert {s}\right\rvert^{1/4})(1+[s+0.1]^{2}) ,$(e)
$(1-s)(s+1)(s^{3}+1)(s^{2}+0.2) ,$
where labels correspond to those in figure 2 and these functional forms are equal to 
$\eta (s)$ subject to normalisation, recalling 
$\eta (s)\in [0,1]$.
Appendix B. Centrelines and radius function of figure 4
The plane curves used as centrelines in figure 4, subject to normalisation, are parameterised in Cartesian 
$xy$-coordinates as
\begin{equation*} \begin{array}{@{}lll@{}} {(a)} & (x(t),y(t))=(t,t^{2}), & t\in[-1,1] , \\ {(b)} & (x(t),y(t))=(\sin{t},\sin{t}\cos{t}), & t\in[{\rm \pi}/2,3{\rm \pi}/2] , \\ {(c)} & (x(t),y(t))=(t^{1/4}\sin{t},t^{1/4}\cos{t}), & t\in[1/10,2{\rm \pi}] , \\ {(d)} & (x(t),y(t))=(t,\tanh{100t}), & t\in[-1/2,1/2] , \\ {(e)} & (x(t),y(t))=\\ &(16\sin^{3}{t},13\cos{t}-5\cos{2t}-2\cos{3t}-\cos{4t}), & t\in[{\rm \pi}/2,3{\rm \pi}/2] , \end{array} \end{equation*}where labels correspond to those in figure 4 and the third Cartesian components are identically zero. The radius function used is given by 
$\sqrt {1-s^{2}}(1.1+\sin {{\rm \pi} s})$, subject to normalisation, as exemplified in figure 5.
