Flocking movement (all particles moving roughly parallelly, figure 1(a)) and vortical movement (all particles moving rotationally around the center-of-mass, figure 1(b)) are two important motion modes of active matter. Flocking movement is widely found in the migratory process of biological flocks, such as locusts and zebras [1]. Thousands of starlings may defend predators by moving in the flocking mode [2]. A swarm of unmanned aerial vehicles (UAV) can search for victims and track targets using flocking motion [3, 4]. The breakdown of ergodicity and the self-averaging of flocking motion is also of great interest [5]. By contrast, the vortical movement is very often a consequence of environmental disturbance or required by biological functions. A colony of bacillus subtilis grows to be the vortical morphotype as a response to a hostile environment, such as a diffusion-limited growth condition or an infertile agar substrate [6]. Ordemann et al [7] studied how a light beam induces Daphnia to move in the vortical mode. Blair et al [8] have confirmed experimentally that vertically vibrated granular rods can move vortically in a circular container with appropriate packing fraction, vibration amplitude, and frequency. Self-propelled microsensors powered by alternating electric fields can form the vortical state to achieve microfluidic mixture in microchannels [9]. A UAV swarm can implement detecting, monitoring, and loitering functions in the vortical motion mode [4, 10]. Besides, the collective behaviors of synthetic active matter are important for biomedicine applications [11].
Figure 1. Spontaneous switch of the system from the flocking mode (a) to the vortical mode (b) is a one-step nucleation process, while the reverse is a two-step nucleation process with a metastable state (c). The particles colored in red in (c) form a metastable flocking nucleus in the vortical background. |
For a finite-size system, the two motion modes may spontaneously switch between each other. The switch between the flocking and vortical swimming modes of golden shiners in a shallow tank was observed to be influenced by boundary and internal factors [12]. The bistable switching behavior was explained by considering the competition between attractive interaction and alignment [13], or between attractive and repulsive interactions [14]. The influences of alignment, attraction, and blind angle on the rate of mode-switching were also investigated [15, 16]. Chen and Hou [17] found that a vortex reversed its chirality through a hierarchical process: particles flipped their moving directions layer by layer from the outside to the inside. For a three-dimensional model of active Brownian particles, it has been found that the switching and hysteresis phenomena cannot exist simultaneously [18].
The switching dynamics between the flocking and vortical motion modes are essential for understanding the physical nature of active matter, which remains elusive. In this work, we modify the Vicsek model to simulate a flock of active matter with 1000 particles under the free boundary condition. With a certain combination of the orientational force and associated noise, we observe that the system spontaneously switches between the flocking mode and the vortical mode. By employing the classical nucleation theory (CNT) [19–21], we identify that the forward mode-switching from flocking motion to vortical motion is a one-step nucleation process with a critical nucleation size of 682 particles, while the backward mode-switching from vortical to flocking is a two-step nucleation process when the system goes over the first effective free energy barrier with a flocking nucleus of 179 particles to reach a metastable state with a flocking nucleus of 433 particles, before it goes over the second barrier with a flocking nucleus of 810 particles to be in the final flocking state.
The two-dimensional active matter is modelled by a Vicsek-like model modified from the Vicsek model [22]. For each particle, an orientational force ${F}^{{\rm{o}}}\left(t\right)$ is applied to mimic the attempt of adjusting its moving direction parallel to its nearest neighbours [23, 24]. A short-range repulsive force ${\mathop{F}\limits^{ \rightharpoonup }}^{{\rm{r}}}\left(t\right)$ is adopted to retain the finite exclusive volume of a particle [25], and the over-damped Langevin equation without translational noise [26] is used to update the particle velocity. With the free boundary condition applied, an attractive force ${F}^{{\rm{a}}}\left(t\right)$ pointing inwards is applied to the particles on the boundary to prevent them from easily fleeing [27]. Each particle has a fixed self-propulsion speed of ${v}_{0}=0.5,$ while the direction of the self-propulsion velocity ${\theta }_{i}\left(t\right)$ at time t is adjusted by ${F}_{i}^{{\rm{o}}}\left(t\right),$ ${F}_{i}^{{\rm{a}}}\left(t\right)$ if the particle is on the boundary, and a Gaussian white noise ${\xi }_{i}\left(t\right)$ mimicking the judging error to moving direction:
$\begin{eqnarray}\begin{array}{c}{\theta }_{i}\left(t\right)={\theta }_{i}\left(t-{\rm{\Delta }}t\right)+{\dot{\theta }}_{i}\left(t\right){\rm{\Delta }}t,\\ \,{\dot{\theta }}_{i}\left(t\right)=\alpha {F}_{i}^{{\rm{o}}}\left(t\right)+\varepsilon {F}_{i}^{{\rm{a}}}\left(t\right)+{\xi }_{i}\left(t\right),\\ {F}_{i}^{{\rm{o}}}\left(t\right)=\displaystyle \displaystyle \sum _{j\in {\rm{VN}}\left(i\right)}\sin \left({\theta }_{j}\left(t-{\rm{\Delta }}t\right)-{\theta }_{i}\left(t-{\rm{\Delta }}t\right)\right),\\ {F}_{i}^{{\rm{a}}}\left(t\right)=\left(1-{u}_{i}^{2}\right)\\ \times \,\sin \left(\theta \left({\mathop{r}\limits^{ \rightharpoonup }}_{i}^{{\rm{c}}}\left(t-{\rm{\Delta }}t\right)-{\mathop{r}\limits^{ \rightharpoonup }}_{i}\left(t-{\rm{\Delta }}t\right)\right)-{\theta }_{i}\left(t-{\rm{\Delta }}t\right)\right),\\ \left\langle {\xi }_{i}\left(t\right)\right\rangle =0,\\ \left\langle {\xi }_{i}\left(t\right){\xi }_{j}\left(t^{\prime} \right)\right\rangle ={\sigma }^{2}{\delta }_{ij}\delta \left(t-t^{\prime} \right),\end{array}\end{eqnarray}$
where ${\mathop{r}\limits^{ \rightharpoonup }}_{i}$ is the position of particle i, ${{\rm{\Sigma }}}_{j\in {\rm{VN}}\left(i\right)}$ means summing over the nearest neighbors of particle i determined by the Voronoi diagram [24], $\alpha $ is the strength of the orientational force, the strength of the attractive force $\varepsilon =0.5$ if the particle is on the boundary and 0 otherwise, the center-of-mass position of the neighbors of particle i is ${\vec{r}}_{i}^{{\rm{c}}}=\tfrac{1}{{n}_{i}}\displaystyle {\sum }_{j\in {\rm{VN}}\left(i\right)}{\mathop{r}\limits^{ \rightharpoonup }}_{j},$ and $\theta \left({\mathop{r}\limits^{ \rightharpoonup }}_{i}^{{\rm{c}}}-{\mathop{r}\limits^{ \rightharpoonup }}_{i}\right)$ is the direction of the vector connecting the particle to the center-of-mass of the neighbors. The reduced local average speed is $\begin{eqnarray}{u}_{i}=\left|\displaystyle \frac{1}{{n}_{i}}\displaystyle \displaystyle \sum _{j\in {\rm{VN}}\left(i\right)}\hat{n}\left({\theta }_{j}\left(t-{\rm{\Delta }}t\right)\right)\right|,\end{eqnarray}$
where ${n}_{i}$ is the number of neighbors of particle i and $\hat{n}\left({\theta }_{j}\right)$ is the normal vector of ${\theta }_{j}.$The velocity of the particle at time t is updated by ${\mathop{v}\limits^{ \rightharpoonup }}_{i}\left(t\right)={v}_{0}\hat{n}\left({\theta }_{i}\left(t\right)\right)+2{v}_{0}{\mathop{F}\limits^{ \rightharpoonup }}_{i}^{{\rm{r}}}\left(t\right),$ and the repulsive force
$\begin{eqnarray}\begin{array}{c}{\mathop{F}\limits^{ \rightharpoonup }}_{i}^{{\rm{r}}}\left(t\right)\\ =\,\displaystyle \displaystyle \sum _{j\in {\rm{VN}}\left(i\right)}{\hat{r}}_{ij}\left(t-{\rm{\Delta }}t\right)/\left(1+\exp \left(\displaystyle \frac{{r}_{ij}\left(t-{\rm{\Delta }}t\right)}{0.5}-1\right)\right),\end{array}\end{eqnarray}$
where ${\hat{r}}_{ij}={\mathop{r}\limits^{ \rightharpoonup }}_{ij}/{r}_{ij},$ ${r}_{ij}=\left|{\mathop{r}\limits^{ \rightharpoonup }}_{ij}\right|=\left|{\mathop{r}\limits^{ \rightharpoonup }}_{i}-{\mathop{r}\limits^{ \rightharpoonup }}_{j}\right|.$ Finally, the position of the particle is updated by ${\mathop{r}\limits^{ \rightharpoonup }}_{i}\left(t\right)={\mathop{r}\limits^{ \rightharpoonup }}_{i}\left(t-{\rm{\Delta }}t\right)+{\mathop{v}\limits^{ \rightharpoonup }}_{i}\left(t\right){\rm{\Delta }}t$ with ${\rm{\Delta }}t=0.05.$To distinguish different states, a motion order parameter (MOP) is defined as
$\begin{eqnarray}\lambda =\varphi -\left|\omega \right|,\end{eqnarray}$
where the reduced average speed $\begin{eqnarray}\varphi =\left|\displaystyle \sum _{i=1}^{N}{\mathop{v}\limits^{ \rightharpoonup }}_{i}\left(t\right)\right|/N{v}_{0},\end{eqnarray}$
and the reduced angular momentum $\begin{eqnarray}\omega =\displaystyle \frac{1}{N}\displaystyle \sum _{i=1}^{N}\displaystyle \frac{\left({\mathop{r}\limits^{ \rightharpoonup }}_{i}\left(t\right)-\displaystyle {\sum }_{j=1}^{N}{\mathop{r}\limits^{ \rightharpoonup }}_{j}\left(t\right)/N\right)\times {\mathop{v}\limits^{ \rightharpoonup }}_{i}\left(t\right)}{\left|{\mathop{r}\limits^{ \rightharpoonup }}_{i}\left(t\right)-\displaystyle {\sum }_{j=1}^{N}{\mathop{r}\limits^{ \rightharpoonup }}_{j}\left(t\right)/N\right|{v}_{0}}.\end{eqnarray}$
The range of $\varphi $ is [0, 1], $\omega $ is [−1, 1], and $\lambda $ is [−1, 1]. A state is identified as flocking if $\lambda \gt 1/3,$ vortical if $\lambda \lt -1/3,$ and disordered otherwise (see the appendix ).
A flock with 1000 particles has been simulated by the above model with $\alpha =0.02$ and $\sigma =1.12$ for $1.6\times {10}^{7}$ steps, when either the flocking mode or the vortical mode can be stable for a period of time (∼105 steps), before it spontaneously switches to the other mode. The time evolution of $\lambda $ during a simulation is plotted in figure 2. This simulation procedure was repeated with different random initial configurations for 24 times, and a sampling interval of 500 steps was used to collect simulation data.
Figure 2. Time evolution of $\lambda $ during one simulation. |
As illustrated in figure 1, the visual examination of the simulation trajectories reveals the following switching dynamics between the two modes. In the globally flocking state (figure 1(a)), some local vortices instantaneously appear and disappear due to thermal fluctuation induced by $\xi $, before a large enough vortex is generated and quickly grows up, turning the whole system into the globally vortical state (figure 1(b)). Reversely, starting from the vortical state (figure 1(b)), some local flocking pieces instantaneously appear and disappear before entering a metastable state with almost half of the particles forming a large flocking nucleus and the others still vortical (figure 1(c)), and after a while the system goes out of the metastable state to form the global flocking state (figure 1(a)).
The physical features and microscopic mechanisms of the above dynamical procedures are studied phenomenologically from the viewpoint of the CNT, which describes the free energy landscape of a nucleus in the one-step nucleation process as7 )) and the bulk free energy (the second term at the RHS of equation (7 )) contributing to the nucleus growth. Since active matter is always far from equilibrium, the term ‘free energy’ is effective when the CNT is applied to study the dynamics of the active matter.
$\begin{eqnarray}{\rm{\Delta }}G\left(n\right)=\gamma A\left(n\right)-{\rm{\Delta }}\mu n,\end{eqnarray}$
where n is the number of particles in the nucleus, $\gamma $ is the surface tension of the nucleus, $A\left(\sim {n}^{1/2}\right)$ is the surface area, and ${\rm{\Delta }}\mu $ is the chemical potential difference between the background and the nucleus. The CNT considers the competition between the surface free energy (the first term at the right hand side (RHS) of equation (A local cluster with two or more particles is identified as a flocking nucleus if the MOP of the cluster ${\lambda }_{{\rm{c}}}\gt 1/3,$ and as a vortical nucleus if ${\lambda }_{{\rm{c}}}\lt -1/3.$ To analyze the forward mode-switching dynamics, the size distribution of vortical nuclei ${p}_{{\rm{v}}}\left(n\right)$ is calculated and plotted in figure 3(a). All simulation data is used in the calculation and the flocking motion is treated as ‘background’ with respect to the vortical nuclei. The corresponding effective free energy landscape ${\rm{\Delta }}{G}_{{\rm{v}}}\sim -\,\mathrm{ln}\left({p}_{{\rm{v}}}\right)$ is plotted in figure 3(b). There is only one maximum, indicating that it is a one-step nucleation process with only one effective free energy barrier to overcome. Under the framework of the CNT, we have found that the effective free energy landscape can be perfectly fitted by
$\begin{eqnarray}\begin{array}{c}{\rm{\Delta }}{G}_{{\rm{v}}}\left(n\right)={\gamma }_{{\rm{v}},0}{n}^{{g}_{{\rm{v}},0}}-{\gamma }_{{\rm{v}},1}{n}^{-{g}_{{\rm{v}},1}}\\ -\,{\rm{\Delta }}{\mu }_{{\rm{v}}}\exp \left({g}_{{\rm{v}},2}n\right),\end{array}\end{eqnarray}$
where ${\gamma }_{{\rm{v}},0}=10.0,$ ${g}_{{\rm{v}},0}=0.0545,$ ${\gamma }_{{\rm{v}},1}=7.71,$ ${g}_{{\rm{v}},1}=0.363,$ ${\rm{\Delta }}{\mu }_{{\rm{v}}}=2.03\times {10}^{-3},$ and ${g}_{{\rm{v}},2}=6.89\times {10}^{-3}.$ The first two terms at the RHS roughly correspond to the effective surface energy while the third term corresponds to the effective bulk energy. By locating the maximal value of the fitted curve, we identify that the critical nucleus size is 682.
Figure 3. (a) Size distributions of the vortical and flocking nuclei. (b) Corresponding effective free energy landscapes along with their fitted curves. |
The backward mode-switching dynamics from vortical to flocking is similarly analysed by treating the vortical motion as the background with respect to the flocking nuclei. The size distribution and the corresponding effective free energy landscape of the flocking nuclei are also plotted in figures 3(a) and (b). The effective free energy landscape is fitted by
$\begin{eqnarray}\begin{array}{c}{\rm{\Delta }}{G}_{{\rm{f}}}\left(n\right)={\gamma }_{{\rm{f}},0}{n}^{{g}_{{\rm{f}},0}}-{\gamma }_{{\rm{f}},1}{n}^{-{g}_{{\rm{f}},1}}-{\rm{\Delta }}{\mu }_{{\rm{f}},0}\exp \left({g}_{{\rm{f}},2}n\right)\\ -\,{\rm{\Delta }}{\mu }_{{\rm{f}},1}\exp \left(\displaystyle \frac{-{\left(n-{n}_{{\rm{m}}}\right)}^{2}}{{g}_{{\rm{f}},3}}\right),\end{array}\end{eqnarray}$
where ${\gamma }_{{\rm{f}},0}=15.3,$ ${g}_{{\rm{f}},0}=0.0914,$ ${\gamma }_{{\rm{f}},1}=12.1,$ ${g}_{{\rm{f}},1}=0.0192,$ ${\rm{\Delta }}{\mu }_{{\rm{f}},0}=1.35,$ ${g}_{{\rm{f}},2}=1.94\times {10}^{-3},$ ${\rm{\Delta }}{\mu }_{{\rm{f}},1}=4.73,$ ${g}_{{\rm{f}},3}\,=1.58\times {10}^{5},$ and ${n}_{{\rm{m}}}=433$ are located at the local minimum of the effective free energy landscape. The first two terms roughly correspond to the effective surface energy, while the third term is to the effective bulk energy. The fourth term reflects the appearance of the metastable state. The two local maxima of the effective free energy landscape locate respectively at 179 and 810 with the minimal value of 433 in between, indicating that the backward mode switch is a two-step nucleation process.Both equations (8 ) and (9 ) have two fitted surface terms, while equation (7 ) has only one surface term for the original CNT. In order to understand this difference, we determine the surface particles for each nucleus by calculating the Voronoi diagram, which are further separated into two categories: those residing on the boundary of the whole system are denoted as boundary particles and the others are denoted as interactive particles. As shown in figure 4(a), by fitting in the whole data range, the number of surface particles can be fitted as ${N}_{{\rm{S}},{\rm{v}}}=1.34{n}^{0.450}$ for the vortical nuclei and ${N}_{{\rm{S}},{\rm{f}}}=1.36{n}^{0.443}$ for the flocking nuclei, respectively, qualitatively consistent with the surface term in equation (7 ) and the first terms at the RHS in equations (8 ) and (9 ). Figure 4(b) shows that, by fitting in the range from beginning to its maximum, the number of interactive particles can be fitted as ${N}_{{\rm{I}},{\rm{v}}}=-11{n}^{-0.248}+11$ for the vortical nuclei and ${N}_{{\rm{I}},{\rm{f}}}=-10{n}^{-0.340}+10$ for the flocking nuclei, respectively, qualitatively consistent with the second terms at the RHS in equations (8 ) and (9 ). Based on the above analysis of the surface particles, the physical interpretation of the appearance of the extra surface term may be that, as the difference of the interactions between flocking and vortical particles is much weaker than that between molecules inside and outside nuclei described by the CNT. The boundary particles applied with the attractive force ${F}^{{\rm{a}}}\left(t\right)$ pointing inwards are at the phenomenological level analogous to the surface particles in the original CNT. Meanwhile, the interactive particles do not experience such an attractive force, so they contribute an extra surface term to be subtracted from the original term.
Figure 4. The number of surface particles (a) and interactive particles (b) of the vortical and flocking nuclei, respectively. |
Through computer simulation with a Vicsek-like model, we have revealed that the dynamics of the forward and backward spontaneous switches between the flocking mode and the vortical mode of a finite-size active matter are asymmetric. The forward mode switch from flocking to vortical is a one-step nucleation process. Small vortical nuclei appear and disappear instantaneously in the flocking state due to thermal fluctuation until one larger than the critical nucleus grows up to turn the whole system into the vortical state by climbing over the single effective free energy barrier. Reversely, the backward mode switch from vortical to flocking is a two-step nucleation process. In the first stage, local pieces of small flocking nuclei generate in the vortical state until the system goes over the first effective free energy barrier to enter the metastable state with almost half of the particles becoming flocking. In the second stage, the system goes over the second barrier to become globally flocking. The physical reason causing the asymmetry of the mode-switching dynamics is that a flocking nucleus is easier to form than a vortical one due to the orientational force, so the forward procedure has to form a very large vortical nucleus simultaneously to fulfil the switch, while the backward one creates flocking structures in two sequential steps to have a lower overall effective switching free energy barrier in comparison with the case of one-step nucleation. This discovery deepens our understanding of motion dynamics in active matter and is anticipated to guide the efficient manipulations of artificial active matter, such as robots and UAVs.
