# Non-equilibrium phases in suspensions of self-propelled colloidal particles controlled by phoretic mobility and hydrodynamics

###### Abstract

We study the collective dynamics of suspensions of self-propelled diffusiophoretic colloids, in (quasi)-2 configurations, by means of a novel mesoscopic method featuring fully resolved hydrodynamics. We show that the system displays a transition when changing the sign of the phoretic mobility (which accounts for the colloid-solute interactions) and, for negative values, it develops a cluster phase. We find that the cluster size distribution follows an exponential behaviour, with a characteristic size growing linearly with the colloid activity, while the density fluctuations grow as a power-law with an exponent depending on the cluster fractal dimension. Our simulations indicate that hydrodynamics plays a relevant role in the aggregation kinetics, cluster morphology and significantly hinders the cluster growth.

Active fluids comprise a variety of systems composed by elements immersed in a
fluid environment which can convert some form of energy into
directed motion MarchettiRMP (); as such
they are intrinsically out-of-equilibrium in the absence of any external forcing .
These systems typically display complex collective phenomena, characterized by long-ranged
correlations and large density fluctuations RamaARCMP (), exhibit
a the ability to form a rich variety
of patterns Budrene (), and display non-equilibrium changes of state,
such as a flocking Vicsek (), clustering Peruani (), or
mobility induced phase transitions Cates (); Buttinoni ().
Fundamental questions arise on how individuals interact and
communicate and what is the
physical origin of their self-organization.
Furthermore, the aim of generating motion in miniaturized devices and the
development of biomimetic environments have recently boosted the
design of micro- and nano-scale self-propelling objects (“microrobots”) artswim ().
In this Letter we present a numerical study of a model system of
particles which move responding to gradients of a self-generated
concentration field; the latter, diffusing, induces a means of
interaction/communication among the active particles.
Motivated by recent experimental studies Bocquet (); Palacci-Science (), we consider the dynamics of
a layer of self-propelled colloids (SPCs) on a flat wall under the
action of gravity embedded in a liquid medium.
We find that the system develops two distinct dynamical regimes, forming either large scale
clusters or a quasi-ordered state with very slow dynamics, depending on the
interaction of the colloidal particles with solute. We characterize the transition between the two
observed non-equilibrium regimes and focus on the morphology and
dynamics of the cluster phase.
With respect to previous studies, we single out quantitatively, for the first time, the
impact of hydrodynamics on the collective dynamics of active Janus colloidal suspensions.

Model and numerical details. The 3 Navier-Stokes equations
for the fluid (solvent solute) velocity field and for
the solute concentration field footnote1 (), which read

(1) | |||||

are integrated by means of a hybrid lattice Boltzmann (LB)/finite difference method Succi (); Ludwig (). Here, is the pressure field, , and are, respectively, the kinematic viscosity, the scalar field diffusivity and a coefficient entering the mixing free energy density (and related to the concentration field compressibility) free-energy (). The source term models the production of solute at the particle surfaces and the sink term has been introduced in order to keep the space average constant in time, thus allowing to reach a steady state. Colloids are described as finite size solid spheres of radius at whose surfaces the proper momentum and torque exchange between particle and fluid is implemented via the bounce-back-on-links scheme for LB probability densities Ladd (), fixing a boundary condition with local slip velocity , where denotes a location on the particle surface. According to the theory of colloidal phoresis AndersonARFM (), in presence of an external concentration field , such effective slip velocity assumes the form

(2) |

where the phoretic mobility contains the details of the colloid-solute interactions. As we will show later on, this parameter plays a crucial role in determining the collective behaviour of SPC suspensions. Setting a constant , the particle moves with a net propulsion velocity , directed towards regions of high concentration of solute (“chemo-attractant”) for negative , and escaping from it (“chemo-repellent”) for positive . Each particle produces solute at a constant rate per unit surface : in the following we will limit ourselves to the case on one hemisphere (defined by , being the particle characteristic unit vector) and on the other (), i.e. Janus-type particles janus (), but in principle the numerical scheme can deal with arbitrarily patchy Sciortino () active colloids. An isolated active Janus colloid of this kind moves with a constant velocity of magnitude GolestanianNJP () (in our simulations the corresponding Péclet number, defined as , is always , for all ’s considered). It is worth stressing that, although we are speaking in terms of a solute/chemical field, the formalism is not specific for self-diffusiophoresis: the only two required ingredients are, in fact, the production of a diffusing scalar field and a regime of linear response of particles to gradients, and as such it enjoys a wider range of applicability, including systems like, for instance, thermophoretic colloids hot () or autochemotactic swimmers Stark () like gliding bacteria glidebact (). We have performed two series of numerical simulations of suspensions with SPCs on lattices of grid points (at a fixed surface area fraction ). Simulations from the two series will be labelled, henceforth, as runs of type A (, ) and B (, ), respectively. Particles are confined on the bottom wall by a gravity field. A hard-core particle-particle and particle-wall repulsion is introduced to prevent overlapping. In every run initially particles are randomly distributed on the surface of the bottom wall.

Non-equilibrium phase transition controlled by the phoretic mobility A number of recent experimental and numerical/theoretical studies of self-propelled particles in (quasi)- have given indication of the emergence of clustering Buttinoni (); Bocquet (); Palacci-Science (); Fily (); Redner (), however the nature of the mechanisms determining the formation of aggregates lacks a consensual agreement and seems to be strongly system-dependent. We deal with spherical particles, which rules out the possibility of alignment-induced collective motion; instead, chemical production and diffusion mediate an effective interaction, analogously to the experimental system studied in Bocquet (). While in the experiments it was surmised that active colloids felt an attractive interaction, here, we can tune the affinity of the particles for the solute via the phoretic mobility , which can be regarded as an effective charge GolestanianSoto (), i.e. positive/negative values induce repulsive/attractive interactions, respectively. Indeed, while for our simulations confirm the formation of clusters, for the system tends to reach a state with long-range order: a non-equilibrium phase transition in the -space takes place. To address the impact of the chemical affinity on the collective dynamics quantitatively, we have performed a Voronoi tessellation Voronoi () analysis of the particle space configurations footnote2 ().

The bottom insets of Fig.1 show the Voronoi diagrams both for the repulsive and cluster-forming regimes; as clearly visible to the naked eye, the geometry of the Voronoi cells for chemo-attractive and -repellent active colloids are distinctively different. The standard deviation of the cell surface distribution (normalized by the square of the mean value ), thus, turns out to be a good indicator to distinguish the two types of dynamics. In the top inset of Fig. 1 we plot , as function of the non-dimensional time (with ), for two cases with different sign of the phoretic mobility. In the attractive case, , cluster formation induces the appearance of very small (and large) cells and, hence, the surface fluctuations grow and eventually saturate at long times. For positive , instead, colloids repel each other and tend to reach an optimal covering of the space; correspondingly attains a (lower) value which remains constant in time.

Therefore, the variation of the asymptotic value of with the phoretic mobility (figure 1) clearly identifies the transition between the two regimes. The phase diagram for chemorepulsive colloids, as recently shown theoretically Benno (), is indeed rather complex and deserves further investigation.

Cluster statistics and morphology. In the following we will focus on chemoattractant colloids, , for which we characterize the cluster size distribution when changing the colloid/solute coupling intensity, . We identify clusters geometrically, based on a distance criterion. Specifically, two particles share a bond if their centers of mass are equal to or less than a cutoff distance apart footnote3 (), and we define as clusters, groups of particles connected to each other through a bond. Fig. 2 plots the cluster size probability density functions (PDFs), which can be in all cases fitted to an exponential, , over a wide range of sizes . The characteristic value and the mean size increase linearly with (see inset of figure 2), hence with the velocity of an isolated particle, in agreement with experimental observations Buttinoni (); Bocquet ().

At increasing , hints are given of fat tails at large
, which are confirmed in the simulations with no hydrodynamic
interactions (discussed later on), where a power law behviour is observed (see left inset of figure 4). A
similar phenomenology was observed in Brownian dynamics simulations of
self-phoretic active particles PohlStark ().

The global attractor for the SPCs dynamics is a set ,
where is the -th cluster with surface of area
and containing SPCs.
Correspondingly, the colloid number density becomes if and zero otherwise.
Since the colloid density fluctuations can be expressed as
(where denotes a surface average), we arrive at

(3) |

where is the measure of the whole plane, and is the probability of having a cluster of area . The number of particles in a cluster is known to scale with the cluster gyration radius as , being the fractal (Hausdorff) dimension fractal (); hence, the density of the -th cluster, , will behave as . The exponential behavior of predicts that

(4) |

where stands for the fluctuations for a inactive particles and is a phenomenological parameter, and where we have used the relation footnote4 () (see inset of figure 2). Fig. 3 shows the quantitative agreement of the theoretical prediction with the numerical observations.

Role of hydrodynamic interactions. Self-propelled colloids interact due to both the chemicals they produce or consume and the hydrodynamic flows they induce. Understanding the relative magnitude and competition between these two sources of dynamic interactions remains challenging. The model put forward allows us to switch off the hydrodynamic interactions (HI), which can be realized setting the fluid velocity to zero at each time step footnote5 (). Interestingly, our numerical study reveals that, although the transition at changing the sign of the phoretic mobility is preserved even without HI (it is determined mainly by the chemical interaction), HI have a profound effect in the kinetics of formation and morphology of the observed aggregates. In the absence of particle induced flows in the solvent, attractive SPCs () show an enhanced tendency to form clusters, as it appears in figure 4, where we compare the time evolution of the mean cluster size , with and without HI (no-HI). For the no-HI case, the growth of follows a power law behaviour (right inset) in the aggregation stage and then saturates, while it is strongly inhibited by hydrodynamics.

This statement is corroborated through the inspection of the radial distribution functions (RDFs) Hansen () (indicated as and , respectively), shown in Fig. 5: without HI the peaks are higher and decay more slowly, associated to the development of clusters larger than those formed when hydrodynamics is switched on. Besides, clusters appear substantially more compact, as appreciated in the snapshots (insets) and quantified by the measurement of a larger fractal dimension ( and , see figure 3).

Hydrodynamics hinders, then, the colloidal aggregation process. Several complex mechanisms can be conjectured to cause this phenomenon: dynamically induced effective repulsion among particles, fluid flow generated disturbances in the chemical field distribution, etc. An important effect, that we could identify quantitatively, is produced by the induced flows close to the wall, leading to hydrodynamic torques that tend to rotate the Janus particles off-plane. The inset in the Fig. 5 (top panel) shows the PDF of the degree of alignment of the particle orientation with the bounding solid wall, . When HI are present, the peak of the PDF around is more pronounced, i.e. there is a larger fraction of colloids pointing out of the plane. Accordingly, in plane particle mobility is reduced, diminisihing their capability to gather and clusterize.

To conclude, we have used a mesoscopic numerical scheme of fully resolved spherical active colloids,
propelled by self-generated gradients of a scalar field (e.g. a
chemical product) where the self-induced hydrodynamic flows can be accounted for. We have identified the role of the phoretic mobility as the key controlling parameter of the non-equilibrium phase
transition between a cluster phase and a quasi-crystal-like state. By means
of a Voronoi tessellation we have characterized the cluster state
finding that the probability distribution of sizes decays
exponentially with a mean size growing linearly with the
particle activity, in agreement with experimental
results Bocquet (); Buttinoni ().
We have quantified the profound effect that hydrodynamic plays
inhibiting clustering for negative phoretic
mobilities. We have identified the interplay between induced flows and
particle reorientation as a key element in the strong slowing down of
cluster coarsening.
This first explorative study shows that our novel
numerical method is powerful and enjoys some unique features, namely
the explicit description of chemical signalling, through the
production and diffusion of a solute concentration field and the
solvent hydrodynamics, to simulate realistic systems.
Moreover, it opens the way to address the dynamics of self-propelled
colloids in general geometries, both in isotropic and unforced situations, where
aggregation can lead to the formation of active colloidal gels, or
under gravity as in the experimental sedimentation setup.

We acknowledge the Spanish MINECO and Generalitat de Catalunya DURSI for financial support under the projects FIS2015-67837-P and 2014SGR-922, respectively. A.S. acknowledges financial support from European Research Council under the EU Seventh Framework Programme (FP7/2007-2013) / ERC Grant Agreement no[279004]. I.P. acknowledges Generalitat de Catalunya under Program Icrea Acadèmia. This work was possible thanks to the access to MareNostrum Supercomputer at Barcelona Supercomputing Center (BSC) and also through the Partnership for Advanced Computing in Europe (PRACE).

## References

- (1) M.C. Marchetti et al, Rev. Mod. Phys. 85, 1143-1189 (2013).
- (2) S. Ramaswamy, Annu. Rev. Cond. Matter Phys. 1, 323-345 (2010).
- (3) E. O. Budrene and H.C. Berg, Nature 349, 630-633 (1991).
- (4) T. Vicsek et al, Phys. Rev. Lett. 75, 1226 (1995).
- (5) F. Peruani, A. Deutsch, and M. Bär, Phys. Rev. E 74, 030904 (2006).
- (6) M.E. Cates et al, Proc. Natl. Acad. Sci. USA 107, 11715â11720 (2010).
- (7) I. Buttinoni et al, Phys. Rev. Lett. 110, 238301 (2013).
- (8) W.F. Paxton et al, J. Am. Chem. Soc. 126, 13424-13431 (2004); R. Dreyfus et al, Nature 437, 862-865 (2005); J. R. Howse et al, Phys. Rev. Lett. 99, 048102 (2007); J. Palacci et al, Phys. Rev. Lett. 105, 088304 (2010); S.J. Ebbens and J.R. Howse Soft Matter 6, 726-738 (2010).
- (9) I. Theurkauff et al, Phys. Rev. Lett. 108, 268303 (2012).
- (10) J. Palacci et al, Science 339, 936 (2013).
- (11) The advection terms play indeed no role, being the Reynolds number always much smaller than unity, typically .
- (12) S. Succi, The lattice Boltzmann equation for fluid dynamics and beyond, Oxford University Press (2001).
- (13) J.-C. Desplat, I. Pagonabarraga and P. Bladon, Comp. Phys. Comm. 134, 273 (2001); K. Stratford et al, J. Stat. Phys. 121, 163-178 (2005); K. Stratford and I. Pagonabarraga, Computers and Mathematics with Applications 55, 1585â1593 (2008).
- (14) M.R. Swift et al, Phys. Rev. Lett. 54, 5041-5052 (1996); V. Kendon et al, J. Fluid Mech. 440, 147-203 (2001).
- (15) A.J.C. Ladd, J. Fluid Mech. 271, 285-309 (1994); A.J.C. Ladd, J. Fluid Mech. 271, 311-339 (1994); C.K. Aidun, Y. Lu and E.-J. Ding, J. Fluid Mech. 373, 287-311 (1998); N.-Q. Nguyen and A.J.C. Ladd, Phys. Rev. E 66, 046708 (2002).
- (16) J.W. Anderson, Ann. Rev. Fluid Mech. 21, 61-99 (1989).
- (17) A. Walther and A.H.E. Müller, Soft Matter 4, 663-668 (2008).
- (18) F. Sciortino, A. Giacometti and G. Pastore, Phys. Rev. Lett. 103, 237801 (2009).
- (19) R. Golestanian, T. Liverpool and A. Ajdari, New J. Phys. 9, 126 (2007).
- (20) R. Golestanian, Phys. Rev. Lett. 108, 038303 (2012); D. Lüsebrink and M. Ripoll, J. Chem. Phys. 136, 084106 (2012); M. Yang, A. Wysocki and M. Ripoll, Soft Matter 10, 6208-6218.
- (21) J. Taktikos, V. Zaburdaev and H. Stark, Phys. Rev. E 84, 041924 (2011); J. Taktikos, V. Zaburdaev and H. Stark, Phys. Rev. E 85, 051901 (2012).
- (22) M.J. McBride, Annu. Rev. Microbiol. 55, 49-75 (2001); F. Peruani et al, Phys. Rev. Lett. 108, 098102 (2012).
- (23) Y. Fily and M.C. Marchetti, Phys. Rev. Lett. 108, 235702 (2012).
- (24) G.S Redner, M.F. Hagan and A. Baskaran, Phys. Rev. Lett. 110, 055701 (2013).
- (25) R. Soto and R. Golestanian, Phys. Rev. Lett. 112, 068301 (2014).
- (26) G. Voronoi, Journal für die Reine und Angewandte Mathematik 133, 97 (1907); R.H. Rycroft et al, Phys. Rev. E 74, 021306 (2006).
- (27) We follow the standard procedure of embedding each particle in a -dimensional cell whose -th edge (face) is set to be equally distant between the reference particle and its -th nearest neighbour.
- (28) B. Liebchen et al, Phys. Rev. Lett. 115, 258301 (2015).
- (29) We set such cutoff to the value of , being the hard-core range of interaction.
- (30) O. Pohl and H. Stark, Phys. Rev. Lett. 112, 238303 (2014).
- (31) T.A. Witten and L.M. Sander,Phys. Rev. Lett. 47, 1400-1403 (1981); P. Meakin Phys. Rev. Lett. 51, 1119-1122 (1983).
- (32) In principle there can be a dependence on also of the fractal dimension ; we assume here, however, that the change in affect only the characteristic cluster size and not its “compactness” (or fractality).
- (33) This procedure removes dynamic interactions mediated by the flow induced in the solvent, without touching the mechanism of exchange of momentum and torque between particle and fluid, in such a way that particles still feel a friction.
- (34) J.-P. Hansen and I.R. McDonald, Theory of simple liquids, Third Edition. Elsevier (2006).