Behavior-dependent critical dynamics in collective states of active particles

We experimentally demonstrate critical collective behavior in a system of active colloidal particles (APs), achieved through variations of their mutual interactions using feedback control. At the transition between a swarm and a swirl we observe an explicit bifurcation dynamics of the rotational order parameter and a critical slowing down, i.e., a growth of the relaxation time by almost one order of magnitude. Additional signatures of critical dynamics, including hysteresis in presence of symmetry-breaking particle interactions, and a maximum of the susceptibility, are measured and characterized in terms of a theoretical model which is based purely on the time-reversal symmetry of the order parameter.

more difficult to change, evidence of critical behavior is essentially based on the properties of groups as a function of their spatial extension and density. Currently it is not clear whether collective states can also be maneuvered to a critical point merely by adjusting the behavior of individuals.
Active colloidal particles (APs) provide an intriguing testbed for investigating collective behavior since they closely mimic many features of living systems [11][12][13][14]. More recently, special attention has been given to systems of APs that are able to replicate interactions similar to those found in living systems, including mutual alignment and attraction but also collision avoidance with neighbors, which results in the spontaneous formation of swarms, flocks and swirls [15][16][17][18]. Notably, interaction rules can be continuously varied in such systems. Here we demonstrate that critical behavior in our system can be achieved in a controlled fashion by altering dynamical interaction rules while keeping other factors such as density and noise strength unchanged. At the bifurcation from swarming to swirling we also observe the emergence of critical slowing down. Additionally, we show that the 64001-p1 rotational order parameter O R follows a hysteresis loop when the symmetry of particle interactions is broken. This allows us to directly measure the system's response to perturbations, i.e., the susceptibility, which has a maximum at the critical point. The observed signatures of critical dynamics are in agreement with a theoretical model which only considers the time-reversal symmetry of O R being obeyed in our system, but is independent of the absolute size and dynamics of single particles. Our results indicate that the collective dynamics in living or responsive systems can be tuned by modifying interactions between constituents. This can be viewed in analogy to critical behavior achieved by changing control parameters (e.g., temperature, pressure, density, etc.) or interaction strengths in various physical systems, which suggests that similar collective behavior may be found in many other systems with comparable symmetries.
Methods. -Light-sensitive APs are fabricated from transparent silica particles with diameter σ = 6.3 μm, and coated on one side with a 80 nm light absorbing carbon cap. They are suspended in a water-lutidine mixture contained in a thin sample cell whose temperature is kept below the fluid's lower demixing point at about 34 • C. Upon light illumination, the caps are selectively heated, which leads to a two-dimensional self-propelled, i.e., active, AP motion due to fluid flows relative to the APs' surfaces [19,20]. Individual translational and orientational interactions between APs are achieved by tuning the intensity and displacement (relative to the cap center) of the laser focus which is scanned across all particles whose positions are determined and continuously updated via image analysis in real time. For details we refer to [17]. This allows particles to stop, move forward or turn to the left and right, depending on a chosen interaction rule with their neighbors [16,17,21]. Because these rules are executed in a real environment, AP motion is not deterministic (similar to living systems) but strongly affected by (thermal) fluctuations, avoidance of particle collisions, and heterogeneities of the motional behavior.
Motivated by perception models that successfully reproduce behaviors in schooling fish [22,23], the following interaction rule is applied: Each AP aims to move towards the center of mass (COM) of the group but with an angular deviation Δ to the left (+) or the right (−). The choice of the sign of Δ depends on the mean swimming direction (polarization) u i of its neighbors within a distance R o = 25 μm ≈ 4σ ( fig. 1(a)). To achieve alignment with neighbors, each AP selects that sign of Δ which is closer to u i . To minimize particle collisions, they turn away from each other when their clearance, i.e., surface-tosurface distance is below 0.25σ. The number of particles was kept constant to N = 48 in all our experiments.

Landau-type description and bifurcation. -
The collective behavior of the APs strongly depends on the value of the deviation angle Δ, which can be freely varied in our experiments (figs. 1(b)-(d)). At small Δ, the motion of APs is essentially towards the group center; this leads to a disordered cohesive swarm. Increasing Δ increases tangential components of the AP motion, which leads to the formation of swirls with a random sense of rotation. To quantify the rotational component of the particle motion, we define the rotational order parameter which is proportional to the angular momentum, and measures the in-plane circular motion of the entire group. Here,r i denotes the unit vector from the group center to the i-th AP, e z the unit vector perpendicular to the sample plane andû i the orientation of the AP, respectively. Although O R does not resolve the motion of individual APs, in the following we show that it correctly describes the distinct dynamical states (swarms and swirls) observed in our experiments.
Because O R is a stochastic quantity which, in the absence of a symmetry-breaking external field, has odd timereversal symmetry (i.e., O R changes its sign under the operation of time-reversal), its equation of motion has the general form where we neglect terms of higher than cubic order. In the spirit of Landau, we have added terms permissible by the symmetries [24]. The noise η has zero mean, and can be considered Gaussian to first approximation: Here the noise strength is set by an effective temperature T eff , and k B is Boltzmann's constant. The minus signs are chosen such that The Langevin equation above is the simplest stochastic model to describe the non-trivial dynamics of the rotational order parameter. As a simple order parameter description, it does not have the complexity of spatially resolved Landau-Ginzburg-type field theories [24,25] which have also been extensively applied to flocking [26]. Figure 2(a) shows the experimentally obtained time derivative of the order parameter dO R /dt for different values of Δ. The grey scale of the data points marks the frequency of occurrence of each value in our experiments. Fitting the data to eq. (1) (red lines), we find excellent agreement with our experiments, which allows us to extract a and b explicitly ( figs. 2(b), (c)). An immediate consequence of eq. (1) is a bifurcation of the average steady state order parameter O R * = ±Re ( −a/b) at the critical point Δ c , where a(Δ c ) ≈ 0. In our system, this condition is met at Δ c ≈ 28 • . Indeed, this bifurcation is observed in the data, as shown in fig. 2(d). This is a first indication of critical slowing down because the timescale τ = a −1 diverges if a → 0. Figures 2(b), (c) also show that for Δ 40 • , a(Δ) ∝ Δ − Δ c , while b(Δ) is essentially constant. This scaling is reminiscent of the behavior of the expansion coefficients in Landau-type potentials near the critical point, where a generically vanishes in linear fashion [24]. For sufficiently large values of Δ, a(Δ) saturates because the measured O R is bounded between −1 and 1; this aspect, which is not captured by the simple Landau-type potential assumed above, would require the inclusion of higher-order terms.
Critical slowing down. -The presence of a bifurcation implies a timescale τ = a −1 which should diverge near the critical point Δ c . Direct evidence for critical slowing down is provided by investigating the dynamics of O R , which should strongly depend on the proximity of Δ to Δ c . Signatures of such dynamical regimes typically appear in autocorrelation functions, as have been observed for natural swarms [27]. Figures 3(a), (b) show the autocorrelation function of O R below the critical point (i.e., Δ < Δ c ), and the corresponding timescales arising from a double-exponential fit to the data, respectively. In general, the decay of the autocorrelation in our experiments is well described by two timescales. At Δ ∼ = 0 • the autocorrelation function of the rotational order parameter O R is mainly governed by a short timescale τ 1 ≈ 10 s. This timescale is observed independently of the value of Δ, and is attributed to the single particle dynamics within the group, which is largely determined by the attraction to the group's center and nearest neighbor repulsion that avoids particle collisions. Upon approaching Δ c , however, an additional time scale τ 2 appears, which eventually dominates the dynamics of O R near Δ c . This collective timescale shows clear evidence for critical slowing down, as it increases from around 200 s up to 1800 s close to the bifurcation. Proximity to the critical point therefore significantly affects relaxation of order parameter correlations, leading to enhanced collective behavior of the swarm. While, at first glance, long relaxation times may appear to be detrimental, they can provide the system with longer collective response and processing times than those available at the single particle level [28]. Above the bifurcation, where the potential has two minima, the collective dynamics is governed by noise-activated transitions between clockwise and counterclockwise rotating swirls. This dynamics, however, is difficult to resolve due to rotational reversals becoming increasingly unlikely as Δ is increased.
Hysteresis. -As shown above, the dynamical collective behavior of APs displays striking similarities with the Landau model of magnetism. This correspondence is further supported by breaking the clockwisecounterclockwise symmetry of swirls in analogy to biasing the spin orientation in a magnet by application of an external magnetic field. Such an asymmetry can be introduced by adding a slight preference to the particles' decision to swim to the left or the right of the center. Formally, this is achieved by making the particle's left/right decision dependent on ũ i = û i −h(r i ×e z ), where h quantifies the strength of the symmetry breaking and h < 0 (h > 0) corresponds to a preference for clockwise (counterclockwise) swirl rotation. Figures 4(a)-(f) show the influence of h on the order parameter for different values of the deviation angle Δ. The data points correspond to the measured mean value of O R in the final state, after h was initially set to +∞ (blue) and −∞ (red) and then suddenly reduced to a finite value of h. While for Δ < Δ c the value of O R is independent of the initial state (figs. 4(a), (b)), a pronounced hysteresis is observed for Δ > Δ c (figs. 4(c)-(f)), as expected. The data has been fitted with a function analogous to the mean-field solution of the Ising model [29], which describes our data well as seen by the solid red and blue curves in figs. 4(a)-(f). The free fit parameters A and B correspond to the saturation value of O R in the limit h → ∞ and the depth of the hysteresis loop. As the noise in our system is expected to be independent of Δ, the scaling of h, which corresponds to the curvature of the function, has been fixed to C = 8, yielding best agreement for all curves. It should be noted that eq. (2) describes the solution in the static T = 0 case and predicts an instantaneous jump at the turning point of the implicit function O R (h). It is known that this description fails to fully describe an Ising model with finite noise, where hysteresis is always a dynamic feature and not present as a static solution. Therefore, only the range from the initialized h to the turning point was taken into account to fit parameters, i.e., only data points where O R (h) > 0 for h initial = ∞ and O R (h) < 0 for h initial = −∞, respectively, were fitted. The data points O data R selected by these criteria are displayed as solid markers in fig. 4. Fitting of the parameters was performed by numerically calculating min A,B |A tanh(BO data R + Ch data )| 2 . The hysteresis curves also enable us to directly calculate the group's susceptibility which is given by χ = dOR dh | h=0 . After numerically solving for O R,0 = O R (A, B, h = 0), we obtain the susceptibility χ = 2AC/(1 − 2AB + cosh(2BO R,0 )). With the parameters obtained from the fit, we find a clear maximum whose position is consistent with the location of the bifurcation point. We remark that the susceptibility here is obtained from a direct response measurement, and not from equilibrium properties via the fluctuation-dissipation theorem; the latter is typically used for the analysis of χ in living systems [5].
Discussion. -Our experiments support a generic relationship between collective dynamical states and a critical point in a system of roughly 50 APs. Despite the finite size of particle groups, the system exhibits a clear dynamical bifurcation and critical slowing at the transition between a swarm and a swirl. We emphasize that the form of collective behavior observed in our system is absent when R o is comparable with the spatial extension of swirls. For the parameters used in this work this corresponds to N ≤ 15. In that regime, the dynamics is complex and cannot be captured by a simple model. Upon increasing the particle number somewhat, swirls and flocks emerge. It is remarkable that at these sizes, this system also exhibits clear bifurcation behavior which is captured very well by a simple coarse-grained model in which a "microscopic" quantity has been averaged over all constituents. Importantly, the system's proximity to the critical point is varied by mere changes of mutual interactions between APs which are similar to interactions typically considered in models of living systems. While there seems to be broad consensus in the literature regarding the advantages for biological systems operating close to a "critical point", most biological systems are of finite size and therefore not necessarily amenable to a proper field-theoretic description in the infinite size limit. Our work demonstrates that several advantages of a "critical point" remain (and can be classified in terms of the more basic concept of a bifurcation of a homogeneous order parameter) in a system consisting of a large but finite number of entities with well-controlled interactions. In general, we expect such critical behavior to persist even when interaction rules are not identical for all APs, but vary slightly within the group, reflecting natural distributions of behaviors in living systems. This notion is motivated by observations in Ising-type systems, where the critical point is not removed by the presence of certain defects or types of disorder [30,31]. Contrary to infinite systems, where the algebraic divergence of the susceptibility and collective relaxation time allows one to precisely estimate the position of a critical point, such dependence does not necessarily hold for finite-sized groups [32]. Even under such conditions, however, the proximity to a critical point may be announced by generic early warning signals as they are suggested to occur in critical systems [33,34]. At larger system sizes, it may be necessary to invoke field theories for a spatially resolved coarse-grained description. * * * The authors acknowledge helpful discussions with M. Krüger, H. Stark and A. Gambassi. CB acknowledges financial support from the ERC Advanced Grant ASCIR (Grant No. 693683). CMR acknowledges financial support from the UCT Emerging Researchers Programme.
Data availability statement: The data generated and/or analysed during the current study are not publicly available for legal/ethical reasons but are available from the corresponding author on reasonable request.