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

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

A. Scagliarini Department of Condensed Matter Physics, University of Barcelona, Carrer de Martí i Franquès 1, Barcelona, Spain Helmholtz-Institute Erlangen-Nürnberg for Renewable Energy Production
Forschungszentrum Jülich GmbH (IEK-11) Fürther Strabe 248, 91429 Nürnberg, Germany
   I. Pagonabarraga Department of Condensed Matter Physics, University of Barcelona, Carrer de Martí i Franquès 1, Barcelona, Spain

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


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


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.

Figure 1: MAIN PANEL: Steady state standard deviation of Voronoi cell surfaces as function of the “coupling constant” (the phoretic mobility) . The two different regimes of the curve, monotically growing for negative and almost constant for positive , indicate the presence of a (non-equilibrium) phase transition between two distinct dynamical regimes. TOP INSET: vs time (, is a characteristic time of self-propulsion) for two cases with positive and negative . BOTTOM INSETS: Snapshots of the colloid distributions and relative Voronoi diagrams in the attractive, (LEFT), and repulsive, (RIGHT), case, respectively. Data come from runs of type A.

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 ().

Figure 2: MAIN PANEL: PDFs of cluster sizes for three values of the phoretic mobility:(), () and (); the dashed lines represent exponential fit which are drawn to guide the eye. INSET: Characteristic () and mean () cluster sizes as function of , i.e. of the intrisic SPC velocity (the dashed line indicates a linear relation). Data come from runs of type B.

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


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


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.

Figure 3: MAIN PANEL: Deviation of the steady state SPC density fluctuations from the value for inactive particles , normalized as , vs the coupling strength from LB (squares) and the phenomenological derivation (4), with fractal dimension , as measured in the simulations, and fitting prefactor . INSET: Mean gyration radius of clusters vs number of particles with () and without () hydrodynamic interactions (HI). The two dashed lines represent the power law : without HI a value of is found, indicating a higher degree of compactness. Data come from runs of type B.

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.

Figure 4: MAIN PANEL: Mean cluster size vs time from simulations (type B) with (HI) and without (no-HI) hydrodynamic interactions. The characteristic time is ; the single particle speed is slightly higher (roughly ) with HI, probably reflecting a finite effect. LEFT INSET: PDF of cluster sizes for and no-HI. The dashed line indicate the power law. RIGHT INSET: Log-log plot of vs , highlighting the scaling in the aggregation stage.

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).

Figure 5: RDFs for at three different times for a type B simulation and corresponding snapshots of the colloid distribution, indicating cluster formation, in a sub-system of size , located in the center of the box. The two panels correspond to simulations with (TOP) and without (BOTTOM) hydrodynamic interactions among colloids. In the top figure also the repulsive case () is reported for comparison. TOP PANEL INSET: PDFs of , the square magnitude projection onto the -plane of the colloid orientations, from simulations with (HI) and without (no-HI) hydrodynamics.

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).


  • (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).
Comments 0
Request Comment
You are adding the first comment!
How to quickly get a good reply:
  • Give credit where it’s due by listing out the positive aspects of a paper before getting into which changes should be made.
  • Be specific in your critique, and provide supporting evidence with appropriate references to substantiate general statements.
  • Your comment should inspire ideas to flow and help the author improves the paper.

The better we are at sharing our knowledge with each other, the faster we move forward.
The feedback must be of minimum 40 characters and the title a minimum of 5 characters
Add comment
Loading ...
This is a comment super asjknd jkasnjk adsnkj
The feedback must be of minumum 40 characters
The feedback must be of minumum 40 characters

You are asking your first question!
How to quickly get a good answer:
  • Keep your question short and to the point
  • Check for grammar or spelling errors.
  • Phrase it like a question
Test description