Issue
Natl Sci Open
Volume 3, Number 4, 2024
Special Topic: Active Matter
Article Number 20230081
Number of page(s) 10
Section Physics
DOI https://doi.org/10.1360/nso/20230081
Published online 22 March 2024

© The Author(s) 2024. Published by Science Press and EDP Sciences.

Licence Creative CommonsThis is an Open Access article distributed under the terms of the Creative Commons Attribution License (https://creativecommons.org/licenses/by/4.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the riginal work is properly cited.

INTRODUCTION

Our understanding of how the microscopic details of the single-particle dynamics lead to different collective behaviors is presently far from satisfactory. For this reason, of late practitioners often resort to non-reciprocal interactions as a tool to model the autonomous clustering of active matter [1-3]. By non-reciprocity we refer to the situation when one particle perceives the presence of a second particle, which, therefore, influences the dynamics of the first particle without modifying its own. Being typically long-ranged, such non-Newtonian interactions involve neither fields of force [4-6], nor (reciprocal) pair-potentials [7], nor steric effects [8, 9], as assumed, for instance, in the early literature on motility-induced phase separation (MIPS) [10]. Indeed, finite-size self-propelled particles may form aggregates by adjusting their velocity according to the direction and the local density of their peers, a mechanism known in biology as quorum sensing (QS) [11, 12].

In its simplest form, MIPS has been shown to be the analogous of a gas-liquid phase separation and can be observed both in biological and synthetic active matter. However, there is still no consensus on how to induce these phases and numerical studies report different behaviors depending on model details [13].

Elementary models of QS active suspensions have been numerically demonstrated to exhibit intriguing clustering capabilities depending on the encoded sensing protocol. Bechinger and coworkers [14] showed that disordered swarms can form when the self-propulsion of individual particles switches off either in regions of high peer concentration, independently of their relative orientation (active-to-passive transitions), or, on reverse, when their sensing cone points away from high density regions (passive-to-active transitions) [15]. Chiral QS protocols can be invoked to generate rotating swarms, or swirls [16]. Moreover, QS particles switching between two active states have been reported to exhibit either social or antisocial behaviors, depending on whether their self-propulsion speed grows smaller or larger for densities above a given concentration threshold (active-to-active transitions) [17].

Inspired by the recent literature on QS active matter, we propose here an even simpler QS protocol, which applies to regular colloids, that is to passive particle suspensions. We assume that the diffusion constant of the overdamped suspension particles can switch from high to low, as the peer density they perceive in their closest surroundings raises above a set threshold. Under certain conditions, such a QS protocol can trigger a clustering process, which eventually leads to a dynamical phase separation.

The content of this paper is organized as follows. In Section “Model”, we detail our QS protocol, characterized by a tunable sensing distance and transition threshold. In Section “Clustering via quorum sensing”, we analyze quantitatively our numerical results for the clustering of a two-dimensional (2D) suspension of hard discs, by introducing appropriate quantifiers. In Section “Model parameter dependence”, we discuss the robustness of the proposed QS clustering mechanism against changes of the model parameters. In Section “Effective diffusivity in a two-phase suspension”, we investigate the effective diffusivity of the suspension particles to stress the dynamical nature of the coexisting phases. Finally, in Section “conclusions”, we briefly mention variations of the present model and possible applications beyond colloids.

MODEL

Single particle dynamics

Any overdamped colloidal particle, labeled here by the index i, executes a regular 2D Brownian motion with Langevin equationr˙i=ξi(t),(1)where ri=(xi, yi) are the coordinates of its center of mass and ξi(t)=[ξix(t),ξiy(t)] is a stationary source of Gaussian noise with ξiq(t)=0 and ξiq(t)ξip(0)=2Diδqpδ(t) for q, p=x, y, modeling the equilibrium thermal fluctuations in the suspension fluid. The particle diffusivity, Di, is allowed to switch between a higher and lower value through the QS protocol detailed below. The stochastic differential Eq. (1) was numerically integrated by means of a standard Euler-Maruyama scheme [18]. To ensure numerical stability, the numerical integrations have been performed using an appropriately short time step, δt=10-3 (See the figure captions for more numerical details).

Colloidal suspension dynamics

We numerically simulated a colloidal suspension by placing N identical, independent colloidal particles of Eq. (1) in a square box of side L and imposing periodic boundary conditions. In a suspension, effects due to the particle finite size cannot be neglected. We modeled steric interactions as either: (i) hard-disk collisions, whereby the particles were modeled as hard discs of radius r0, or (ii) short-range pair-repulsions with truncated Lennard-Jones potential [19],Vij={4ϵ[(σrij)12(σrij)6],  if  rijrm,0,otherwise,where rij is the distance between particles i and j, rm=21/6σlocates the potential minimum, and σr=2r0 is the particle effective diameter. Further reciprocal interactions and hydrodynamical interactions [20, 21] have been neglected.

Quorum sensing protocol

We then assumed that the diffusivity of each particle depends on the spatial distribution of its neighbors. In biological systems this process is mediated by some form of inter-particle communication (mostly chemical in bacteria colonies [11, 12]). On the other hand, the diffusivity of colloids is suppressed with increasing density [10, 8]. Without entering the details of the specific sensing mechanisms, we can define the sensing function of particle i as [14]Pi(dc)=jVidc12πrij,(2)where Vidc denotes its perception area of radius dc. This means that each suspension particle senses the presence of its peers in all directions within a restricted horizon, rijdc (see Figure 1A).

thumbnail Figure 1

Coexisting phase patterns in a suspension of N=470 hard discs of radius r0=1 in a square box of side L=100. (A) Schematics of the QS protocol for the hot-cold transition with DM=0.1 (red dots) and Dm=0.01 (blue dots); (B)–(H) suspension snapshots at t=105 for Lp=10 and increasing QS distances, respectively, dc=5, 8, 10, 12, 15, 17, and 19.

The particle diffusivity is then governed by the following simple quorum sensing protocol:Di={DM,P(dc)iPth,Dm,Pi(dc)>Pth,(3)where DM>Dm and the transition threshold is defined as Pth=ρ0Lp, with ρ0=N/L2 denoting the suspension density in the case of uniform spatial distribution. Note that for a uniform suspension the particles’ sensing function of Eq. (2) is approximatively P¯(dc)=ρ0dc. Clearly, this form of particle interaction is non-reciprocal, since j may be perceived by i without being itself affected by the presence of i.

Tunable model parameters

A simple rescaling of the space and time variables, xx/L, yy/L, and tt/τ with τ=L2/Dm, shows that the free dynamical parameters for our suspension are Lp/L, dc/L, DM/Dm, and σ/L. Accordingly, the quantities defining the QS protocol scale like ρ0L2ρ0=N, PthLPth=N(Lp/L), and P¯(dc)LP¯(dc)=N(dc/L). The choice of the diffusion constants, Dm and DM, affects the stationary properties of the suspension only through the ratio DM/Dm. The particle diameter, σ, enters the simulation output mostly through the definition of the actual sensing function Pi(dc) in Eq. (2), where it plays the role of a natural cut-off when rij0. On the other hand, for sensing distances, dc, larger than the mean interparticle distance, lN=L2/N, QS is expected to control the suspension diffusivity irrespective of any short-ranged reciprocal interactions. For this reason, in our simulations we fix the suspension parameters, N and L, the particle size, σ=2r0 and, therefore, the packing fraction ϕ=N(σ/L)2, and vary the remaining tunable parameters Lp, dc, and DM/Dm.

CLUSTERING VIA QUORUM SENSING

We start presenting our numerical results for a suspension of N hard discs of radius r0 confined to a square box of side L. As illustrated in Figure 1, on increasing the ratio dc/Lp the suspension separates into two phases, a hot and a cold one respectively with Di=DM and Di=DM. For dc/Lp≪1(dc/Lp≫1), all particles turn hot (cold). As dc/Lp approaches 1, we observe the appearance of several small cold clusters, which eventually merge into one large cold cluster. Analogously, upon lowering dc/Lp toward 1, small cavities open up in the uniform cold phase, which host sparse hot particles. As dc/Lp approaches 1, these cavities also merge into one large cavity. One can regard the cavity formation in the cold phase as a sort of symmetric counterpart of the cluster formation in the hot phase. It is worth noticing that for dc/Lp~1 the transition between cluster and cavity phase separation is marked by the emergence of band structures that cut across the entire simulation box. This is a symmetry breaking effect induced by the periodic conditions at the box boundaries. The orientation of the emerging bands can be either vertical or horizontal, depending on the system initial conditions. The coexistence of extended compact hot and cold phases would require an infinitely large simulation box.

The corresponding phase diagram in the (dc, Lp) plane is shown in Figure 2A for a fixed ratio DM/Dm. The pattern sequence of Figure 1 repeats itself for any value of Lp, except LplN and LpL/2, where the QS protocol becomes respectively to discreteness (finite interparticle separation) and periodicity effects (sensing area encroaching the box boundaries).

thumbnail Figure 2

(A) Phase diagram in the space parameter (dc, Lp) for the hard-discs suspension of Figure 1. Seven distinct phase topologies are distinguishable, say, for fixed Lp and increasing dc: all hot, small clusters, one cluster, bands, one cavity, small cavities, all cold. Configurations with many small clusters or cavities are hard to detect at too low or too large values of Lp, i.e., of the QS threshold in Eq. (2) (see text). (B) Fraction of cold particles, Nm/N, vs. dc for Lp=10. All remaining parameters are as in Figure 1. Vertical colored bands in (B) locate the different phase topologies from (A).

Under stationary conditions, that is for sufficiently long simulation runs (typically t>105), the number of cold particles fluctuates around a stable average value, Nm, also displayed in Figure 2Bvs. dc at fixed Lp. The corresponding topology of the coexisting phases is marked by the same color code as in Figure 2A for reader’s convenience. As already apparent in Figure 2A, the transition from a predominantly hot to a predominantly cold suspension occurs for dc~Lp. This result comes as no surprise in view of the definition of the QS threshold in Eq. (3), Pth=ρ0Lp, and recalling that the particles’ average sensing function for a uniform distribution would be P¯(dc)=ρ0dc.

To further characterize the QS induced hot-to-cold transition, we considered the suspension-averaged radial distribution function g(r) [22]. A few g(r) curves for the hard disc suspension of Figures 1 and 图2B are displayed in Figure 3A for different values of the QS distance, dc. The strength of the QS effect is quantified by the height, gmax, of the first g(r) peak around rσ, reported in Figure 3B as a function of dc. The clustering effect turns out to be the strongest as the one cold cluster grows into a periodic band, namely, for the two-phase topology that separates the predominantly hot from the predominantly cold configurations. Recalling that for a uniform particle distribution g(∞)=1, we are not surprised to observe that in the presence of clustering, dc<Lp, the tails of g(r) approach asymptotic values g(∞)<1. Vice versa for dc>Lp, the formation of hot cavities increases the density of the predominant cold phase, so that g(∞)>1.

thumbnail Figure 3

(A) Radial distribution function g(r) for the hard-disc suspension of Figure 2B with different values of dc (see legend in (C)). (B) Maximum amplitude of the g(r) curves in (A), gmaxvs. dc. The color bands denote the different phase topologies as in Figure 2B. (C) Distribution of the local suspension density for the same model parameters as in (A). All curves reported in (A)–(C) have been averaged over 1000 stationary configuration taken 100 time units apart, starting from t=105. The local densities, ρ in (C) have been computed over square bins of side δL=10.

Finally, we also considered the distribution of the local suspension density ρ. We divided the simulation box into square bins of size δL×δL and counted the number of particles, δN, in each bin under stationary conditions, i.e., for t>105. We then averaged the bin density, ρ=δN/δL2 over time for better statistics. Some ρ distributions, P(ρ), are plotted in Figure 3C. The coexistence of two phases is signaled by bimodal ρ distributions, with peaks to the left (right) of ρ0=N/L2 (uniform particle distribution) respectively for the hot (cold) phase. The two peaks equilibrate in the intermediate band phase topology. For a more sophisticated description of the clustered suspension we implemented Voronoi’s tessellation to consistently partition the simulation box into non-overlapping local polyhedra without introducing the arbitrary mesh length δL [23]. The area distribution of the Voronoi cells (not shown) proved consistent with the curves of Figure 3C.

MODEL PARAMETER DEPENDENCE

To establish the robustness of the QS induced phase transition we discuss here the dependence of the quantifiers introduced above on the model parameters. In Figure 4 we focus on the dependence of the average cold disc fraction, Nm/N, on the QS threshold, Pth, i.e., on the parameter Lp. As anticipated in Section “Clustering via quorum sensing”, the simulation output is insensitive to the actual value of Lp as long as we keep the ratio dc/Lp fixed. Vice versa, this scaling property fails for exceedingly small (LplN) or large Lp values (LpL/2), where discreteness and finite box size effects come into play.

thumbnail Figure 4

Cold particle fraction, Nm/N, for the hard-disc suspension of Figure 1, but different Lp (solid symbols). Crosses represent the results obtained for Lp=10, i.e., same parameters as in Figure 2B, except here particles interact via the potential of Eq. (2) with σ=2 and ϵ=1/6.

We then replaced the disc hard-core collision with particle pair interactions mediated by the truncated Lennard-Jones potential of Eq. (2) of equal (effective) diameter σ=2r0: no appreciable changes in the two phase fractions were detected. In the following we will keep reporting data for suspensions of particles interacting via the potential of Figure 2 to further substantiate this conclusion.

As anticipated in Section “Model”, we expect that under stationary conditions QS induced diffusivity phase transitions only depend on the ratio DM/Dm, one diffusivity constant, say Dm, entering the time scaling factor τ. Our expectations are confirmed by the simulation data for Nm/N and gmaxvs. dc plotted in Figure 5. Curves for the same ratios DM/Dm overlap, no matter what the choice of DM and Dm. It is worth pointing out that on increasing DM/Dm the hot-to-cold transition is slightly anticipated (Figure 5A), and the cluster packing enhanced (Figure 5B). This is an effect of the diffusivity difference between coexisting phases, as also detected in hard-disc suspensions.

thumbnail Figure 5

Dependence of (A) Nm/N and (B) gmax on the diffusivity parameters, DM and Dm. Other suspension parameters: N=382, Lp=10, and particles interacting via the potential of Eq. (2) with σ=2 and ϵ=1/6.

Finally, we looked at the dependence of the curves Nm/N versus dc on the suspension parameters N and σ (not shown). Due to local density fluctuations, the onset of the QS transition happens for sensing distances, dc, appreciably lower than Lp (see Figure 4). However, we noticed that the transition onset shifts to higher dc values, dcLp, with increasing N. Moreover, the transition onset distance also decreases (almost proportionally) with σ. This effect is due to the fact that, as anticipated in Section “Model”, the particle size enters the QS protocol as a natural cut-off of the sensing function of Eq. (2). These remarks concern the transition onset alone, the overall hot-to-cold transition being localized by the condition dc~Lp, as discussed above.

EFFECTIVE DIFFUSIVITY IN A TWO-PHASE SUSPENSION

We conclude our analysis by investigating the effective diffusivity of a particle in a two-phase suspension. Of course, hot and cold particles have instantaneous diffusivity DM and Dm, respectively. The question then rises of how they diffuse within the suspension. Our numerical simulations show that all particles obey the same normal diffusion lawlimt[x(t)x(t)]2=2Defft,(4)where ... denotes a time average taken over the trajectory of a single particle. For a better statistics and without loss of generality, in Figure 6 the single-particle effective diffusivity has been ensemble averaged over all N suspension particles. When plotted versus dc, the single-particle Deff drops from DM* at dc=0 down to Dm* for dc>>Lp. These two limiting situations refer to the hot and cold uniform one-phase suspensions, respectively. Subject to QS transitions, a single particle thus switches from hot to cold, and vice versa, without belonging permanently to either phase. Moreover, one notices by inspection that DM* and Dm* are smaller than the relevant free-particle diffusivity, DM and Dm, but their ratio remains the same, DM*/Dm*=DM/Dm. Indeed, we checked that the effective diffusivity in a uniform one-phase suspension decays with the packing fraction, φ, according to the well known power law [10, 8, 24], DM(m)*=DM(m)(1λϕ)2 with λ0.9, which holds for φ≪1.

thumbnail Figure 6

Effective particle diffusivity, Deff, vs. dc in a two-phase suspension with N=382, Lp=10, DM=0.1, Dm=0.01, and pair potential of Eq. (2) with σ=2 and ϵ=1/6. Single-particle Deff data (crosses) are compared with better statistics ensemble averaged data (dots) and the fitting law of Eq. (5) (solid curve), obtained from the corresponding Nm/N vs. dc curve in Figure 5. Inset: MSD vs. t for different phase topologies.

In the intermediate regime, dc~Lp, for a stationary two-phase suspension a simple two-state argument suggests that Deff=(NmDm*+NMDM*)/N with NM=N-Nm; henceDeff=DM*[1NmN(1DmDM)].(5)This equation relates the dc dependence of Deff and Nm/N, in good agreement with the simulation data of Figure 6 (even for a single particle).

CONCLUSIONS

We have analyzed a simple model of QS induced phase transition, whereby the particles of a colloidal suspension switch from high to low diffusivity as the concentration of their peers in their immediate vicinity grows above a certain threshold. We have numerically proven that such mechanism suffices to trigger the formation of cold clusters (hot cavities) in a hot (cold) suspension. The coexistence of two phases under stationary conditions has been characterized by means of several quantifiers, like particle phase fractions, radial distributions, an effective diffusivity. In the current literature, like in the present paper, QS synthetic matter has been investigated in 2D, only, because no actual 3D experiments have been performed so far. However, experiments of living swimmers suggest that QS clustering also occurs in 3D [11, 12].

While in our presentation we refer for simplicity to the simulated system as an ideal colloidal suspension, the proposed model can be extended to describe assemblies of less elementary units, where sensing-at-distance is much easier to implement. For instance, it is everyday’s experience that individuals in a social setting tend to slow down the pace as they perceive, for instance visually, the presence of an even casual gathering, thus triggering the formation of large, persisting crowds. Stated otherwise, the proposed clustering mechanism relies on non-reciprocal interactions, which can operate at much longer distances than any pair interaction, steric, electrostatic and hydrodynamic, alike.

By the same token, one can think of reversing the QS mechanism of Eq. (3), asDi={DM,Pi(dc)Pth,Dm,Pi(dc)<Pth.(6)

In this case we can expect an anti-associative behavior with all individuals trying to distance themselves from one another so as to avoid local gatherings. This variation of the present model is the subject of ongoing investigations.

Funding

Y.L. was supported by the National Natural Science Foundation of China (12375037 and 11935010).

Author contributions

Y.Z. collected the data; Y.L. and F. M. conceived and analyzed the work; F.M. wrote the original draft. All authors have read and are agreed to publish the manuscript.

Conflict of interest

The authors declare no conflict of interest.

References

All Figures

thumbnail Figure 1

Coexisting phase patterns in a suspension of N=470 hard discs of radius r0=1 in a square box of side L=100. (A) Schematics of the QS protocol for the hot-cold transition with DM=0.1 (red dots) and Dm=0.01 (blue dots); (B)–(H) suspension snapshots at t=105 for Lp=10 and increasing QS distances, respectively, dc=5, 8, 10, 12, 15, 17, and 19.

In the text
thumbnail Figure 2

(A) Phase diagram in the space parameter (dc, Lp) for the hard-discs suspension of Figure 1. Seven distinct phase topologies are distinguishable, say, for fixed Lp and increasing dc: all hot, small clusters, one cluster, bands, one cavity, small cavities, all cold. Configurations with many small clusters or cavities are hard to detect at too low or too large values of Lp, i.e., of the QS threshold in Eq. (2) (see text). (B) Fraction of cold particles, Nm/N, vs. dc for Lp=10. All remaining parameters are as in Figure 1. Vertical colored bands in (B) locate the different phase topologies from (A).

In the text
thumbnail Figure 3

(A) Radial distribution function g(r) for the hard-disc suspension of Figure 2B with different values of dc (see legend in (C)). (B) Maximum amplitude of the g(r) curves in (A), gmaxvs. dc. The color bands denote the different phase topologies as in Figure 2B. (C) Distribution of the local suspension density for the same model parameters as in (A). All curves reported in (A)–(C) have been averaged over 1000 stationary configuration taken 100 time units apart, starting from t=105. The local densities, ρ in (C) have been computed over square bins of side δL=10.

In the text
thumbnail Figure 4

Cold particle fraction, Nm/N, for the hard-disc suspension of Figure 1, but different Lp (solid symbols). Crosses represent the results obtained for Lp=10, i.e., same parameters as in Figure 2B, except here particles interact via the potential of Eq. (2) with σ=2 and ϵ=1/6.

In the text
thumbnail Figure 5

Dependence of (A) Nm/N and (B) gmax on the diffusivity parameters, DM and Dm. Other suspension parameters: N=382, Lp=10, and particles interacting via the potential of Eq. (2) with σ=2 and ϵ=1/6.

In the text
thumbnail Figure 6

Effective particle diffusivity, Deff, vs. dc in a two-phase suspension with N=382, Lp=10, DM=0.1, Dm=0.01, and pair potential of Eq. (2) with σ=2 and ϵ=1/6. Single-particle Deff data (crosses) are compared with better statistics ensemble averaged data (dots) and the fitting law of Eq. (5) (solid curve), obtained from the corresponding Nm/N vs. dc curve in Figure 5. Inset: MSD vs. t for different phase topologies.

In the text

Current usage metrics show cumulative count of Article Views (full-text article views including HTML views, PDF and ePub downloads, according to the available data) and Abstracts Views on Vision4Press platform.

Data correspond to usage on the plateform after 2015. The current usage metrics is available 48-96 hours after online publication and is updated daily on week days.

Initial download of the metrics may take a while.