Self-organization of active particles by quorum sensing rules

Many microorganisms regulate their behaviour according to the density of neighbours. Such quorum sensing is important for the communication and organisation within bacterial populations. In contrast to living systems, where quorum sensing is determined by biochemical processes, the behaviour of synthetic active particles can be controlled by external fields. Accordingly they allow to investigate how variations of a density-dependent particle response affect their self-organisation. Here we experimentally and numerically demonstrate this concept using a suspension of light-activated active particles whose motility is individually controlled by an external feedback-loop, realised by a particle detection algorithm and a scanning laser system. Depending on how the particles’ motility varies with the density of neighbours, the system self-organises into aggregates with different size, density and shape. Since the individual particles’ response to their environment is almost freely programmable, this allows for detailed insights on how communication between motile particles affects their collective properties.

O ne of the most intriguing properties of living matter is its ability to spontaneously organise from random into complex structures on different length scales such as biofilms 1 , swarms 2,3 and flocks [4][5][6][7] . This requires communication between individual group members, which is typically realised by complex internal signal pathways. In the case of cells or bacterial colonies, communication can be achieved by extracellular signalling molecules which are produced, released and sensed by all group members 8,9 . Such biochemical communication, socalled quorum sensing, enables organisms to measure their local population density and to regulate their response accordingly (quorum sensing should be distinguished from chemotaxis describing the response of organisms in concentration gradients). The first example of quorum sensing was observed in the bioluminescent bacteria Aliivibrio fischeri which start to luminesce once their population density exceeds a certain density threshold 10 . By now, many other examples of quorum sensing have been found and it is considered to be a generic cell-to-cell communication mechanism, which is relevant, e.g., for the secretion of virulence factors 11 , biofilm formation 12 and motility control 10,13 .
In contrast to living systems, where the organism's response to a molecular concentration is determined by internal signal pathways, the motility of synthetic self-propelling active particles (APs) can be externally adjusted, e.g. by optical 14,15 , electrical 16,17 or thermal 18,19 fields. Experiments with such systems show a wealth of dynamical states ranging from living crystals 20 to phase separation [20][21][22] and swarming 16 that can be controlled by geometry 23 or boundary conditions 24 . In addition to suspensions with homogeneous motility, theoretical studies also considered APs whose motility and orientation changes upon variations of their local density [25][26][27] or their own chemical concentration gradients 28 . Under such conditions, not only motility-induced phase separation but also the occurrence of moving clumps, lanes and asters has been observed.
Here we present an experimental realisation of an active suspension, whose individual particle motion varies depending on its neighbouring density. This is achieved using APs whose individual motility, i.e. magnitude of propulsion, is controlled by the intensity of an incident focused laser beam. In contrast, the propulsion direction, which is given by the particle orientation, remains unaffected by the laser illumination and undergoes free Brownian rotational diffusion. With an external feedback-loop which consists of a real-time optical particle detection algorithm which controls the position of a scanned laser beam, we are able to explore how a specific choice of a density-dependent particle motility affects the cooperative behaviour. With experiments, numerical simulations and theory we demonstrate that small variations on how particles change their motility in response to their environment can strongly affect their self-organisation.

Results
Experimental realisation. Active particles are made from silica spheres with diameter σ = 4.4 μm which are half-coated by a 30 nm carbon film and suspended in a critical mixture of water-lutidine several degrees below its lower demixing temperature T c . Upon laser illumination, the carbon caps are selectively heated above T c , which leads to local demixing of the solvent near the caps 14 . As a result of compositional flows within the solvent, the particles self-propel opposite to the cap with their speed determined by the incident laser intensity (Methods section). Because particles are separately illuminated by a scanned focused laser beam, we are able to dynamically assign each AP an individual propulsion velocity. It should be mentioned that the direction of the propulsion velocity is opposite to the cap orientation which is undergoing free Brownian rotational diffusion. Accordingly, in our experiments, only the magnitude but not the direction of motion is controlled externally. The entire suspension is contained in a thin sample cell with height 200 μm, where particles form a two-dimensional system due to gravitational forces. To avoid variations of the total particle density within our field of view, we have applied a lateral circular confinement with reflective boundary conditions and radius R = 65 μm (Fig. 1a). However, cluster formation by quorum-sensing rules, as reported here, is also observed using periodic boundary conditions (Methods section). In the following, we investigate a suspension with constant particle density ρ 0 = 0.0092 μm −2 (corresponding ' 15% of close packing). To obtain steady-state conditions, we have allowed the system to evolve in time for at least 30 min before taking data. For further details of the experimental setup, (see Methods section).
Quorum sensing in living systems requires communication between individuals, which is achieved by the release of signalling molecules with production rate γ 29 . Because such molecules have a finite lifetime, each organism senses the molecular concentration Here, r ij ¼ jr i À r j j is the distance to the neighbour labelled j, σ the linear size of the individuals, and λ ¼ ffiffiffiffiffiffiffi ffi D c τ p (with D c the diffusion coefficient of the signalling molecules and τ their lifetime) a decay length, defining the range of the concentration profile and thus the distance over which particles communicate. The prefactor is given byc ¼ γ= 4πD c σ ð Þ. To externally introduce particle communication in a suspension of synthetic APs that are lacking the ability to release or detect such molecules, we first determine the APs' positions in our sample at periodic time intervals. Next, we use Eq. (1) and calculate the hypothetical molecule concentration c i (t) at each AP To trigger a density-dependent motility response, we apply the following rule: when the concentration c i (t) 'sensed' by an AP exceeds a threshold, i.e. c i > c th , it becomes non-motile (i.e. the laser illumination is set to zero and thus no self-propulsion occurs, v = 0), otherwise it is motile and propels with velocity v = v 0 (Fig. 1c, d). Such a sharp threshold is in agreement with the conditions of many living systems exhibiting quorum sensing 10,13 . In the following, the propulsion velocity was set to v 0 = 0.2 μm s −1 (Methods section). Due to diffusive and active motion, particles undergo configurational variations over time which lead to continual changes of the APs' behaviour from motile to non-motile. Particle motilities are updated every 500 ms (Methods section). During this time interval, particles typically move by a distance <5% of their diameter. Note that the reduction of the motility update time by a factor of 100 yields identical results as confirmed by numerical simulations (Supplementary Figure 1).
For c th = 0 (i.e. an entirely Brownian suspension), particles are homogeneously distributed. Increasing c th leads to an inhomogeneous particle distribution and formation of clusters with growing density (Fig. 2a-c). In the following, clusters are defined by regions where the particle density ρ is 20% larger than ρ 0 . When c th = ∞ (i.e. permanently motile APs), cluster formation again disappears. This is in contrast to APs with constant motility, where cluster formation is observed at much higher velocities and densities of APs 22,31 (Supplementary Figure 2). This clearly demonstrates that here the organisation into densely packed regions is entirely due to the presence of quorum-sensing rules. The formation of clusters under such conditions is qualitatively understood as follows: isolated and motile particles approaching dense regions slow down; as a result, they become entirely diffusive since they sense super-threshold concentrations. This facilitates the aggregation of particles leading to cluster growth. Because particle self-propulsion is turned off when joining a cluster (in contrast to permanently motile particles), the density within clusters is rather loose as confirmed by the radial density profiles ρ(r), where r is measured relative to the clusters' centre of mass (Fig. 2e). Such loose packing is in strong contrast to the closely packed (even crystalline) aggregates which are observed in dense suspensions of APs with constant motility 32 . Our experimentally measured density profiles are in excellent agreement with those obtained by numerical simulations (solid lines in Fig. 2e) of a simple model neglecting hydrodynamic interactions (Methods section).
Numerical and analytical analysis. To gain further insights into the influence of the quorum-sensing rule on the collective particle behaviour, we performed numerical simulations where the concentration threshold c th and the decay length λ were systematically changed. In addition, we developed a mean-field theory from which analytical results can be obtained. It describes a passive, circular homogeneous cluster in coexistence with a motile gas due to the balance of a diffusive current away and an active current towards the cluster (Methods section). Our findings are summarised in Fig. 3a. In agreement with Fig. 2a-d, clustering occurs only within a well-defined range of concentration thresholds c th . Our numerical simulations predict cluster formation between the solid red and blue line (the latter also in agreement with the onset of cluster formation within the meanfield calculations). With increasing decay length λ, each particle senses a larger number of neighbours, i.e. a higher concentration (Eq. (1)). Accordingly, the maximal threshold c th for which only non-motile particles are observed, becomes larger. The same trend also applies to the conditions where only motile particles exist. As can be seen, cluster formation occurs only within these limiting cases. With increasing c th (λ = const.) the clusters are compressed, i.e. they become denser and smaller, as seen by the open symbols in Fig. 3b, c. This behaviour is in good agreement with our experiments (closed symbols) and qualitatively reproduced by the mean-field calculations (solid line) in which fluctuations and the excluded volume are neglected. As expected, the agreement with mean-field generally becomes better with increasing decay length since particles sense over larger distances and are thus less sensitive to fluctuations of the particle  Figure 3). The upperboundary for cluster formation suggested by mean-field theory is shown as a dashed red line and largely overestimates the limit where stable clusters are formed according to our experimental and numerical data. The reason is that clusters shrink as we increase the threshold, and clusters composed of a few particles become unstable with respect to fluctuations (which are neglected in mean-field). A small density fluctuation (e.g. because one particle leaves the cluster) might lead to a drop of the chemical concentration below the threshold for other particles in the cluster, which then become motile and also leave the cluster. For small clusters, this positive feedback between density and motility fluctuations is sufficient to dissolve the entire cluster while larger clusters remain stable. In the simulations, we observe that clusters with less than N p ≈ 65 passive particles spontaneously dissolve (and also re-form) independent of λ (Supplementary Figure 4). This effect lowers the upper boundary in Fig. 3a to much smaller concentration thresholds. Larger global densities stabilise clusters at larger thresholds and move the upper boundary closer to the theoretical limit of closely packed clusters (Supplementary Figure 5). The actual number below which clusters dissolve depends on the system size (Supplementary Figure 6). Our experiments and simulations suggest that cluster formation resulting from quorum-sensing rules requires not only the coexistence of motile and non-motile particles (Fig. 3a), but also the ability to change their motility (Fig. 4a, b). To investigate the importance of such motility changes in more detail, we performed simulations with mixtures of motile and non-motile particles but without the possibility of motility changes. Independent of the mixing ratio, no clustering is observed for our packing fraction (data not shown).
In order to quantify the motility changes, we introduce the motility change density _ n p$a ðrÞ, which is the number of motilepassive change events within a circular ring at distance r from the centre of the cluster per time and area. Figure 4c compares how _ n p$a ðrÞ and the particle density profile ρ(r) depend on the radial position within a cluster for λ = 10σ and c th = 8.6c. The rate _ n p$a ðrÞ becomes largest at the interface between the dense cluster and the (dilute) gas because configurational changes of the APs are most pronounced near the interface. Accordingly, the distributions are shifted to smaller radii (Fig. 4d) when the cluster size becomes smaller by increasing the sensing threshold c th (cf. Figure 2a-d). In Fig. 4e, we have plotted the total rate of motility changes _ N p$a (density _ n p$a ðrÞ integrated over the whole area of the system) vs. c th . Obviously, induced cluster formation requires a minimal rate of motility changes. Starting from a homogeneous particle distribution (c th = 0), a growing threshold c th leads to clusters with increasing density (cf. Fig. 3a). Accordingly, the deviation of the particle distribution compared to thermal equilibrium increases. Sustaining such nonequilibrium conditions requires an increasing amount of motility control, i.e. switching rate _ N p$a . As the clusters become unstable with respect to fluctuations, the switching rate drops and becomes again zero when all particles remain active. The behaviour shown in Fig. 4e is robust with respect to changing of λ (Supplementary Figure 7) and thus demonstrates that the dynamics of motility changes plays an important role for cluster formation induced by quorum-sensing rules.
We have also studied how the collective behaviour changes upon further variations of the particles' response to their environment. Figure 5a shows a cluster that is formed by introducing a second concentration threshold c th,2 > c th above which particles recover their motility 27 . This leads to enhanced active particle motion near the cluster centre and preferential escape from this region. As a result, the particle density near the centre decreases which leads to a 'ring'-like structure.
We also considered the situation where the concentration profile of signalling molecules becomes non-isotropic, e.g. due to thermophoretic forces induced by external temperature gradients. To account for the angle-dependence of the concentration profile around each particle j, we have modified Eq. (1) by an angle- where Θ ij denotes the angle between the connecting vector of particles i and j and the x-axis. For f(Θ) = cos 4 (Θ) we find an elliptically elongated shape, whose axis is rotated by 90°for f(Θ) = sin 4 (Θ) (Fig. 5b, c). Choosing f(Θ) = cos 4 (2⋅Θ) we obtain an almost quadratically shaped cluster as shown in Fig. 5d.

Discussion
We have demonstrated the collective behaviour of synthetic particles which interact via quorum sensing rules. Contrary to living systems, where the interpretation of signalling molecules can be rather complex 33 , in our experiments well-defined perception-response relations are imposed externally by a feedbackloop. Because our approach does not rely on specific particle interactions, this enables us to freely vary not only the type of stimuli but also how particles respond to them. In addition to a density-dependent incentive which controls the particles' motility, other types of particle responses including active torques 34 but also a time-delay between a stimulus and the particles' response and non-reciprocal interaction rules can be realised.
Finally, since our experiments are performed at low Reynolds number akin to the surroundings of bacteria, comparison of the collective behaviour of synthetic with living systems may allow to unveil what type and over what distances information must be exchanged to initiate collective behaviour.

Methods
Experimental details. In our experiments we use silica particles with diameter σ = 4.4 μm and a 30-nm-thick carbon film on one hemisphere, which are suspended in a critical mixture of water-lutidine at temperature T = 25°C. Particles are individually illuminated with a scanned laser beam (beam waist w = 5 μm) aiming at the centre of each particle. Under such conditions, particles self-propel with velocity v opposite to the orientation of the capped hemisphere which is undergoing rotational Brownian diffusion. It is important to notice that this reorientation is not altered by the illumination. The translational and rotational diffusion coefficients have been determined to D 0,exp = (0.0208 ± 0.0012) μm 2 s −1 and D R,exp = (112.6 ± 26.1) −1 s −1 by evaluating the translational mean-square displacement in a dilute system 35 . While D R,exp agrees well with the theoretical value, D 0,exp is about 50% below the bulk Stokes-Einstein value. The reduction is due to the enhanced viscous friction near a surface 36 and in agreement with previous studies 14 .
To avoid changes in the particle density during experiments, we employ reflective boundary conditions at the edge of a circular confinement with radius R = 65 μm by the application of an active torque which leads to a particle reorientation 14,34 . This is accomplished by displacing the illuminating laser beam relative to the particle centre by ≈2.6 μm resulting in a local intensity gradient. This causes an effective torque and thus particle reorientation. Such torques are applied to particles when leaving the circular confinement until their swimming direction Cluster formation by motility switching. a Schematic illustration of the influence of motility switches. An active particle (red) is approaching a cluster of three passive particles (blue) and becomes passive which favours joining the cluster. On the opposite, a passive particle slightly diffusing away from a cluster becomes active which facilitates leaving it. b Same situation without motility changes. c Numerically obtained radial density profile ρ(r) and motility change density _ points towards the confinement centre. To reduce variations of the calculated concentration by particles leaving and re-entering the confinement, we take all particles into account which leave the confinement by <10 μm. The confinement contains on average 122 particles, corresponding to a density of ρ 0 = 0.0092 μm −2 .
Experimental realisation of feedback-controlled particle motility. The propulsion velocity of the particles is controlled independently from each other by a feedback loop. A video camera acquires images with a repetition rate of 2 Hz which are then evaluated on a computer by a real-time particle detection algorithm to determine the particle positions. With an acousto-optical deflector (AOD) a laser beam, with beam waist w = 5 μm, is consecutively directed to the prior determined particle positions and each particle is illuminated for a period of 8 μs which is repeated every 4 ms. Since the remixing timescale of the binary mixture is on the order of 100 ms 14 , the repetition is fast enough to produce stable particle selfpropulsion conditions. With this approach, the motilities of up to 400 particles can be controlled independently.
The time-averaged illumination intensity of each particle is set to either I = 0.2 W mm −2 (resulting in a propulsion velocity v 0 = 0.2 μm s −1 ) or to I = 0 (diffusive motion) (Fig. 6). Particle configurations are updated every 500 ms, typical positional changes during this time interval are below 5% of the particle diameter. Therefore, the feedback-loop can be considered as quasi-instantaneous.
Illumination intensity corrections. To avoid velocity changes due to the overlap of the illuminating Gaussian beams of neighbouring particles, we have additionally adjusted the laser intensity depending on the relative particle positions. The intensity profile of the laser beam illuminating particle i is given by where w = 5 μm is the beam waist, I 0,i its intensity and r the radial distance. Accordingly, when particles (diameter 4.4 μm) get close to each other, they will also receive light from the beams centred on their neighbours. This leads to an increased intensity at the centre of particle i by with r ij the distance between particles i and j. To avoid such configurationdependent variations in the effective illuminating intensity (and thus of the propulsion velocity), we corrected the illuminating intensity by reducing the intensity of beam i toÎ For computational reasons, we only consider particles j with r ij < 2w in the calculation of ΔI i . Using this empirical relationship, the illumination integrated over each particle becomes independent of the positional configuration. To demonstrate the validity of Eq. (5), we have numerically tested this for a huge number of arbitrary particle configurations including situations with and without applied quorum sensing interaction. Figure 7 shows an example (corresponding to the particle configuration in Fig. 1a-c) with motile and non-motile APs. We have considered an illumination intensity of non-motile and motile particles of I = 0 and I = 0.2 W mm −2 , respectively. Without the correction discussed, we obtain a bimodal illumination intensity distribution of the particles as shown in Fig. 7b. Obviously, some of the motile particles receive up to 0.3 W mm −2 , i.e. 50% more than the nominal illumination intensity. After application of the correction procedure, this unwanted effect is almost completely suppressed (Fig. 7c). Note, that the correction is less effective for passive particles. However, due to the presence of a minimal intensity to initiate active motion (dotted line, I = 0.1 W mm −2 ), this will not change the behaviour of non-motile particles.
Simulations. In our numerical simulations, we integrate the coupled equations of motion (assuming overdamped dynamics and neglecting hydrodynamic interactions) for N particles at positions r i . A particle is propelled along its orientation vector e i with velocity v 0 = 0.2 μm s −1 if it senses a concentration c i < c th (with c i determined by Eq. (1)), for c i > c th the propulsion is set to zero and the particle is non-motile. The particle orientations undergo rotational diffusion with rotational diffusion coefficient D R = (1/120) s −1 . Translational diffusion is modelled by the random force ξ i with zero mean and variance <ξ i (t)ξ j (t′) > = 2D 0 δ ij δ(t − t′) with translational diffusion coefficient D 0 = 0.02 μm 2 s −1 . We model steric particle interactions via the repulsive Weeks-Chandler-Andersen potential cut off at r cut = 2 1/6 σ wca . We set ϵ ¼ 100 k B T and σ wca = 3.98 μm, which implies an effective (Barker-Henderson) particle diameter σ = 4.4 μm 37 . The particles are randomly initialised and the equations of motion are integrated with time step Δt = 40 ms. The concentrations c i are updated every time step. To compare the rate of motility changes to the experiment, motility changes are recorded every 480 ms. We have performed two types of simulations: employing periodic boundary conditions (shown in Fig. 8) and modelling the experimental system through a circular confinement with N = 132. For the latter, instead of applying a torque to particles reaching the boundary, in the simulations their orientation vectors are instantaneously reoriented towards the centre of the confinement with R = (65 + 10) μm. As shown in Supplementary Figure 1 Analytical theory. Neglecting the excluded volume of particles, some insights can be obtained from a simplified mean-field theory of our model extending a previous approach 27 . The evolution of an ensemble of APs with joint probability ψ(r, φ, t) is governed by with scalar speed v(c) and orientation e = (cos φ, sin φ) T . We strongly simplify the experimental situation and assume that APs interact only through the chemical concentration profile c(r) generated by the APs, which we assume to adapt instantaneously to a change of particle positions (exploiting the huge difference between colloidal, D 0 , and molecular, D c , diffusion coefficients). This is the standard model of active Brownian particles (ABPs) extended by interactions through an additional scalar field c(r).
It is sufficient to consider only the first two moments We consider stationary profiles with rotational symmetry and vanishing angular polarisation. Switching to polar coordinates with distance r from the origin (the centre of the cluster), we obtain with density profile ρ(r) and radial polarisation p(r). The first equation expresses the balance between active and diffusive particle currents. Eliminating the density, we obtain the Bessel differential equation for the radial polarisation, the solution of which is a function of r/ξ with the length We have two solutions, one for the inner passive region with speed v = 0 and radius r * , the other for the active gas with speed v = v 0 . In the inner region, from Eq. (10) we find that the density gradient is zero and thus the density ρ = ρ c is constant with vanishing polarisation. In the outer region r > r * , the polarisation decays with p(r) = −bK 1 (r/ξ) with integration constant b and modified Bessel function of the second kind K n (x). Through integration, we obtain the density profile for r > r * . Equation (10) also implies that the density is continuous at r * and thus the polarisation has to jump. The jump condition and conservation of the total density allow to determine ρ c and b.
The remaining unknown r * is determined by the concentration threshold. The concentration profile reads c r ð Þ ¼ Z dr′ρ r′ ð Þu r À r′ j j ð Þ ; uðrÞ ¼c e Àr=λ r=σ ð13Þ with condition c(r * ) = c th , which yields r * . We solve the resulting system of equations iteratively to obtain the cluster size r * and density ρ c for given threshold c th and interaction range λ (Fig. 3).
Data availability. The experimental and numerical data that support the findings of this study are available from the corresponding author upon reasonable request. Cluster formation with periodic boundary conditions. a Snapshot of a simulation with N = 1000 particles with density ρ 0 = 0.0075 μm −2 (λ = 5σ, c th = 5.5c). As in the confined system, we observe the formation of a cluster of non-motile particles (blue) surrounded by a dilute gas of motile particles (red). b Blue: Corresponding radial density profile ρ(r)/ρ 0 (with respect to the particles' centre of mass) with R half the box length. The cluster density is the same as for a circularly confined system with N = 1000 particles at the same density (dashed line). For the smaller system with N = 132, the cluster density is slightly higher (dotted line). The reason is that due to the smaller confinement, quorum concentrations are lower and thus the particles have to form a denser cluster to overcome the threshold