Fluctuating hydrodynamics and the Brownian motion of an active colloid near a wall

TL;DR

This paper systematically derives the fluctuating hydrodynamics of an active colloid near a wall, revealing activity-induced breakdown of detailed balance and providing tools for Brownian dynamics simulations.

Contribution

It introduces a systematic derivation of passive, active, and Brownian traction contributions for active colloids, including fluctuation-dissipation relations and a regularization method for simulations.

Findings

Derived exact relations between Brownian traction and thermal forces.

Showed the stationary distribution is non-Gibbsian near a wall.

Provided a regularization for friction tensors ensuring positive variances.

Abstract

The traction on the surface of a spherical active colloid in a thermally fluctuating Stokesian fluid contains passive, active, and Brownian contributions. Here we derive these three parts systematically, by "projecting out" the fluid using the boundary-domain integral representation of slow viscous flow. We find an exact relation between the statistics of the Brownian traction and the thermal forces in the fluid and derive, thereby, fluctuation-dissipation relations for every irreducible tensorial harmonic traction mode. The first two modes give the Brownian force and torque, from which we construct the Langevin and Smoluchowski equations for the position and orientation of the colloid. We emphasize the activity-induced breakdown of detailed balance and provide a prescription for computing the configuration-dependent variances of the Brownian force and torque. We apply these general results to an active colloid near a plane wall, the simplest geometry with configuration-dependent variances, and show that the stationary distribution is non-Gibbsian. We derive a regularization of the translational and rotational friction tensors, necessary for Brownian dynamics simulations, that ensures positive variances at all distances from the wall. The many-body generalization of these results is indicated.

Tables

Table 1. Table 1: Expressions for the 3 × 3 3 3 3\times 3 friction tensors at a height h ℎ h from a plane wall, using the result of the first Jacobi iteration Singh and Adhikari ( 2016a ) for h > b ℎ 𝑏 h>b and its regularized form for h ≤ b ℎ 𝑏 h\leq b . Here ℱ 1 = 1 + b 2 6 ​ ∇ 𝐑 2 superscript ℱ 1 1 superscript 𝑏 2 6 superscript subscript ∇ 𝐑 2 \mathcal{F}^{1}=1+\tfrac{b^{2}}{6}\nabla_{{\scriptscriptstyle\mathbf{R}}}^{2} is an operator encoding the finite size of the colloid.

$h>b$ (Jacobi iteration)	$h\le b$ (Regularized)
${𝜸}^{TT}={\gamma }_0^T\left(𝑰-{\gamma }_0^T{ℱ}^1{ℱ}^1{𝐆}^{\ast }\right)$	${𝜸}^{TT}={\gamma }_0^T\left[\left(1+\frac{9r^{\ast }}{32b}\right)𝑰-\frac{3{𝐫}^{\ast }{𝐫}^{\ast }}{32b{r}^{\ast }}\right]$
${𝜸}^{RR}={\gamma }_0^R\left(𝑰-\frac{{\gamma }_0^R}{4}{\nabla }_{𝐑}\times {\nabla }_{𝐑}\times {𝐆}^{\ast }\right)$	${𝜸}^{RR}={\gamma }_0^R\left[\left(1+\frac{27r^{\ast }}{32b}-\frac{27r^{{\ast }^3}}{64b^3}\right)𝑰-\left(\frac{9r}{32b}-\frac{3r^{{\ast }^3}}{64b^3}\right)\frac{{𝐫}^{\ast }{𝐫}^{\ast }}{r^{\ast 2}}\right]$
${𝜸}^{RT}=-\frac{{\gamma }_0^R{\gamma }_0^T}{2}{\nabla }_{𝐑}\times {𝐆}^{\ast }$	${𝜸}^{RT}=-{\gamma }_0^R\left(\frac{2}{b^2}-\frac{3r^{\ast }}{4b^3}\right)𝜺\cdot {𝐫}^{\ast }$

Equations

$$\int \mathbf{f}dS+{\mathbf{F}}^P=0,$$

$$\nabla \cdot \mathbf{u}=0,\nabla \cdot \mathbf{\sigma }+\mathbf{\xi }=0,$$

$$\mathbf{u}=\text{v},$$

$$⟨\int {\mathbf{u}}^{\mathcal{D}}(\mathbf{r})\cdot \mathbf{\xi }(\mathbf{r},t)dV\int {\mathbf{u}}^{\mathcal{D}}({\mathbf{r}}^{\prime })\cdot \mathbf{\xi }({\mathbf{r}}^{\prime },t^{\prime })d{V}^{\prime }⟩=2k_{B}T\dot{\mathcal{E}}({\mathbf{u}}^{\mathcal{D}})\delta (t-t^{\prime }),$$

$$\dot{\mathcal{E}}({\mathbf{u}}^{\mathcal{D}})=\eta \int \left[\nabla {\mathbf{u}}^{\mathcal{D}}+\left(\nabla {\mathbf{u}}^{\mathcal{D}}\right){}^T\right]{}^{2}dV.$$

$$\frac{1}{2}{v}_{\alpha }(\mathbf{r})=$$

$${\nabla }_{\alpha }{G}_{\alpha \beta }(\mathbf{r},{\mathbf{r}}^{\prime })=0,-{\nabla }_{\alpha }{P}_{\beta }(\mathbf{r},{\mathbf{r}}^{\prime })+\eta {\nabla }^2{G}_{\alpha \beta }(\mathbf{r},{\mathbf{r}}^{\prime })=-\delta \left(\mathbf{r}-{\mathbf{r}}^{\prime }\right){\delta }_{\alpha \beta },$$

$$K_{\alpha \beta \gamma }(\mathbf{r},{\mathbf{r}}^{\prime })=-{\delta }_{\alpha \gamma }{P}_{\beta }(\mathbf{r},{\mathbf{r}}^{\prime })+\eta \left({\nabla }_{\gamma }{G}_{\alpha \beta }(\mathbf{r},{\mathbf{r}}^{\prime })+{\nabla }_{\alpha }{G}_{\beta \gamma }(\mathbf{r},{\mathbf{r}}^{\prime })\right),$$

$$\mathbf{G}\cdot \mathbf{f}=\int \mathbf{G}(\mathbf{r},{\mathbf{r}}^{\prime })\cdot \mathbf{f}({\mathbf{r}}^{\prime })d{\mathrm{S}}^{\prime },\mathbf{K}\cdot \mathbf{v}=\int {\hat{\mathbf{\rho }}}^{\prime }\cdot \mathbf{K}(\mathbf{r},{\mathbf{r}}^{\prime })\cdot \mathbf{v}({\mathbf{r}}^{\prime })d{\mathrm{S}}^{\prime },$$

$$\mathbf{w}=\int \mathbf{G}(\mathbf{r},{\mathbf{r}}^{\prime })\cdot \mathbf{\xi }({\mathbf{r}}^{\prime })d{V}^{\prime },$$

$$\mathbf{f}=$$

$${\mathbf{v}}^{\mathcal{A}}=\sum _{l=1}^{\infty }\frac{1}{(l-1)!(2l-3)!!}{\mathbf{V}}^{(l)}\cdot {\mathbf{Y}}^{(l-1)}(\hat{\mathbf{\rho }}).$$

$$\mathbf{f}=\sum _{l=1}^{\infty }\frac{2l-1}{4\pi b^2}{\mathbf{F}}^{(l)}\cdot {\mathbf{Y}}^{(l-1)}(\hat{\mathbf{\rho }}).$$

$$\text{∫fd}S={\mathbf{F}}^{(1s)},\int \mathbf{\rho }\times \mathbf{f}dS=b{\mathbf{F}}^{(2a)}.$$

$$\mathbf{f}={\mathbf{f}}^{\mathcal{D}}+{\mathbf{f}}^{\mathcal{A}}+\hat{\mathbf{f}},$$

$$V_{\alpha }+{ϵ}_{\alpha \beta \gamma }{\Omega }_{\beta }{\rho }_{\gamma }$$

$$\frac{1}{2}{v}_{\alpha }^{\mathcal{A}}(\mathbf{r})=-\int G_{\alpha \beta }(\mathbf{r},{\mathbf{r}}^{\prime })f_{\beta }^{\mathcal{A}}({\mathbf{r}}^{\prime })d{\mathrm{S}}^{\prime }$$

$$0=\int G_{\alpha \beta }(\mathbf{r},{\mathbf{r}}^{\prime }){\xi }_{\beta }({\mathbf{r}}^{\prime })d{V}^{\prime }$$

$$-{\mathbf{G}}^{(l,l^{\prime })}(\mathbf{R})\cdot {\mathbf{F}}^{\mathcal{D}(l^{\prime })}={\mathbf{V}}^{\mathcal{D}(l)}.$$

$${\mathbf{G}}^{(l,l^{\prime })}(\mathbf{R})=\frac{2l-1}{4\pi b^2}\frac{2l^{\prime }-1}{4\pi b^2}\int {\mathbf{Y}}^{(l-1)}(\hat{\mathbf{\rho }}){\mathbf{G}}^{}(\mathbf{R}+\mathbf{\rho },\mathbf{R}+{\mathbf{\rho }}^{\prime }){\mathbf{Y}}^{(l^{\prime }-1)}({\hat{\mathbf{\rho }}}^{\prime })d\mathrm{S}d{\mathrm{S}}^{\prime },$$

$${\mathbf{F}}^{\mathcal{D}(l\sigma )}=$$

$${\mathbf{\gamma }}^{(l\sigma ,1s)}={\mathbf{P}}^{(l\sigma )}\cdot {\left[{\mathbf{G}}^{-1}(\mathbf{R})\right]}^{(l,1)}\cdot {\mathbf{P}}^{(1s)},{\mathbf{\gamma }}^{(l\sigma ,2a)}={\mathbf{P}}^{(l\sigma )}\cdot {\left[{\mathbf{G}}^{-1}(\mathbf{R})\right]}^{(l,2)}\cdot {\mathbf{P}}^{(2a)}.$$

$${\mathbf{F}}^{\mathcal{D}}=-{\mathbf{\gamma }}^{TT}\cdot \mathbf{V}-{\mathbf{\gamma }}^{TR}\cdot \mathbf{\Omega },$$

$${\mathbf{T}}^{\mathcal{D}}=-{\mathbf{\gamma }}^{RT}\cdot \mathbf{V}-{\mathbf{\gamma }}^{RR}\cdot \mathbf{\Omega }.$$

$${\mathbf{G}}^{(l)}(\mathbf{r},\mathbf{R})=\frac{2l-1}{4\pi b^2}\int {\mathbf{G}}^{}(\mathbf{r},\mathbf{R}+\mathbf{\rho }){\mathbf{Y}}^{(l-1)}(\hat{\mathbf{\rho }})d\mathrm{S},$$

$$-{\mathbf{G}}^{(l,l^{\prime })}(\mathbf{R})\cdot {\mathbf{F}}^{\mathcal{A}(l^{\prime })}+{\mathbf{K}}^{(l,l^{\prime })}(\mathbf{R})\cdot {\mathbf{V}}^{(l^{\prime })}=\frac{1}{2}{\mathbf{V}}^{(l)},$$

$${\mathbf{K}}^{(l,l^{\prime })}(\mathbf{R})=\frac{2l-1}{4\pi b^2}\frac{1}{(l^{\prime }-1)!(2l^{\prime }-3)!!}\int {\mathbf{Y}}^{(l-1)}(\hat{\mathbf{\rho }})\mathbf{K}(\mathbf{R}+\mathbf{\rho },\mathbf{R}+{\mathbf{\rho }}^{\prime })\cdot {\hat{\mathbf{\rho }}}^{\prime }{\mathbf{Y}}^{(l^{\prime }-1)}({\hat{\mathbf{\rho }}}^{\prime })d\mathrm{S}d{\mathrm{S}}^{\prime }.$$

$$4\pi a^2{\mathbf{V}}^{\mathcal{A}}$$

$${\mathbf{F}}^{{}^{\mathcal{A}}(l\sigma )}=-$$

$${\mathbf{\gamma }}^{(l\sigma ,l^{\prime }{\sigma }^{\prime })}={\mathbf{P}}^{(l\sigma )}\cdot {\left[{\mathbf{G}}^{-1}\cdot \left(\frac{1}{2}\mathbf{I}-\mathbf{K}\right)\right]}^{(l,l^{\prime })}\cdot {\mathbf{P}}^{(l^{\prime }{\sigma }^{\prime })}.$$

$${\mathbf{F}}^{\mathcal{A}}=-\sum _{l^{\prime }{\sigma }^{\prime }=1s}^{\infty }{\mathbf{\gamma }}^{(T,l^{\prime }{\sigma }^{\prime })}\cdot {\mathbf{V}}^{(l^{\prime }{\sigma }^{\prime })},$$

$${\mathbf{T}}^{\mathcal{A}}=-\sum _{l^{\prime }{\sigma }^{\prime }=1s}^{\infty }{\mathbf{\gamma }}^{(R,l^{\prime }{\sigma }^{\prime })}\cdot {\mathbf{V}}^{(l^{\prime }{\sigma }^{\prime })}.$$

$$\mathbf{G}{}^{(l,l^{\prime })}(\mathbf{R})\cdot {\hat{\mathbf{F}}}^{(l^{\prime })}(t)$$

Full text

DOI: 10.1080/17797179.2017.1294829

I Introduction

There has been a renewal of interest in the study of colloids with “active” boundaries, on which the usual no-slip boundary condition is replaced by one involving a “slip velocity”. This slip is the macroscopic manifestation of microscopic non-equilibrium processes in a thin boundary layer surrounding the colloid. Classic examples of slip-driven motion are the multitude of phoretic phenomena including electro-, thermo- and diffusio-phoresis and the motion of ciliated microorganisms Anderson (1989); Ebbens and Howse (2010). More recently, the slip model has been adapted to provide an effective description of flagellated microorganisms Ghose and Adhikari (2014). It provides a very general framework for the dynamics of phenomena where colloidal motion occurs without external influence.

The traction, that is the force per unit area, on the surface of an active colloid has three components: the Stokes drag proportional to the rigid body motion, the active thrust proportional to the slip velocity and the Brownian stress proportional to the temperature and the viscosity of the fluid. When the inertia of the colloid is negligible, its rigid body motion is obtained by setting the net force and net torque due to these tractions to zero. A central problem, then, is to derive explicit expressions for each part of the traction and to obtain, thereby, the force, the torque and its remaining moments.

At low Reynolds number and at finite temperature the Cauchy stress in the fluid obeys the fluctuating Stokes equation, a linear stochastic partial differential equation containing zero-mean Gaussian random fluxes with variances determined by the fluctuation-dissipation relation. These represent thermal fluctuations in the fluid. The solenoidal fluid velocity obeys the slip boundary condition at the colloid-fluid boundary. For a given rigid body motion and slip, the solution of the fluctuating Stokes equation provides the Cauchy stress in the fluid and, hence, the traction on the boundary. This solution for the traction only contains the boundary condition and the random fluxes; the fluid is, therefore, “projected out”. The Brownian forces and torques on the colloid are the first and second antisymmetric moments of the stochastic part of the traction. Their variance, by linearity of the governing equations, is proportional to the variance of the random fluxes. The Langevin equations for the position and orientation of the colloid follow straightforwardly by inserting the expressions for the net force and torque in the corresponding Newton’s equation.

Fox and Uhlenbeck were the first to derive the Langevin equation for the position of a passive spherical colloid from the fluctuating hydrodynamic equations for the fluid Fox and U. (1970). The fluid was taken to satisfy the no-slip boundary condition on the surface of the colloid and to be quiescent at the remote boundaries. The Lorentz reciprocal identity was used to relate the deterministic (Stokes) and stochastic (Einstein) parts of the force and to derive, thereby, the fluctuation-dissipation relation for the Brownian force from that of the random fluxes in the fluid. This approach was extended by several authors to include fluid inertia, particle inertia, Brownian fluxes at the colloid-fluid boundary and to many colloidal particles Hauge and Martin-Löf (1973); Bedeaux and Mazur (1974); Beenakker and Mazur (1983); Noetinger (1990); Roux (1992). Zwanzig, in an earlier contribution, derived the variance of the Brownian force on a spherical colloid in an unbounded fluid using the Faxén relation. The use of this variance in the Green-Kubo relation for the transport coefficient recovered the correct value of the Stokes friction Zwanzig (1964).

In this contribution, we derive the traction on the surface of an $active$* *colloidal particle near a plane wall by projecting out the fluid degrees of freedom. Our derivation differs in three important ways from previous work. First, we derive the complete distribution of the Brownian traction, and not just its first two moments as has been customary. This provides, in addition to the Brownian force and torque, the Brownian stresslet, a quantity important in suspension rheology. Second, the no-slip boundary condition at the colloid-fluid boundary is replaced by the slip boundary condition. This introduces an additional uncompensated source of dissipation and is the source, as we shall see, of the breakdown of detailed balance in the Langevin equations. Third, a no-slip boundary condition is introduced at the plane wall. This results in Stokes drags and active thrusts that depend on the distance of the colloid from the wall. The configuration-dependent friction requires, by the fluctuation-dissipation relation, configuration-dependent Brownian forces and torques and necessitates a “prescription” to render the Langevin equations mathematically meaningful. The solution to this old chestnut, the so-called Itô-Stratonovich dilemma, has been provided by multiple authors on multiple occasions but tends to be forgotten Lau and Lubensky (2007); Volpe et al. (2010); Volpe and Wehr (2016); Mannella and McClintock (2011). The procedure is to adiabatically eliminate the momentum, considered as a fast variable, from underdamped Langevin equations (where no such dilemma exists). The interpretation of the Brownian forces and torques in the overdamped Langevin equations is then unambiguous, though not necessarily conforming to either the Itô or Stratonovich prescriptions van Kampen (1981); Gardiner (1985); van Kampen (1992); Klimontovich (1990, 1994). Due to the linearity of the governing equations, this method of imputing meaning to configuration-dependent Brownian forces and torques remains valid in the presence of activity, as we show below. To the best of our knowledge, ours is the first systematic derivation of the Langevin equations for an active (and, as a special case, of a passive) colloid when the Stokes friction is configuration-dependent.

The remainder of the paper is organized as follows. In section II, we transform the fluctuating Stokes equation to its boundary-domain integral representation and hence obtain a linear integral equation for the traction on the colloid-fluid boundary. A formal solution, expressed in terms of the inverse of the single-layer operator of the integral equation (see below), clearly identifies the passive, active and Brownian parts of the traction. In section III, we consider a spherical colloid and provide an explicit solution to the boundary-domain integral equation. We derive, through a Galerkin discretization, an equivalent linear algebraic system for the coefficients of the expansion of the traction in a complete orthogonal basis of tensorial spherical harmonics. The solution of the linear system shows that the each tensorial coefficient of the Brownian traction is a zero-mean Gaussian random variable and provides their variances in terms of the variance of the stochastic term in the fluctuating Stokes equation. Variances of the Brownian force, torque and stresslet follow immediately. In section IV, we use the previous results to construct the overdamped Langevin equations for the position and orientation of the colloid. The Smoluchowski equations, corresponding to the prescription for the Brownian forces and torques provided by the adiabatic elimination procedure, are then presented and the activity-induced breakdown of the fluctuation-dissipation relation is pointed out. In section V, we specialize to the case of Stokes flow bounded by a plane wall and provide a leading order solution in terms of the Green’s function first identified by Lorentz and subsequently derived systematically by Blake Blake (1971). We show that the usual barometric distribution of height is no longer the stationary solution, a consequence of the breakdown of detailed balance. We also provide a regularization of the friction tensors, for what would correspond to heights in which the sphere overlaps with the wall. This ensures positive variances and is a necessary ingredient in Brownian dynamics simulations of spheres without hard steric potentials. Finally, the extension of the above results to many active colloids is indicated.

II Boundary-domain integral representation of fluctuating Stokes flow

We consider, in this section, the motion of an active colloid of arbitrary shape in an incompressible fluid of viscosity $\eta$ at a temperature $k_{B}T$. The boundary condition induces a local force per unit area $\boldsymbol{f}=\mathbf{\hat{\boldsymbol{\rho}}}\cdot\boldsymbol{\sigma}$ on the colloid boundary, where $\boldsymbol{\sigma}$ is the Cauchy stress in the fluid and $\mathbf{\boldsymbol{\hat{\rho}}}$ is the unit surface normal Landau and Lifshitz (1959). We shall henceforth refer to $\boldsymbol{f}$ as the traction. In addition, the colloid may also be acted upon by body forces $\mathbf{F}^{P}$ and body torques $\mathbf{T}^{P}$. In the absence of inertia of both particle and fluid, Newton’s equations for the colloid reduces to instantaneous balance of the surface forces (and torques) with the body forces (and torques):

[TABLE]

These are overdamped Langevin equations when $\boldsymbol{f}$ contains Brownian contributions.

To obtain the traction $\boldsymbol{f}$ it is necessary to know the flow field $\mathbf{u}$. At low-Reynolds number this satisfies the Stokes equation Happel and Brenner (1965)

[TABLE]

where $\bm{\sigma}=-p\boldsymbol{I}+\eta(\bm{\nabla}\mathbf{u}+(\bm{\nabla}\mathbf{u})^{T})$ is the Cauchy stress, $p$ is the fluid pressure, $\boldsymbol{\xi}$ is the thermal force acting on the fluid, $\boldsymbol{v}$ is the boundary velocity and $S$ is the surface of the colloid and $V$ is the domain of flow. The thermal force is a zero-mean Gaussian random field whose variance is given by the fluctuation-dissipation relation

[TABLE]

where $\mathbf{u}^{\mathcal{D}}$ is any flow field that satisfies the rigid body boundary condition on $S$ and $\dot{\mathcal{E}}$, the rate of dissipation of the fluid kinetic energy due to rigid body motion, is given by

[TABLE]

This manner of describing the thermally fluctuating fluid, due to Hauge and Martin-Löf Hauge and Martin-Löf (1973), is specially suited for flows with boundaries. The addition of a random flux, the more conventional manner of description first introduced by Landau and Lifshitz Landau and Lifshitz (1959), contains ambiguities in the presence of boundaries and is best used, therefore, when the fluid is unbounded in all directions Hauge and Martin-Löf (1973).

The key property of the above problem, which makes it possible to eliminate the fluid degrees of freedom exactly, is linearity. It is most clearly expressed in terms of the boundary-domain integral representation of slow viscous flow which provides the fluid flow in the bulk in terms of the boundary condition and the thermal force. The traction $\boldsymbol{f}$ is obtained as a solution of the boundary-domain integral equation

[TABLE]

where $G_{\alpha\beta}(\mathbf{r},\,\mathbf{r}^{\prime})$ is a Green’s function of the Stokes equation and $K_{\alpha\beta\gamma}(\mathbf{r},\,\mathbf{r}^{\prime})$ is the associated stress tensor. These kernels satisfy the Stokes system

[TABLE]

where $P_{\alpha}$($\mathbf{r},\,\mathbf{r}^{\prime})$ is the pressure vector. Implicit in the above is a choice of Green’s function that satisfies no-slip boundary conditions at any boundary of the fluid that is not part of $S$.

Defining the single-layer and double-layer integral operators $\boldsymbol{G}$ and $\boldsymbol{K}$, which act, respectively, on tractions and velocities, as

[TABLE]

and a Brownian velocity field $\boldsymbol{w}$ as

[TABLE]

the solution of the boundary-domain integral equation can be expressed formally in terms of the inverse of the single-layer integral operator as

[TABLE]

This formal solution shows that: (i) the traction is a sum of a Brownian contribution $\boldsymbol{\hat{f}}=\,\boldsymbol{G}^{-1}\cdot\boldsymbol{w}$ from the fluctuations in the fluid and a deterministic contribution from the boundary condition, containing both the rigid body motion $\boldsymbol{f}^{\mathcal{D}}=-\boldsymbol{G}^{-1}\cdot\boldsymbol{v}^{\mathcal{D}}$ and the active slip $\boldsymbol{f}^{\mathcal{A}}=\boldsymbol{G}^{-1}\cdot\left(-\tfrac{1}{2}\boldsymbol{I}+\boldsymbol{K}\right)\cdot\boldsymbol{v}^{\mathcal{A}}$ (ii) the Brownian traction is a zero-mean Gaussian random variable whose variance is linearly related to the variance of the thermal force $\boldsymbol{\xi}$ and (iii) the variance of the Brownian traction can be determined from the inverse of the single-layer operator and the fluctuation-dissipation relation for the thermal force. In the next section, we provide a solution for the boundary-domain integral equation in a basis adapted for a spherical active colloid and, thereby, derive the explicit form of the traction in terms of generalized friction tensors.

III Traction on a spherical active colloid

We now consider a spherical colloid of radius $b$ whose center is at $\mathbf{R}$ and whose orientation is specified by the unit vector $\mathbf{p}$. A point on the boundary is $\mathbf{r}=\mathbf{R}+\boldsymbol{\rho}$, where $\boldsymbol{\rho}$ is the radius vector. The boundary velocity, \boldsymbol{v}=\text{\boldsymbol{v}}^{\mathcal{D}}+\boldsymbol{v}^{\mathcal{A}}, is the sum of its rigid body motion \text{\boldsymbol{v}}^{\mathcal{D}}=\mathbf{V}+\boldsymbol{\Omega}\times\boldsymbol{\rho}, specified by the velocity $\mathbf{V}$ and angular velocity $\mathbf{\Omega}$, and the active slip velocity $\boldsymbol{v}^{\mathcal{A}}$. The only restriction on the active slip is that its integral over the surface of the sphere is zero. This ensures conservation of mass in the fluid.

We choose the tensorial spherical harmonics $\mathbf{Y}^{(l)}(\bm{\hat{\rho}})=(-1)^{l}\rho^{l+1}\bm{\nabla}^{(l)}\rho^{-1}$, where $\bm{\nabla}^{(l)}=\bm{\nabla}_{\alpha_{1}}\dots\bm{\nabla}_{\alpha_{l}}$, as complete orthogonal basis functions on the sphere. In this basis, the active slip is expanded as**

[TABLE]

The coefficients $\mathbf{V}^{(l)}$ are reducible Cartesian tensors of rank $l$, with three irreducible parts of ranks $l,$ $l-1,$ and $l-2$, corresponding to symmetric traceless, antisymmetric and pure trace combinations of the reducible indices. These irreducible parts are $\mathbf{V}^{(l\sigma)}=\mathbf{P}^{(l\sigma)}\cdot\mathbf{V}^{(l)}$, where the index $\sigma=s,\,a\text{ and }t$, labels the symmetric irreducible, antisymmetric and pure trace parts of the reducible tensors Brunn (1976); Schmitz (1980); Ghose and Adhikari (2014); Singh et al. (2015); Singh and Adhikari (2016a, b). The operator $\mathbf{P}^{(ls)}=\boldsymbol{\Delta}^{(l)}$ extracts the symmetric irreducible part, $\mathbf{P}^{(la)}=\boldsymbol{\Delta}^{(l-1)}\boldsymbol{\boldsymbol{\varepsilon}}$ the antisymmetric part and $\mathbf{P}^{(lt)}=\boldsymbol{\delta}$ the trace of the operand. Here $\mathbf{\Delta}^{(l)}$ is a tensor of rank $2l$, projecting any $l$-th order tensor to its symmetric irreducible form Hess (2015), $\boldsymbol{\varepsilon}$ is the Levi-Civita tensor and $\boldsymbol{\mathbf{\delta}}$ is the Kronecker delta.

The traction can be similarly expanded in the tensorial harmonic basis as**

[TABLE]

and, as with the velocity coefficients, each traction coefficient at order $l$ has three irreducible parts indexed by $\sigma$. The net force and the torque are given by first two irreducible coefficients

[TABLE]

Note that each $\mathbf{V}^{(l)}$ and $\mathbf{F}^{(l)}$ have the dimension of velocity and force respectively.

It is convenient to express the traction as a sum of rigid body, active and Brownian contributions,

[TABLE]

with corresponding expansion coefficients $\mathbf{F}^{\mathcal{D}(l)}$, $\mathbf{F}^{\mathcal{A}(l)}$ and $\hat{\mathbf{F}}^{(l)}$. By linearity, the three parts of the traction satisfy* independent *boundary integral equations. Recalling that rigid body motion is an eigenvector of the double-layer integral operator, these are

[TABLE]

Addition of the above three equations recovers the boundary integral equation, Eq. (5), for the net traction.

The first integral equation for the Stokes drag has been well-studied in the literature on suspension mechanics. The second integral equation for the active thrust has been studied recently in the context of active colloids in an athermal fluid Singh et al. (2015); Singh and Adhikari (2016b, a). The third integral equation for the fluctuating traction is studied here for the first time. Physically, the fluctuating traction corresponds to the distribution of surface forces necessary to keep the sphere stationary in the incident Brownian velocity field $\boldsymbol{w}(\mathbf{r)}$. From this, it is particularly clear that the Brownian traction is a zero-mean Gaussian random variable whose variance is related to that of $\boldsymbol{\xi}$. We now present explicit solutions for each of the integral equations using Galerkin’s method of discretization. Linear algebraic systems are derived by inserting the basis expansions on each side of the integral equations, weighting the result by another basis function and integrating over the boundary. We refer the reader to Singh et al. (2015); Singh and Adhikari (2016b, a) where the procedure is explained in detail.

III.1 Rigid body traction

The linear algebraic system for the rigid body contribution to the traction, with the summation convention for repeated indices, is

[TABLE]

Here the matrix elements of the single-layer operator $\boldsymbol{G}(\mathbf{R)}$ in the tensorial harmonic basis are

[TABLE]

and $l$-th tensorial harmonic coefficients of the rigid body motion and the traction are, respectively, $\mathbf{V}^{\mathcal{D}(l)}$ and $\mathbf{F}^{\mathcal{D}(l)}$. Clearly, the only non-zero coefficients of $\mathbf{V}^{\mathcal{D}(l)}$ are $\mathbf{V}$ and $\boldsymbol{\Omega},$ corresponding to $l\sigma=1s$ and $l\sigma=2a$ respectively. The solution of the linear system is

[TABLE]

and the friction tensors $\boldsymbol{\gamma}^{(l\sigma,\,1s)}$ and $\boldsymbol{\gamma}^{(l\sigma,\,2a)}$ are given by

[TABLE]

The friction tensors above give the contribution to the traction from rigid body motion Brenner (1963); Felderhof (1976). The traction modes $l\sigma=1s$ and $l\sigma=2a$ correspond to the forces and torques

[TABLE]

where we have introduced the indices $T,R=1s,2a$ to make contact with the usual notation. The inverse of the single-layer operator can be computed in several ways both analytically and numerically. The Jacobi iteration used in Singh and Adhikari (2016b, a) is convenient for analytical expressions.

The fluid flow $\mathbf{u}^{\mathcal{D}}$ for a rigid body motion of the sphere has the integral representation $\mathbf{\mathbf{u}^{\mathcal{D}}(\mathbf{r})}=-\,\boldsymbol{G}{}^{(l)}(\mathbf{R},\,\mathbf{r})\cdot\mathbf{F^{\mathcal{D}}}^{(l)}$ where the boundary integral

[TABLE]

is the contribution to the external flow from the $l$-th tensorial coefficient of the traction. The double layer has no contribution from rigid body motion to the external flow. This result will be used below in deriving a key identity necessary for deriving the variance of the Brownian traction.

III.2 Active traction

The linear algebraic system for the active contribution to the traction is

[TABLE]

where the matrix elements of the double-layer operator $\boldsymbol{K}(\mathbf{R)}$ in the tensorial harmonic basis are

[TABLE]

The irreducible parts of the slip coefficients are $\mathbf{V}^{(l\sigma)}$. The first two modes, $\mathbf{V}^{(1s)}\equiv-\mathbf{V}^{\mathcal{A}}$ and $\mathbf{V}^{(2a)}\equiv-b\mathbf{\Omega}^{\mathcal{A}}$ are given by the integrals

[TABLE]

and are equal to the self-propulsion velocity and self-rotation angular velocity of an isolated active colloid in unbounded flow Anderson (1989); Ghose and Adhikari (2014). The solution of the linear system is

[TABLE]

The generalized friction tensors $\boldsymbol{\gamma}^{(l\sigma,\,l^{\prime}\sigma^{\prime})}$ give the active contribution to the traction. These tensors were first introduced in Singh and Adhikari (2016a) and are given by

[TABLE]

They can be interpreted as Onsager coefficients of the linear response of the traction to the surface slip. The above expression is identical to Eq. (18) for $l^{\prime}\sigma^{\prime}=1s,2a$ since rigid body motion is an eigenvector of the double layer operator. Therefore, both active *and *passive friction tensors can be recovered from the general expression above.

The active force and torque on the colloid are the $l\sigma=1s$ and $2a$ coefficients of the traction. These are given as

[TABLE]

The active forces and torques depend on all modes of the slip. In contrast to passive colloids, where only four friction tensors determine the force and torque, there are, in general, infinitely many friction tensors for active colloids, corresponding to the infinitely many modes of the slip. These infinitely many friction tensors account for the diversity of phenomenon seen in active suspensions.

III.3 Brownian traction

The linear algebraic system for the Brownian contribution to the traction is

[TABLE]

where $\mathbf{\boldsymbol{W}}^{(l)}$ are coefficients of the irreducible tensor expansion

[TABLE]

of the Brownian velocity field incident on the surface of the colloid. From the definition of the Brownian velocity field, Eq.(8), these coefficients are given by

[TABLE]

The solution of the linear system for the Brownian traction is

[TABLE]

The coefficients of the Brownian traction are proportional to the coefficients of the Brownian velocity field incident on the surface of the colloid and, by Eq. (8), to the thermal force in the fluctuating Stokes equation. It is clear, then, that the traction coefficients are zero-mean Gaussian random variables and to fully specify their distribution it is necessary, then, to only determine their variance. By the previous equation, their variance is related to that of the Brownian velocity coefficients as

[TABLE]

To determine the variance of Brownian velocity coefficients we use the boundary integral representation of $\boldsymbol{u}^{\mathcal{D}}$ given above (see also Singh et al. (2015); Singh and Adhikari (2016a)). Inserting this on the left of the fluctuation-dissipation relation for the thermal force, Eq. (3), gives

[TABLE]

On the other hand, the power dissipated by the rigid body motion, on the right of the fluctuation-dissipation relation can be expressed as

[TABLE]

The first equality is obtained by using the divergence theorem to reduce the volume integral of the quadratic form in Eq. (4) to the boundary of the colloid and then using the constitutive relation between the stress and the strain rate in Stokes flow Landau and Lifshitz (1959). The second equality is an elementary consequence of the orthogonality of the tensorial spherical harmonics Singh et al. (2015); Singh and Adhikari (2016a). The third equality is obtained by eliminating the velocity using the linear algebraic system, Eq.(15) for rigid body motion. Comparing the above two equations, we obtain the key identity for the variance of the Brownian velocity coefficients,

[TABLE]

This expression when used in Eq. (31) yields the variance of Brownian traction coefficients

[TABLE]

These are an infinite number of fluctuation-dissipation relations between the variance of the tensorial harmonic modes of the fluctuating traction and the matrix elements, in the irreducible tensorial harmonic basis, of the *inverse *of the single-layer operator. To the best of our knowledge, the complete statistics of the Brownian traction is derived here for the first time and is the central result of this paper.

The variance of the irreducible parts of the fluctuating traction follow straightforwardly as

[TABLE]

The first two coefficients of the fluctuating traction are the force and torque, and choosing $l\sigma=1s,\,2a$ we obtain

[TABLE]

where $\boldsymbol{\gamma}^{\alpha\beta}$, with $\alpha,\beta=T,R$, are the one-particle friction tensor and are ($l\sigma=1s,\,2a$) elements of \mathbf{P}^{(l\sigma)}\cdot\big{[}\mathbf{\boldsymbol{G}}^{-1}(\mathbf{R})\big{]}^{(l,\,l^{\prime})}\cdot\mathbf{P}^{(l^{\prime}\sigma^{\prime})}. The fluctuation-dissipation relation for the Brownian force and torque are, thus, derived from the fluctuation-dissipation relation for the thermal force on the fluid.

We make the following remarks about the above derivation. First, the explicit form of the inverse of the single-layer operator is is not necessary to obtain Eq. (35); it is sufficient to know that the inverse exists. Therefore, the fluctuation-dissipation relation for the irreducible coefficients, Eq. (36), is valid for any geometry bounding the fluid, provided the flow vanishes there. In particular, it holds for a colloid near a plane wall. Second, our derivation provides the fluctuation-dissipation relation for *all *modes of the Brownian traction. Earlier derivations have focused on only the force and torque. Therefore, our derivation provides the fluctuating stresslet ($l\sigma=2s)$ which is needed to compute the Brownian contribution to the suspension stress. Third, the configuration-dependence of both the fluctuation, Eq. (27), and the dissipation, Eq. (36), is made explicit in our derivation. The configuration-dependent “noise” variance follows from this automatically. The interpretation of the resulting multiplicative noise in the Langevin equation that we derive below is obtained by recalling that the momentum and angular momentum of the colloid are fast variables that have, implicitly, been adiabatically eliminated Gardiner (1984). The form of the Smoluchowski equation for this problem is well-known and is used below to consistently interpret the multiplicative noise in the Langevin equation Chandrasekhar (1949); Murphy and Aguirre (1972); Wilemski (1976).

IV Langevin equations for a spherical active colloid

In this section we derive the Langevin and Smoluchowski equations for the Brownian motion of an active colloid. We use the results derived above for three kinds of forces and torques acting on the colloid. With these, the balance conditions, Eq. (1), become

[TABLE]

The above can be inverted to obtain the rigid body motion of the colloid in explicit form. This gives the Langevin equations for a Brownian active colloid with hydrodynamic interactions, first derived heuristically in Laskar and Adhikari (2015),

[TABLE]

Here $\boldsymbol{\eta}^{\alpha}$, and $\boldsymbol{\zeta}^{\alpha}$ are Gaussian white noises with zero-mean and variances $1$ and $1/b$ respectively. The matrix square roots are to be interpreted as Cholesky factors. The mobility matrices $\boldsymbol{\mu}^{\alpha\beta}$ are inverses of the friction matrices $\boldsymbol{\gamma}^{\alpha\beta}$ Happel and Brenner (1965); Felderhof (1977); Mazur and Saarloos (1982); Schmitz and Felderhof (1982); Nunan and Keller (1984); Ladd (1988); Durlofsky et al. (1987); Brady et al. (1988); Kim and Karrila (1992); Cichocki et al. (1994). The propulsion tensors $\boldsymbol{\pi}^{(\alpha,\,l\sigma)}$, first introduced in Singh et al. (2015), relate the rigid body motion to modes of the active velocity. They are related to the generalized friction tensors, introduced in Singh and Adhikari (2016a), by

[TABLE]

The translational propulsion tensors $\bm{\pi}^{(\text{T},\,l\sigma)}$ are dimensionless while the rotational propulsion tensors $\bm{\pi}^{(R,\,l\sigma)}$ have dimensions of inverse length. Stochastic trajectories of motion can be obtained from the kinematic equations

[TABLE]

using the standard Ermak-McCammon integrator Ermak and McCammon (1978). In this integrator, the noise variances are computed using mobilities in the configuration at time $t$ but a “spurious” drift, proportional to the configurational divergence of the mobilities, is added to to arrive at the configuration at time $t+\Delta t$. There is nothing particularly spurious about this drift; it is simply the residual effect of the adiabatically eliminated degrees of freedom.

The Smoluchowski equation for the distribution function $\Psi(\mathbf{R};\,\mathbf{p})$ of positions and orientations follows immediately from the Langevin equations. We write it in the form of a conservation law in configuration space

[TABLE]

where the “velocities” $\boldsymbol{\mathcal{V}}_{{\scriptscriptstyle\mathbf{R}}}$ and $\boldsymbol{\mathcal{V}}_{{\scriptscriptstyle\mathbf{p}}}$ are

[TABLE]

Here $\mathbf{F}^{P}=-\boldsymbol{\nabla}_{{\scriptscriptstyle\mathbf{R}}}U$, $\mathbf{T}^{P}=-\mathbf{p}\times\boldsymbol{\nabla}_{{\scriptscriptstyle\mathbf{p}}}U$, and $U$ is a potential that contains both positional and orientational interactions. Note the position of the mobility in the second derivative terms: the $\boldsymbol{\nabla\mu\nabla}$ order (in contrast with two other inequivalent permutations) is provided unambiguously by the adiabatic elimination of momenta.

In the absence of activity, the drift and diffusion coefficients in the Smoluchowski equation obey the fluctuation-dissipation relation and the Gibbs distribution $\Psi\sim\exp(-U/k_{B}T$) is the stationary solution. However, the Gibbsian form is not a stationary solution in the presence of the active terms, as can easily be verified by substitution. This applies, a fortiori, to an active colloid near a plane wall discussed in the next section, where the barometric distribution of a passive suspension is no longer the stationary distribution.

V Brownian active colloid near a plane wall

We now apply the preceding general results to the specific case of an active colloid near a plane wall. The Green’s function for the problem, denoted by $\mathbf{G}^{\text{w}}$, is taken to vanish at the location of the wall, $z=0$. The form of the Green’s function, first derived by Lorentz Lorentz (1896), can be written in the following form due to Blake Blake (1971)

[TABLE]

where $G_{\alpha\beta}^{0}(\mathbf{r})=(\nabla^{2}\delta_{\alpha\beta}-\nabla_{\alpha}\nabla_{\beta})(r/8\pi\eta$) is the Oseen tensor and the correction necessary to satisfy the boundary condition is

[TABLE]

The correction is interpreted as a sum of three images, a Stokeslet, a dipole, and a degenerate quadrupole, located at $\mathbf{\mathbf{R}}^{*}=\boldsymbol{\mathcal{M}}\cdot\mathbf{\mathbf{R}}$, where $\boldsymbol{\mathcal{M}}=\boldsymbol{I}-2\mathbf{\hat{z}}\mathbf{\hat{z}}$ is the mirror operator with respect to the wall Brennen and Winet (1977); Berke et al. (2008); Lauga and Powers (2009); Goldstein (2015); Mathijssen et al. (2015); Singh and Adhikari (2016b). Here $h$ is the height of the colloid from the wall and $\mathbf{r}^{*}=\mathbf{\mathbf{R}}^{\prime}-\mathbf{\mathbf{R}}^{*}$. The correction has no singularities in the domain of flow and the only singular contribution there is from the Oseen tensor. The flows produced by a sphere translating or rotating near the wall are shown in the top panels of Fig. (1).

The leading terms of the inverse of the single-layer operator are obtained here using Jacobi’s iterative method Singh and Adhikari (2016a). The results for the friction tensors $\boldsymbol{\gamma}^{\alpha\beta}$ are given in the left column of Table (1). The grand resistance tensor, formed out each of the $3\times 3$ blocks, is manifestly symmetric and, for $h>b$, positive-definite. This is established by an explicit computation of the eigenvalues, shown in the first figure of the bottom panel of Fig. (1). Thus positive-definiteness is ensured for all configurations in which the colloid does not overlap with the wall.

In Brownian dynamics simulations, however, it is convenient to avoid hard-sphere potentials as these require special integrators. In the absence of such potentials, it is no longer possible to maintain the constraint $h>b$ during integration, and the colloid may substantially overlap with the wall. For such configurations, a naive continuation of the friction tensors to the domain $h<b$ is untenable: as the dashed lines in the eigenvalue plot in Fig. (1) show, the grand resistance tensor has negative eigenvalues in this region. This implies negative entropy production and is clearly unphysical. The cure, first proposed by Rotne and Prager in the context of bead-spring models of polymers Rotne and Prager (1969), is to regularize the matrix elements of the single-layer, by computing the surface integrals in their definition, over the union of overlapping surfaces. In this case, the integral is to be computed over the union of the surface of the colloid and its image. Such overlap integrals have been computed recently Wajnryb et al. (2013) and we use those results to obtain the regularized forms of the friction in the right column of Table (1). The use of the regularized form provides friction tensors that are positive-definite for $all$ heights, as can be seen from the solid curves for $h<b$ in the eigenvalue plot in Fig. (1). Their use results in positive-definite variances for the Brownian force and torque at $all$ heights above the wall.

In the second and third figures of the bottom panel, we plot the condition number of the grand resistance tensor as a function of height and its pseudo-spectrum at $h=1.3b$. The condition number remains small for all heights and the pseudo-spectrum shows no sign of non-normality. Thus, computing the “square-root” Cholesky factors of the grand friction tensor (or, its inverse, the grand mobility tensor) poses no problem and iterative methods are expected to converge rapidly. Brownian dynamics simulations of active colloids near a plane wall can then be efficiently performed using the above results.

VI Discussion

We have shown how to obtain the three parts of the traction, due to rigid body motion, activity and thermal fluctuations, on a spherical colloid in a fluctuating Stokesian fluid. We then applied the general results to the specific case of an active colloid near a plane wall. We provided a regularization of the friction tensors to enable Brownian dynamics simulations in which the spheres can overlap with the walls. We conclude with three remarks about our results.

The first is that expressions in the left column of Table (1) are valid not only for the Lorentz-Blake Green’s function, but for any Green’s function that can be expressed as $\mathbf{G}^{0}+\mathbf{G}^{\ast}$, where $\mathbf{G}^{0}$ is the Oseen tensor and $\mathbf{G}^{*}$ is the correction necessary to satisfy the boundary conditions. A variety of Green’s function can be expressed in this form, including those for flow between parallel walls and in periodic domains. As mentioned before, the correction term does not contain singularities in the domain of flow and, therefore, the boundary integrals in the definition of the matrix elements can be expressed in terms of derivatives of the correction. Singular terms from the Oseen tensor can be calculated explicitly using well-known results for integrals of Bessel functions. Therefore, our results for $h>b$ are of broader validity than might have been anticipated.

The second is that active colloids show fascinating behaviour in the proximity of a wall Palacci et al. (2013); Petroff et al. (2015). While Brownian motion appears to be negligible in comparison to activity-induced motion for a large class of these colloids, there is considerable theoretical interest in understanding the interplay between passive friction, thermal fluctuations, and activity, especially when the friction is configuration-dependent. Brownian dynamics, with the regularized frictions provided here, will be a powerful tool to study such interplays.

Finally, the extension of the method presented here to determine the tractions on the surface of many active colloids is straightforward in principle but tedious in practice. The result for the active traction has been obtained recently Singh and Adhikari (2016a) and a heuristic argument has been provided to determine the first two moments of the Brownian traction. It will be instructive to obtain all moments and, in particular, the symmetric second moment to determine their contribution to the stress in a non-dilute active Brownian suspension.