Unified mobility expressions for externally driven and self-phoretic propulsion of particles

Abstract

The mobility of externally driven phoretic propulsion of particles is evaluated by simultaneously solving the solute conservation equation, interaction potential equation and the modified Stokes equation. While accurate, this approach is cumbersome, especially when the interaction potential decays slowly compared with the particle size. In contrast to external phoresis, the motion of self-phoretic particles is typically estimated by relating the translation and rotation velocities with the local slip velocity. While this approach is convenient and thus widely used, it is only valid when the interaction decay length is significantly smaller than the particle size. Here, by taking inspiration from Brady (J. Fluid Mech., vol. 922, 2021, A10), which combines the benefits of two approaches, we reproduce their unified mobility expressions with arbitrary interaction potentials and show that these expressions can conveniently recover the well-known mobility relationships of external electrophoresis and diffusiophoresis for arbitrary double-layer thickness. Additionally, we show that for a spherical microswimmer, the derived expressions relax to the slip velocity calculations in the limit of the thin interaction length scales. We also employ the derived mobility expressions to calculate the velocities of an autophoretic Janus particle. We find that there is significant dampening in the translation velocity even when the interaction length is an order of magnitude larger than the particle size. Finally, we study the motion of a catalytically self-propelled particle, while it also propels due to external concentration gradients, and demonstrate how the two propulsion modes compete with each other.

1. Introduction

Phoretic phenomena, i.e. the movement of particles in response to an external field (Anderson 1989; Velegol et al. 2016; Khair 2022), is pivotal for a range of applications such as separation of biomacromolecules (Heller 2001; Lee et al. 2012), purification of nucleic acids from whole blood (Persat, Marshall & Santiago 2009), measurement of zeta potential (Doane et al. 2012; Shin et al. 2017a), banding of colloidal particles (Abécassis et al. 2008; Banerjee & Squires 2019; Raj, Shields & Gupta 2023b), membraneless water filtration (Shin et al. 2017b) and understanding biological pattern formation (Alessio & Gupta 2023), among others. In contrast, self-phoretic particles, also known as microswimmers, respond to a self-generated field gradient (Paxton et al. 2004; Howse et al. 2007; Ebbens & Howse 2010; Popescu, Uspal & Dietrich 2016; Ganguly, Alessio & Gupta 2023). Microswimmers are extensively studied for applications in targeted drug delivery (Xuan et al. 2014; Luo et al. 2018), environmental remediation (Gao et al. 2013; Wang et al. 2019), remote sensing of toxic chemicals (Esteban-Fernández de Ávila et al. 2016), the autonomous motion of microbots (Zarei & Zarei 2018; Hu, Liu & Sun 2022) and collective behaviour of active colloids (Palacci et al. 2013; Takatori & Brady 2016; Illien, Golestanian & Sen 2017). The mobility expressions of phoretic and self-phoretic processes are identical, with the key distinction being that the origin of field gradients in the two processes is different. In phoretic processes, this field gradient is externally imposed on colloidal particles, while in self-phoretic particles they are locally generated by the particles themselves typically through surface reactions or other mechanisms.

Studies on external electrophoretic motion have focused on the dependence of electrophoretic mobility on the effect of particle shape (Yoon & Kim 1989; Solomentsev & Anderson 1994), surface heterogeneity (Fair & Anderson 1992; Velegol, Anderson & Garoff 1996), finite double-layer thickness (Henry 1931; O'Brien & White 1978) and, more recently, strong deformation of double-layers (Khair 2018, 2022) and charge reversal (Kubíčková et al. 2012; Gupta et al. 2020a). Similarly, researchers have predicted the dependence of diffusiophoretic mobility (Anderson 1989; Brady 2011) on finite double-layer thickness (Prieve et al. 1984; Keh & Wei 2000), surface chemistry (Gupta, Shim & Stone 2020b) and multiple electrolytes (Gupta, Rallabandi & Stone 2019; Alessio et al. 2021). Studies on self-phoretic systems (Ramaswamy 2010; Moran & Posner 2017) focus on the impact of particle shape (Shklyaev, Brady & Córdova-Figueroa 2014; Nourhani & Lammert 2016; Poehnl, Popescu & Uspal 2020; Daddi-Moussa-Ider et al. 2021; Ganguly & Gupta 2023; Lee et al. 2023; Raj et al. 2023a), active patch shape (Lisicki, Reigh & Lauga 2018; Lee et al. 2021), surface interaction (Sharifi-Mood, Koplik & Maldarelli 2013) and finite Péclet number (Michelin & Lauga 2014).

Broadly speaking, there are two approaches for predicting the mobilities described above. The first approach solves the coupled solute conservation equations and the modified Stokes equation (Henry 1931; O'Brien & White 1978; Prieve et al. 1984; Prieve & Roman 1987; Anderson 1989; Keh & Wei 2000; Sharifi-Mood et al. 2013; Khair 2018, 2022; Gupta et al. 2019) and employs a force-free and torque-free condition to arrive at the translation velocity, $\boldsymbol {U}$, and rotational velocity, $\boldsymbol {\varOmega }$, of the particle. Thus, the above approach requires resolving the interaction potential simultaneously with the hydrodynamic equations. While exact and powerful, the methodology described above is cumbersome for analytical results when the particle–fluid interaction potential decays at much larger length scales than the particle size. Further, the solution strategy needs to be revised whenever the interaction potential changes, making it less convenient to be integrated into other analyses.

The second approach employs the reciprocal theorem in the thin interaction length limit (Stone & Samuel 1996; Brady 2011; Michelin & Lauga 2014; Lisicki et al. 2018; Masoud & Stone 2019; Poehnl et al. 2020; Ganguly & Gupta 2023; Raj et al. 2023a). In this limit, it is assumed that there exists a slip velocity, $\boldsymbol {u}_s$, at the interface of the inner region where the interaction potential is non-zero, and the outer region where the interaction potential is zero. This allows one to treat the outer problem as a classical Stokes flow problem with a slip boundary condition. Consequently, $\boldsymbol {U}$ and $\boldsymbol {\varOmega }$ can be represented as surface integrals of appropriate functions of $\boldsymbol {u}_s$. This approach based on the reciprocal theorem was first utilized by Stone & Samuel (1996) to study the impact of distortions in spherical microswimmers. This methodology is particularly powerful because, unlike the first method, computing $\boldsymbol {U}$ and $\boldsymbol {\varOmega }$ is relatively straightforward and agnostic to the mechanistic origin of $\boldsymbol {u}_s$. However, this approach is valid only when the interaction potential length is significantly smaller than the particle size, restricting its applicability. Additionally, it requires the knowledge of $\boldsymbol {u}_s$ a priori, and most studies have to rely on a lumped mobility parameter to estimate the value of $\boldsymbol {u}_s$ (Michelin & Lauga 2014; Lisicki et al. 2018; Poehnl & Uspal 2021; Ganguly & Gupta 2023; Raj et al. 2023a).

In this work, we seek to unify the benefits realized through the two approaches by employing the results of Brady (2021) and demonstrating how it can reconcile a large volume of mobility results for externally driven and self-phoretic propulsion of particles, and using these results for additional analyses. This approach is similar to the prior literature on the inertial correction to Stokes flow (Brenner & Cox 1963; Hinch 1991; Leal 2007), and swimming through non-Newtonian fluids (Datt et al. 2015; Elfring & Goyal 2016; Datt et al. 2017), where the first-order corrections include a body force term from the leading order and reciprocal theorem is employed to find the resulting motion. In § 2, we obtain a general mobility expression for an arbitrary particle shape, subjected to an osmophoretic (combination of osmotic and phoretic force) body force $\boldsymbol {b}$, identical to the results in Brady (2021). Subsequently, we take our expression to the thin interaction length scale limit and retrieve the mobility expressions in Stone & Samuel (1996), see § 3. In § 3, we also retrieve the expression for the electrophoretic mobility of translation of a charged spherical particle in an external electric field, obtained by Henry (1931). Additionally, we derive the diffusiophoretic mobility of a charged spherical particle in an externally imposed solute gradient at finite interaction lengths, as first obtained by Keh & Wei (2000). Finally in § 4, we apply our result to study the autophoretic motion of spherical microparticles with catalytic caps. We study how translation velocity depends on the cap size, the surface interaction potential and the interaction length relative to particle size. Since our methodology works for both externally driven and self-propelling particles, we also study particle propulsion by both modes simultaneously, see § 5. These model problems demonstrate the wide applicability of the expressions derived by Brady (2021) and reproduced in this manuscript. Finally, in § 6 we summarize the key findings of our work and outline future ideas.

2. Derivation of the unified mobility expression

In this section, we derive the translation velocity ($\boldsymbol {U}$) and rotational velocity ($\boldsymbol {\varOmega }$) of an arbitrary particle with surface $S_p$ immersed in a fluid of volume $V$ due to an arbitrary osmophoretic body force $\boldsymbol {b}$, see figure 1(a). The particle surface is defined via the vector $\boldsymbol {x}_S$ relative to the centre of mass of the particle. We define $\boldsymbol {e}_r$ as the outward unit normal to the particle surface, $r$ is the distance from the centre of mass of the particle and $\boldsymbol {r}$ is the position vector defined from the centre of mass of the particle. The fluid velocity around the particle can be resolved through the modified Stokes equation, defined as

(2.1)\begin{equation} \boldsymbol{\nabla}\boldsymbol{\cdot}\boldsymbol{\sigma} + \boldsymbol{b} = \boldsymbol{0}, \end{equation}

where $\boldsymbol {\sigma }$ is the hydrodynamic stress tensor. The velocity field $\boldsymbol {u}$ is assumed to decay to zero in the far-field, $\boldsymbol{u}|_{S_\infty} \to {\boldsymbol{0}}$. At the particle surface, the fluid obeys a no-slip, rigid body boundary condition, $\boldsymbol {u}|_{S_p} = \boldsymbol {U} + \boldsymbol {\varOmega } \times \boldsymbol {x}_{S}$. To obtain $\boldsymbol {U}$ and $\boldsymbol {\varOmega }$ using Lorentz reciprocal theorem (Masoud & Stone 2019), we define an auxiliary Stokes flow ($\hat {\boldsymbol {U}},\hat {\boldsymbol {\varOmega }}$) while preserving particle geometry with the same no-slip rigid surface, $\hat {\boldsymbol {u}}|_{S_p} = \hat {\boldsymbol {U}} + \hat {\boldsymbol {\varOmega }} \times \boldsymbol {x}_{S}$, and far-field decay, $\hat{\boldsymbol{u}}|_{S_\infty} \to {\boldsymbol{0}}$, boundary conditions.

Figure 1. Two approaches to finding the velocity of a particle of characteristic length scale a by resolving the fluid velocity at the particle surface. The particle surface is defined by a vector $\boldsymbol{x}_s$ relative to the centre of mass (COM) of the particle. (a) Obtain the fluid velocity near the particle surface by resolving the modified Stokes equation with an arbitrary body force, $\boldsymbol {b}$. The body force, $\boldsymbol {b}$, depends on charge, $\rho$, salt, $s$, and interaction potential, $\phi$. (b) When the interaction length is small $\lambda /a \ll 1$, the velocity near the fluid surface, at the outer edge of the interaction layer, $\boldsymbol {u}_{s}$, is taken to be the velocity at the particle surface. The slip velocity, $\boldsymbol {u}_{s}$, depends on the lumped mobility, $\mathcal {M}$, which depends on the interaction between the surface and solute, and the solute concentration at the vicinity of the surface, $c$.

Using the Lorentz reciprocal theorem, we can relate the phoretic problem ($\boldsymbol {U}$,$\boldsymbol {\varOmega }$,$\boldsymbol {b}$) and the auxiliary problem ($\hat {\boldsymbol {U}}$,$\hat {\boldsymbol {\varOmega }}$,$\hat {\boldsymbol {b}}$) to be

(2.2)\begin{equation} \int_{S_p} \boldsymbol{e}_r\boldsymbol{\cdot}\boldsymbol{\sigma}\boldsymbol{\cdot}\hat{\boldsymbol{u}}\,{\rm d}S - \int_{S_p} \boldsymbol{e}_r\boldsymbol{\cdot}\hat{\boldsymbol{\sigma}}\boldsymbol{\cdot}\boldsymbol{u}\,{\rm d}S = \int_V \hat{\boldsymbol{u}}\boldsymbol{\cdot}\boldsymbol{b}\,{\rm d}V - \int_V \boldsymbol{u}\boldsymbol{\cdot}\hat{\boldsymbol{b}}\,{\rm d}V. \end{equation}

Substituting in the expressions of the fluid velocities, $\boldsymbol {u}|_{S_p}$ and $\hat {\boldsymbol {u}}|_{S_p}$, at the particle surface we can simplify (2.2). Additionally, we assume that there is no body force in the auxiliary problem, ${\boldsymbol {\hat b}} = \boldsymbol {0}$. Thus, we can rewrite (2.2) to be

(2.3)\begin{equation} \int_{S_p} \boldsymbol{e}_r\boldsymbol{\cdot}\boldsymbol{\sigma}\boldsymbol{\cdot} \left(\hat{\boldsymbol{U}}+\hat{\boldsymbol{\varOmega}}\times \boldsymbol{x}_S\right) {\rm d} S - \int_{S_p} \boldsymbol{e}_r\boldsymbol{\cdot}\hat{\boldsymbol{\sigma}}\boldsymbol{\cdot} \left(\boldsymbol{U}+\boldsymbol{\varOmega}\times\boldsymbol{x}_S\right) {\rm d}S = \int_V \hat{\boldsymbol{u}}\boldsymbol{\cdot}\boldsymbol{b}\,{\rm d}V. \end{equation}

Since the inertia of the particle is negligible, for both the phoresis and the auxiliary problem, the particle is force-free and torque-free. For the phoretic propulsion, the hydrodynamic force and torque is balanced by the osmophoretic force and torque, or

(2.4)$$\begin{gather} \underbrace{\displaystyle \int_{S_p} \boldsymbol{e}_r\boldsymbol{\cdot}\boldsymbol{\sigma}\, {\rm d}S}_\textit{hydrodynamic} \ \underbrace{- \displaystyle \int_V \boldsymbol{b}\, {\rm d}V}_\textit{osmophoretic} = \boldsymbol{0}, \end{gather}$$

(2.5)$$\begin{gather}\underbrace{\displaystyle \int_{S_p} \boldsymbol{x}_s \times \boldsymbol{e}_r \boldsymbol{\cdot} \boldsymbol{\sigma} \,{\rm d}S}_\textit{hydrodynamic} \ \underbrace{- \displaystyle \int_V \boldsymbol{r}\times\boldsymbol{b}\, {\rm d}V}_\textit{osmophoretic} = \boldsymbol{0}, \end{gather}$$

where the negative sign in front of the osmophoretic term appears because the osmophoretic force on the particle is equal in magnitude to the osmophoretic force on the fluid but opposite in sign (Brady 2011). For the auxiliary system, we balance the hydrodynamic and external forces and torques, or

(2.6)$$\begin{gather} \underbrace{\displaystyle \int_{S_p} \hat{\boldsymbol{\sigma}}\boldsymbol{\cdot}\boldsymbol{e}_r\,{\rm d}S }_\textit{hydrodynamic} \ + \underbrace{\hat{\boldsymbol{F}}_{ext}}_\textit{external}= \boldsymbol{0}, \end{gather}$$

(2.7)$$\begin{gather}\underbrace{\displaystyle \int_{S_p} \boldsymbol{x}_s \times \hat{\boldsymbol{\sigma}}\boldsymbol{\cdot}\boldsymbol{e}_r\,{\rm d}S }_\textit{hydrodynamic} \ + \underbrace{\hat{\boldsymbol{L}}_{ext}}_\textit{external}= \boldsymbol{0}, \end{gather}$$

where $\hat {\boldsymbol {F}}_{ext}$ and $\hat {\boldsymbol {L}}_{ext}$ are the external force and torque required to move the particle in the auxiliary problem, which is to be determined. We note that (2.4)–(2.5) are different from (2.6)–(2.7) since the phoretic propulsion is induced by $\boldsymbol {b}$ whereas in the auxiliary problem motion is caused by $\hat {\boldsymbol {F}}_{ext}$ and $\hat {\boldsymbol {L}}_{ext}$.

To calculate $\boldsymbol {U}$ and $\boldsymbol {\varOmega }$ for a given $\boldsymbol {b}$ from (2.3)–(2.7), we need to express, $\hat {\boldsymbol {F}}_{ext}$, $\hat {\boldsymbol {L}}_{ext}$ and $\boldsymbol {\hat {u}}$ in the auxiliary problem as functions of $\boldsymbol {\hat {U}}$ and $\boldsymbol {\hat {\varOmega }}$. To do so, we use a resistance formulation to write

(2.8) \begin{equation} \begin{bmatrix} \hat{\boldsymbol{F}}_{ext}\\ \hat{\boldsymbol{L}}_{ext} \end{bmatrix} = \begin{bmatrix} \boldsymbol{\mathsf{R}}_{FU} & \boldsymbol{\mathsf{R}}_{F\varOmega} \\ \boldsymbol{\mathsf{R}}_{LU} & \boldsymbol{\mathsf{R}}_{L\varOmega} \end{bmatrix}\boldsymbol{\cdot}\begin{bmatrix} \hat{\boldsymbol{U}} \\ \hat{\boldsymbol{\varOmega}} \end{bmatrix}, \end{equation}

where the resistance matrices $\boldsymbol{\mathsf{R}}_{FU}$, $\boldsymbol{\mathsf{R}}_{F\varOmega }$, $\boldsymbol{\mathsf{R}}_{LU}$ and $\boldsymbol{\mathsf{R}}_{L\varOmega }$ relate the driving force ($\hat {\boldsymbol {F}}_{ext}$) and torque ($\hat {\boldsymbol {L}}_{ext}$) to the translational ($\boldsymbol {\hat {U}}$) and rotational velocity ($\boldsymbol {\hat {\varOmega }}$). Further, we describe $\hat {\boldsymbol {u}}$ as

(2.9)\begin{equation} \hat{\boldsymbol{u}}=\boldsymbol{\mathsf{D}}_T\boldsymbol{\cdot}\hat{\boldsymbol{U}}+ \boldsymbol{\mathsf{D}}_R\boldsymbol{\cdot}\hat{\boldsymbol{\varOmega}}\times\boldsymbol{r}, \end{equation}

where $\boldsymbol{\mathsf{D}}_T$ is the translation disturbance tensor and $\boldsymbol{\mathsf{D}}_R$ is the rotation disturbance tensor. Equations (2.8)–(2.9) combined with (2.6)–(2.7) provide necessary information to simplify (2.3) as a function of $\boldsymbol {\hat {U}}$ and $\boldsymbol {\hat {\varOmega }}$.

Next, we choose convenient values of $\boldsymbol {\hat {U}}$ and $\boldsymbol {\hat {\varOmega }}$ to simplify (2.3). Specifically, we use six auxiliary flow problems, pure translation ($\hat {\boldsymbol {\varOmega }} = \boldsymbol {0}$ and $\hat {\boldsymbol {U}}$ = $U_0 \boldsymbol {e}_1$, $U_0 \boldsymbol {e}_2$, $U_0 \boldsymbol {e}_3$) and pure rotation ($\hat {\boldsymbol {U}}=\boldsymbol {0}$ and $\hat {\boldsymbol {\varOmega }}=U_0/a \boldsymbol {e}_1$, $U_0/a \boldsymbol {e}_2$, $U_0/a \boldsymbol {e}_3$) with $U_0$ being the characteristic velocity scale and $a$ being the characteristic particle length, to obtain

(2.10) \begin{equation} \begin{bmatrix} \boldsymbol{\mathsf{R}}_{FU} & \boldsymbol{\mathsf{R}}_{F\varOmega} \\ \boldsymbol{\mathsf{R}}_{LU} & \boldsymbol{\mathsf{R}}_{L\varOmega} \end{bmatrix}\boldsymbol{\cdot}\begin{bmatrix} {\boldsymbol{U}} \\ {\boldsymbol{\varOmega}} \end{bmatrix} = \begin{bmatrix} \displaystyle \int_V \left(\boldsymbol{\mathsf{D}}_T - \boldsymbol{I}\right)\boldsymbol{\cdot}\boldsymbol{b}\,{\rm d}{V} \\ \displaystyle \int_V \left(\boldsymbol{\mathsf{D}}_R - \boldsymbol{I}\right)\boldsymbol{\cdot}\boldsymbol{r}\times\boldsymbol{b}\,{\rm d}{V} \end{bmatrix} . \end{equation}

By inverting the resistance tensor, we obtain a mobility formulation that resolves $\boldsymbol {U}$ and $\boldsymbol {\varOmega }$ in terms of volume integrals of $\boldsymbol {b}$,

(2.11) \begin{equation} \begin{bmatrix} \boldsymbol{U} \\ \boldsymbol{\varOmega} \end{bmatrix} = \begin{bmatrix} \boldsymbol{\mathsf{M}}_{UF} & \boldsymbol{\mathsf{M}}_{UL} \\ \boldsymbol{\mathsf{M}}_{\varOmega F} & \boldsymbol{\mathsf{M}}_{\varOmega L} \end{bmatrix}\boldsymbol{\cdot}\begin{bmatrix} \displaystyle \int_V \left(\boldsymbol{\mathsf{D}}_T - \boldsymbol{I}\right)\boldsymbol{\cdot}\boldsymbol{b}\,{\rm d}{V} \\ \displaystyle \int_V \left(\boldsymbol{\mathsf{D}}_R - \boldsymbol{I}\right)\boldsymbol{\cdot} \boldsymbol{r}\times\boldsymbol{b}\,{\rm d}{V} \end{bmatrix} \end{equation}

where the matrices $\boldsymbol{\mathsf{M}}_{UF}$, $\boldsymbol{\mathsf{M}}_{\varOmega F}$, $\boldsymbol{\mathsf{M}}_{UL}$ and $\boldsymbol{\mathsf{M}}_{\varOmega L }$ are the corresponding mobility tensors. For an in-depth mathematical analysis and mechanistic discussion regarding the various forms of $\boldsymbol {b}$, we redirect the reader to Brady (2021).

Physically, (2.11) is insightful as it helps parse apart the difference between the phoretic problem and the auxiliary problem. The rightmost term is the effective force and torque on the particle due to phoretic interactions and has two contributions: (i) the term associated with the identity tensor ($\boldsymbol{\mathsf{I}}$) is the osmophoretic force and torque acting on the particle, and (ii) the term associated with the disturbance tensors ($\boldsymbol{\mathsf{D}}_T$, $\boldsymbol{\mathsf{D}}_R$) is the hydrodynamic correction to the distribution of body forces around the particle. This correction arises because the phoretic interactions near the particle surface lead to an additional compensating fluid motion (Brady 2011) causing a long-range hydrodynamic disturbance. This effect is not captured in the definition of the hydrodynamic mobility tensor and thus manifests separately. If the terms associated with disturbance tensors were not present, (2.11) is essentially identical to (2.8) with osmophoretic force on the particle as the external force.

We note that (2.11) takes an explicit form only when $\boldsymbol {b}$ is independent of $\boldsymbol {U}$ and $\boldsymbol {\varOmega }$ and is thus most convenient for systems with small Péclet number (${Pe} \ll 1$), where ${Pe} = {U_0 a}/{D}$, and $D$ is the diffusivity of the solute. The distinction from prior work, such as Stone & Samuel (1996), Michelin & Lauga (2014), Lisicki et al. (2018), Poehnl et al. (2020), Poehnl & Uspal (2021) and Ganguly & Gupta (2023), that utilize Lorentz reciprocal theorem is that they invoke the thin interaction length limit and apply the analysis in the outer region where $\boldsymbol {b}=\boldsymbol {0}$, see figure 1(b). Consequently, they do not arrive at (2.11) but rather represent $\boldsymbol {U}$ and $\boldsymbol {\varOmega }$ in terms of a slip velocity at the particle surface $\boldsymbol {u}_s$.

We acknowledge that similar results have been presented in Khair (2018) and Brady (2021). However, in Khair (2018), $\boldsymbol {b}$ only focused on the electrophoretic contributions. In contrast, Brady (2021) argued that $\boldsymbol {b}$ should include both osmotic and phoretic contributions, and we thus refer to $\boldsymbol {b}$ as an osmophoretic body force. Care should be taken that the osmotic contribution only includes an excess osmotic effect since a particle cannot move without a phoretic interaction; interested readers are referred to Brady (2021). The phoretic contribution arises from the interaction of the particle with a macroscopically established potential field. The nature of this field depends on the specific model problem under consideration. We refer the readers to Brady (2021) for an in-depth mathematical analysis and a general discussion on the mechanistic origin of $\boldsymbol {b}$. Building on the work by Brady (2021), we systematically illustrate how both phoretic and osmotic contributions to the body force term are required to reconcile a broad range of results in the literature and arrive at universal mobility relationships. Additionally, through this framework, we quantify the impact of interaction length on microswimmer motion in electrolytic solutions, elaborating on the suggestion established in Brady (2021).

For a spherical particle (see Duprat & Stone (2016) for derivation), the relevant hydrodynamic parameters are $\boldsymbol{\mathsf{D}}_T = {3a}/{4r}(\boldsymbol{\mathsf{I}}+\boldsymbol {e}_r\boldsymbol {e}_r)+{a^3}/{4r^3}(\boldsymbol{\mathsf{I}}-3\boldsymbol {e}_r \boldsymbol {e}_r)$, $\boldsymbol{\mathsf{D}}_R = ({a^3}/{r^3})\boldsymbol{\mathsf{I}}$, $\boldsymbol{\mathsf{M}}_{UF} = ({1}/{6{\rm \pi} \mu a})\boldsymbol{\mathsf{I}}$, $\boldsymbol{\mathsf{M}}_{UL}=0$, $\boldsymbol{\mathsf{M}}_{\varOmega F}=0$ and $\boldsymbol{\mathsf{M}}_{\varOmega L} = ({1}/{8{\rm \pi} \mu a^3}) \boldsymbol{\mathsf{I}}$; where $a$ is the radius of the sphere, $\mu$ is the fluid viscosity, $r$ is the radial distance from the centre of the sphere and $\boldsymbol {e}_r$ is the radial vector pointing away from the centre. Substituting, these definitions of the hydrodynamic disturbance and mobility in (2.11) we obtain

(2.12)$$\begin{gather} \boldsymbol{U} = \frac{1}{6{\rm \pi}\mu a} \displaystyle \int_V \left[\left(\frac{3a}{2r}-\frac{a^3}{2r^3}-1\right)\boldsymbol{b}_\perp{+} \left(\frac{3a}{4r}+\frac{a^3}{4r^3}-1\right)\boldsymbol{b}_\parallel\right]{\rm d}{V}, \end{gather}$$

(2.13)$$\begin{gather}\boldsymbol{\varOmega} = \frac{1}{8{\rm \pi}\mu a^3} \displaystyle \int_V r\left(\frac{a^3}{r^3}-1\right)\boldsymbol{e}_r \times \boldsymbol{b}_\parallel\,{\rm d}{V}, \end{gather}$$

where the body force is decomposed into $\boldsymbol {b}=\boldsymbol {b}_\perp +\boldsymbol {b}_\parallel$. The perpendicular subscript denotes the component normal to the sphere and the parallel subscript denotes the component parallel to the surface. Equations (2.12)–(2.13) were also reported in the prior literature for phoretic systems (Brady 2021) as well as for different physical systems (Brenner & Cox 1963; Hinch 1991; Leal 2007; Datt et al. 2015; Elfring & Goyal 2016; Datt et al. 2017). We extensively validate this result in the next section and show it relaxes to the various well-known expressions present in the literature, for both microswimmers and externally driven particles. Further, we employ this expression to study a microswimmer in the arbitrary interaction layer limit and a microswimmer driven by an external gradient in addition to its self-propelling mode of swimming.

3. Validation

3.1. Simplification at the limit of the thin interaction length scale

In this subsection, we aim to simplify (2.12)–(2.13) in the limit of the thin interaction length scale for a spherical microswimmer and recover the equations discussed in Stone & Samuel (1996).

We proceed to simplify (2.12)–(2.13) at the thin interaction length limit, ${\lambda }/{a} \ll 1$, where $\lambda$ is the interaction length scale. To this end, we define a stretched radial coordinate $\rho = {(r-a)}/{\lambda }$. Next, we expand and rewrite (2.12)–(2.13) in orders of ${\lambda }/{a}$. Subsequently, the leading-order contribution to the translation and rotation velocities are obtained to be

(3.1)$$\begin{gather} \boldsymbol{U} ={-}\frac{\lambda}{4{\rm \pi}\mu a^2} \displaystyle \int_V \rho \boldsymbol{b}_\parallel\,{\rm d}{V}, \end{gather}$$

(3.2)$$\begin{gather}\boldsymbol{\varOmega} ={-}\frac{3\lambda}{8{\rm \pi}\mu a^3} \displaystyle \int_V \rho \boldsymbol{e}_r \times \boldsymbol{b}_\parallel\,{\rm d}{V} . \end{gather}$$

Since the volume of interest at the thin interaction limit is a spherical shell of thickness $\lambda$ surrounding the particle, we can rewrite the differential volume element to be ${\rm d}{V} = \lambda \,{\rm d}\rho \,{\rm d}{S}$ and the volume integrals as

(3.3)$$\begin{gather} \boldsymbol{U} ={-}\frac{\lambda^2}{4{\rm \pi}\mu a^2} \displaystyle \int_{S_p}\left[\displaystyle \int_0^\infty \rho \boldsymbol{b}_\parallel\,{\rm d}\rho\right]{\rm d}{S}, \end{gather}$$

(3.4)$$\begin{gather}\varOmega ={-}\frac{3\lambda^2}{8{\rm \pi}\mu a^3} \displaystyle \int_{S_P} \boldsymbol{e}_r \times \left[\displaystyle \int_0^\infty \rho \boldsymbol{b}_\parallel\,{\rm d}\rho\right]{\rm d}S. \end{gather}$$

In the thin interaction limit, the shear force is balanced by the parallel body force, or

(3.5)\begin{equation} \frac{\mu}{\lambda^2}\frac{\partial^2 \boldsymbol{u}_\parallel}{\partial \rho^2} + \boldsymbol{b}_\parallel= 0 . \end{equation}

We note that in (3.5), since $\boldsymbol {b}$ is the osmophoretic force, $\boldsymbol {b}_\parallel$ also includes the excess osmotic term. Multiplying (3.5) by $\rho$ and employing integration by parts, we arrive at

(3.6)\begin{equation} \displaystyle \int_0^\infty \rho \boldsymbol{b}_\parallel\,{\rm d}{\rho} = \frac{\mu}{\lambda^2} \boldsymbol{u}_\parallel |_0^\infty = \frac{\mu}{\lambda^2} \boldsymbol{u}_s, \end{equation}

where $\boldsymbol {u}_s = \boldsymbol {u}_{\parallel,\infty }-\boldsymbol {u}_{\parallel,0}$, is the phoretic slip velocity. Substituting (3.6) into (3.3)–(3.4), we get the widely used result derived in Stone & Samuel (1996),

(3.7)$$\begin{gather} \boldsymbol{U} ={-}\frac{1}{4{\rm \pi} a^2} \displaystyle \int_{S_P} \boldsymbol{u}_s\,{\rm d}{S}, \end{gather}$$

(3.8)$$\begin{gather}\boldsymbol{\varOmega} ={-}\frac{3}{8{\rm \pi} a^3} \displaystyle \int_{S_P} \boldsymbol{e}_r \times \boldsymbol{u}_s\,{\rm d}{S}, \end{gather}$$

where the integral is over the surface of the sphere.

3.2. Electrophoretic mobility at arbitrary interaction length scales

To further validate (2.12)–(2.13) by the determination of the electrophoretic mobility of a sphere in the Debye–Hückel limit for an arbitrary Debye length (Henry 1931; Teubner 1982; Kim & Karrila 2013), we assume a homogeneous sphere of radius, $a$, immersed in a binary monovalent electrolytic solution such that the electrical permittivity of the solution is denoted as $\varepsilon$. Our objective is to analyse the motion of the particle with a given surface zeta potential driven by an external electric field. First, we assume that the surface zeta potential, $\zeta$, falls in the Debye–Hückel limit, ${e\zeta }/{k_B T} \ll 1$, where $e$ is the charge of an electron, $k_B$ is the Boltzmann constant and $T$ is the absolute temperature. We also assume an electric field disturbance of $\boldsymbol {E}_\infty$ far away from the particle such that $\boldsymbol {E}_\infty = \epsilon E_0 \boldsymbol {e}_z$, where $\epsilon$ is a small parameter and physically indicates that the length scale of the far-field potential decay is much larger than the particle size, see figure 2(a). Note that $\varepsilon$ is the electrical permittivity and should not be confused with $\epsilon$, which is a small parameter in our analysis. Finally, the total osmophoretic body force ($\boldsymbol {b}$) driving the particle arises through a combination of the electrostatic interaction and net excess osmotic pressure in the fluid, or

(3.9)\begin{equation} \boldsymbol{b} ={-}e\left(c_{+} - c_{-}\right) \boldsymbol{\nabla} \phi - k_B T \boldsymbol{\nabla}\left(c_{+} + c_{-}\right), \end{equation}

where $c_{+}$ and $c_{-}$ are the concentrations of the positive and negative electrolytic species, respectively, and $\phi$ is the electric potential. As a convenient choice, we can represent the solute concentrations in terms of net charge, $\rho = e ( c_{+} - c_{-} )$, and salt, $s = c_{+} + c_{-}$, and rewrite the body force to be $\boldsymbol {b} = - \rho \boldsymbol {\nabla } \phi - k_B T \boldsymbol {\nabla } s$. It should be noted that care should be exercised in choosing the appropriate expression for $\boldsymbol {b}$. Specifically, the osmotic contribution $-k_B T \boldsymbol {\nabla }s$ refers to the excess osmotic contribution arising out of an interaction that locally drives the solute out of equilibrium. This effectively implies that in the absence of such interactions, an external salt gradient on its own cannot induce net particle motion, as demonstrated in Brady (2021). The equivalence of (2.12) and the results of Brady (2021) can be seen by defining an additional surface stress contribution, $\boldsymbol {\sigma }_p = -k_B T s \boldsymbol{\mathsf{I}}$, as per (2.17) in Brady (2021), due to the excess osmotic pressure. The divergence of $\boldsymbol {\sigma }_p$ leads to the second term in (3.9), $-k_B T \boldsymbol {\nabla } s$. We refer the readers to our discussion on the mechanistic origin of $\boldsymbol {b}$ in § 2 and redirect them to Brady (2021) for more details.

Figure 2. Methodology to validate proposed mobility expressions for a charged particle with a zeta potential ($\zeta$) in the Debye–Hückel limit for (a) electrophoresis with an external field $\boldsymbol {E}_\infty =\epsilon E_0 \boldsymbol {e}_z$ and (b) diffusiophoresis with externally imposed solute gradient $\boldsymbol {\nabla }s_\infty = 2\epsilon c_0 \boldsymbol {e}_z$. The expressions of dimensionless osmophoretic force $\tilde {\boldsymbol {b}}$ are provided. Substituting the appropriate $\boldsymbol {b}$ in (2.12) enables us to recover mobility relationships that otherwise require cumbersome calculations.

To appropriately derive the particle motion, we are required to obtain the solutions to $\rho$, $s$ and $\phi$ for a given $\boldsymbol {E}_\infty$ and $\zeta$. As mentioned earlier, we assume that the species are monovalent, $z_\pm =\pm 1$ and have diffusivities $D_\pm$, the species balance is given by the steady Nernst–Planck equations,

(3.10)$$\begin{gather} \boldsymbol{u}\boldsymbol{\cdot}\boldsymbol{\nabla}c_+=D_+ \left[\nabla^2 c_{+} + \frac{e}{k_B T} \boldsymbol{\nabla}\boldsymbol{\cdot}\left(c_{+}\boldsymbol{\nabla}\phi\right)\right] , \end{gather}$$

(3.11)$$\begin{gather}\boldsymbol{u}\boldsymbol{\cdot}\boldsymbol{\nabla}c_-= D_- \left[\nabla^2 c_{-} - \frac{e}{k_B T} \boldsymbol{\nabla}\boldsymbol{\cdot}\left(c_{-}\boldsymbol{\nabla}\phi\right)\right]. \end{gather}$$

For ${Pe} = aU/D \ll 1$ ($U$ is the velocity scale for the particle, $D = {2 D_+ D_-}/({D_+ + D_-})$ is the ambipolar diffusivity), we ignore the convective effects. Consequently, (3.10) and (3.11) can be rewritten in terms of $\rho$ and $s$ as

(3.12)$$\begin{gather} \nabla^2 s + \frac{1}{k_B T} \boldsymbol{\nabla}\boldsymbol{\cdot}\left(\rho\boldsymbol{\nabla}\phi\right) = 0, \end{gather}$$

(3.13)$$\begin{gather}\nabla^2\rho + \frac{e^2}{k_B T} \boldsymbol{\nabla}\boldsymbol{\cdot}\left(s\boldsymbol{\nabla}\phi\right) = 0. \end{gather}$$

Finally, the system of equations is closed by using Poisson's equation to resolve the electric potential,

(3.14)\begin{equation} - \varepsilon \nabla^2 \phi = \rho. \end{equation}

In the far-field, at $r\to \infty$, the potential gradient is the externally imposed electric field and the fluid is electroneutral, or

(3.15)$$\begin{gather} \left. -\boldsymbol{\nabla}\phi \right| _{r \rightarrow \infty}= \epsilon E_0 \boldsymbol{e}_z. \end{gather}$$

(3.16)$$\begin{gather}\left. \rho \right|_{r \rightarrow \infty} = 0. \end{gather}$$

Moreover, in the far-field $c_{+} = c_{-} = c_0$, where $c_0$ is a characteristic solute concentration, we can write $s$ to follow

(3.17)\begin{equation} s |_{r \rightarrow \infty} = 2c_0. \end{equation}

At the particle surface, at $r=a$, the electrostatic potential is equal to the zeta potential at the surface, or

(3.18)\begin{equation} \phi |_{r = a} = \zeta. \end{equation}

Additionally, there is no salt or charge flux normal to the particle surface,

(3.19)$$\begin{gather} \boldsymbol{e}_r\boldsymbol{\cdot}\left[\boldsymbol{\nabla}s + \frac{\rho}{k_B T} \boldsymbol{\nabla}\phi\right]_{r = a} = 0, \end{gather}$$

(3.20)$$\begin{gather}\boldsymbol{e}_r\boldsymbol{\cdot}\left[\boldsymbol{\nabla}\rho + \frac{s e^2}{k_B T} \boldsymbol{\nabla}\phi\right]_{r = a} = 0. \end{gather}$$

Equations (3.12)–(3.20) are non-dimensionalized using the following appropriate scales:

(3.21a–f)\begin{equation} \tilde{\boldsymbol{\nabla}} = a \boldsymbol{\nabla},\quad \tilde{\nabla}^2 = a^2 \nabla^2,\quad \tilde{\phi} = \frac{e\phi}{k_B T},\quad\tilde{\rho} = \frac{\rho}{e c_0},\quad \tilde{s} = \frac{s}{c_0}, \quad\tilde{r} = \frac{r}{a}. \end{equation}

Thus, the non-dimensional Poisson–Nernst–Planck equations are given as

(3.22)$$\begin{gather} \tilde{\nabla}^2\tilde{s} + {\tilde{\boldsymbol{\nabla}}}\boldsymbol{\cdot} \left(\tilde{\rho}\tilde{\boldsymbol{\nabla}}\tilde{\phi}\right) = 0, \end{gather}$$

(3.23)$$\begin{gather}\tilde{\nabla}^2\tilde{\rho} + \tilde{\boldsymbol{\nabla}}\boldsymbol{\cdot} \left(\tilde{s}\tilde{\boldsymbol{\nabla}}\tilde{\phi}\right) = 0, \end{gather}$$

(3.24)$$\begin{gather}\tilde{\nabla}^2 \tilde{\phi} ={-}\frac{\kappa^2}{2} \tilde{\rho}, \end{gather}$$

where $\kappa = ({2 a^2 e^2 c_0 }/{\varepsilon k_B T})^{1/2}$ is the dimensionless inverse Debye length. In the far-field, thus

(3.25)$$\begin{gather} -\left. \tilde{\boldsymbol{\nabla}}\tilde{\phi} \right|_{\tilde{r} \rightarrow \infty} = \epsilon \tilde{E}_0 \boldsymbol{e}_z, \end{gather}$$

(3.26)$$\begin{gather}\left. \tilde{\rho} \right|_{\tilde{r} \rightarrow \infty} = 0 , \end{gather}$$

(3.27)$$\begin{gather}\left. \tilde{s}\right|_{\tilde{r} \rightarrow \infty} = 2 , \end{gather}$$

where $\tilde{E}_0 = a e E_0/(k_B T )$ is the non-dimensional electric-field. Similarly, at the particle surface, the non-dimensional boundary conditions read

(3.28)$$\begin{gather} \boldsymbol{e}_r \boldsymbol{\cdot} \left[\tilde{\boldsymbol{\nabla}}\tilde{s} + \tilde{\rho} \tilde{\boldsymbol{\nabla}}\tilde{\phi}\right]_{\tilde{r} = 1} = 0, \end{gather}$$

(3.29)$$\begin{gather}\boldsymbol{e}_r \boldsymbol{\cdot} \left[\tilde{\boldsymbol{\nabla}}\tilde{\rho} + \tilde{s} \tilde{\boldsymbol{\nabla}} \tilde{\phi}\right]_{\tilde{r} = 1} = 0, \end{gather}$$

(3.30)$$\begin{gather}\left. \tilde{\phi} \right|_{\tilde{r} = 1} = \frac{e \zeta}{k_B T} = \tilde{\zeta}. \end{gather}$$

As we discuss later, it is more appropriate to write (3.30) as a constant charge boundary condition, which renders the gradient of the potential to be constant instead. However, for the weak-field, these surface boundary conditions are equivalent and hence we retain the constant potential boundary condition for simplicity.

For the remainder of the calculation until (3.67), we will drop the tilde superscript in (3.22)–(3.30) for convenience and restore the dimensions once the non-dimensional calculations are complete. We expand $\phi$, $\rho$ and $s$ in the small parameters $\zeta$ and $\epsilon$ as

(3.31)$$\begin{gather} \phi = \phi_{00} + \zeta \phi_{01} + \epsilon \left(\phi_{10} + \zeta\phi_{11}\right), \end{gather}$$

(3.32)$$\begin{gather}\rho = \rho_{00} + \zeta \rho_{01} + \epsilon\left(\rho_{10} + \zeta \rho_{11}\right), \end{gather}$$

(3.33)$$\begin{gather}s = s_{00} + \zeta s_{01} + \epsilon \left(s_{10} + \zeta s_{11}\right). \end{gather}$$

The asymptotic expansions in (3.31)–(3.33) are substituted into (3.22)–(3.30), and the corresponding equations are solved at each asymptotic order.

Order $O(1)$: an uncharged particle without any electric field. The governing equations and boundary conditions are obtained to be

(3.34)$$\begin{gather} \nabla^2 s_{00} + \boldsymbol{\nabla}\boldsymbol{\cdot}\left(\rho_{00} \boldsymbol{\nabla}\phi_{00}\right) = 0, \end{gather}$$

(3.35)$$\begin{gather}\nabla^2 \rho_{00} + \boldsymbol{\nabla}\boldsymbol{\cdot}\left(s_{00}\boldsymbol{\nabla}\phi_{00}\right) = 0, \end{gather}$$

(3.36)$$\begin{gather}\nabla^2 \phi_{00} ={-}\frac{\kappa^2}{2} \rho_{00}, \end{gather}$$

(3.37)$$\begin{gather}\boldsymbol{e}_r\boldsymbol{\cdot}\left[\boldsymbol{\nabla}s_{00}+\rho_{00}\boldsymbol{\nabla}\phi_{00}\right]=0,\quad\mathrm{at}\ r=1, \end{gather}$$

(3.38)$$\begin{gather}\boldsymbol{e}_r \boldsymbol{\cdot} \left[\boldsymbol{\nabla} \rho_{00} + s_{00} \boldsymbol{\nabla}\phi_{00}\right]=0,\quad \mathrm{at}\ r=1, \end{gather}$$

(3.39)$$\begin{gather}\phi_{00} = 0,\quad \mathrm{at}\ r=1, \end{gather}$$

(3.40)$$\begin{gather}s_{00}=2,\quad \mathrm{at}\ r\to\infty, \end{gather}$$

(3.41)$$\begin{gather}\rho_{00}=0,\quad \mathrm{at}\ r\to\infty, \end{gather}$$

(3.42)$$\begin{gather}\boldsymbol{\nabla}\phi_{00}=0,\quad\mathrm{at}\ r\to\infty. \end{gather}$$

The system of (3.34)–(3.40) have a trivial solution, i.e. $\phi _{00}=0$, $\rho _{00}=0$, $s_{00}=2$. Physically, the solution simply implies that the ion concentration is uniform because the particle is uncharged and there is no electric field.

Order $O(\epsilon )$: perturbation due to the electric field for an uncharged particle. The governing equations for the salt and charge dynamics, and electrostatic potential after substituting the expressions of $\rho _{00}$, $s_{00}$, and $\phi _{00}$ are given as

(3.43)$$\begin{gather} \nabla^2 s_{10}=0, \end{gather}$$

(3.44)$$\begin{gather}\nabla^2 \rho_{10} + 2 \nabla^2 \phi_{10}=0, \end{gather}$$

(3.45)$$\begin{gather}\nabla^2 \phi_{10} ={-}\frac{\kappa^2}{2} \rho_{10}. \end{gather}$$

The corresponding reduced boundary conditions when $r \to \infty$ are

(3.46)$$\begin{gather} s_{10} = 0, \end{gather}$$

(3.47)$$\begin{gather}\rho_{10}=0, \end{gather}$$

(3.48)$$\begin{gather}-\boldsymbol{\nabla}\phi_{10} = E_0 \boldsymbol{e}_z, \end{gather}$$

and when $r=1$, they read

(3.49)$$\begin{gather} \boldsymbol{e}_r\boldsymbol{\cdot}\boldsymbol{\nabla}s_{10}=0, \end{gather}$$

(3.50)$$\begin{gather}\boldsymbol{e}_r\boldsymbol{\cdot}\left[\boldsymbol{\nabla}\rho_{10}+2\boldsymbol{\nabla}\phi_{10}\right]=0. \end{gather}$$

(3.51)$$\begin{gather}\boldsymbol{e}_r \boldsymbol{\cdot} \boldsymbol{\nabla} \phi_{10}=0, \end{gather}$$

where the derivative of the potential is set to be zero to ensure that there is no excess charge on the surface, see discussion below (3.30). Mathematically, at $O(\epsilon )$, the charge on the particle surface is zero. Through Gauss's law, the no surface charge boundary condition necessitates that $\boldsymbol {e}_r \boldsymbol{\cdot}\boldsymbol {\nabla }\phi _{10} = 0$. Hence, (3.50) implies that $\boldsymbol {e}_r\boldsymbol{\cdot}\boldsymbol {\nabla }\rho _{10}=0$ at the particle surface. Along with (3.44) and (3.47), we obtain $\rho _{10}=0$. Similarly, (3.43), (3.46) and (3.49) reveal $s_{10}=0$.

Since $\rho _{10}=0$, $\phi _{10}$ is governed by the Laplace equation, and the solution with appropriate boundary conditions reads (Griffiths 2005)

(3.52)\begin{equation} \phi_{10}\left(r,\theta\right) ={-}E_0 r \left(1+\frac{1}{2r^3}\right) \cos \theta. \end{equation}

Physically, this order implies that the perturbed potential and corresponding electric field lines get modified due to the geometry of the particle but there is no charge and salt accumulation.

Order $O(\zeta )$: perturbation of a charged particle without an external electric field. The equations governing the charge, salt and potential at $O(\zeta )$ are analogous to (3.43)–(3.45) at $O(\epsilon )$ but have different boundary conditions. The governing equations are

(3.53)$$\begin{gather} \nabla^2 s_{01} = 0, \end{gather}$$

(3.54)$$\begin{gather}\nabla^2 \rho_{01} + 2 \nabla^2 \phi_{01} =0, \end{gather}$$

(3.55)$$\begin{gather}\nabla^2 \phi_{01} ={-}\frac{\kappa^2}{2}\rho_{01}. \end{gather}$$

As $r \to \infty$ the net charge, salt and electric potential gradient all decay to zero, or

(3.56)$$\begin{gather} s_{01} = 0, \end{gather}$$

(3.57)$$\begin{gather}\rho_{01} = 0, \end{gather}$$

(3.58)$$\begin{gather}\boldsymbol{\nabla} \phi_{01} = 0. \end{gather}$$

At the particle surface, $r=1$, we obtain

(3.59)$$\begin{gather} \boldsymbol{e}_r\boldsymbol{\cdot}\boldsymbol{\nabla}s_{01} = 0, \end{gather}$$

(3.60)$$\begin{gather}\boldsymbol{e}_r\boldsymbol{\cdot}\left[\boldsymbol{\nabla}\rho_{01} + 2\boldsymbol{\nabla}\phi_{01}\right]=0, \end{gather}$$

(3.61)$$\begin{gather}\phi_{01}=1. \end{gather}$$

Equations (3.53), (3.56) and (3.59) yield $s_{01} = 0$. The remainder of the equations reveal

(3.62)$$\begin{gather} \phi_{01}(r) = \frac{1}{r} \exp({-\kappa(r-1)}) , \end{gather}$$

(3.63)$$\begin{gather}\rho_{01}(r) ={-}\left(\frac{2}{r}\right)\exp({-\kappa(r-1)}) . \end{gather}$$

Physically, the results indicate the distribution of potential and charge arising due to a charged particle.

Order $O(\epsilon \zeta )$: perturbation of both imposed electric field and surface zeta potential. To fully resolve the body force to the order of $O(\epsilon \zeta )$, we have to obtain $s_{11}$, this can be observed by the expansion of (3.9) and collecting the respective orders. The governing equations and boundary conditions for salt, $s_{11}$, are given by

(3.64)$$\begin{gather} \nabla^2 s_{11} + \boldsymbol{\nabla}\boldsymbol{\cdot}\left(\rho_{01}\boldsymbol{\nabla} \phi_{10}\right) = 0, \end{gather}$$

(3.65)$$\begin{gather}s_{11} \to 0,\quad \mathrm{at}\ r \to \infty, \end{gather}$$

(3.66)$$\begin{gather}\boldsymbol{e}_r \boldsymbol{\cdot} \left[\boldsymbol{\nabla}s_{11} + \rho_{01}\boldsymbol{\nabla}\phi_{10}\right]=0\quad \mathrm{at}\ r=1. \end{gather}$$

Substituting in the expression of $\rho _{01}$ from (3.63) and $\phi _{10}$ from (3.52) we can solve for $s_{11}(r,\theta )$. The salt dynamics are in the form of $s_{11}=E_0\,f(r)\cos \theta$ where $f(r)$ is

(3.67) \begin{align} f(r) &= \frac{1}{r^2}\left(\frac{\kappa}{6}-\frac{5}{3}-\frac{10}{\kappa}-\frac{2}{\kappa^2}\right)+ \frac{\kappa^2}{3}\exp\left({\kappa}\right)\,{Ei} \left(-\kappa r\right)+\frac{2}{\kappa^2 r^2} \exp({-\kappa(r-1)})\nonumber\\ &\quad-\frac{4}{3\kappa r^2}\exp({-\kappa(r-1)})+\frac{ \kappa}{3} \exp({-\kappa(r-1)})\nonumber\\ &\quad -\frac{1}{3r}\exp({-\kappa(r-1)})+\frac{2}{\kappa r} \exp({-\kappa(r-1)}), \end{align}

where ${Ei}()$ is the elliptic integral.

Reintroducing dimensions: at this stage, we restore the dimensions and reintroduce the tilde for dimensionless variables. The total osmophoretic body force $\boldsymbol {b}$ is made dimensionless by writing $\boldsymbol {b} = ({k_B T c_0}/{a} )\tilde { \boldsymbol {b} } = ({ \varepsilon (k_B T )^2 \kappa ^2 }/{2 e^2 a^3} )\tilde {\boldsymbol {b}}$. The first relevant order of $\tilde {\boldsymbol {b}}$ for electrophoresis is $\epsilon \tilde {\zeta }$ because it is the order at which a charged particle is being driven by an electric field. To this end, we write

(3.68)\begin{equation} \tilde{\boldsymbol{b}} ={-} \epsilon \tilde{\zeta} \left[ \tilde{\rho}_{00} \tilde{\boldsymbol{\nabla}} \tilde{\phi}_{11} + \tilde{\rho}_{11} \tilde{\boldsymbol{\nabla}} \tilde{\phi}_{00} + \tilde{\rho}_{01} \tilde{\boldsymbol{\nabla}} \tilde{\phi}_{10} + \tilde{\rho}_{10} \tilde{\boldsymbol{\nabla}} \tilde{\phi}_{01} + \tilde{\boldsymbol{\nabla}} \tilde{s}_{11} \right]. \end{equation}

Based on the solutions at different orders, it is straightforward to see that $\boldsymbol {\tilde {b}}$ reduces to

(3.69)\begin{equation} \boldsymbol{\tilde{b}} ={-}\epsilon \tilde{\zeta} \left[ \tilde{\rho}_{01} \tilde{\boldsymbol{\nabla}} \tilde{\phi}_{10} + \tilde{\boldsymbol{\nabla}} \tilde{s}_{11} \right]. \end{equation}

We note that the body force term of $\tilde {s}_{11}$ integrates out to zero in the calculation of $\boldsymbol {U}$ and $\boldsymbol {\varOmega }$ for electrophoresis and is generally not included in prior analyses. However, we retain this term for consistency as it does become crucial for diffusiophoretic phenomena, as we detail in § 3.3.

After substituting the values of $\tilde {\rho }_{01}$, $\tilde {\phi }_{10}$ and $\tilde {s}_{11}$, the resultant equation in dimensional form is

(3.70)\begin{align} \boldsymbol{b}&={-}\frac{\epsilon \varepsilon E_0 \zeta \kappa^2} {2 a^2} \left\{\left[\frac{2}{\tilde{r}}\exp({-\kappa(\tilde{r}-1)})\left(1-\frac{1}{\tilde{r}^3} \right)\cos\theta+\frac{{\rm d} f(\tilde{r})}{{\rm d} \tilde{r}}\cos\theta\right]\boldsymbol{e}_r \right. \nonumber\\ &\quad \left.-\left[\frac{2}{\tilde{r}}\exp({-\kappa(\tilde{r}-1)}) \left(1+\frac{1}{2 \tilde{r}^3}\right)\sin\theta+\frac{f(\tilde{r})}{\tilde{r}}\sin \theta\right]\boldsymbol{e}_\theta \right\}. \end{align}

Equation (3.70) is substituted into (2.12) to obtain the translational velocity to be

(3.71)\begin{equation} \boldsymbol{U}=\mathcal{M} \boldsymbol{E}_{\infty}, \end{equation}

where the mobility $\mathcal {M}$ after restoring dimensions is

(3.72)\begin{equation} \mathcal{M}=\frac{\varepsilon\zeta}{6\mu}\left[(1+\kappa)+(12-\kappa^2)\displaystyle \int_1^\infty \frac{e^{\kappa (1-t)}}{t^5} {\rm d}t\right], \end{equation}

which is the seminal result of Henry (1931) for arbitrary double layer thickness, and has also been reported by Teubner (1982) and Kim & Karrila (2013). For a homogeneous sphere, our analysis reveals $\boldsymbol {\varOmega }=\boldsymbol {0}$, as expected. We emphasize that it is straightforward to extend the calculations to heterogeneous spheres (Velegol et al. 1996; Teubner 1982) and obtain results for electrorotation in arbitrary double-layer thicknesses, which otherwise requires considerable effort.

3.3. Electrolytic diffusiophoretic mobility at arbitrary interaction length scales

Next, we focus on the process of electrolytic diffusiophoresis in the Debye–Hückel limit and for arbitrary double-layer thickness. We assume that the external concentration gradient of a binary monovalent electrolyte is given as $\boldsymbol {\nabla } s_{\infty } = 2 \epsilon \boldsymbol {\nabla } c_0$, where $\epsilon$ is a small parameter, much like in § 3.2, see figure 2(b). Here, $\boldsymbol {b}$ is required to be expanded to an additional higher order of $\epsilon \tilde {\zeta }^2$. The term of order $O ( \epsilon \tilde {\zeta } )$ is identical to electrophoresis and represents the electrophoretic component of the diffusiophoretic mobility. The second term of order $O ( \epsilon \tilde {\zeta }^2 )$ denotes the chemiphoretic component. We employ the expression of $\boldsymbol {b}$, derived in this subsection, to (2.12)–(2.13). This allows us to retrieve the expression of the translation velocity of a charged spherical particle in an unbounded solution of a symmetrically charged electrolyte for an arbitrary double-layer thickness, which otherwise requires considerable efforts, see Keh & Wei (2000).

We acknowledge the electrokinetic equations used to describe such diffusiophoretic systems are analogous to our treatment of the motion of electrophoretically propelled particles in § 3.2. However, the key mechanistic difference is the presence of an external gradient of solute $\boldsymbol {\nabla } s_{\infty }$ instead of an imposed electric field $\boldsymbol {E}_\infty$; this results in a change of boundary conditions and subsequently the solutions at different asymptotic orders. To preserve the pedagogical nature of our manuscript, we will rederive the electrophoretic contribution and subsequently solve for the chemiphoretic contribution to the osmophoretic body force term and attempt to emphasize key physical and mathematical differences between the derivations laid out in §§ 3.2 and 3.3.

Consider a colloidal particle with a surface zeta potential, $\zeta$, in an external solute gradient of a symmetric binary electrolyte. We assume that the electrolytes are monovalent such that $z_\pm = \pm 1$. The ions are assumed to have different diffusivities $D_+ \neq D_-$. The governing equations of the concentration of the ionic species, $c_\pm$, are identical to (3.10)–(3.11) and the interaction potential is governed by Poisson's equation, as given in (3.14). However, the far-field boundary conditions are different. Specifically, as $r \to \infty$, we assume that the concentration of the ionic species is linear with position $z$, or

(3.73)\begin{equation} c_\pm= c_0\left(1+ \frac{\epsilon z}{a}\right). \end{equation}

Further, it is assumed that the electric current in the far-field is zero, which yields (Prieve et al. 1984; Velegol et al. 2016; Gupta et al. 2019)

(3.74)\begin{equation} - {\boldsymbol{\nabla}} \phi |_{r\to\infty} = \epsilon \beta \frac{k_B T}{a e} = \epsilon E_0 \boldsymbol{e}_z, \end{equation}

where $\beta = ({D_+ - D_-})/({D_+ + D_-})$ and $E_0 = \beta ({k_B T}/{a e})$. The electric field is thus induced due to unequal diffusivities and a non-zero salt gradient.

Similar to electrophoresis, modification of the governing equations in terms of charge, $\rho = e(c_+ - c_-)$, salt, $s=c_+ + c_-$, and potential, $\phi$, result in (3.12)–(3.14). In the far-field, the boundary conditions for $\phi$ and $\rho$ are identical to (3.15) and (3.16) with the aforementioned definition of $E_0$. However, the boundary condition of $s$ is modified to

(3.75)\begin{equation} s|_{r \to\infty } = 2c_0 \left( 1+\frac{\epsilon z}{a} \right). \end{equation}

Note that $\epsilon /a = \boldsymbol {\nabla } \log s_\infty$. The boundary conditions are identical at the particle surface, i.e. (3.18)–(3.20). The objective is to solve $\rho$, $s$ and $\phi$ with the modified boundary conditions above and subsequently evaluate the total osmophoretic body force, following the procedure used to obtain (3.9). We non-dimensionalize the equations using the same scales as (3.21a–f) and also define $\boldsymbol {b} = ({k_B T c_0}/{a} )\tilde {\boldsymbol {b}} = ({\varepsilon (k_B T)^2 \kappa ^2}/{2 e^2 a^3}) \tilde {\boldsymbol {b}}$, where the definition of $\kappa$ is also identical.

For simplicity, we drop the tilde from our analysis until (3.92) and reintroduce them afterward. Thus, the non-dimensional osmophoretic body force is given as $\boldsymbol {b} = -\rho \boldsymbol {\nabla }\phi -\boldsymbol {\nabla }s$. We expand $\rho$, $s$ and $\phi$ until $O(\epsilon \zeta ^2)$, and solve the equations at each order.

Order $O(1)$: an uncharged particle without any external salt gradient. The results at this order are identical to electrophoresis and thus yield $s_{00}=2$, $\rho _{00}=0$, and $\phi _{00}=0$, indicating a uniform concentration of ion with no charge and potential.

Order $O(\epsilon )$: perturbation of the external salt concentration to an uncharged particle. This order is distinct compared with electrophoresis since the far-field boundary condition for salt is different, while the remainder of the equations and boundary conditions are identical. We note that the boundary condition for the electric field is similar to electrophoresis since we have defined $E_0$. The solution simply reduces to zero and uniform charge density $\rho _{10}=0$, while both salt and potential follow the Laplace equation. The results read

(3.76)$$\begin{gather} s_{10}\left(r,\theta\right)=2r\left(1+\frac{1}{2r^3}\right)\cos\theta, \end{gather}$$

(3.77)$$\begin{gather}\phi_{10} ={-}E_0 r \left(1+\frac{1}{2r^3}\right)\cos\theta. \end{gather}$$

Physically, at this order, a gradient in the salt concentration far away perturbs the salt field and induces a potential field if diffusivity asymmetry is present ($E_0 \neq 0$ only when $\beta \neq 0$). However, since the surface is uncharged, $\rho _{10}=0$.

Order $O(\zeta )$: perturbation in the surface charge of the particle without an external field. Since there is no external field at this order, the solution is identical to electrophoresis with $s_{01}=0$ and

(3.78)$$\begin{gather} \phi_{01}(r) = \frac{1}{r}\exp({-\kappa(r-1)}), \end{gather}$$

(3.79)$$\begin{gather}\rho_{01}(r) ={-}\frac{2}{r} \exp({-\kappa(r-1)}). \end{gather}$$

This order represents the potential and charge profiles due to the surface charge of the particle. However, there is no salt accumulation at this order since the reduction in the coion concentration is balanced by the increase in the counter-ion concentration.

Order $O(\epsilon \zeta )$: perturbation in both the imposed salt concentration and surface charge. The governing equations for charge ($\rho _{11}$), salt ($s_{11}$), and potential ($\phi _{11}$) are

(3.80)$$\begin{gather} \nabla^2 s_{11} + \boldsymbol{\nabla}\boldsymbol{\cdot}\left[\rho_{01}\boldsymbol{\nabla}\phi_{10} \right]=0, \end{gather}$$

(3.81)$$\begin{gather}\nabla^2 \rho_{11}+\boldsymbol{\nabla}\boldsymbol{\cdot}\left[s_{10}\boldsymbol{\nabla}\phi_{01}+ 2 \boldsymbol{\nabla}\phi_{11}\right]=0, \end{gather}$$

(3.82)$$\begin{gather}\nabla^2 \phi_{11} ={-}\frac{\kappa^2}{2}\rho_{11}. \end{gather}$$

The boundary conditions at the particle surface, $r=1$, are

(3.83)$$\begin{gather} \boldsymbol{e}_r\boldsymbol{\cdot}\left[\boldsymbol{\nabla}s_{11}+\rho_{01}\boldsymbol{\nabla}\phi_{10}\right]=0, \end{gather}$$

(3.84)$$\begin{gather}\boldsymbol{e}_r\boldsymbol{\cdot}\left[\boldsymbol{\nabla}\rho_{11}+s_{10}\boldsymbol{\nabla}\phi_{01}+2\boldsymbol{\nabla}\phi_{11}\right]=0, \end{gather}$$

(3.85)$$\begin{gather}\boldsymbol{e}_r\boldsymbol{\cdot}\boldsymbol{\nabla} \phi_{11}=0. \end{gather}$$

Again, there is no external field as $r\to \infty$. Substituting in the solutions obtained in $O(\epsilon )$, and $O(\zeta )$ we begin to solve $s_{11}$, $\rho _{11}$ and $\phi _{11}$. Further, we can separate the $r$ and $\theta$ contributions by redefining $s_{11}(r, \theta ) = f_{s_{11}}(r)\cos \theta$, $\rho _{11}(r,\theta )=f_{\rho _{11}}(r)\cos \theta$, and $\phi _{11}(r, \theta )=f_{\phi _{11}}(r)\cos \theta$ and solve the equations numerically, see Appendix A. It is possible to find analytical solutions to (A1)–(A6), similar in form to (3.67). In the scope of our current paper, we choose to resolve the dynamics at $O(\epsilon \zeta )$ and $O(\epsilon \zeta ^2)$ numerically. For an analytical derivation of such higher-order effects the reader is directed to Keh & Wei (2000).

Order $O(\epsilon \zeta ^2)$: first-order perturbation in salt field and second-order perturbation in surface charge. We only seek to solve $s_{12}$ at this order since it is the only quantity required to resolve the body force up to $O(\epsilon \zeta ^2)$, see (3.91). The equation governing the dynamics of $s_{12}$ is

(3.86)\begin{equation} \nabla^2 s_{12} + \boldsymbol{\nabla}\boldsymbol{\cdot}\left[\rho_{01}\boldsymbol{\nabla}\phi_{11}+\rho_{11}\boldsymbol{\nabla}\phi_{01}\right]=0, \end{equation}

with

(3.87)\begin{equation} \boldsymbol{e}_r\boldsymbol{\cdot}\left[\boldsymbol{\nabla}s_{12}+\rho_{01}\boldsymbol{\nabla}\phi_{11}+\rho_{11} \boldsymbol{\nabla}\phi_{01}\right]_{r=1}=0, \end{equation}

and $s_{12}(r \rightarrow \infty, \theta ) =0$. We write $s_{12}$ as $s_{12}(r,\theta )=f_{s_{12}}(r)\cos \theta$ and solve the equations numerically, see Appendix A. As discussed previously, we solve (A7), (A8), and the far-field constraint numerically.

Restoring dimensions: we now restore dimensions and reintroduce tilde to describe dimensionless variables. Therefore, we write the body force $\boldsymbol {b} = ({\varepsilon (k_B T)^2 \kappa ^2}/{2 e^2 a^3} )\tilde {\boldsymbol {b}}$ such that

(3.88)$$\begin{gather} \tilde{\boldsymbol{b}} = \epsilon \tilde{\zeta} \boldsymbol{\tilde{b}}_{11} + \epsilon \tilde{\zeta}^2 \boldsymbol{\tilde{b}}_{12}, \end{gather}$$

(3.89)$$\begin{gather}\tilde{\boldsymbol{b}}_{11} ={-} \left( \tilde{\rho}_{01} \boldsymbol{\tilde{\nabla}}\tilde{\phi}_{10}+\boldsymbol{\tilde{\nabla}}s_{11} \right), \end{gather}$$

(3.90)$$\begin{gather}\tilde{\boldsymbol{b}}_{12} ={-} \left( \tilde{\rho}_{01}\boldsymbol{\tilde{\nabla}} \tilde{\phi}_{11}+\tilde{\rho}_{11}\boldsymbol{\tilde{\nabla}}\tilde{\phi}_{01}+\boldsymbol{\tilde{\nabla}}\tilde{s}_{12} \right). \end{gather}$$

At this point, some comments are in order. We note that $\tilde { \boldsymbol {b}}$ could also include a term at $O(\epsilon )$ since $s_{10} \neq 0$. However, $\tilde {\boldsymbol {b}}$ only includes excess osmotic pressure and not the osmotic pressure itself. This is because the osmotic pressure contribution due to $\boldsymbol {\nabla }s_\infty$ would lead to particle motion even in the absence of phoretic interactions. Only the terms at subsequent orders are included to ignore this effect. Further, we highlight that both $\tilde {\phi }_{10}$ and $\tilde {s}_{11}$ are proportional to $\tilde {E}_0$. Therefore, the $O (\epsilon \tilde {\zeta } )$ term only depends on $E_0$ and is referred to as the electrophoretic contribution. Since $\tilde {\rho }_{01}$, $\tilde {\phi }_{10}$ and $\tilde {s}_{11}$ are identical to electrophoretic solution, the $O(\epsilon \tilde {\zeta })$ is equal to the one described earlier in (3.72).

However, in contrast, for the $O(\epsilon \tilde {\zeta }^2)$ contribution, $\tilde {\rho }_{01}$, $\tilde {\phi }_{11}$, $\tilde {\rho }_{11}$, $\tilde {\phi }_{01}$ and $\tilde {s}_{12}$ are independent of $E_0$. Furthermore, $\tilde {\phi }_{11}$, $\tilde {\rho }_{11}$ and $\tilde {s}_{12}$ are all proportional to $\tilde {s}_{10}$, which is consistent with the prior literature (Anderson 1989; Keh & Wei 2000; Gupta et al. 2019).

Since $\tilde {\boldsymbol {b}}_{11}$ is identical to the electrophoretic motion, we focus our attention on $\tilde {\boldsymbol {b}}_{12}$, which reads

(3.91)\begin{equation} \tilde{\boldsymbol{b}}_{12}={-} \left(\rho_{01}\frac{{\rm d}f_{\phi_{11}}}{{\rm d}r}+f_{\rho_{11}}\frac{{\rm d}\phi_{01}}{{\rm d}r}+\frac{{\rm d}f_{s_{12}}}{{\rm d}r} \right) \cos\theta \boldsymbol{e}_r - \left(\frac{\rho_{01} f_{\phi_{11}}}{r}+\frac{f_{s_{12}}}{r}\right)\sin\theta \boldsymbol{e}_\theta . \end{equation}

Upon substituting the value of $\boldsymbol {b}$ in (2.12), the translation velocity of the particle could be simplified to read

(3.92)\begin{equation} \boldsymbol{U} = \mathcal{M} \boldsymbol{\nabla} \log s_{\infty}, \end{equation}

where

(3.93)\begin{equation} \mathcal{M} = \frac{\varepsilon}{\mu}\left[\frac{k_B T}{e}\beta \zeta \varTheta_1\left(\kappa\right) + \frac{\zeta^2}{8}\varTheta_2\left(\kappa\right)\right], \end{equation}

where $\varTheta _1$ and $\varTheta _2$ are evaluated numerically, see Appendix A. Figure 3 demonstrates good quantitative agreement between the values obtained in Keh & Wei (2000) and our results.

Figure 3. Comparison of proposed mobility expressions of diffusiophoretic mobility in (3.93) with the mobility reported in Keh & Wei (2000). Quantitative agreement of both (a) the electrophoretic component $\varTheta_1 (\kappa)$ and (b) the chemiphoretic component $\varTheta_2 (\kappa)$ is observed.

This section highlighted the generality of (2.12) and (2.13). We showed they are able to recover the mobilities for microswimmers in thin interaction layer limit, electrophoresis for arbitrary double-layer thickness and electrolytic diffusiophoresis for arbitrary double-layer thickness.

4. Autophoretic motion of microswimmers

In this section, we use the formula obtained in (2.12)–(2.13) to study the translation of Janus-like particles with a spherical cap, see figure 4(a). The key novelty of our analysis is that (2.12)–(2.13) do not impose the restriction on interaction length scale. As we show later, if the interaction length is comparable to particle size, the particle velocity is significantly impacted.

Figure 4. (a) Self-phoretic Janus particle where propulsion is controlled by the size of the spherical cap $\theta _0$, the reactive flux $J$ and the interaction length scale $\kappa ^{-1}$. (b) Here, $U$ versus $\theta _0$ for different $\kappa$ values demonstrates a maximum velocity for $\theta _0={{\rm \pi} }/{2}$ irrespective of $\kappa$. (c) Here, ${U\kappa ^2}/{U_0}$ versus $\kappa$ asymptotically approaches the thin limit calculations as $\kappa \to 100$. However, considerable dampening is observed even for $\kappa = O(10)$. The values are reported for $\theta _0={{\rm \pi} }/{2}$. (d) The dampening of ${U\kappa ^2}/{U_0}$ with $\kappa$ is observed for all $\theta _0$. The value $U_0$ is the asymptotic limit of $U \kappa ^2$ from the thin interaction layer calculations.

We define catalytic surface activity through a non-dimensional outward surface flux of strength, $J$ (scaled by a characteristic flux $D c_0 / a$, where $D$ is the diffusivity of the solute, $c_0$ is a reference concentration of the solute and $a$ is the particle radius). The non-dimensional interaction length is characterized by $\kappa ^{-1}$ (scaled by $a$). Lastly, the size of the catalytic cap is controlled by the polar angle $\theta _0$, such that $\theta _0=0$ indicates no catalytic cap on the particle and $\theta _0={{\rm \pi} }/{2}$ represents a hemispherical cap.

We assume a Helmholtz-like equation governs the interaction potential ($\phi$, scaled by $k_B T$) with a constant surface potential ($\phi _0$) and a far-field decay. We take this opportunity to highlight the choice of Helmholtz-like potential. While the potential is not representative of different surface interactions possible, it provides a convenient choice to explore the impact of $\kappa$ and thus has been chosen for this analysis. We note that our analysis can be easily extended to other interaction potentials provided that the integrals in (2.11) are convergent. For a detailed analysis of phoretic motion due to a general particle–solute interaction the reader is directed to Brady (2021).

To resolve particle translation for a given surface activity and interaction, we write

(4.1)$$\begin{gather} \nabla^2 \phi = \kappa^2 \phi, \end{gather}$$

(4.2)$$\begin{gather}\phi = \phi_0,\quad {\rm at}\ r=1, \end{gather}$$

(4.3)$$\begin{gather}\phi \to 0,\quad {\rm at}\ r \to \infty. \end{gather}$$

Solute transport is governed by diffusion and phoretic interactions with the following boundary conditions:

(4.4)$$\begin{gather} \boldsymbol{\nabla} \boldsymbol{\cdot} \left(\boldsymbol{\nabla} c + c \boldsymbol{\nabla} \phi \right) = 0, \end{gather}$$

(4.5)$$\begin{gather}-\boldsymbol{n}\boldsymbol{\cdot}\left(\boldsymbol{\nabla} c + c\boldsymbol{\nabla}\phi\right) = J,\quad {\rm at}\ r=1, \end{gather}$$

(4.6)$$\begin{gather}c\to0,\quad {\rm at}\ r\to\infty. \end{gather}$$

Our model problem is illustrated in figure 4(a). We solve the coupled equations (4.1)–(4.6) numerically. A scaled model geometry was constructed with the particle radius given to be $a=1$ and an outer radius of $r|_\infty =20$, representing the far-field. The interaction potential ($\phi$) and solute concentration ($c$) were defined by (4.1)–(4.6). To obtain the translation velocity through (2.12), we define the body force $\boldsymbol {b} = -c\boldsymbol {\nabla }\phi -\boldsymbol {\nabla } c$. Note that the arbitrary body force has both a phoretic and an osmotic contribution. Details of the computational method are provided in Appendix B.

We first analyse the effects of the spherical cap size ($\theta _0$) and interaction length scale ($\kappa ^{-1}$) for a fixed surface flux ($J=1$) and surface potential ($\phi _0=-1$), see figure 4. First, we note that the particle moves in the direction of the catalytic cap, since when $\phi _0=-1$ the particle is attracted towards regions of higher solute concentration. From figure 4(b), we observe that the propulsion velocity is maximum for $\theta _0={{\rm \pi} }/{2}$. This is in agreement with prior observations in the literature (Golestanian, Liverpool & Ajdari (2007), Michelin & Lauga (2014), Popescu et al. (2018) and others). Additionally, the dependence of translation velocity on the cap size is symmetric about $\theta _0={{\rm \pi} }/{2}$. As described in Michelin & Lauga (2014), when $\theta _0={{\rm \pi} }/{2}$, for small Péclet numbers, the sharpest concentration gradients are located near the equator, consequently leading to a larger slip velocity over an extended inert surface. In contrast, for smaller or larger catalytic caps, the aforementioned solute front is closer to the pole and thus involves a smaller share of the particle surface in generating slip velocities and therefore a smaller swimming speed.

Thin interaction layer: we also study the dependence of propulsion velocity on the interaction length scale ($\kappa ^{-1}$), see figure 4(c). To compare our calculations with the thin interaction layer limit, we perform analogous calculations following the approach of Anderson (1989) and Derjaguin et al. (1947). We consider diffusive transport of solute $\nabla ^2 c = 0$ through the fluid volume. The surface flux condition is given by $-\boldsymbol {n}\boldsymbol{\cdot}\boldsymbol {\nabla }c=J$ with a far-field decay condition. The phoretic slip at the particle surface is defined as

(4.7)\begin{equation} \boldsymbol{u}_{slip} ={-} \boldsymbol{\nabla}_s c \displaystyle \int_1^\infty (r-1)\left[\exp\left(-\phi\right)-1\right] {\rm d}r, \end{equation}

where the non-dimensional $\boldsymbol {u}_{slip}$ is scaled with $k T c_0 a / \mu$. After that, we refer to Stone & Samuel (1996), also derived in (3.7), to obtain the translation velocity $\boldsymbol {U}_{thin}$ to be

(4.8)\begin{equation} \boldsymbol{U}_{thin} ={-}\frac{1}{4{\rm \pi}} \displaystyle \int_S \boldsymbol{u}_{slip} \,{\rm d}S. \end{equation}

It is well known in the diffusiophoretic literature (Golestanian 2019) that for $\kappa \gg 1$, $|\boldsymbol {U}_{thin}| \propto 1/ \kappa _{thin}^2$. Consequently, we introduce a scaled velocity expression, $\boldsymbol {U}_0 = \boldsymbol {U}_{thin} \kappa ^2_{thin}$ which becomes constant as $\kappa \to \infty$. As shown in figure 4(c), we observe that the velocity ratio $U\kappa ^2/U_0 = |\boldsymbol {U}|\kappa ^2/|\boldsymbol {U}_0|$ significantly decreases when the interaction limit becomes comparable to particle size, i.e. $\kappa \rightarrow 1$. We find that even for $\kappa = O(10)$, the velocity is reduced by almost a factor of two. We only reach the thin interaction limit for $\kappa = O(10^2)$. This observation is consistent with the passive diffusiophoretic literature where even moderately thin double layers can significantly reduce the diffusiophoretic velocity (Prieve et al. 1984), see figure 3(b). Our analysis highlights that even for autophoretic swimmers, this effect can be observed, and using a thin interaction limit could overestimate the velocity for moderately thin interaction thickness such as $\kappa \rightarrow 50$. Over the past decade, there has been an increasing interest in nanoparticles or very dilute systems where relative interaction length scales are large (Leunissen et al. 2007; Shin et al. 2016; Gupta et al. 2020b; Wu, Greydanus & Schwartz 2021; Shi, Wu & Schwartz 2023). Thus, these results could be crucial for future experimental studies. To ensure that our trends of reductions in velocity are consistent for other conditions, in figure 4(d), we observe that reduction in velocity for smaller $\kappa$ values is consistent for all $\theta _0$ values.

5. Autophoretic swimmers with external solute gradients

To further demonstrate the applicability of (2.12) in scenarios where it is difficult to use a slip velocity approach, we modify our model problem to include an external solute flux in the fluid bulk. The model geometry is preserved as shown in figure 5(a). The concentration in the far-field is $\boldsymbol {\nabla } c_\infty = J_{ext} \boldsymbol {e}_z$, where $\boldsymbol {e}_z$ is the $z$-coordinate basis vector in the universal Cartesian frame of reference and $J_{ext}$ is a free parameter used to control the direction and strength of this external field. The magnitude of the external gradient is thus $\left |{\boldsymbol{\nabla}}c_{\infty} \right | = J_{ext}$. We note that all additional parameters have been appropriately non-dimensionalized as per the discussion in § 4. We analyse the variations of the local surface flux ($J$), external solute flux ($J_{ext}$) and interaction length scale $\kappa$. For a given surface flux $J=1$, it is observed that the translation velocity is linear with the imposed solute flux, $\boldsymbol {\nabla } c_\infty$. The particle moves with a velocity $\boldsymbol {U} = U_z \boldsymbol {e}_z$. The direction of the propulsion is governed by the relative magnitudes of different contributions of the body force terms. In the absence of an external flux, the body force contributions due to the phoretic activity cause the particle to move in the positive $z$-direction. In our analysis, we have ignored the excess osmotic contribution ($\boldsymbol {\nabla } c$) from the external solute gradient and only considered the phoretic contribution. This is done to ignore the motion of the particle only due to $\boldsymbol {\nabla } c_\infty$ when $\phi =0$. For $J_{ext}>0$, one can observe that force terms arising from the gradients in active and passive concentration fields are in the same direction, hence $U_z >0$ as in figure 5(b). Alternatively when $J_{ext}\,{<}\,0$, the phoretic force terms arising out of the passive solute concentration field compete with forcing arising due to activity. This leads to a direction reversal below a threshold concentration gradient, $(\boldsymbol {\nabla } c_\infty )_{crit}$, where $U_z < 0$. Here, $(\boldsymbol {\nabla } c_\infty )_{crit}$ is the external solute gradient necessary to arrest particle motion due to activity. The linearity of the results in figure 5(b) is due to the computations being performed in the weak field limit. To inspect the effects of interaction potential length, we obtained the $(\boldsymbol {\nabla } c_\infty )_{crit}$ values at $J=1$ while varying over $\kappa$. A qualitative agreement is observed between the results in figure 5(c) and the $\kappa$-dependence of velocity ($U$). Further, we see that a stronger external flux is necessary to arrest motion as we approach the thin interaction length limit. This effect arises because, for $\kappa \to \infty$, the catalytically ejected solute decays more sharply, $\boldsymbol {\nabla } c$ increases locally near the particle surface. This necessitates the need for a larger $(\boldsymbol {\nabla } c_\infty )_{crit}$ to counteract activity-induced gradients to arrest propulsion.

Figure 5. (a) A self-phoretic particle also driven by external $\boldsymbol {\nabla } c_\infty$. (b) Here, $U_z$ versus $\boldsymbol {\nabla } c_\infty$ shows that both the particle speed and direction depend on the competing effects of external and self-propulsion modes. We define $(\boldsymbol {\nabla } c_\infty )_{crit}$ as the value when the particle motion was arrested despite the presence of concentration gradients due to surface activities. (c) Here, $(\boldsymbol {\nabla } c_\infty )_{crit}$ versus $\kappa$. We observe an increase in the magnitude of the external flux needed to arrest motion as we approach thin interaction length limits.

6. Conclusion

The main result of our work is that (2.12)–(2.13), reported in the prior literature for phoretic systems (Brady 2021) and different physical set-ups (Brenner & Cox 1963; Hinch 1991; Leal 2007; Datt et al. 2015; Elfring & Goyal 2016; Datt et al. 2017), can retrieve the mobility expressions for both electrophoretic and electrolytic diffusiophoretic motion of colloidal spherical particles at arbitrary interaction lengths at the low potential limit. Additionally, the asymptotic limit of (2.12)–(2.13) for thin interaction length recovers net translation and rotation of the particle in terms of a phoretic slip velocity as in Stone & Samuel (1996). Finally, we employ these mobility relationships for self-propulsion of spherical microswimmers where we observe a peak in translational velocity as we approach a hemispherical catalytic coverage with a velocity reduction for lower or higher coverages. At moderate interaction lengths, $\kappa = O(10)$, a dampening in the translation velocity is observed. However, for higher $\kappa$ values our calculations retrieve the thin interaction length limit results arising from equivalent slip velocity calculations. Further, we add an external passive solute concentration gradient to our problem to understand the competing effects of surface-generated and externally imposed solute concentration gradients. We find configurations where the external flux arrest motion induced by surface activity $(\boldsymbol {\nabla } c_\infty )_{crit}$ and regions where the propulsion induced by the passive external solute gradient aids with or competes against the propulsion induced due to surface activity. A nonlinear decay of $(\boldsymbol {\nabla } c_\infty )_{crit}$ with $\kappa$ is observed from our analysis which is qualitatively similar to the $\kappa$-dependence of the translation velocity, $|\boldsymbol {U}|$.

Beyond the result described in this paper, our work can be utilized to predict $\boldsymbol {U}$ and $\boldsymbol {\varOmega }$ for an arbitrarily shaped particle. For such a calculation, one would need the appropriate expressions of the mobility ($\boldsymbol{\mathsf{M}}$) and the disturbance tensors ($\boldsymbol{\mathsf{D}}_T$, $\boldsymbol{\mathsf{D}}_R$), which might be possible to obtain analytically or numerically. While this manuscript focuses on the mobility of a single particle, a similar analysis could be extended to multiple particles and particles under confinement.

The derived expression will be particularly useful for multiphysics propulsion. For instance, it might be possible to induce propulsion of particles using multiple modes such as a combination of electric fields and other fields, such as electrodiffusiophoresis (Wang, Behdani & Silvera Batista 2022; Jarvey, Henrique & Gupta 2023). Another possibility is the inclusion of chemical kinetics at the particle surface (Davis & Yariv 2022) which will modify the solute problem and thus consequently change $\boldsymbol {b}$.