Interplay between Brownian and hydrodynamic tracer diffusion in suspensions of swimming microorganisms

Abstract

The general problem of tracer diffusion in non-equilibrium baths is important in a wide range of systems, from the cellular level to geographical length scales. In this paper, we revisit the archetypical example of such a system: a collection of small passive particles immersed in a dilute suspension of non-interacting dipolar microswimmers, representing bacteria or algae. In particular, we consider the interplay between thermal (Brownian) diffusion and hydrodynamic (active) diffusion due to the persistent advection of tracers by microswimmer flow fields. Previously, it has been argued that even a moderate amount of Brownian diffusion is sufficient to significantly reduce the persistence time of tracer advection, leading to a significantly reduced value of the effective active diffusion coefficient $D_A$ compared to the non-Brownian case. Here, we show by large-scale simulations and kinetic theory that this effect is in fact practically relevant only for microswimmers that effectively remain stationary while still stirring up the surrounding fluid – so-called shakers. In contrast, for moderate and high values of the swimming speed, relevant for biological microswimmer suspensions, the effect of Brownian motion on $D_A$ is negligible, leading to the effects of advection by microswimmers and Brownian motion being additive. This conclusion contrasts with previous results from the literature, and encourages a reinterpretation of recent experimental measurements of $D_A$ for tracer particles of varying size in bacterial suspensions.

1. Introduction

Understanding the mass transport of colloidal and molecular species in non-equilibrium environments is crucial for various processes, ranging from active intracellular transport (Mogre, Brown & Koslover 2020) to the dispersion of nutrients in world oceans (Katija 2012). Apart from their practical importance, the transport properties of tracer particles in generic ‘active baths’ have attracted much interest from a statistical physics perspective, where they can be viewed as a minimal example of particles driven by external, non-equilibrium noise (Argun et al. 2016; Park et al. 2020). Beyond the level of tracer particles driven by generic non-equilibrium noise, the archetypical example of a tracer particle in an active bath is a collection of point-like tracers being advected by a set of microswimmers such as bacteria or algae (Lauga & Powers 2009). When swimming through a viscous fluid, these swimmers create long-range flow fields that advect the tracers, leading to tracer dynamics that is ballistic at short times and diffusive over time scales longer than the autocorrelation time of the local flow field (Lin, Thiffeault & Childress 2011). Realisations of this system have been studied extensively both experimentally, typically in suspensions of E. coli bacteria (Wu & Libchaber 2000; Kim & Breuer 2004; Drescher et al. 2011; Jepson et al. 2013; Koumakis et al. 2013; Miño et al. 2013, 2011; Patteson et al. 2016; Peng et al. 2016; Semeraro, Devos & Narayanan 2018) or Chlamydomonas algae (Leptos et al. 2009; Yang et al. 2016; Ortlieb et al. 2019; von Rüling, Kolley & Eremin 2021), and theoretically, with microswimmers typically being modelled either as force dipoles acting on the surrounding fluid (Pushkin, Shum & Yeomans 2013; Pushkin & Yeomans 2013; Morozov & Marenduzzo 2014; Nordanger, Morozov & Stenhammar 2022), as spherical ‘squirmers’ with an imposed slip velocity along their body (Thiffeault & Childress 2010; Lin et al. 2011; Thiffeault 2015), or as needle-shaped ‘slender swimmers’ with imposed stresses along their body lengths (Saintillan & Shelley 2012; Krishnamurthy & Subramanian 2015). While the details of these three microswimmer models differ, the results regarding enhanced tracer diffusion are largely generic and consistent with experimental results, which have shown the swimmer-induced, hydrodynamic diffusivity $D_A$ to scale linearly with microswimmer density $n$ in the dilute limit where swimmer–swimmer correlations can be neglected (Thiffeault & Childress 2010; Lin et al. 2011; Miño et al. 2013). In this limit, a fruitful way of calculating $D_A$ is to consider the net displacement due to binary swimmer–tracer scattering events (Pushkin & Yeomans 2013; Morozov & Marenduzzo 2014); two examples of resulting (deterministic) tracer trajectories for scattering events are shown in figure 1. For a non-tumbling swimmer, starting and ending at $x = \pm \infty$, the resulting loop is closed, leading to a vanishing tracer net displacement $\varDelta$ (figure 1a). For tumbling swimmers with a finite persistence length, the trajectory is, however, punctuated midway through the tracer loop, leading to significantly larger values of $\varDelta$ (figure 1b). The resulting value of $D_A$ due to a large set of such scattering events can then be obtained by summing over all possible sets of scattering parameters.

Figure 1. Tracer trajectories for infinite and finite swimmer paths. (a) Typical trajectory for a non-diffusing tracer advected by a non-tumbling, point-dipole swimmer following an effectively infinite, straight path. (b) Corresponding trajectory terminated due to a tumbling event. Note that, per (1.1), the effective tracer diffusion $D_A$ is independent of tumbling rate $\lambda$ for sufficiently high swimming speeds when averaged over all possible swimmer–tracer configurations, even though the net displacement $\varDelta$ is much larger for the tumbling swimmer. Tracer trajectories were obtained through direct numerical integration using a non-regularised dipolar flow field, as described by Morozov & Marenduzzo (2014).

In spite of the dependence of $\varDelta$ on the microswimmer tumbling rate $\lambda$ for single scattering events such as that in figure 1, Pushkin & Yeomans (2013) showed that, in the limit of large swimming speeds $v_s$, $D_A$ is in fact independent of $\lambda$ when summed over all possible swimmer trajectories. This result was later generalised by Škultéty et al. (2020) to arbitrary swimming speeds, leading to the following approximate expression for $D_A$:

(1.1)\begin{equation} D_A \approx \frac{7\kappa^2 n}{2048 \lambda \varepsilon + 336{\rm \pi} v_s}. \end{equation}

Here, $\kappa$ is the magnitude of the microswimmer dipole (in units of volume over time), and $\varepsilon$ is a characteristic linear size of the microswimmer, which we take to be equal to the short-range regularisation length of the dipolar flow field; unless stated otherwise, in the following we will use $\kappa$ and $\varepsilon$ to non-dimensionalise the numerical data.

Rather than using the scattering approach outlined above, (1.1) was derived by formulating a kinetic theory for the spatiotemporal correlations of the disturbance velocity field ${\boldsymbol U}$ created by the swimmers. Due to the linearity of Stokes flow, we can write ${\boldsymbol U}$ as a superposition of the individual swimmer flow fields ${\boldsymbol u}_s$:

(1.2)\begin{equation} {\boldsymbol U}({\boldsymbol r},t) = \sum_{i=1}^{N} {\boldsymbol u}_s({\boldsymbol r};{\boldsymbol r}_i,{\boldsymbol p}_i), \end{equation}

where ${\boldsymbol r}_i$ and ${\boldsymbol p}_i$ are, respectively, the position and orientation of swimmer $i$. Knowing the statistical properties of ${\boldsymbol U}$, $D_A$ can be calculated via the Green–Kubo relation

(1.3)\begin{align} D_A = \frac{1}{3} \int_{0}^{\infty}\langle \dot{{\boldsymbol r}}_T(t) \boldsymbol{\cdot} \dot{{\boldsymbol r}}_T(0)\rangle \, {\rm d} t = \frac{1}{3} \int_{0}^{\infty}\langle {\boldsymbol U}({\boldsymbol r}_T,t) \boldsymbol{\cdot} {\boldsymbol U}({\boldsymbol r}_T,0) \rangle \, {\rm d} t \equiv \frac{1}{3} \int_{0}^{\infty}C_T(t) \, {\rm d} t, \end{align}

where, in the second equality, we have assumed point-like tracers advected by the disturbance flow, so that $\dot {{\boldsymbol r}}_T = {\boldsymbol U}({\boldsymbol r}_T)$, and the third equality defines the velocity autocorrelation function $C_T(t)$ in the co-moving tracer frame. While yielding results for dilute suspensions identical to those for the scattering approach discussed above, kinetic theories are more readily extended to accommodate the effect of swimmer–swimmer correlations due to the mutual advection and reorientation of swimmers (Škultéty et al. 2020). Importantly, these interactions break the symmetry between rear-actuated (pusher) swimmers, such as most bacteria, and front-actuated (puller) ones, such as Chlamydomonas, leading to a super-linear scaling of $D_A$ with $n$ for pushers, and a corresponding sub-linear scaling for pullers (Stenhammar et al. 2017).

Equation (1.1) shows two qualitatively different regimes for high and low $v_s$: for $v_s \!\rightarrow\! 0$ – the so-called shaker limit – the dominant mechanism controlling the decorrelation of $C_T(t)$ is tumbling of the swimmer. For fast swimmers, with $v_s \gg \lambda \varepsilon$, the decorrelation of the fluid velocity is instead dominated by the swimmer self-propulsion and is thus independent of $\lambda$, and (1.1) reduces to the expression derived by Pushkin & Yeomans (2013). A third, somewhat less explored, mechanism affecting $C_T(t)$ is Brownian translational diffusion of the tracer: even though Brownian diffusion does not affect the statistics of the flow field ${\boldsymbol U}({\boldsymbol r},t)$ as measured in the lab frame, the positional noise of the tracer particle will cause it to cross the streamlines of the disturbance flow, thus perturbing its trajectory compared to the athermal case shown in figure 1, leading to a lower $D_A$. Since the Brownian diffusion constant $D_0$ depends inversely on the tracer radius per the Stokes–Einstein relation, the magnitude of this effect is expected to be significant primarily for small tracer particles, and it has been hypothesised to explain the non-monotonic size dependence of enhanced tracer diffusion seen in experiments with colloids in E. coli suspensions, where Patteson et al. (2016) observed a maximum in $D_A$ for a tracer radius of approximately 5 $\mathrm {\mu }$m. Furthermore, the effect of Brownian motion on enhanced tracer diffusion was analysed theoretically by Kasyap, Koch & Wu (2014) for a model of slender-body swimmers, showing that $D_A$ is a non-monotonic function of $D_0$, with $D_A$ first showing a small increase for intermediate $D_0$, before falling below the athermal value as $D_0$ grows larger. In a more recent study of finite-size spherical tracers in microswimmer suspensions, Dyer & Ball (2021) analysed numerically the combined effect of thermal fluctuations and near-field flows on the size-dependent tracer dynamics, finding non-monotonic behaviour similar to that observed experimentally. In this paper, we will revisit the problem of the interplay between Brownian and hydrodynamic diffusion for the simple case of point-like tracers immersed in a dilute suspension of microswimmers described via a regularised dipolar flow field. Using kinetic theory and large-scale lattice Boltzmann (LB) simulations of E. coli-like suspensions, we show that the effect of Brownian diffusion on active diffusion is practically relevant only whenever $v_s < \lambda \varepsilon$, which corresponds to extremely slow (or frequently tumbling) swimmers. For biologically relevant values of $v_s$ and $\lambda$, swimming is instead the dominant decorrelation mechanism, so that $D_A$ becomes independent of both $D_0$ and $\lambda$. In contrast to previous studies, our results thus indicate that the effect of Brownian motion on the enhanced diffusion is in fact negligible for most microswimmer realisations, and thus that the experimentally observed non-monotonic size dependence of $D_A$ on tracer size reported by Patteson et al. (2016) has other explanations.

2. Model and method

We consider a collection of $N$ non-interacting microswimmers at number density $n = N/V$ moving through a three-dimensional viscous fluid of viscosity $\mu$. Each microswimmer is composed of two equal and opposite point forces of magnitude $F$ separated by a length $\ell$, and swims with a constant speed $v_s$. The resulting reduced hydrodynamic dipole strength is $\kappa = F\ell /\mu$. Furthermore, the swimming direction ${\boldsymbol p}_i$ relaxes through Poisson-distributed random tumbles with uncorrelated directions, occurring with average frequency $\lambda$.

The position ${\boldsymbol r}_T$ of a point-like tracer obeys the equation of motion

(2.1)\begin{equation} \dot{{\boldsymbol r}}_T = {\boldsymbol U}({\boldsymbol r}_T) + \sqrt{2D_0}\,\boldsymbol{\eta}, \end{equation}

where $\boldsymbol {\eta }$ is a unit-variance white noise, $\delta$-correlated in time, and $D_0$ is the Brownian diffusion constant. Thus, the effect of Brownian motion is fully contained in the tracer dynamics, while we assume the effects of thermal fluctuations on the fluid and on the pairwise swimmer–tracer dynamics to be subdominant. The fluid disturbance velocity ${\boldsymbol U}({\boldsymbol r}_T)$ due to the presence of all swimmers is given by (1.2) and can, in principle, be summed up explicitly on each time step. However, to avoid the (prohibitively costly) pairwise summation over all swimmers and tracers, we instead solve numerically for the flow field using an efficient point-force implementation of the LB method described previously (Nash, Adhikari & Cates 2008; Bárdfalvy et al. 2020). Our simulations comprise a system with periodic boundaries and size $100^3$ lattice units. In LB units, set by the LB lattice spacing ${\rm \Delta} l$ and time step ${\rm \Delta} t$, the microswimmer density was kept constant at $n = 0.01$ unless stated otherwise, corresponding to $N = 10^4$ microswimmers. Furthermore, $N_T = 5 \times 10^4$ point tracers were included for statistical averaging. In addition to the LB simulations, in § 3 we will furthermore extend the kinetic theory developed previously by Škultéty et al. (2020) to the case of microswimmers undergoing Brownian diffusion, and show that in the limit of non-interacting swimmers, the derived expression for $D_A$ is equivalent to that for Brownian tracers in a suspension of non-Brownian microswimmers.

Inserting (2.1) into the Green–Kubo relation (1.3) yields

(2.2)\begin{equation} D_{Tot} = D_0 + \frac{1}{3} \int_{0}^{\infty} C_T(t) \, {\rm d} t = D_0 + D_A, \end{equation}

where $C_T$ was defined in (1.3). Thus, to obtain $D_A$, we evaluate numerically the time correlation of the disturbance velocity measured in the co-moving tracer frame. Since the tracer position ${\boldsymbol r}_T$ depends on $D_0$, $C_T(t)$, and thus $D_A$, will depend implicitly on $D_0$.

To characterise the system, we will use three dimensionless quantities. First, we define the Péclet number, which measures the relative importance of active and thermal forces, as

(2.3)\begin{equation} {Pe} \equiv \frac{D_A(D_0=0)}{D_0}, \end{equation}

where $D_A(D_0=0)$ is the active diffusivity of a tracer immersed in an equivalent microswimmer suspension but in the absence of Brownian tracer motion. It should be noted that our definition of $Pe$ is qualitatively different from that of Kasyap et al. (2014), who, instead of $D_A(D_0=0)$, use the swimming speed $v_s$ to characterise the active forces. We, however, argue that the activity experienced by the tracers depends on the magnitude of the velocity fields generated by the swimmers, and is thus dependent on $\kappa$ and encoded in $D_A(D_0 = 0)$, while $v_s$ is instead a measure of the swimmers’ self-propulsion. For experimental realisations of microswimmers, $\kappa$ and $v_s$ are directly proportional to each other; however, the specific relation between them will nevertheless be specific to each type (or species) of swimmer, and decoupling them enables us to study separately the effects of self-propulsion and fluid forcing, as we demonstrate further below.

Second, we measure the change in active diffusion due to Brownian motion through the quantity

(2.4)\begin{equation} \xi \equiv \frac{D_A(D_0)}{D_A(D_0 = 0)}. \end{equation}

In the limit ${Pe} \rightarrow \infty$, where Brownian motion becomes negligible, we thus expect that $\xi \rightarrow 1$ as $D_A$ approaches its non-Brownian value. Finally, in accordance with Škultéty et al. (2020), we account for the effect of microswimmer self-propulsion using the reduced swimmer persistence length $L$, defined by

(2.5)\begin{equation} L \equiv \frac{v_s}{\varepsilon \lambda}. \end{equation}

3. Kinetic theory

In this section, we will outline the main steps in the derivation of $D_A$ for a suspension of Brownian tracer particles immersed in a dilute microswimmer suspension, whose dynamics are governed by (2.1). Just as in our previous works (Stenhammar et al. 2017; Škultéty et al. 2020), we describe the flow field measured at ${\boldsymbol r}$ due to a swimmer with position ${\boldsymbol r}_i$ and orientation ${\boldsymbol p}_i$ by a regularised dipolar flow field ${\boldsymbol u}_s({\boldsymbol r})$:

(3.1)\begin{equation} {\boldsymbol u}_s({\boldsymbol r};{\boldsymbol r}_i,{\boldsymbol p}_i) = \frac{\kappa}{8{\rm \pi}} \left[ 3\,\frac{({\boldsymbol p}_i \boldsymbol{\cdot} {\boldsymbol r}')^2 {\boldsymbol r}' + \varepsilon^2 ({\boldsymbol p}_i \boldsymbol{\cdot} {\boldsymbol r}'){\boldsymbol p}_i }{(r'^2 + \varepsilon^2)^{5/2}} - \frac{{\boldsymbol r}'}{(r'^2 + \varepsilon^2)^{3/2}} \right], \end{equation}

where ${\boldsymbol r}' = {\boldsymbol r}-{\boldsymbol r}_i$, $r' = |{\boldsymbol r}'|$, and $\varepsilon$ is the regularisation length. Our starting point is the derivation of Škultéty et al. (2020), where we formulated and solved a kinetic theory describing the fluctuations of the velocity field ${\boldsymbol U}({\boldsymbol r},t)$ due to a superposition of single-swimmer flow fields. In the limit of non-interacting swimmers, which is the case that we consider here, the temporal correlations of the steady-state velocity field ${\boldsymbol U}$ measured in the lab frame, $C_U(t) \equiv \langle {\boldsymbol U}({\boldsymbol r},t) \boldsymbol {\cdot } {\boldsymbol U}({\boldsymbol r},0) \rangle$, are given by

(3.2)\begin{equation} C_U(t) = \frac{\kappa^2 n}{15 {\rm \pi}^2 \varepsilon} \int_0^{\infty} A^2(\zeta) \, {\rm e}^{-\tau}\, \mathcal{F}(L\zeta\tau) \, {\rm d} \zeta, \end{equation}

where

(3.3)\begin{equation} A(x) = \tfrac{1}{2} x^2\,K_2(x), \end{equation}

with $K_2$ being the modified Bessel function of the second kind, and

(3.4)\begin{equation} \mathcal{F}(x) = 15\,\frac{(5x^2-12)\sin x - x(x^2-12)\cos x}{x^5}, \end{equation}

defined such that $\mathcal {F}(0) = 1$. Furthermore, we used the dimensionless variables ${L = v_s/(\lambda \varepsilon )}$, $\tau = t\lambda$ and $\zeta = k\varepsilon$, where $k = |{\boldsymbol k}|$ is the wavevector magnitude. Equation (3.2) can be expressed equivalently in closed form in terms of elliptic integrals, as given by (72) of Škultéty et al. (2020). Equation (3.2) contains decorrelation of the flow field due to two separate mechanisms: exponential decay of $C_U(t)$ due to tumbling, and a more complex, oscillatory behaviour due to swimming, encoded in the function $\mathcal {F}$.

To obtain the hydrodynamic diffusivity of a passive tracer, Škultéty et al. (2020) used a stationary tracer approximation, implying that the tracer advection by the swimmer flow field is negligible compared to the self-propulsion of the swimmer. This implies that ${\boldsymbol r}_T$ remains effectively constant over the time it takes for ${\boldsymbol U}$ to relax, such that $C_T(t) = \langle {\boldsymbol U}({\boldsymbol r}_T[t],t) \boldsymbol {\cdot } {\boldsymbol U}({\boldsymbol r}_T[t=0],0) \rangle \approx \langle$ ${\boldsymbol U}({\boldsymbol r}_T[t=0],t) \boldsymbol {\cdot } {\boldsymbol U}({\boldsymbol r}_T[t=0],0) \rangle = C_U(t)$, where, in the last equality, we have made the additional assumption that the tracers are distributed homogeneously in space so that they sample an unweighted spatial average of the flow field. Thus if we can replace the correlation function $C_T$ in the co-moving tracer frame with that in the stationary lab frame, $C_U$, then we can insert (3.2) into the Green–Kubo relation (1.3) and integrate over time to yield the following expression for $D_A$, identical to (85) of Škultéty et al. (2020):

(3.5)\begin{equation} D_A = \frac{\kappa^2 n}{45 {\rm \pi}^2 \lambda \varepsilon} \int_0^{\infty} A^2(\zeta)\, \mathcal{G} ( L\zeta ) \, {\rm d} \zeta, \end{equation}

where

(3.6)\begin{equation} \mathcal{G}(x) = \frac{5}{2}\,\frac{3x + 2x^3 -3(1+x^2)\arctan x}{x^5}, \end{equation}

defined such that $\mathcal {G}(0) = 1$. Furthermore, by matching the asymptotic behaviours for $L\rightarrow 0$ and $L \rightarrow \infty$, (3.5) can be approximated by the simple expression given in (1.1)

The generalisation of (3.5) to the case of Brownian tracers might seem straightforward, but unfortunately is not. Since the approximation $C_T(t) = C_U(t)$ amounts to the tracer remaining effectively stationary during a swimmer–tracer scattering event, this approximation will, by construction, not capture any effects on $D_A$ coming from Brownian diffusion across streamlines. This can be realised easily by noticing that $C_U(t)$ is solely a property of the swimmer suspension, and will be strictly unaffected by the tracer dynamics; thus $D_A$ in (3.5) remains unaffected by the inclusion of tracer diffusion. Instead of going beyond the stationary tracer approximation, we circumvent this problem by noticing that in a suspension of non-interacting microswimmers, the dynamics of a Brownian tracer will be statistically identical to that of a non-Brownian tracer in a suspension of Brownian swimmers with the same translational diffusivity $D_0$. This is because the single-swimmer flow field in (3.1) depends solely on the separation vector ${\boldsymbol r}_T - {\boldsymbol r}_i$, implying that the flow field experienced by a diffusing tracer (i.e. noise acting on ${\boldsymbol r}_T$) is identical to that experienced by a non-Brownian tracer sampling the flow field from a swimmer with the same noise instead applied to ${\boldsymbol r}_i$. In a non-interacting microswimmer suspension, this equivalence is exact as long as the noise has zero mean and identical spectral properties, and we verify it numerically in figure 3. However, it breaks down as soon as swimmer–swimmer correlations become significant, since swimmer diffusion will then affect the magnitude of such correlations, which tracer diffusion will not.

Thus, as outlined in Appendix A, we instead calculate $C_U(t)$ for the case of a suspension of diffusing swimmers, described by the dynamics

(3.7)\begin{equation} \dot{{\boldsymbol r}}_i = v_s {\boldsymbol p}_i + \sqrt{2D_0}\,\boldsymbol{\eta},\end{equation}

including the same tumbling dynamics as before. This yields the following generalisation of (3.2):

(3.8)\begin{equation} C_U(t;D_0) = \frac{\kappa^2 n}{15 {\rm \pi}^2 \varepsilon} \int_0^{\infty} A^2(\zeta) \, {\rm e}^{-(1+\tilde{D}\zeta^2) \tau}\, \mathcal{F}(L\zeta\tau) \, {\rm d} \zeta, \end{equation}

where we have additionally defined the non-dimensional diffusivity $\tilde {D} = D_0/(\lambda \varepsilon ^2)$. Since the effect of Brownian diffusion is now fully incorporated into the properties of $C_U$, we again use the stationary-tracer approximation and insert this expression into the Green–Kubo relation (1.3), leading to the following expression for $D_A$ in the presence of Brownian diffusion:

(3.9)\begin{equation} D_A(D_0) = \frac{\kappa^2 n}{45 {\rm \pi}^2 \lambda \varepsilon} \int_0^{\infty} \frac{A^2(\zeta)}{1 + \tilde{D}\zeta^2}\,\mathcal{G}\left( \frac{L\zeta}{1 + \tilde{D}\zeta^2} \right) {\rm d} \zeta. \end{equation}

In § 4, we evaluate (3.8) and (3.9) numerically, and compare the results with direct numerical simulations of microswimmer suspensions.

4. Results and discussion

In figure 2, we begin by verifying (3.5) for the hydrodynamic diffusion coefficient $D_A$ in the limit $D_0 = 0$. In figure 2(a), we demonstrate the $\lambda ^{-1}$ dependence of $D_A$ in the shaker limit $v_s = 0$, while figure 2(b) shows its more complex dependence on $v_s$ for constant $\lambda$. Apart from numerically verifying the kinetic theory expression (3.5), these results illustrate how $D_A$ decreases abruptly due to the temporal decorrelation of the flow field ${\boldsymbol U}({\boldsymbol r}_T,t)$ induced respectively by tumbling and swimmer self-propulsion. The slight deviation between the theoretical curve and simulation results at small $v_s$ in figure 2(b) is likely due to the specific form for the short-range regularisation, which becomes important as $v_s \rightarrow 0$. In our derivation of (3.5), we use the regularised flow field (3.1), based on the regularisation first introduced by Cortez, Fauci & Medovikov (2005). In the LB simulations, we instead use a numerical interpolation scheme based on a regularisation of the $\delta$ function (Peskin 2002) acting separately on the two point forces that make up each microswimmer. Unlike the expression in (3.1), this numerical regularisation does not allow a direct mapping (or adjustment) of the regularisation length $\varepsilon$. We thus do not expect perfect agreement between kinetic theory and simulation in the low-$v_s$ regime where the short-range regularisation becomes important, and therefore treat $\varepsilon$ as a fitting parameter when comparing data from LB simulations with kinetic theory predictions. However, we find that the fitted value of $\varepsilon$ varies only slightly ($\varepsilon \in [1.9\,{\rm \Delta} l, 2.5\,{\rm \Delta} l]$) for the values of $v_s$ used throughout this work, in good accordance with the fact that the regularised $\delta$ function is interpolated over a support of $2\,{\rm \Delta} l$ in each Cartesian direction; for a more in-depth discussion of the effect of the interpolation scheme on the tracer dynamics, see de Graaf & Stenhammar (2017).

Figure 2. Effective tracer diffusion in the absence of Brownian motion. Panel (a) shows the $\lambda ^{-1}$ dependence of ${D_A(D_0=0)}$ for shakers with $v_s = 0$, and panel (b) its $v_s$ dependence at constant $\lambda = 10^{-4}$. Symbols denote simulation results, and solid lines show results from (3.5) using $\varepsilon$ as a fitting parameter. Error bars represent one standard deviation as obtained from averaging over four separate runs with different initial conditions. The results are non-dimensionalised in terms of $\kappa$ and $\varepsilon$.

In figure 3, we study the additional effect of varying the Brownian diffusion coefficient $D_0$, as encoded in (3.9). From the data in figure 3(a), it is clear that for large enough $D_0$, the active diffusivity $D_A$ decreases compared to its non-Brownian value. To enable an easier analysis of the effect of varying swimming speed, in figure 3(b) we present the same data instead plotted as a function of the reduced variables $\xi$ and $Pe$. For shakers with $L = 0$ (blue curve in figure 3b), $D_A$ is reduced compared to its non-Brownian value ($\xi < 1$) as soon as ${Pe} < 1$, reaching a value as low as $\xi = 0.2$ for ${Pe} \approx 10^{-3}$. However, for finite values of $v_s$, this effect on $D_A$ occurs for gradually lower values of $Pe$; for the fastest swimmers considered here, with $L = 4.0$, no significant reduction of $D_A$ is observed even for $Pe$ as low as $10^{-3}$. Instead, we observe a small but significant increase in the active diffusion compared to its non-Brownian value, in accordance with what was observed previously for slender swimmers by Kasyap et al. (2014); we discuss this effect further below. Crucially, a reduced persistence length of $L = 4$ nevertheless corresponds to relatively slow swimming from a biological perspective: According to the approximate calculation in Škultéty et al. (2020), $L$ for E. coli bacteria lies somewhere in the range between 5 and 20, indicating that the effect of Brownian motion on active diffusion is likely negligible in suspensions of swimming bacteria due to their fast self-propulsion. In our LB simulations, studying values higher than $L \approx 4$ is challenging, as these large swimming speeds both require very large systems to avoid significant finite-size effects and yield artefacts due to the effect of finite Reynolds number in the swimmer–tracer scattering dynamics (de Graaf & Stenhammar 2017). We nevertheless studied $\xi ({Pe})$ numerically using (3.9) for larger values of $L$, verifying that both the peak and the subsequent decrease in $\xi$ continue to move to even lower values of $Pe$ as $L$ is increased.

Figure 3. Brownian motion suppresses active diffusion for slow swimming speeds. (a) Values of $D_A$ measured from LB simulations (diamonds) and calculated from (3.9) (solid lines), both expressed in LB units. (b) The same data but expressed in the dimensionless quantities $\xi$ and $Pe$. For very slow swimmers with $L \lesssim 1$, $D_A$ is reduced compared to the non-Brownian value ($\xi = 1$) whenever ${Pe} < 1$, while for faster swimmers, significantly lower values of $Pe$ are necessary to affect $D_A$. The circles for $L = 0$ and 2.4 correspond to the hydrodynamic diffusion of non-Brownian tracers measured in a suspension of Brownian swimmers of the same $D_0$, verifying the statistical equivalence between tracer and swimmer diffusion in the non-interacting limit. Error bars represent one standard deviation as obtained from averaging over four separate runs with different initial conditions. The results in (a) are non-dimensionalised using $\kappa$ and $\varepsilon$.

To understand further the mechanism behind the reduction in $D_A$ with $D_0$, we consider the two autocorrelation functions $C_T(t)$ and $C_U(t)$, which measure the fluid autocorrelation in the co-moving tracer frame and in the lab frame, respectively. Figure 4 shows these correlation functions for $L = 0$ and $L = 2.4$, with figures 4(a,b) corresponding to LB results for $C_T$, and figures 4(c,d) to kinetic theory results from (3.8) for $C_U$. First, we notice that the two sets of curves are very similar, implying that the stationary tracer approximation $C_T \approx C_U$ is indeed accurate. Second, we notice that the decay of the correlation function is significantly faster for swimmers than for shakers, again illustrating that self-propulsion acts an efficient decorrelation mechanism for ${\boldsymbol U}$. The effect of finite $D_0$ for shakers (figures 4a,c) is simply to decrease the relaxation time of the exponential decay, in accordance with the $L = 0$ limit of (3.8). For swimmers, the situation is, however, more complex: For short times, the flow field decays faster with decreasing $Pe$, while the long-time tail of $C_T$ and $C_U$ instead becomes somewhat more extended with decreasing $Pe$. For fast enough swimmers, the latter effect leads to the local maximum at $\xi > 1$ for intermediate $Pe$ observed in figure 3(b) for the two highest $L$. Finally, we note that the equal-time fluid velocity variance $\langle U^2({\boldsymbol r}_T) \rangle$, corresponding to the $t = 0$ values of $C_T$ and $C_U$, is independent of $Pe$. This means that, regardless of the ratio between diffusive and active motion, the tracer particles sample the overall flow field homogeneously. This fact is non-trivial, since dry active particles that move autonomously on a solid substrate are known to preferentially sample regions where they move slowly (Stenhammar et al. 2016). Our results thus highlight that this generic mechanism is absent for inertialess point tracers advected by an incompressible fluid. However, for microswimmer systems where entrainment due to tracers being captured by the near-field flows of passing swimmers (Jeanneret et al. 2016), we would expect the fluid flow sampled by tracers to be significantly different from the average flow field in the system.

Figure 4. Brownian motion decorrelates tracer trajectories. (a,b) Fluid velocity autocorrelation $C_T(t)$ in the co-moving tracer frame measured from LB simulations for (a) shakers ($L=0$), and (b) swimmers with $L = 2.4$ at indicated values of $Pe$. The dotted line shows the correlation function $C_U(t)$ of the fluid velocity in the lab frame, demonstrating that the stationary-tracer approximation $C_T(t) \approx C_U(t)$ is excellent in the absence of Brownian tracer diffusion ($Pe \rightarrow \infty$). (c,d) Corresponding lab-frame correlation functions $C_U(t)$ obtained from kinetic theory (3.9) for a suspension of diffusing swimmers, as described in § 3. Throughout, $C_T$ and $C_U$ are non-dimensionalised using $\kappa$ and $\varepsilon$.

In figure 5, we study the dependence of the suppression of active diffusion on the microswimmer density $n$ in the shaker limit $L=0$. At first sight, this dependence might appear trivial, since $D_A$ is well known to be linearly dependent on $n$ (Jepson et al. 2013; Miño et al. 2013) in the limit of non-interacting swimmers, a fact that is unaffected by Brownian motion as shown by (3.9). Since $\xi$ measures the ratio between the Brownian and non-Brownian values of $D_A$, one would naively expect $\xi$ to be independent of $n$. However, since the Péclet number itself, as defined in (2.3), increases with $n$ for constant $D_0$, the relative effect of Brownian motion on $D_A$ is in fact a complex function of $n$ even for non-interacting swimmers, as shown in figure 5. More specifically, according to figure 5(b), the suppression of active diffusion becomes more significant with increasing microswimmer density. In physical units, the highest concentration considered here (${n\varepsilon ^3 \approx 0.08}$) corresponds approximately to a bacterial concentration of $10^9$ ml$^{-1}$, which is somewhat higher than the highest E. coli concentration considered by Jepson et al. (2013), but still within the range of concentrations where swimmer–swimmer correlations are reasonably small (Stenhammar et al. 2017). In summary, our results thus show that for Brownian diffusion to have any measurable influence on hydrodynamic diffusion, it is necessary to create a system with a relatively high density of very slow microswimmers; as we discuss in § 5, this set of parameters is likely not achievable for suspensions of biological microswimmers.

Figure 5. Reduction of $D_A$ varies with microswimmer density. (a) Plots of $\xi$ as a function of $Pe$ for shakers ($L=0$) at various densities $n$, as indicated. (b) Plot of $\xi$ as a function of $n$, at a fixed value of ${Pe} = 0.013$. Simulation data are given by the symbols, with error bars obtained from averaging over four separate runs with different initial conditions, while solid lines are computed from (3.9).

5. Summary and conclusions

In this study, we have demonstrated a number of theoretical and computational results regarding the effect of Brownian diffusion on the swimmer-induced hydrodynamic diffusion of tracer particles in a suspension of dipolar microswimmers. Our key finding is that the effect of Brownian diffusivity $D_0$ on the activity-induced, hydrodynamic diffusivity $D_A$ is significant only when the Péclet number as defined by (2.3) is below unity, meaning that Brownian diffusivity needs to dominate hydrodynamic diffusivity. However, the necessary requirement ${Pe} < 1$ is sufficient only in the shaker limit $v_s \rightarrow 0$: for swimmers with persistence lengths larger than the organism size, significantly lower values of $Pe$ are required to perturb the tracer trajectories sufficiently to affect $D_A$. This conclusion is analogous to the independence of $D_A$ from the tumbling rate $\lambda$ for large $v_s$ illustrated in (1.1). Whenever $v_s$ is large, the decorrelation of ${\boldsymbol U}$ by swimming will dominate the decorrelation due to tumbling and translational diffusion, and the dependence on $\lambda$ and $D_0$ will thus vanish in the limit $v_s \rightarrow \infty$. While this effect is expected, what is perhaps surprising is the rather moderate values of $L = v_s / (\varepsilon \lambda )$ necessary to render the coupling between $D_0$ and $D_A$ negligible, as illustrated in figure 3(b). To put these values into perspective, we use the conservative estimate $L = 5$ for E. coli. By virtue of figure 3(b), we require that ${Pe} \leq 10^{-3}$ for Brownian motion to have an effect of ${\sim }5$ % on $D_A$. Using as an example the minimum value $D_A \approx 10^{-2}\ \mathrm {\mu }{\rm m}^2\ {\rm s}^{-1}$ measured by Jepson et al. (2013) in a three-dimensional E. coli suspension, this requirement thus implies that $D_0 \geq 10\ \mathrm {\mu }{\rm m}^2\ {\rm s}^{-1}$, which by virtue of the Stokes–Einstein relation corresponds to a tracer radius $R_0 \sim 20$ nm. While significantly smaller than the particles used in typical measurements on colloidal tracers (Leptos et al. 2009; Miño et al. 2013; Patteson et al. 2016), this diffusion coefficient is close to the value of $D_0$ measured for dextran in E. coli suspensions by Kim & Breuer (2004). It is also fully feasible to realise such low Péclet numbers for micron-sized spheres by instead decreasing the bacterial density to very low values; however, measuring the correction to $D_A$ for ${Pe} \sim 10^{-3}$ represents a major difficulty, since it amounts to measuring a $\sim$5 % deviation of an effective diffusivity that is itself a thousand times smaller than the Brownian diffusion. It would thus require an extremely accurate determination of $D_0$, which then needs to be subtracted from the total measured diffusion constant to determine $D_A$. Obtaining this accuracy in a colloidal suspension would be very challenging due to particle polydispersity, interactions with boundaries, temperature gradients, and other system-specific complications. We thus conclude that, for typical $L$ values relevant for biological microswimmers, our results imply that the effect of Brownian motion on $D_A$ is likely negligible for all practical purposes. To observe experimentally the reduction in $D_A$, one would instead need to study a system of dipolar shakers, which stir up the surrounding fluid without self-propelling. While this is a somewhat exotic type of system, it could potentially be realised by anchoring molecular motors or biological microswimmers to a surface. In a biological setting, the shaker limit furthermore resembles previously developed models of enzymes anchored to lipid bilayers that induce dipolar flows through cyclical conformation changes (Hosaka, Komura & Mikhailov 2020).

Our results are qualitatively consistent with the previous theoretical results obtained by Kasyap et al. (2014), including the non-monotonic behaviour of $D_A$ with respect to $D_0$ at high swimming speed. Their results, however, differ in three important ways. First, they consider a ‘slender swimmer’ model, where the fluid is forced via a stress applied along a line representing the swimmer body, rather than by two point forces. While this model should lead to a dipolar flow in the far field, the near-field differences are significant, as illustrated in figure 8 of Kasyap et al. (2014). Furthermore, they consider only the fast-swimming limit where, according to (1.1), $D_A$ for dipolar microswimmers is independent of $\lambda$. Finally, they parametrise their model in a qualitatively different way to us: as discussed above, they adopt a definition of ${Pe}$ based on $v_s$ rather than on the tracer diffusivity as in (2.3). In this description, $v_s$ is furthermore directly coupled to $\kappa$, so that changing $Pe$ simultaneously changes the activity of the bath (via $\kappa$) and the swimming speed $v_s$. These differences make it difficult to compare directly with our results, as we consider the effects of fluid advection and self-propulsion separately via $Pe$ and $L$, respectively.

Thus, even though Brownian motion is unlikely to provide a significant dependence of $D_A$ on $R_0$ for tracers immersed in suspensions of biological microswimmers, there are several other mechanisms that need to be studied to explain the non-monotonic dependence observed experimentally (Patteson et al. 2016) and computationally (Dyer & Ball 2021). First, the effect of tracer entrainment by the near field flows of the swimmer is strongly dependent on the size ratio between the swimmer and the tracer (Jeanneret et al. 2016), although we expect this term to be small for micron-sized tracer particles in E. coli suspensions. Second, the finite size of the tracer will change the equation of motion (2.1) into the Faxén equation that takes into account the nonlinearity of the flow field (Kim & Karrila 1991), an effect that was included implicitly in the wavelet Monte Carlo simulations by Dyer & Ball (2021), and, together with tracer entrainment, is a significant explanation of their observed $R_0$ dependence of $D_A$. Finally, non-hydrodynamic interactions such as direct collisions, electrostatic interactions and artefacts due to container walls are likely to depend in a non-trivial manner on the tracer size for each system in question. We thus conclude by noting that significant further experimental and theoretical work is necessary to disentangle the system-specific properties from the generic properties of tracer dynamics in microswimmer suspensions.