Direct-forcing fictitious domain method for simulating non-Brownian active particles

We present a direct-forcing fictitious domain method for simulating non-Brownian squirmer particles with both the hydrodynamic interactions and collisions being fully resolved. In this method, we solve the particle motion by distributing collocation points inside the particle interior domain that overlay upon a fixed Eulerian mesh. The fluid motions, including those of the “fictitious fluids” being extended into the particle, are solved on the entire computation domain. Pseudo-body forces are used to enforce the fictitious fluids to follow the particle movement. A direct-forcing approach is employed to map physical variables between the overlaid meshes, which does not require additional iterations to achieve convergence. We perform a series of numerical studies at both small and finite Reynolds numbers. First, accuracy of the algorithm is examined in studying benchmark problems of a free-swimming squirmer and two side-by-side squirmers. Then we investigate statistic properties of the quasi-two-dimensional collective dynamics for a monolayer of squirmer particles that are confined on a surface immersed in a bulk flow. Finally, we explore the physical mechanisms of how a freely moving short cylinder interacts with a monolayer of active particles, and find out that the cylinder movement is dominated by collision. We demonstrate that a more directional migration of cylinder can be resultant from an inhomogeneous distribution of active particles around the cylinder that has an anisotropic shape.

Figure 1

Free-swimming motion of a single squirmer particle at $\mathrm{Re}=1.0$. (a) Translational velocity $U_p$ as a function of time for pusher ($\beta =-3$), puller ($\beta =3$), and neutral ($\beta =0$) particle, respectively. Solid lines and open circles represent the simulation results obtained by choosing the Eulerian grid size and time step as ($\mathrm{\Delta }x=\frac{a}{8},\mathrm{\Delta }t=0.002$) and ($\mathrm{\Delta }x=\frac{a}{16},\mathrm{\Delta }t=0.001$), respectively. The black dashed lines represent the analytical predictions [49]. Inset: Lagrangian collocation points distribution. (b–d) The in-plane velocity field for (b) pusher, (c) puller, and (d) neutral particles superimposed on the $u$ velocity component (colorfield) when choosing ($\mathrm{\Delta }x=\frac{a}{8},\mathrm{\Delta }t=0.002$).

Figure 2

(a) Comparison of the trajectories of two colliding squirmer particles ($\beta =\pm 5,a=1.0$). (b) Trajectories of two side-by-side pushers ($\beta =-3,a=1.0$) at various $\mathrm{Re}$ numbers. Inset: Variation of the separation distance between the two particles as a function of time by changing $\delta y$. (c) Trajectories of two side-by-side pullers ($\beta =3,a=1.0$) at various $\mathrm{Re}$ numbers. Inset: Net force exerted on the puller particle in the $y$ direction as a function of time.

Figure 3

The normalized velocity-velocity correlation function $C_{uu}$ (a) and pair-distribution function $g\left(r\right)$ (b) computed when using two different domain sizes with $L_z=51.2,102.4$ at $\mathrm{Re}=0.01$ for a concentrated puller suspension when choosing $\beta =3,a=1.0,\varphi =24\%$.

Figure 4

Velocity vector fields for a monolayer of pushers ($\beta =-3$, top row) and pullers ($\beta =3$, bottom row) at three concentrations (area fraction) at $\mathrm{Re}=0.01$. (d, h) In-plane velocity fields on the $(x,z)$ plane extracted at $y=25.6$ [marked by the dashed line in panels (c) and (g)] at $\varphi =48\%$ to reveal the 3D flow structures. Inset in panel (h): Comparison of velocity decay in the $z-\text{directionl}$ for the pusher and puller cases. For panels (a)–(c) and (e)–(g), the particle color shows the bond parameter $|q_6^{\left(k\right)}{|}^2\in [0,1]$ [60].

Figure 5

Normalized probability density function $PDF\times {\sigma }^{\prime }$ as a function of in-plane fluctuation quantity $u^{\prime }/{\sigma }^{\prime }$ at the different concentrations and $\mathrm{Re}=0.01,1.0$: (a) puller ($\beta =3$), (b) pusher ($\beta =-3$).

Figure 6

Normalized velocity-velocity correlation function $C_{uu}$ at the different concentrations and $\mathrm{Re}=0.01,1.0$: (a) puller ($\beta =3$), (b) pusher ($\beta =-3$).

Figure 7

Particle velocity correlation function $I_p$ at the different concentrations and $\mathrm{Re}=0.01,1.0$: (a) puller ($\beta =3$), (b) pusher ($\beta =-3$).

Figure 8

Pair distribution function $g\left(r\right)$ at the different concentrations and $\mathrm{Re}=0.01,1.0$: (a) puller ($\beta =3$), (b) pusher ($\beta =-3$).

Figure 9

Probability density functions (PDF) of average speed of the squirmers at the different concentrations and $\mathrm{Re}=0.01,1.0$: (a) puller ($\beta =3$), (b) pusher ($\beta =-3$).

Figure 10

Average speed of the squirmers $\langle \overline{|{\mathbf{U}}_p|}\rangle$ and the in-plane fluid flow $\langle \overline{\left|\mathbf{u}\right|}\rangle$ computed at the different concentrations and $\mathrm{Re}=0.01,1.0$: (a) puller ($\beta =3$), (b) pusher ($\beta =-3$).

Figure 11

Comparison of transport of a short circular cylinder (CC) and wedgelike cylinder (WC, $\gamma =\pi /3$) in a pusher monolayer ($\beta =-3,a=1.0,\mathrm{Re}=0.01,\varphi =12\%,36\%$): (a) Trajectory and (b) MSD. Inset in panel (a): Schematic for WC as well as the in-plane Lagrangian points distribution. Inset in panel (b): Mean migration velocity (i.e., C.O.M. velocity) magnitude of WC as a function of $\gamma$ at two different concentrations.

Figure 12

Orientation probability distribution $p_{\theta }$ for the in-plane migration velocity ${\mathbf{U}}_{\mathrm{cyl}}$, the applied viscous (${\mathbf{F}}_H$) and collision (${\mathbf{F}}_C$) force, as well as the estimated contact force through a surface integral ${\mathbf{F}}_S$ [in Eq. (33)] for the cases in Fig. 11: (a, b) CC at $\varphi =12\%,36\%$; (c, d) WC with $\gamma =\pi /3$ at $\beta =-3,a=1.0,\mathrm{Re}=0.01,\varphi =12\%,36\%$.

Figure 13

The time-averaged number density field of pusher particles ($\beta =-3,a=1.0,\mathrm{Re}=0.01$): (a, e) CC at $\varphi =12\%,36\%$; (b, f) WC with $\gamma =\pi /3, \varphi =12\%,36\%$; (c, g) WC with $\gamma =\pi /2, \varphi =12\%,36\%$; (d, h) WC with $\gamma =2\pi /3, \varphi =12\%,36\%$.