Anisotropic diffusion of ellipsoidal tracers in microswimmer suspensions

Tracer particles immersed in suspensions of biological microswimmers such as E. coli or C. reinhardtii display phenomena unseen in conventional equilibrium systems, including strongly enhanced diffusivity relative to the Brownian value and non-Gaussian displacement statistics. In dilute, 3-dimensional suspensions, these phenomena have typically been explained by the hydrodynamic advection of point tracers by isolated microswimmers, while, at higher concentrations, correlations between pusher microswimmers such as E. coli can increase the effective diffusivity even further. Anisotropic tracers in active suspensions can be expected to exhibit even more complex behaviour than spherical ones, due to the presence of a nontrivial translation-rotation coupling. Using large-scale lattice Boltzmann simulations of model microswimmers described by extended force dipoles, we study the motion of ellipsoidal point tracers immersed in 3-dimensional microswimmer suspensions. We find that the rotational diffusivity of tracers is much less affected by swimmer-swimmer correlations than the translational diffusivity. We furthermore study the anisotropic translational diffusion in the particle frame and find that, in pusher suspensions, the diffusivity along the ellipsoid major axis is higher than in the direction perpendicular to it, albeit with a smaller ratio than for Brownian diffusion. Thus, we find that far field hydrodynamics cannot account for the anomalous coupling between translation and rotation observed in experiments, as was recently proposed. Finally, we study the probability distributions (PDFs) of translational and rotational displacements. In accordance with experimental observations, for short observation times we observe strongly non-Gaussian PDFs that collapse when rescaled with their variance, which we attribute to the ballistic nature of tracer motion at short times.


I. INTRODUCTION
Active transport of particles is of importance in many biological contexts, such as intracellular transport [1], absorption of nutrients in intestines [2] and by microorganisms [3], and possibly mass transport in oceans [4]. It also has a more fundamental relevance, as tracers can be used to probe the nature of non-equilibrium fluctuations in active systems, in particular in relation to (equilibrium) Brownian motion [5]. Arguably, the simplest model system for studying tracer diffusion in such "active baths" is a spherical colloidal particle immersed in a dilute suspension of biological microswimmers, such as E. coli [6][7][8][9][10][11][12][13][14][15] or C. reinhardtii [16][17][18][19]. A generic feature of tracer motion in these systems is that the displacement is characterized by ballistic motion at short times (approximately 0.02 − 2 seconds [13]) and diffusive motion at long times, although with a significantly increased translational diffusion coefficient D T compared to the corresponding Brownian one. Furthermore, in the concentration regime where swimmer-swimmer correlations are sufficiently weak, D T displays a linear scaling with swimmer density [8,11,16]. For dilute, 3dimensional microswimmer suspensions, this enhanced tracer diffusion has been rationalised from the superposition of independent hydrodynamic swimmer-tracer scattering events due to the swimmers' long-ranged dipolar * joakim.stenhammar@fkem1.lu.se flow fields [8,[20][21][22][23][24][25], while, at much higher densities, for finite-size particles, or in confined geometries, direct swimmer-tracer collisions [26][27][28] and near-field hydrodynamic effects [13,29] also become significant. Furthermore, for the case where microswimmers are significantly larger than the tracer particles, so-called tracer entrainment can occur, leading to isolated, very large displacements of individual tracers and a qualitative change in the displacement statistics [25,30,31]. Another wellstudied non-equilibrium feature of tracer dynamics in dilute microswimmer suspensions is the probability distribution function (PDF) of tracer displacements within a fixed time window. For Brownian diffusion, these PDFs are Gaussian, while in biological microswimmer suspensions they become strongly non-Gaussian at low swimmer densities and short observation times [16,17,32]. This effect is explained by the fact that tracer displacements result from just a small number of swimmer-tracer scattering events, so that the central limit theorem does not apply [33][34][35], and has also been reproduced in computational models of microswimmers [21,23,35].
l v s p s −F p s F p s Figure 1. Schematic illustration of a pusher-type model swimmer and its parameters.
Even though the enhanced tracer diffusion at sufficiently low microswimmer densities can be approximated as arising from independent swimmer-tracer scattering events, swimmer-swimmer correlations will strictly affect the enhanced diffusion for any finite microswimmer concentration [36,37]. When such correlations are taken into account, the pusher-puller symmetry and the linear scaling of D T with swimmer density are broken: In puller suspensions, such as C. reinhardtii, the scaling progresses slower than linearly, while in pusher suspensions, such as E. coli, the increase becomes superlinear [36][37][38][39]. At a well-defined critical density of pushers, D T diverges, corresponding to the onset of active turbulence, a phenomenon in which the system exhibits large-scale vortices and jets, with fluid velocities surpassing the swimming speed of an individual swimmer. [36,[40][41][42].
Going beyond the case of spherical tracers, a few recent experimental studies have focused on the dynamics of anisotropic passive particles in the form of ellipsoids [14,18] or dumbbells [19] immersed in active suspensions consisting of swimming bacteria or algae confined to liquid films. In addition to the enhanced overall diffusion seen for spherical tracers, these studies found significantly increased rotational diffusion coefficients D R . More surprisingly, they found that ellipsoidal tracers display a qualitatively anomalous anisotropic diffusion compared to the Brownian case, in that the ratio D /D ⊥ of diffusion coefficients parallel and perpendicular to the particle major axis is below unity in the case of pusher-type swimmers (E. coli ) at high densities [14]. This contrasts with the Brownian diffusion of ellipsoids in bulk suspensions, where D /D ⊥ → 2 for large aspect ratios due to their anisotropic friction [43][44][45]. For puller-type swimmers (C. reinhardtii ), the same authors instead found that the corresponding ratio remains above unity [18]. These results were partially rationalised as due to the different symmetries of pusher and puller flow fields in quasi-2D geometries [14,18], but their precise origin is still unclear.
To build a better theoretical understanding of anisotropic tracer dynamics in active suspensions, we will here consider a simple computational model of ellipsoidal tracers immersed in 3-dimensional microswimmer suspensions. In terms of their translational dynamics, tracers are described as point particles, while the orientational tracer dynamics are governed by Jeffery's equation for extended ellipsoids. The swimmers are modelled as extended force dipoles, which affect the tracers only via their flow fields, which we implement using a lattice Boltzmann framework. The simplicity of the model allows for large-scale particle resolved simulations with N > 10 5 swimmers and N = 10 5 tracers, while accurately capturing the far-field contribution to the tracer dynamics.
Our main finding is that translational diffusion in the particle frame is anisotropic only for the case of pusher suspensions, and furthermore only at densities around and above the transition to active turbulence. In con-trast to experimental results for ellipsoidal tracers in E. coli suspensions [14], we find that the ratio D /D ⊥ is above unity, indicating that far-field hydrodynamics cannot account for the anomalous coupling between translation and rotation seen in experiments. 2 Even though the enhanced tracer diffusion at sufficiently low microswimmer densities can be approximated as arising from independent swimmer-tracer scattering events, swimmer-swimmer correlations will strictly affect the enhanced diffusion for any finite microswimmer concentration [35,36]. When such correlations are taken into account, the pusher-puller symmetry and the linear scaling of D T with swimmer density are broken: In puller suspensions, such as C. reinhardtii, the scaling progresses slower than linearly, while in pusher suspensions, such as E. coli, the increase becomes superlinear [35][36][37][38]. At a well-defined critical density of pushers, D T diverges, corresponding to the onset of active turbulence, a phenomenon in which the system exhibits large-scale vortices and jets, with fluid velocities surpassing the swimming speed of individual swimmer [35,[39][40][41].
Going beyond the case of spherical tracers, a few recent experimental studies have focused on the dynamics of anisotropic passive particles in the form of ellipsoids [14,18] or dumbbells [19] immersed in active suspensions. In addition to the enhanced overall diffusion seen for spherical tracers, these studies found significantly increased rotational diffusion coefficients D R . More surprisingly, they found that ellipsoidal tracers display a qualitatively anomalous anisotropic diffusion compared to the Brownian case, in that the ratio D /D ⊥ of diffusion coefficients parallel and perpendicular to the particle major axis is below unity in the case of pushertype swimmers (E. coli ) at high densities [14]. This contrasts with the Brownian diffusion of ellipsoids in bulk suspensions, where D /D ⊥ → 2 for large aspect ratios due to their anisotropic friction [42][43][44]. For pullertype swimmers (C. reinhardtii ), the same authors instead found that the corresponding ratio remains above unity [18]. These results were partially rationalised as due to the different symmetries of pusher and puller flow fields in quasi-2D geometries [14,18], but their precise origin is still unclear.
To build a better theoretical understanding of anisotropic tracer dynamics in active suspensions, we will here consider a simple computational model of ellipsoidal tracers immersed in 3-dimensional microswimmer suspensions. In terms of their translational dynamics, tracers are described as point particles, while the orientational tracer dynamics are governed by Jeffrey's equation for extended ellipsoids. The swimmers are modelled as extended force dipoles, which affect the tracers only via their flow fields, which we implement using a lattice Boltzmann framework. The simplicity of the model allows for large-scale particle resolved simulations with N > 10 5 swimmers and N = 10 5 tracers, while accurately capturing the far-field contribution to the tracer dynamics.
Our main finding is that translational diffusion in the particle frame is anisotropic only for the case of pusher suspensions, and furthermore only at densities around and above the transition to active turbulence. In contrast to experimental results for ellipsoidal tracers in E.
coli suspensions [14], we find that the ratio D /D ⊥ is above unity, indicating that far-field hydrodynamics cannot account for the anomalous coupling between translation and rotation seen in experiments.

II. MODEL AND METHOD
We consider a three-dimensional suspension containing N s swimmers and N t = 10 5 tracer particles, moving in a box with periodic boundaries. Swimmers are represented by extended force dipoles, where the body and flagella exert two equal and oppositely directed forces ±F p s separated by a distance l on the fluid, where the unit vector p s represents the orientation of the swimmer (see Fig. 1). Each swimmer is characterized by its dipole strength κ = ±F l/µ, where µ is the dynamic viscosity

II. MODEL AND METHOD
We consider a three-dimensional suspension containing N s swimmers and N t = 10 5 tracer particles, moving in a box with periodic boundaries. Swimmers are represented by extended force dipoles, where the body and flagella exert two equal and oppositely directed forces ±F p s separated by a distance l on the fluid, where the unit vector p s represents the orientation of the swimmer (see Fig. 1). Each swimmer is characterized by its dipole strength κ = ±F l/µ, where µ is the dynamic viscosity of the fluid, and the sign of κ distinguishes pushers from pullers, represented by κ > 0 and κ < 0, respectively.
The time evolution of each swimmer's position r s and orientation p s is governed by the following equations of motion (EOMs) [36,40]: where U(r s ) is the velocity of the fluid at the position of the swimmer's body, v s is the swimming speed of an individual swimmer, and I is the unit tensor. As discussed in [36], the orientational EOM (2) describes the dynamics of a flow-aligning swimmer composed of a spherical body of radius a connected by a thin rod of length l to the flagellum, in the limit a l (see Fig. 1). In addition to the reorientation caused by hydrodynamic interactions, swimmers undergo run-and-tumble motion, implemented by imposing Poisson distributed random reorientations with an average frequency λ. In the absence of other microswimmers, these dynamics yield a persistent random walk with the persistence length v s /λ. Since the effective diffusion of tracer particles in microswimmer suspensions is independent of λ for persistence lengths v s /λ l [23,46], we for simplicity choose to keep λ fixed at a value close to that of wild-type E. coli, corresponding to v s /λl = 5.
x y z p t θ y θ z Figure 3. Illustration of the rotation scheme used for evaluating particle displacements in the particle coordinate frame.
The EOMs for the positions r t and orientations p t of the ellipsoidal point tracers are given bẏ where Eq. (4) is Jeffery's equation, describing the reorientation of an ellipsoidal particle in a shear flow. Here, β = (q 2 − 1)/(q 2 + 1) with q denoting the ellipsoid aspect ratio, and E = (∇U + ∇U † )/2 and W = (∇U − ∇U † )/2 are, respectively, the rate-of-strain and vorticity tensors evaluated at r t . For simplicity we set β = 1, corresponding to infinite aspect ratio, although we also checked that our results only depend weakly on the value of β for aspect ratios q ≥ 2. In this limit, Eq. (4) simplifies tȯ where the gradient (∇U) · p t along the tracer orientation was numerically evaluated as the first-order central difference evaluated at the tracer position, similar to the swimmer EOM (2). Notably, the translational EOM (3) is identical to that of a point tracer, so that the anisotropic behavior is fully encoded in the orientational EOM (5). This is in contrast to Brownian motion of elongated particles, which is intrinsically anisotropic due to the different friction coefficients along the minor and major axes. As Brownian motion is not explicitly included in our simulations, all observed anisotropic behaviour is thus solely due to this intricate translation-rotation coupling. The translational and rotational displacements of tracers were tracked both in the lab frame and by considering a separate, co-moving frame of reference for each individual particle. Our approach is equivalent to that of Ref. [44], but extended from two to three dimensions. We start from the tracer position and orientation in the lab frame at time step t n , denoted by r t (t n ) and p t (t n ), respectively. During a time interval δt = t n − t n−1 , the particle undergoes a translation δr t,n = r t (t n )−r t (t n−1 ), which is transformed into its body frame counterpart δr t,n by applying two subsequent rotations (Fig. 3): with R y =   cos θ y,n 0 sin θ y,n 0 1 0 − sin θ y,n 0 cos θ y,n   , Here, θ i,n = [θ i (t n−1 ) + θ i (t n )]/2, θ y is the angle between p t and the lab xy plane and θ z is the azimuthal angle between the x-axis and a projection of p t onto the xy plane (see Fig. 3). The total body frame displacement of the tracer during a macroscopic time interval t n is obtained by summing over all displacements δr t,n : Using Eq. (9), we consider body frame displacements for trajectories of duration ∆t via ∆r t (∆t) =r t (t 0 + ∆t) − r t (t 0 ). The first component of this quantity represents the displacement along the tracer major axis, with the other two components representing the displacement in the plane perpendicular to this axis.  For the evaluation of translational diffusion coefficients, mean-square displacements (MSDs) in the lab frame and in the particle frame were obtained by averaging over all particle trajectories with a given duration ∆t. Due to the persistent character of the fluid flows, the MSDs can be fitted to a persistent random walk with a ballistic regime at short times and a diffusive regime beyond a crossover time τ T , resulting in [47] where the MSD was obtained by averaging over all tracer particles and all time origins, after removing initial transients due to relaxation of the system to steady state. The rotational diffusion coefficient D R was calculated by fitting the orientational autocorrelation function C(∆t) ≡ p(t 0 ) · p(t 0 + ∆t) to the approximate expression suggested by Wilkinson and Pumir [48], which interpolates between the correct ballistic (C = exp −(D R /τ R )∆t 2 ) and diffusive (C = exp [−2D R ∆t]) behaviours: The numerical evaluation of fluid flows were implemented through the D3Q15 BGK lattice Boltzmann (LB) method as developed by Nash et al. [49,50], and further discussed in [36,40]. We employed a cubic box of size (100 3 ) lattice points, which is large enough to minimise finite-size effects even in the active turbulence regime [40], and a simulation length of up to 5 × 10 5 time steps. In terms of LB units, defined by the lattice spacing ∆L and time step ∆t, the swimmer parameters were set to v s = 10 −3 , F = 1.57 · 10 −3 , l = 1, λ = 2 · 10 −4 , and µ = 1/6, where the latter value corresponds to the fluid relaxing to local equilibrium on each timestep. The resulting (non-dimensionalised) value of κ corresponds closely to that measured for E. coli [7,40]. In the following, rather than using LB units, we will present all results in terms of the swimmer length l and the swimming time scale l/v s . These can in turn be related to physical units by rescaling with the corresponding dipole lengths and swimming speeds of the experimental system in question: for E. coli, l ≈ 2 µm, v s ≈ 20 µm/s, and F ≈ 0.4 pN [7]. The model introduced in the previous Section has been extensively used to study the onset of collective motion in dilute suspensions of pusher-like microswimmers [36,40,51], and its phenomenology is wellunderstood. Below the threshold value n c of the number density n = N/V of microswimmers, the suspension remains homogeneous and isotropic while exhibiting significant correlations between microswimmers. As the number density approaches its critical value, the interparticle correlations diverge, culminating in a collective state. The associated velocity fields are correlated across the whole domain and comprise large-scale jets and vortices [40,42,52,53]. In puller-like suspensions, on the other hand, long-range hydrodynamic interactions suppress collective motion and yield suspensions less correlated than their non-interacting counterparts.
In Fig. 2a, we show representative trajectories of a single ellipsoidal tracer immersed in a suspension of pusher microswimmers below (n = 0.01) and above (n = 0.2) the onset of collective motion (see also movies in [54]). For both densities, the tracer exhibits what looks like a "diffuse-and-tumble" motion: its center of mass and its orientation perform random walks punctuated by rare significant re-orientation events. The origin of these dynamics lies in the decay law of the velocity fields gener-ated by individual microswimmers. The tracer advection due to a superposition of dipolar fields (∼ r −2 ) generated by the whole suspension is a truly long-range effect in 3D. The tracer rotation, on the other hand, is only marginally long-range since the superposition of the associated velocity gradients (∼ r −3 ) diverges as a logarithm of the system size.
To further visualise this difference, in Fig. 2b, we plot U and its gradient in the dilute regime, n = 0.001. The fluid velocity U exhibits patches of significant magnitude, though much smaller than v s , embedded in a background of small but finite amplitude. The velocity gradient, on the other hand, has significant values only in the direct vicinity of individual microswimmers. A tracer sampling such fields, performs, effectively, a translational random walk with a constant orientation until a close encounter with a microswimmer (rare in dilute suspensions) changes its orientation significantly. Such sudden reorientations due to close encounters with individual microswimmers can be seen as rotational analogues of entrainment events disscussed above [30,46].
Above the onset of collective motion (n = 0.2, Fig. 2c) the fluid velocity is significantly larger than v s and is correlated across the whole domain. The resulting translational motion is again a random walk with persistence length set by the magnitude of the fluid velocity and the correlation time of the fluid flow, which is signficantly longer than the tumble time λ −1 of individual microswimmers [37,40], yielding a significantly enhanced translational diffusivity D T . On the other hand, the velocity gradient is still correlated over short distances only and its maximum values are only somewhat bigger than their dilute counterparts in Fig. 2b. This again implies rare reorientation events interspaced with periods of weak rotational random walk.
To further quantify these observations, in Fig. 4 we show the tracer diffusion coefficients D T and D R and corresponding crossover times τ T and τ R as functions of the swimmer concentration n, all measured in the lab coordinate system. In all cases, results are presented for suspensions of pushers, pullers, and non-interacting swimmers, where the latter refers to simulations in which hydrodynamic interactions between swimmers have been disabled by setting all terms containing the fluid velocity U in their EOMs (1)-(2) to zero, while keeping the tracer EOMs (3)-(5) unchanged; note that, in this limit, pushers and pullers are statistically equivalent [40]). Thus, at the lowest swimmer concentrations (n < 10 −2 ), where swimmer-swimmer correlations are small, D T and D R are both equivalent for the three swimmer types, and all show a linear increase with n; for D T , this is in accordance with previous results from experiments [6,8,11,16,17], simulations [21,36,39], and theory [23,35,37,46].
For intermediate concentrations (0.01 < n < 0.2), significant deviations from the linear dependence develop due to swimmer-swimmer correlations, with pushers showing a more steep increase, and pullers a slower increase than the noninteracting swimmers [36,38,39]; for pusher suspensions, this corresponds to the buildup of correlations leading up to the transition to active turbulence. Above the transition (n ≥ 0.2), the increase in both D T and D R becomes less steep, indicating a different scaling behavior inside the turbulent regime. While this behaviour for D T is similar to the behaviour of the RMS fluid velocity U 2 studied in [40], the three corresponding regimes for D R are somewhat less well-defined. It is furthermore clear from Fig. 4ab that the relative effect of correlations is approximately an order of magnitude smaller for D R than for D T , in accordance with the qualitative observations in Fig. 2.
Furthermore, the corresponding persistence times τ T and τ R (Fig. 4cd) for pushers show sharp peaks at concentrations corresponding to the transition to active turbulence (n ≈ 0.15), with the plot for τ T highly similar to the one for the persistence time of the fluid velocity in [40]. Notably, however, the peak in τ R is significantly flatter and shifted towards higher densities than for τ T , highlighting the complex nature of the flow fields near the transition from disordered to collective motion. The latter peak furthermore shows values of the persistence time approximately ten times smaller than for translation, indicating that velocity gradients decorrelate faster than the velocity field itself.
To quantify the anisotropy of the translational diffusion, in Fig. 5 we show the ratio D /D ⊥ of diffusion coefficients along the tracers' major and minor axes, as described in Section II. First, we observe that suspensions of pullers and non-interacting swimmers display no measurable anisotropic diffusion, since the ratio D /D ⊥ remains close to unity across the whole range of swimmer concentrations. Pusher suspensions, on the other hand, exhibit a steep increase in this ratio around n = 0.15, before stabilizing at a value of ∼ 1.1. This sharp increase coincides with the transition to active turbulence (c.f. Fig. 4ab), and thus appears to be a signature of this transition. Notably, the abrupt change in the ratio D /D ⊥ appears significantly sharper than the more gradually growing deviation from linear behaviour seen in Fig. 4ab. Furthermore, although the ratio D /D ⊥ exceeds unity, corresponding to a higher diffusion coefficient parallel to the major axis than perpendicular to it, it is still significantly smaller than the value D /D ⊥ = 2 expected for Brownian diffusion [43,45]. This is not suprising per se, as the origin of the latter value is the anisotropic drag of an elongated particle, which does not influence the hydrodynamic advection underlying the increased active diffusion studied here. The behaviour is however qualitatively different from what was experimentally observed in Peng et al. [14] for ellipsoidal particles immersed in a quasi-2D E. coli suspension, where D /D ⊥ < 1 at densities corresponding to collective motion; we discuss this apparent discrepancy further below.
In Fig. 6, we show the rescaled probability distributions (PDFs) of translational and rotational tracer displacements ∆r α t and ∆θ z as a function of the observation time ∆t. For short observation times and low densities, the PDFs are strongly non-Gaussian, in accordance with experimental results [16,17]. Theoretically, this effect can be attributed to the fact that significant displacements are driven by only a small number of swimmer-tracer encounters within the observation window, so that the central limit theorem does not apply [33][34][35]. Notably, the PDFs of rotational displacements deviate significantly more from the Gaussian form than the translational counterpart at the same density, again showing that rotational diffusion is driven by fewer swimmer-tracer scattering events than translational diffusion. This is in line with both the observed tracer trajectories (Fig. 2a) and the relatively small enhancement of the rotational diffusion observed in Fig. 4. For short observation windows ∆t, the non-Gaussian PDFs furthermore collapse when rescaled with their standard deviations σ, while, as ∆t is increased, both the translational and rotational distributions start approaching Gaussian distributions due to the growing number of swimmertracer scattering events within the time window. This data collapse is qualitatively in accordance with previous experimental [16,17] and theoretical [21,23,35] results, where it has been observed in both the diffusive (large ∆t) and the ballistic (short ∆t) regimes. In Fig. 6, the collapse generally appears to occur in the regime ∆t < τ T /R , i.e., as long as the length of the observation window is within the ballistic regime of the translational or rotational motion. This observation can be simply understood as follows: For purely ballistic translational motion, i.e., for ∆t τ T , the displacement ∆r t is given by ∆r t ≈ U(r t )∆t. When rescaled by ∆t, the PDF of the tracer displacement in the ballistic regime therefore has to collapse onto that of the fluid velocity. For the individual components of ∆r t and U, we thus expect that P (∆r α t /∆t) = P (U α ), with α denoting Cartesian components. For rotational displacements, the same scaling argument holds for the short-time PDFs of the orientation displacement, P (∆p α t /∆t). As per Eq. (5), these should furthermore coincide with the PDF of (δ αβ −p α p β )∂ γ U β n γ , i.e., the components of the velocity gradient tensor projected onto randomly oriented unit vectors n with random positions in the simulation box. From the above reasoning, it obviously follows that the tracer displacement PDFs in the ballistic regime should collapse also when rescaled by their standard deviation, as in Fig. 6. This argument shows that, for small enough ∆t, the "diffusive" collapse of the PDFs observed in Fig. 6 is in fact a consequence of ballistic particle motion at short times. This picture is verified in Fig. 7: For very short times, the collapse of P (∆r α t /∆t) and P (∆p α t /∆t) for different ∆t is essentially perfect, while it gradually worsens as ∆t approaches τ T or τ R , as expected.
It should finally be noted that the above collapse mechanism is different from the one studied theoretically in [35], which focusses on time windows within the diffusive regime of the MSD, i.e., for ∆t > τ T . While experimental studies on tracers in C. reinhardtii suspensions [16,17] have focussed on the diffusive collapse, the  squirmer simulations by Lin et al. [23] showed a corresponding data collapse also in the short-time, ballistic regime. Short-time ballistic tracer MSDs have furthermore been observed in most experimental studies of tracer dynamics in E. coli suspensions [6,11,13,26], with measured values of τ T between 0.02 and 2 seconds. These values are comparable to our E. coli -like system (Fig. 4c) and well within the accessible range of particle tracking experiments. Thus, verification of this "ballistic collapse" in E. coli suspensions would be an interesting topic for further experimental investigations.

IV. SUMMARY AND CONCLUSIONS
In this study, we have computationally investigated the dynamics of ellipsoidal tracer particles in 3-dimensional microswimmer suspensions, using large-scale lattice Boltzmann simulations. Our model accurately incorporates far-field hydrodynamics due to dipolar flow fields which have previously been shown to dominate the enhanced diffusion of tracer particles in 3-dimensional E. coli suspensions [8].
In line with previous theoretical and computational results [36][37][38], we found positive deviations from a linear dependence of D T and D R on density in pusher suspensions, coinciding with the onset of significant swimmerswimmer correlations. Notably, D T deviates from linear behaviour by an order of magnitude more than D R , indicating that the fluid velocity, which advects the particles, is significantly more affected by these correlations than the velocity gradient, which rotates them. This is in accordance with the visual nature of the tracer dynamics (see Fig. 2a and movies in [54]), which show particles being advected over long distances without reorienting, punctuated by short reorientation events when the tracer moves through regions with sufficiently large velocity gradients. This qualitative picture of the particle dynamics is further quantified by the fact that the rotational displacement PDFs in Figs. 6bd are more strongly non-   Figure 7. PDF of translational (a) and rotational (b) displacements in the short-time, ballistic regime. In (a), the displacement data is the same as the short-time data in Fig. 6a, but displacements have instead been rescaled by ∆t. Note the collapse of P (∆r α t /∆t) to the PDF of the components of the instantaneous fluid velocity U (solid black line) for short ∆t. Panel (b) shows the corresponding data for the rescaled orientational displacement ∆p α t /∆t for ∆t < τR. The solid black line shows the PDF of the corresponding velocity gradient, as described in the text.
Gaussian than the translational equivalents in Fig. 6ac.
By separately evaluating diffusion coefficients parallel (D ) and perpendicular (D ⊥ ) to the tracer major axes, we were furthermore able to quantify the anisotropy of translational diffusion through the ratio D /D ⊥ . In the case of pushers, this ratio displays a sharp increase coinciding with the onset of active turbulence. In fact, this increase appears sharper than what is seen for D T or D R , indicating that D /D ⊥ might be used as a robust fingerprint for locating the transition to active turbulence, similar to the previous proposal by Krishnamurthy and Subramanian [39] who suggested to use D T as an observable sensitive to the transition to collective motion.
As highlighted above, this behaviour is furthermore in qualitative contrast with experimental measurements on ellipsoidal tracer particles in a thin film of E. coli suspension [14], which showed a monotonic decrease in D /D ⊥ with density. Our results show that, in unbounded 3D systems, this type of anisotropy cannot be attributed to the generic far-field hydrodynamic advection by dipolar swimmers. The hydrodynamic argument for the observed anisotropy given in [14] for pushers and [18] for pullers is instead based on the advection by a single swimmer moving in the same 2-dimensional plane as the tracers, while advecting and rotating them with a 3-dimensional dipolar flow field. This specific geometry breaks the 3D pusher-puller equivalence at the single-swimmer level and therefore yields different translation-rotation couplings for tracers in pusher and puller flow fields, as was shown in [14]. Nevertheless, given that we observe no measurable anisotropic diffusion in 3D until near the onset of active turbulence, we argue that the far-field, singleswimmer contribution to the experimentally observed anisotropic diffusion is generally small, although it would be interesting to investigate how 2D tracer confinement affects this finding, since such confinement is known to strongly affect the hydrodynamic interactions between dipolar microswimmers [55][56][57]. Furthermore, the fact that the anomalous translation-rotation coupling in [14] is only observed for high bacterial densities where collective motion occurs indicates that the single-swimmer picture is not sufficient to describe the phenomenology. Instead, we hypothesize that the anomalous anisotropic diffusion in [14] is due to specific effects such as nearfield hydrodynamic interactions or direct swimmer-tracer collisions that are enhanced by the quasi-2D experimental geometry, where the latter have indeed been shown to yield to D /D ⊥ < 1 [28]. Another indication that the translation-rotation coupling is strongly affected by details of the system geometry and other experimental parameters is the contrasting results obtained for a confined, stiff colloidal chain in an E. coli suspension, where D /D ⊥ was observed to exceed unity, and even exceed the maximal Brownian value of D /D ⊥ = 2 for large enough microswimmer concentrations [58].
In addition, we studied the probability distributions of translational and rotational displacements of tracers as functions of the length ∆t of the observation window. In accordance with previous works [16,17,21,23,[33][34][35], these distributions show a universal and strongly non-Gaussian behaviour, although they here only occur for short observation times. In the limit ∆t → 0, we show that this collapse is a generic consequence of the nature of ballistic particle motion: as long as particles move in straight paths, the displacement statistics must be identical when renormalized with the average length of these paths. Notably, this "ballistic collapse" is different from that discussed in [35], where an intermittent collapse of the displacement PDFs in the diffusive regime was rationalised, showing that there are in fact two separate mechanisms leading to universal scaling of tracer PDFs in the two time regimes.
While this study forms a natural starting point for the computational investigation of anisotropic tracer motion in active suspensions, more work is needed to quantitatively connect with experimental results, which are likely to be significantly influenced by both near-field flows and non-hydrodynamic interactions between swimmers and tracers, neither of which are included in our model. An especially interesting aspect of the anisotropic tracer dynamics to be investigated in the future is the effects from changing the system geometry, which changes both the pusher-puller symmetry at low densities and affects the nature of swimmer-swimmer correlations and collective behaviour, and might therefore shift both the relative importance of rotational and translational diffusion and the nature of the translation-rotation coupling.
This work was financed through the Knut and Alice Wallenberg Foundation (project grant KAW 2014.0052). J.S. gratefully acknowledges financial support from the Swedish Research Council (Project grant 2019-03718). All simulations were performed on resources provided by the Swedish National Infrastructure for Computing (SNIC) at LUNARC.