# Consequence of anisotropy on flocking: the discretized Vicsek model

Mintu Karmakar⋆<EMAIL_ADDRESS>School of Mathematical & Computational
Sciences, Indian Association for the Cultivation of Science, Kolkata – 700032,
India.  Swarnajit Chatterjee⋆<EMAIL_ADDRESS>Center for
Biophysics & Department for Theoretical Physics, Saarland University, 66123
Saarbrücken, Germany.  Raja Paul<EMAIL_ADDRESS>School of Mathematical
& Computational Sciences, Indian Association for the Cultivation of Science,
Kolkata – 700032, India.  Heiko Rieger<EMAIL_ADDRESS>Center for
Biophysics & Department for Theoretical Physics, Saarland University, 66123
Saarbrücken, Germany.  INM – Leibniz Institute for New Materials, Campus D2 2,
66123 Saarbrücken, Germany.

###### Abstract

We numerically study a discretized Vicsek model (DVM) with particles orienting
in $q$ possible orientations in two dimensions. The study probes the
significance of anisotropic orientation and microscopic interaction on the
macroscopic behavior. The DVM is an off-lattice flocking model like the active
clock model [ACM; EPL 138, 41001 (2022)] but the dynamical rules of particle
alignment and movement are inspired by the prototypical Vicsek model (VM). The
DVM shows qualitatively similar properties as the ACM for intermediate noise
strength where a transition from macrophase to microphase separation of the
coexistence region is observed as $q$ is increased. But for small $q$ and
noise strength, the liquid phase appearing in the ACM at low temperatures is
replaced in the DVM by a configuration of multiple clusters with different
polarization which does not exhibit any long-range order. We find that the
dynamical rules have a profound influence on the overarching features of the
flocking phase. We further identify the metastability of the ordered liquid
phase subjected to a perturbation.

Flocking; Vicsek Model discretization; Metastability

$\star$$\star$footnotetext: These authors have contributed equally to this
work

##  I Introduction

Active matter systems represent a fascinating class of materials composed of
self-propelled entities that convert energy into mechanical motion, giving
rise to complex and often out-of-equilibrium behaviors Ramaswamy (2010);
Marchetti  _et al._ (2013); Gompper  _et al._ (2020); Needleman and Dogic
(2017). The emergent phenomena in active matter systems has gained significant
attention in recent years due to their potential applications in various
physical, biological, and engineering systems Marchetti  _et al._ (2013);
Vernerey  _et al._ (2019); Ghosh  _et al._ (2021). Active matter exhibits
dynamic behaviors like collective motion Marchetti  _et al._ (2013), pattern
formation Bär  _et al._ (2020); Shaebani  _et al._ (2020), and even the
ability to exhibit controlled transport Sanchez  _et al._ (2012). These
systems encompass a wide range of physical, chemical, and biological entities,
from swimming bacteria Peruani  _et al._ (2012), mammalian herds Garcimartín
_et al._ (2015), fish schools Becco  _et al._ (2006); Calovi  _et al._ (2014),
and sterling flocks to amoeba and bacteria colonies Steager  _et al._ (2008),
to the cooperative behavior of cytoskeletal filaments and molecular motors in
living cells Schaller  _et al._ (2010); Sumino  _et al._ (2012); Sanchez  _et
al._ (2012) or in vitro environments to synthetic colloidal particles equipped
with motors Veigel and Schmidt (2011); Wong  _et al._ (2016). To understand
and unravel the fundamental principles governing active matter systems, new
models Shaebani  _et al._ (2020) have been developed in the last two decades.

The Vicsek model (VM), introduced by Vicsek and collaborators in 1995 Vicsek
_et al._ (1995), provides a fundamental framework for studying the collective
behavior of particles under aligning interactions. In this model, particles
adjust their velocities to align with the average velocities of neighboring
particles, leading to the emergence of coherent motion and ordered patterns.
The VM has played a crucial role in advancing our understanding of flocking
dynamics Toner and Tu (1995, 1998); Toner (2012); Solon  _et al._ (2015a). At
low particle density or high noise, particles move in random directions, and
no long-range order is observed. The transition from the gas phase at high
noise and low density to the polar ordered Toner-Tu phase at low noise and
high density, displaying long-range order (LRO) by a coherent motion of all
particles, is first order Chaté  _et al._ (2008). But, in contrast to
conventional first-order phase transition scenarios, the coexistence phase of
the VM can manifest as either multiple bands of particles moving collectively,
a phenomenon known as microphase separation Solon  _et al._ (2015a, b), or a
polar-ordered cross sea phase Kürsten and Ihle (2020), primarily driven by
giant number fluctuations (GNF) Solon  _et al._ (2015a).

Nonetheless, it is important to note that the VM posits a continuous range of
possible directions for the motion of particles. However, when considering a
scenario in which particles are constrained to discrete, equidistant angular
orientations within a two-dimensional plane, such as in the active clock model
(ACM) Solon  _et al._ (2022); Chatterjee  _et al._ (2022), the VM-inspired
dynamical principles governing particle alignment and movement remain
uncharted territory. In a recent study on the ACM Solon  _et al._ (2022), it
was revealed that in large systems, any values of discrete orientations result
in significant changes in phenomenology when compared to the VM. These changes
include the loss of long-range correlations, the pinning of global order, and
the transformation of coexistence bands into a single moving domain.
Additionally, another study on the ACM Chatterjee  _et al._ (2022) with
different dynamical rules shows that for a small number of directions, the ACM
mirrors the active Potts model (APM) Chatterjee  _et al._ (2020); Mangeat  _et
al._ (2020), exhibiting macrophase separation and reorientation transition.
Conversely, with more directions, the ACM transitions towards the VM,
displaying microphase separation and transverse bands without reorientation.
Remarkably, the transition in the $q\to\infty$ limit of the ACM Chatterjee
_et al._ (2022), known as the active XY model, shares the same universality
class as the VM. Motivated by these findings, in this paper, we undertake an
extensive computational investigation that examines in detail a $q$-state
discretized version of the Vicsek model (DVM) where the rules of particle
alignment and movement follow the Vicsek protocol.

We ask several intriguing questions that persist within the context of the
DVM, e.g., (a) How does discretizing the directions of particles in the VM,
affect the overall diverse collective dynamics and steady-state phases? (b)
What is the impact of $q$ and system size on the coexistence region (micro- or
macro-phase separation)? (c) How do the behavior of the density fluctuations,
the direction of the system’s global order, and the behavior of correlation
functions in the liquid phase correspond with the self-organized patterns in
the phase-coexistence region? (d) What is the nature of the DVM liquid phase
as a function of $q$? To answer these questions, we study the DVM in an off-
lattice domain, focusing on the three key factors: the anisotropy and
degeneracy parameter $q$, noise level $\eta$, and system size.

The paper is organized as follows: after introducing the model in Sec. II, we
present our numerical results in Sec. III. Finally, we summarize and discuss
the implications of our findings in Sec. IV.

##  II Model

We consider $N$ self-propelled particles within a two-dimensional off-lattice
domain of size $L_{x}\times L_{y}$ ($L_{x}>L_{y}$ for rectangular domain and
$L_{x}=L_{y}=L$ for square domain) with periodic boundary conditions. Akin to
the two-dimensional VM, each point particle $i$ is endowed with an off-lattice
position vector $\bm{r}_{i}=(x_{i},y_{i})$ and moves with a constant speed
$v_{0}$ in individual directions given by a unit orientation vector
${\bm{\sigma}}_{i}=(\cos\theta_{i},\sin\theta_{i})$ with an orientation angle
$\theta_{i}\in(0,2\pi)$ where

$\theta_{i}=\frac{2\pi n_{i}}{q}\,,$ (1)

and $n_{i}=\\{0,1,2,\cdots,(q-1)\\}$ denote discrete orientations of the
particles. $q$ denotes the ground state degeneracy where each particle can
only have discrete orientations allowed by the $q$ value and therefore, the
continuous $U(1)$ symmetry of the VM is replaced by the discrete $Z_{q}$
symmetry.

At each discrete time step $\Delta t=1$, a particle $i$ moves a fixed distance
$v_{0}\Delta t$ and interacts with $\mathcal{N}_{i}$ neighboring particles
within a circle of unit radius around it. The position evolves in the
following way:

$\bm{r}_{i}^{t+\Delta t}=\bm{r}_{i}^{t}+v_{0}\bm{\sigma}_{i}^{t+\Delta
t}\Delta t\,,$ (2)

while the new orientation is determined by a projection of the updated
orientation proposed by the Vicsek rule onto one of the $q$ allowed
directions:

$\theta_{i}^{t+\Delta t}=\mathbb{P}(\bar{\theta}_{i}^{t}+\eta\xi_{i}^{t})\,,$
(3)

where $\mathbb{P}$ is the projection and $\bar{\theta}_{i}^{t}$ is the
orientation angle of a spin-weighted sum

$\bar{\bm{\sigma}}_{i}^{t}=\sum_{j\in\mathcal{N}_{i}}\bm{\sigma}_{j}^{t},$ (4)

of orientation vectors of neighboring particles. $\xi_{i}^{t}\in[-\pi,\pi]$ is
a scalar noise uniformly distributed and uncorrelated for all sites and times.
Such noise is often called $white$ since it has a flat Fourier spectrum.
$\eta$ is the noise amplitude.

We define the projection onto the allowed directions probabilistically by

$\mathbb{P}(\theta)=\begin{cases}\theta_{1}&\text{with probability}\ 1-p\,,\\\
\theta_{2}&\text{with probability}\ p\,,\end{cases}$ (5)

where $\theta_{1}$ and $\theta_{2}$ are the two allowed directions which are
closest to $\theta$, such that $\theta_{1}=2\pi n/q<\theta$ and
$\theta_{2}=2\pi(n+1)/q>\theta$ for some $n$. The probability $p\in[0,1]$ is
given by $p=\frac{q}{2\pi}(\theta-\theta_{1})$, i.e. minimal (0) for $\theta$
close to $\theta_{1}$ and maximal (1) for $\theta$ close to $\theta_{2}$ (c.f.
Fig. 1). Note that then for all particles going into the direction, say,
$\theta_{1}$, the probability to turn stochastically into another direction,
say $\theta_{2}$, is of the order $\sim$ $\frac{q}{2}\eta$, i.e. small for
small $q$ and noise $\eta$.

It should be mentioned here that the dynamical rules governing DVM are
different from the dynamical rules of the $q$-state ACM Chatterjee  _et al._
(2022). In the DVM, as per Eq. (2), the hopping probability of a particle with
orientation angle $\theta$ along another discrete direction $\phi$ is zero
whereas the probability of hopping along non-preferred directions is nonzero
in the ACM and depends on the self-propulsion velocity of the particle.
Therefore, the transverse fluctuations, say in $q=4$ ACM, are stronger than
the $q=4$ DVM. In the ACM, transverse fluctuations mainly originate from the
nonzero hopping probability of a particle along its non-preferred directions.
Here thermal fluctuation, through the inverse temperature $\beta$, also plays
an important role as it controls the flipping dynamics. However, in the large
$\beta$ limit of the ACM, similar to the small $\eta$ limit of the DVM, the
probability that all particles moving in a particular direction will flip
their orientation to another direction is also very small. But, unlike DVM, an
ACM particle can move along a direction different than its orientation angle.
This crucial difference in the hopping dynamics, as we will see, plays an
essential role in the steady-state pattern formation of the DVM at low $\eta$
and $q$.

Figure 1: (Color online) Schematic of the DVM for $q=4$. (a) Four allowed
orientations for the particle are $0,\pi/2,\pi,$ and $3\pi/2$ with
$\theta_{i}^{t}=0$. The circular neighborhood (of radius $R$) of particle $i$
contains neighboring particles to calculate $\theta_{i}^{t+\Delta t}$. (b) The
new orientation proposed by the Vicsek rule (blue dotted arrow) will be
projected either along $\theta_{1}=0$ or $\theta_{2}=\pi/2$ probabilistically.

DVM control parameters are the average particle density
$\rho_{0}=N/L_{x}L_{y}$, noise strength $\eta$, particle velocity $v_{0}=0.5$
(unless mentioned otherwise), and the measure of anisotropy $q$. According to
Eq. (1) a large $q$ signifies weak anisotropy while a small $q$ signifies
strong anisotropy.

We performed numerical simulations of the stochastic process with parallel
updates. The initial configuration is prepared homogeneously by assigning
random orientations and positions to the particles as defined in Eq. (1) and
Eq. (2), respectively. After the initialization, we let the system evolve
under various control parameters for $t_{\rm eq}\sim 10^{5}$ to reach the
steady state followed by measurements of various quantities until the maximum
simulation time $t_{\rm max}\sim 10^{6}$.

##  III Numerical Results

Collective motion & phase diagram. We present the typical non-equilibrium
steady-state configurations of the DVM in Fig. 2 for $q=9$ and density
$\rho_{0}=2$. The system exhibits phases of polar ordered liquid (a), liquid-
gas coexistence (b–c), and disordered gas (d) as noise strength $\eta$ is
increased from 0.1 to 0.5. The phase-coexistence region is characterized by a
low-noise cross-sea phase (b) and a high-noise band phase (c). The cross-sea
phase has particle density much higher at the crossing points than anywhere
else and has recently been reported as the fourth phase of the VM Kürsten and
Ihle (2020). We notice that such patterns assemble spontaneously without an
external drive at certain parameter values. Conversely, the band phase is a
collection of high-density bands moving parallelly along a specific direction
at a constant speed. Polar flocks [the homogeneous ordered liquid phase shown
in Fig. 2(a)] can be observed in a large class of active matter systems and
have been considered robust to fluctuations. But recent studies have argued
that liquid polar flocks are metastable to the presence of a small obstacle
Codina  _et al._ (2022) or to the nucleation of an opposite-phase droplet
Benvegnen  _et al._ (2023) and triggers counter-propagating dense clusters
leading to the reversal of the liquid flow. In light of these observations, in
the subsequent part of this paper, we will investigate the stability of the
DVM liquid phase. Fig. 2 suggests that the system exists in distinct phases,
which we will characterize next.

Figure 2: (Color online) Snapshots on a square domain of size $1024^{2}$
exhibiting the four phases of the DVM for $q=9$. (a) Polar ordered liquid
phase ($\eta=0.1$). (b) Cross-sea state ($\eta=0.3$). (c) Band state
($\eta=0.35$). (d) Disordered gas phase ($\eta=0.5$). Dark color represents
high particle density and red arrows indicate the average direction of motion.
Parameters: $\rho_{0}=2$, $v_{0}=0.5$. Figure 3: (Color online) $\eta-q$ phase
diagram of the DVM illustrated by snapshots on a $1024^{2}$ domain at time
$t=10^{5}$. As a function of $\eta$ and $q$, we observe six distinct self-
organized patterns: cluster ($\eta=0.1$, $q=4$), macrophase ($\eta=0.3$,
$q=4$), microphase ($\eta=0.4$, $q=6\to 16$), cross-sea ($\eta=0.3$, $q=8\to
10$), ordered liquid ($\eta=0.2$, $q=6\to 16$), and disordered gas
($\eta=0.5$, $q=4\to 16$). Parameters: $\rho_{0}=2$, $v_{0}=0.5$.

The non-equilibrium steady-state behavior of the DVM is illustrated by
representative late-stage snapshots as a function of noise strength $\eta$ and
anisotropy parameter $q$ (see Fig. 3). We observe six distinct self-organized
patterns in the ($q$, $\eta$) plane which completely describe the DVM. At low
noise and $q$, we observe a locally ordered cluster phase. Although each
cluster is highly dense and polar, the system as a whole does not possess any
net polarization (see Appendix A). This cluster phase, which has not been
observed earlier in any flocking models, appears only for small $q$ and $\eta$
when the system is strongly discretized and the impact of fluctuation is
insignificant.

The appearance of a cluster phase in the $q=4$ DVM for small $(q,\eta)$ and
high density instead of a polar ordered phase similar to the 4-state APM
Chatterjee  _et al._ (2020); Mangeat  _et al._ (2020) and ACM Solon  _et al._
(2022); Chatterjee  _et al._ (2022) can be attributed to the absence of
transverse fluctuations through hopping. In DVM, similar to the VM, particle
movement is solely controlled by the orientation $\theta$ [see Eq. (2)].
Therefore, $q=2$ DVM only manifests one-dimensional movement (along the
$x$-axis) of high-density clusters of self-propelled particles having
orientations $\theta=0$ and $\theta=\pi$. We observe that these clusters never
coalesce due to the lack of transverse fluctuations and thus never form a band
or polar liquid phase. Similar observations are made when the constant
transverse diffusion is switched off in the two-dimensional active Ising model
(AIM) Solon and Tailleur (2015), although the one-dimensional AIM Benvegnen
_et al._ (2022) exhibits a flocking state consisting of a single dense ordered
aggregate. This signifies that diffusion along the non-preferred hopping
directions plays a crucial role in the formation of large liquid domains. For
APM and ACM at large $\beta$ and small $q$ and DVM with small $\eta$ and $q$,
the probability of transverse flipping is very small as fluctuations are weak.
However, the nonzero finite hopping rates along the unbiased directions
facilitate the formation of large liquid domains in the APM and ACM whereas an
absence of this forms a cluster phase in the DVM. Finally, it is the interplay
of $q$ and $\eta$ that determines the fate of the cluster phase. If $q$ is
small but fluctuation is large, the rigid cluster phase can relax and form a
large ordered domain (see the snapshot for $q=6$ and $\eta=0.2$). On the
contrary, if $\eta$ is small but $q$ is large, the weak anisotropy helps the
cluster phase merge into a large ordered domain (see the snapshot for $q=10$
and $\eta=0.1$). This also explains why cluster size increases with $q$ for a
fixed $\eta$ (see the snapshots for $\eta=0.1$).

Beyond the cluster phase, for $\eta=0.2$, we observe the emergence of the
polar ordered liquid phase ($q\geqslant 6$). For intermediate noise strength
($\eta=0.3-0.4$), a transition from macrophase separation ($q=4$, a single
liquid band moving through the gaseous background) to microphase separation
(multiple bands moving parallelly or in a cross-sea pattern through the gas
phase) in the coexistence region is observed where the number of bands
increases with $q$. This is a consequence of having more allowed particle
orientations through $q$. The cross-sea phase, where the interactions become
more intensive for the characteristics of the band structure, appears between
the polar liquid phase and the band state Xue  _et al._ (2023) for a fixed
$q$. This phase is not simply a superposition of two band waves, but an
independent self-organized complex pattern with an inherently selected
crossing angle. One should note here that the cross-sea phase can not be
formed with the single band macrophase separation, a single cross-sea pattern
at least requires two bands where they cross each other approximately at an
angle $\sim\pi/4$. It can be noticed that for $\eta=0.3$, as we increase $q$
from $q=7$ to $q=16$, we observe that the self-organized patterns change from
band structure $(q=7)$ to cross-sea pattern $(q=8$, 9, 10) to a pattern
showing the outset of the polar liquid phase ($q=16$). For the same values of
the control parameters, the VM $(q\to\infty)$ exhibits a similar phase to the
$q=16$ DVM where the phase point on the $(\eta,\rho_{0})$ diagram almost lies
on the liquid binodal Solon  _et al._ (2015a). This happens because anisotropy
becomes weaker with $q$ and, for $q\geqslant 16$, the characteristic of the
system becomes similar to the VM. The VM does exhibit a cross-sea phase but at
a different parameter regime Kürsten and Ihle (2020). For very large noise,
e.g., $\eta=0.5$, we observe a disordered gas phase irrespective of $q$.

Since the emergence of different phases in the DVM depends on the spatial
anisotropy, in Appendix B, we show that the cluster phase is ubiquitous for
small ($q,\eta$) and present a cluster size analysis of the $q=4$ DVM for
varying noise. See also Appendix C and Appendix D for more discussions on the
DVM phases with spatial anisotropy (rectangular domain).

Based on the above observations, we can quantify the DVM phase diagram. Notice
that from Fig. 3, we have visually distinguished four different phases so far
Kürsten and Ihle (2020): the ordered phase, the cross-sea phase, the band
phase, and the disordered gas phase. It is, however, challenging to define the
phase boundaries by visual inspection.

Figure 4: (Color online) The structural order parameter $\langle C_{2}\rangle$
versus the noise amplitude $\eta$ for several values of $q$. The dotted curves
between two successive points are approximated using spline interpolation.
Parameters: $L=1024$, $v_{0}=0.5$, and $\rho_{0}=2$.

Earlier, phase-separated density profiles were used to identify the binodals
that delimit the gas and liquid phases from the coexistence region Solon  _et
al._ (2015a). Yet, this technique can not differentiate between the distinct
self-organized patterns we observe in the coexistence region (i.e. macrophase
separation, microphase separation, and cross-sea) of the DVM. One might also
want to distinguish different structures in terms of the global polar order
parameter or the global magnetization defined as:

$\displaystyle m=|{\bf m}|$
$\displaystyle=\frac{1}{N}\left|\sum_{i}{\bm{\sigma}}_{i}\right|\,.$ (6)

The polar order parameter can be used to study the transition between the
disorder gas phase and the band states but can not differentiate between the
band and the cross-sea states Kürsten and Ihle (2020). The Binder cumulant
constructed from this order parameter shows only the transition between the
disordered phase and the band state Kürsten and Ihle (2020). The polar order
parameter $m$ is also insensitive around the cross-sea state. Thus for precise
quantification of different phases and their boundaries, we compute the
structural order parameter $C_{2}$ Kürsten and Ihle (2020); Xue  _et al._
(2023); Kürsten  _et al._ (2020) as follows:

$\displaystyle C_{2}$ $\displaystyle=N^{2}\int G_{2}({\bf r}_{1},{\bf
r}_{2})d{\bf r}_{1}d{\bf r}_{2}\Theta(R-|{\bf r}_{1}|)\Theta(R-|{\bf r}_{2}|)$
$\displaystyle=\left(\frac{N}{L^{2}}\right)^{2}\int_{\mathbb{R}^{2}}\left[g(|{\bf
r}_{1}-{\bf r}_{2}|)-1\right]$ $\displaystyle\times\Theta(R-|{\bf
r}_{1}|)\Theta(R-|{\bf r}_{2}|)d{\bf r}_{1}d{\bf r}_{2}\,,$ (7)

where $G_{2}({\bf r}_{1},{\bf r}_{2})=P_{2}({\bf r}_{1},{\bf
r}_{2})-P_{1}({\bf r}_{1})P_{2}({\bf r}_{2})$ for one- and two-particle
probability density functions $P_{1}$ and $P_{2}$, $\Theta$ is the Heaviside
step function, and $R$ is the distance around an arbitrary fixed point in
space. For macroscopically isotropic systems, $G_{2}$ can be expressed in
terms of the pair correlation function $g(r)$ which physically signifies the
probability of finding a particle at a distance $r$ relative to that of a
given reference particle and provides a statistical description of the local
packing and particle density of the system. It has been shown Kürsten and Ihle
(2020); Kürsten  _et al._ (2020) that $C_{2}$ performs better than $m$ in
capturing the structural change and the Binder cumulant of the $C_{2}$ is also
more efficient in distinguishing different features than the Binder cumulant
of the order parameter $m$. In practice, to compute $C_{2}$, we take all pairs
of particles, draw circles of radius $R$ around them, and calculate the
overlap area of the two circles. The overlap area is given by $A_{\rm
overlap}=2R^{2}\cos^{-1}\left(\frac{d}{2R}\right)-\frac{d}{2}\sqrt{4R^{2}-d^{2}}$
Kürsten  _et al._ (2020), where $d$ is the distance between the centers of the
circles.

Fig. 4 shows the dependence of the time-averaged structural order parameter
$\langle C_{2}\rangle$ with respect to $\eta$ for different $q$. The value of
$\langle C_{2}\rangle$ is the lowest for disorder gas and becomes maximum when
particles are clustered. For the band phase, $\langle C_{2}\rangle$ values for
macrophase separation are larger than the microphase separation and cross-sea
phase. Between microphase separation and the cross-sea phase, $\langle
C_{2}\rangle_{\rm cross}>\langle C_{2}\rangle_{\rm micro}$. Although a change
in $\langle C_{2}\rangle$ is not very significant in these two phases, it is
still a better candidate for distinguishing the cross-sea from the microphase
separation than the traditional polar order parameter. In Fig. 5, we plot the
$\eta-q$ phase diagram by computing $\langle C_{2}\rangle$ for the six
different phases where in the coexistence region, $\langle C_{2}\rangle_{\rm
macro}>\langle C_{2}\rangle_{\rm cross}>\langle C_{2}\rangle_{\rm micro}$.
This phase diagram complements Fig. 3, which has been constructed using the
density field and depicts the nature of phase separation in the DVM as a
function of $q$.

Now, the nature of the coexistence region for any finite, large $q$ has been a
subject of discussion in the context of $q$-state ACM Solon  _et al._ (2022);
Chatterjee  _et al._ (2022). It was argued in Ref. Solon  _et al._ (2022) that
spatial anisotropy plays a crucial role in determining the macro/micro-phase
separation of the coexistence region in the ACM and one should observe a
macrophase separation of the coexistence region for a finite $q$ beyond a
characteristic length scale which diverges for large $q$. ACM with a different
set of dynamical rules than Ref. Solon  _et al._ (2022) was studied in Ref.
Chatterjee  _et al._ (2022) where the flocking transition in the ACM was
argued to be equivalent to the VM at large $q$. The $q$-state DVM is governed
by a completely different set of microscopic rules than both the ACM Solon
_et al._ (2022); Chatterjee  _et al._ (2022) and we will therefore investigate
next the impact of microscopic rules on the DVM steady-state as a function of
$q$. We plan to do this through the analysis of number fluctuations, the
pinning-unpinning property of the system’s global order, and the structure
factor manifesting the correlation of polarization.

Figure 5: (Color online) $\eta-q$ phase diagram of the DVM by computing
$\langle C_{2}\rangle$. The colorbar represents the range of $\langle
C_{2}\rangle$ values for different phases. Parameters: $L=1024$, $v_{0}=0.5$,
and $\rho_{0}=2$.

Number fluctuations.—It is known that fluctuations play an essential role in
selecting the phase-separated patterns in flocking models. In the AIM (Solon
and Tailleur, 2015), the number fluctuation $\Delta n^{2}=\langle
n^{2}\rangle-\langle n\rangle^{2}$ in the liquid phase was found to be normal
($\Delta n^{2}=n$) and the coexistence region shows a macrophase separation.
In the VM Solon  _et al._ (2015a), on the contrary, giant density fluctuations
($\Delta n^{2}=n^{1.6}$) are observed which break large liquid domains and
prevent the bands from coarsening further. This results in a microphase-
separated coexistence region. In the $q$-state ACM Solon  _et al._ (2022);
Chatterjee  _et al._ (2022), which can be thought of as a bridge between the
discrete AIM and continuous VM, density fluctuations show a transition from
normal to giant fluctuation as $q$ increases and can be explained as a
transition from AIM physics to VM physics.

In Fig. 6, therefore, we show the number fluctuations $\Delta n^{2}$ versus
average particle number $\langle n\rangle$ computed in the ordered liquid
phase of the DVM for various $q$. $n(\ell)$ is the number of particles in
boxes of different sizes $\ell$ included in a $300^{2}$ domain (with
$\ell\leqslant 150$), with $\langle n\rangle=\rho_{0}\ell^{2}$. As shown in
Table 1, the number fluctuation behaves like $\langle n\rangle^{\xi}$ with the
fluctuation exponent $\xi$ increasing with $q$, from $\xi\simeq 1.22$ for
$q=4$ to $\xi\simeq 1.64$ for large $q$. This transition from normal
fluctuations for small $q$ to giant fluctuations for larger $q$ was also
observed in the ACM Solon  _et al._ (2022); Chatterjee  _et al._ (2022)
although $\xi$ for $q=4$ and 5 are moderately larger in the DVM than the ACM.
The increase of the density fluctuations with $q$ can be attributed to the
fact that for large $q$, particles have more rotational degrees of freedom due
to the weak anisotropy and therefore more directional freedom to propel. The
existence of giant number fluctuations (GNF) and its connection with the
microphase separation in the VM was hypothesized in Ref. Solon  _et al._
(2015a). It was argued that GNF ($\xi\simeq 1.6$) breaks bulk liquid domains
and produces a smectic-like microphase separation in the coexistence regime
whereas the system stabilizes in the bulk phase when the density fluctuations
are normal ($\xi\simeq 1$) Solon and Tailleur (2015). Exploiting the same
logic for ACM Chatterjee  _et al._ (2022), GNF was argued to be responsible
for the microphase separation in the coexistence regime for $q\geqslant 8$,
although a direct relation between the existence of GNF in the ordered phase
and the microphase separation in the coexistence phase is still ambiguous.
Nonetheless, if we compare $\xi$ in Table 1 and the snapshots in Fig. 3, we
observe a correspondence between the fluctuation exponent and the pattern
formation in the coexistence region of DVM.

One should also consider the finite size effect on the fluctuation exponents
as discussed in Ref. Chatterjee  _et al._ (2022). In Fig. 6, the data can be
fitted to two different power-law regimes (the extracted exponents depend on
the interval along the $x$-axis to which the fits are restricted): (a) $\xi$
tabulated in Table 1 in the interval $[10^{3},5\times 10^{4}]$ and (b) a
smaller $\xi$ in the interval $[5\times 10^{4},10^{5}]$. Around the second
interval, the plot shows a “saturation” because $\xi$ must decrease with
increasing $\langle n\rangle$ due to the finite-size cut-off at $\langle
n\rangle=N=\rho_{0}L^{2}$, where $\Delta n^{2}$ vanishes.

$q$ | 4 | 5 | 6 | 7 | 8 | 9 | 10 | 16
---|---|---|---|---|---|---|---|---
$\xi$ | 1.22 | 1.28 | 1.41 | 1.48 | 1.67 | 1.64 | 1.64 | 1.64
Table 1: (Color online) Number fluctuation exponents $\xi$ for several values
of $q$, reported from Fig. 6. The typical error on the fluctuation exponents
is 0.03.  Figure 6: (Color online) Number fluctuations $\Delta n^{2}=\langle
n^{2}\rangle-\langle n\rangle^{2}$ versus average particle number $\langle
n\rangle$ for several $q$ values in a $300^{2}$ domain. Parameters:
$\eta=0.3$, $v_{0}=0.5$, and $\rho_{0}=6$.

Pinned property of the order parameter.—The global order parameter defined in
Eq. (6) quantifies the overall ordering of the particles in the system. In the
ACM (Solon  _et al._ , 2022), at finite size, the direction of the global
polar order $\Phi\equiv\rm{arg}\langle{\bf m}\rangle$ exhibits distinct
behaviors in the liquid phase depending on the value of $q$: it displays AIM-
like properties (pinned along an angle) for small $q$ and VM-like behavior
(unpinned over time) for large $q$.

In Fig. 7(a), we show the time evolution of $\Phi(t)$ in the DVM liquid phase
for varying $q$. Akin to the observation made in the ACM Solon  _et al._
(2022), $\Phi(t)$ starts wandering slowly with $q$ and becomes an unpinned
variable of $t$ for large $q$. In other words, for weak anisotropy, the global
ordering does not remain constrained to a specific orientation. While, for
small $q$, $\Phi(t)$ remains pinned and tends to exhibit a stable global
ordering. Microscopically, this refers to a picture that at large $q$, a
proportional number of degrees of freedom allows the particles to choose
between adjacent directions facilitated by fluctuation while it is not the
obvious choice for particles at small $q$ requiring a significantly larger
jump to switch directions. Translating this to the global polar order
parameter and comparing Fig. 7(a) with Fig. 6 we propose that GNF corresponds
to the unpinned behavior of $\Phi(t)$. Likewise, the unpinning nature of the
direction of the global polar order is a characteristic of microphase
separation.

In addition to $q$, the finite system size also affects the evolution of
$\Phi(t)$, which is shown in Fig. 7(b) for $q=9$. Similar to the ACM Solon
_et al._ (2022), we observe a transition from unpinned behavior to pinned
behavior in $\Phi(t)$ as the system size is increased. In larger systems with
polar order, a particle interacts with more particles in the neighborhood and
correlates over longer distances. Higher connectivity promotes stronger
alignment and cooperative motion among the particles. As a result, the
direction of the global order becomes more pronounced and persistent in larger
systems, leading to the pinned state. Fig. 7(b) further signifies that if
$\Phi(t)$ is pinned for $L=300$, it must also be pinned for $L=1024$. This is
inconsistent with the microphase separation and cross-sea patterns observed in
the coexistence region of the $q=9$ DVM (Fig. 3). Evidently, the correlation
proposed earlier between the pinned property of the system’s ordered liquid
phase and system morphology observed in the coexistence region
(macro/micro/cross-sea) is not conclusive in the DVM (see Appendix E for more
details).

Figure 7: (Color online) (a) Time series of the orientation of global polar
order $\Phi$ in the liquid phase showing a transition from unpinned to pinned
as $q$ increases ($\eta=0.3$, $\rho_{0}=6$, $v_{0}=0.5$, and $L=200$). (b)
$\Phi$ shows an unpinning to pinning transition as a function of system size
$L$ for $q=9$ ($\eta=0.3$, $v_{0}=0.5$ and $\rho_{0}=6$).

Structure factor.— To explore the correlation among the particle
polarizations, we consider the transverse magnetization structure factor
$S_{\perp}({\bf k})=\langle m_{\perp}({\bf k})m_{\perp}(\bf{-k})\rangle$
against wavelength ${\bf k}$ and plot it in Fig. 8. The structure factor has
been calculated in the liquid phase on a $300^{2}$ domain for (a) various $q$
values (the same behavior is observed for the structure factor of the density
field) and (b) for a fixed $q=9$ but for various system sizes. The results
presented in Fig. 8(a) show that for small $q$, the structure factor $S(k)$
converges to finite values as the wave vector ${\bf k}\to 0$. This convergence
indicates an AIM-like behavior or a macrophase separation of the coexistence
region Solon  _et al._ (2022). However, one can notice that this convergence
is achieved only beyond certain length scales and these length scales are
functions of $q$. The structure factor captures the correlations between
particle orientations and therefore the magnitude of ordering at different
length scales. A saturation of $S(k)$ for small $q$ when ${\bf k}$ approaches
zero thus signifies a strongly correlated liquid domain whereas for large $q$,
due to weak anisotropy or more allowed orientations for ordering, particles
inside the liquid domain are not as strongly correlated as for small $q$. Fig.
8(b) shows $S(k)$ for several system sizes and manifests that with larger
system sizes $(L\geq 300)$, $S(k)$ tends to converge to a finite value when
${\bf k}\to 0$. Our earlier argument that stronger interactions between
particles (with larger $L$) promote robust ordering also applies here.
Comparing Fig. 8 with Fig. 7, we conclude that the pinning behavior (unpinning
behavior) of $\Phi(t)$ and the saturation of $S(k)$ (algebraic scaling of
$S(k)$) compliment each other and convey the same physics.

It was argued in Ref. (Solon  _et al._ , 2022) considering the pinning
properties of the order parameter and behavior of the structure factor in the
liquid phase that for the large $q$ ACM, VM behavior (microphase separation)
will be observed only up to large finite sizes. But the asymptotic large
length scale behavior will be AIM-like (macrophase separation) where the
length scale diverges with $q$ as $\exp(q^{2})$. Ref. (Solon  _et al._ , 2022)
also showed that in the phase-coexistence region of the ACM, a microphase to
macrophase transition happens upon increasing the linear system size along the
transverse direction at fixed $q$ which we do not observe (in the DVM,
multiple bands do not merge to a single band upon increasing $L_{y}$ for a
fixed $L_{x}$). We observe a similar behavior of $\Phi(t)$ and $S(k)$ as Ref.
(Solon  _et al._ , 2022) but can not draw a conclusive correspondence between
the large length-scale liquid phase behavior of $\Phi(t)$ and $S(k)$ to the
phase-coexistence behavior of the DVM as shown in Fig. 3.

Figure 8: (Color online) Structure factor $S_{\perp}({\bf k})$ vs ${\bf
k}=(k_{\parallel},0)$ in the ordered liquid phase for (a) several $q$ values
($L=300$) and (b) varying system size ($q=9$). Parameters: $\eta=0.3$,
$\rho_{0}=6$, and $v_{0}=0.5$.

In Fig. 3, the snapshots on a large square domain (without any spatial
anisotropy) show the existence of a microphase separation and cross-sea
patterns (for which one needs at least two bands) of the phase-coexistence
region for large $q$ values. The number fluctuation plotted in Fig. 6
corroborates this observation by exhibiting GNF for those large $q$ values.
The large length scale asymptotic behavior (for large $q$) of the direction of
global order $\Phi(t)$ (Fig. 7) and the structure factor $S(k)$ (Fig. 8) in
the ordered liquid phase respectively shows a pinned behavior and saturation
for $k_{\parallel}=0$ which as argued in Ref. (Solon  _et al._ , 2022)
signifies an AIM phenomenology. However, our numerical investigation of the
DVM does not show a cross-over from micro- to macrophase separation for higher
$q$-values. Therefore, we argue that the impact of dynamical rules governing
flipping and hopping in a flocking model has a significant influence over the
system dynamics.

Figure 9: (Color online) Order parameter distributions of the $q=9$ DVM for
the liquid phase. Parameters: $\eta=0.3$, $v_{0}=0.5$ and $\rho_{0}=6$. Ring-
like distributions for smaller system sizes ($L=50$ and $L=100$) in (a) and
(b) change to distinct isolated spots for larger system sizes ($L=200$ and
$L=300$) in (c) and (d), which correspond to the ordered liquid phase’s
$q$-fold degeneracy.

Order parameter distribution.—In Vicsek-like models, where particles are
active and move with a constant velocity, the ordered state exhibits a true
long-range order (LRO) in two dimensions because of a spontaneous symmetry
breaking due to the out-of-equilibrium nature of the models. In Fig. 9, we
show the time- and ensemble-averaged distribution of the order parameter
$\textbf{m}=(m_{x},m_{y})$ for increasing system sizes, where
$m_{x}=\frac{1}{N}\sum_{i=1}^{N}\cos\theta_{i}$ and
$m_{y}=\frac{1}{N}\sum_{i=1}^{N}\sin\theta_{i}$. In Fig. 9(a) and Fig. 9(b),
ring-like distributions (unpinned orientations) characteristic of the quasi
long–range ordered (QLRO) phase is observed. But, this is only due to the
finite-size effect and similar to the impact of finite-size on $\Phi(t)$ and
$S(k)$. We recover the LRO for larger system sizes ($L=200$ and $L=300$) where
the distributions display nine distinct isolated spots (pinned orientations)
that correspond to the 9-fold degeneracy of the ordered liquid phase, each
spot having equal probability. One can expect that the finite size effect will
be much weaker for $L=1024$ and the spread of the distribution in the LRO
phase around the allowed ordering angles will also be more precise.

The DVM for large $q$ exhibits a QLRO phase when $v_{0}=0$ (see Appendix F).
We argue that DVM with immobile particles reduces to the two-dimensional
$q$-state clock model (with a quenched bond-disorder as only particles within
a fixed distance interact) which approaches the XY model for large $q$ with
vanishing LRO regime Chatterjee  _et al._ (2018). For $v_{0}>0$ and a fixed
$L$, as flocking directions increase with $q$, we again observe ring-like
distributions for large $q$ (see Appendix F) but beyond a length scale which
is proportional to $q$, the order parameter distributions for large $q$ show
$q$ isolated spots characteristic of the LRO phase. This is similar to the
unpinned to the pinned transition of $\Phi(t)$ and convergence of $S(k)$ to a
finite value at ${\bf k}\to 0$ for large $q$ values as the system size is
increased. In the VM ($q\to\infty$), even for a large $L$, the order parameter
distribution shows a ring-like structure because of the continuous symmetry.

Stability of the ordered liquid phase.—

Figure 10: (Color online) Time evolution snapshots of the orientation field
showing the reversal of the initial ordered liquid phase in (a) $q=6$ and (b)
$q=16$-state DVM following the introduction of a tiny counter-propagating
high-density droplet (motion of direction is denoted by black arrows, $t=1$)
of radius $r_{d}=5$ and density $\rho_{0}^{d}=10\rho_{0}$. Colorbar represents
particle orientation field. Parameters: $L=1024$, $\eta=0.4$, $v_{0}=1$ and
$\rho_{0}=2$.

As discussed in the context of Fig. 2, here we present a brief analysis of the
stability of the DVM ordered liquid phase by inserting a small high-density
counter-propagating liquid droplet. Initially, the average direction of the
particles in the polar liquid phase is aligned along the direction $\Phi_{\rm
liq}=\pi$, and the average particle orientation of the liquid droplet is
$\Phi_{d}=0$. The radius of the droplet is $r_{d}$ and density is
$\rho_{0}^{d}$ and it is inserted in an ordered phase of density $\rho_{0}$
($\rho_{0}^{d}\gg\rho_{0}$). We take several $q$ values and calculate the
probability ($P_{\rm rev}$) that the droplet grows against the main order and
reverses the ordered phase as a function of $\eta$, $r_{d}$, and
$\rho_{0}^{d}$. The simulation protocol follows Ref. Benvegnen  _et al._
(2023), where the initial configuration is prepared by particles with
$\theta=\pi$ and we let the system evolve up to a certain time $t$ to reach
the steady state. The system retains the average global polarization in the
direction $\theta=\pi$. Then, a circular region of radius $r_{d}$ centered at
$(L/2,L/2)$ is selected and an additional $\Delta N=(\rho_{0}^{d}-\rho_{0})\pi
r_{d}^{2}$ number of particles are added to make a high-density circular
droplet. Finally, the orientation of all the particles within the droplet is
changed to $\theta=0$.

In Fig. 10, we study the fate of polar flocks in $q=6$ and 16-state DVM by
introducing a small high-density counter-propagating droplet against the
initial polar ordered liquid phase of the main flow and observe the subsequent
time evolution. One should note here that the perturbation through the droplet
is very small i.e. the ratio of droplet diameter to the linear length of the
simulation box is $\sim 10^{-2}$ ($r_{d}=5$, $L=1024$). We observe, similar to
Ref. Codina  _et al._ (2022), that the droplet grows with time leaving behind
a dilute region ($t=10^{2}$) and adds more and more particles as it moves
along forming a principal dense, curved band (followed by several other curved
bands) that invades the whole system ballistically ($t=10^{3}$ and
$t=10^{4}$). In the final stage, this principal curved band connects itself
over the system boundaries and widens until a steady state liquid phase of a
different $\Phi_{\rm liq}$ emerges (at $t=10^{5}$), signifying the
metastability of the DVM liquid phase. The time-evolution is similar for both
small and large $q$, which signifies that both discrete and continuous-
symmetry flocks are metastable.

The growth pattern of the DVM droplet is similar to the VM Codina  _et al._
(2022) but distinct from the AIM Benvegnen  _et al._ (2023). In the DVM, after
its introduction, the droplet front interacts with the liquid particles
outside and creates a curved band of particles having several different
orientations (impact of $q$). In this process, the droplet seizes to exist and
it is this high-density curved band that destroys the initial flow. In AIM
Benvegnen  _et al._ (2023), the droplet grows along the direction transverse
to the propagation (due to the constant transverse diffusion) creating a
comet-like trail of disordered particles that can not be observed in the DVM.

Fig. 11 quantifies $P_{\rm rev}$, the probability of reversing the main flow
upon the introduction of a given droplet, for several control parameters. For
each set of control parameters, we have taken 20 independent realizations to
calculate $P_{\rm rev}$. Akin to Ref. Codina  _et al._ (2022), we observe that
the noise strength $\eta$ has a strong influence on the reversal dynamics and
$P_{\rm rev}$ increases from 0 to 1 as $\eta$ is increased [Fig. 11(a)]. This
is because for small $\eta$, the ordered phase is very stable, and thus, the
counter-propagating dense bands find it difficult to reverse its flow. For
large $\eta$, fluctuations are stronger, and therefore, the probability of
reversal increases. Fig. 11(a) also exhibits that the transitional value of
$\eta$ ($P_{\rm rev}=\frac{1}{2}$) above which a droplet triggers a reversal
is a decreasing function of $r_{d}$ although there is a critical $\eta$
($\eta\sim 0.35$), below which no droplet can trigger a reversal irrespective
of its density. The $r_{d}-\rho_{0}^{d}$ phase diagram in Fig. 11(b) has been
constructed by calculating $P_{\rm rev}$ for several $(r_{d},\rho_{0}^{d})$.
Unsurprisingly, large $r_{d}$ combined with large $\rho_{0}^{d}$ are shown to
facilitate the reversal of the initial liquid phase. The droplet-induced
reversal of the liquid flow is found independent of $q$ where $P_{\rm rev}$ is
found to behave similarly for each $q$ under certain values of $\eta$ [Fig.
11(c)].

Figure 11: (Color online) (a) Reversal probability $P_{\rm rev}$ versus noise
$\eta$ for $q=16$ and two different droplet radius $r_{d}=5$ and 10. (b)
$r_{d}-\rho_{0}^{d}/\rho_{0}$ phase diagram for $\eta=0.45$ and $q=16$. The
colorbar represents reversal probability $P_{\rm rev}$. (c) $P_{\rm rev}$
versus $q$ for various noise, $r_{d}=5$. Parameters: $L=400$,
$\rho_{0}^{d}=10\rho_{0}$ (a & c), $v_{0}=1$ and $\rho_{0}=2$.

In addition, a study of the metastability of the DVM liquid phase on
rectangular geometry ($L_{x}=800,L_{y}=100$), produces a similar outcome. One
can also consider droplet movement in other directions than opposite to the
global flow, e.g. transversely propagating droplets with $\Phi_{d}=\pi/2$ for
$q=8$ or $\Phi_{d}=2\pi/3$ for $q=6$ and perform a similar analysis to check
whether the liquid phase is susceptible to droplets irrespective of their
propagation direction.

##  IV Discussion

Our study motivated by the active clock model (Chatterjee  _et al._ , 2022;
Solon  _et al._ , 2022), considers a true $q$-state discrete version of the
Vicsek model where $q$ defines the strength of orientation anisotropy. At
small $q$, the system is highly anisotropic which however vanishes in the
limit $q\rightarrow\infty$ when we recover the Vicsek model. The DVM shows
qualitatively similar features as the ACM Chatterjee  _et al._ (2022) for
intermediate noise strength $\eta$ where a transition from macrophase to
microphase separation is observed in the coexistence region as $q$ is
increased. But for small $q$ and $\eta$, the liquid phase appearing in the ACM
at low temperatures is replaced in the DVM by a cluster phase. The cluster
phase consists of multiple clusters with different polarization (see Fig. 12)
which does not exhibit a long-range order. The clusters grow and merge with
increasing $q$ leading to a homogeneous ordered phase at large $q$. For small
$q$, a long-range ordered phase can be achieved by increasing the noise
strength. For low noise and small $q$, the probability of transverse flipping
is very small as fluctuations are weak. In addition to that, transverse
fluctuations through hopping are also absent. Consequently, the combined
influence of these factors results in clusters failing to grow continuously
for low $\eta$ and small $q$, preventing the system from reaching a
homogeneous liquid state.

The self-organized patterns in the coexistence region of the discretized VM
indicate a transition from AIM-like patterns to VM-like patterns as anisotropy
becomes weaker. This observation is corroborated by the giant density
fluctuations for large $q$. However, the large length scale behavior of the
direction of global order $\Phi(t)$, the structure factor $S(k)$, and the
order parameter distribution in the liquid phase do not correspond with the
phase-coexistence patterns of the large $q$ DVM. It shows microphase or cross-
sea pattern for large length scales without any spatial anisotropy as $q$
increases.

We also find that the DVM liquid phase is susceptible to perturbation applied
through a counter-propagating droplet. The liquid phase reorients and
propagates along the direction of the droplet. The reversal dynamics is
significantly impacted by the noise strength $\eta$ Codina  _et al._ (2022)
but remains independent of $q$. The stability of the high-density flocking
ordered phase at low noise is still an open problem and will be addressed in a
subsequent study Chatterjee  _et al._ (unpublished).

As a final remark, we add that the rotational flexibility of the particles and
microscopic details of the dynamical rules can significantly impact the
macroscopic properties of the ordered phase. It would be interesting to
compare the model predictions of the DVM with suitable experiments where
anisotropy in the particle orientation may be controlled by a finite number of
motility directions.

##  V Acknowledgments

MK acknowledges financial support in the form of a research fellowship from
CSIR, Govt. of India (Award Number: 09/080(1106)/2019-EMR-I). RP, MK thanks
the Indian Association for the Cultivation of Science (IACS) for the
computational facility. SC and HR are financially supported by the German
Research Foundation (DFG) within the Collaborative Research Center SFB 1027.
SC wants to thank Dr. Matthieu Mangeat for many valuable discussions.

##  Appendix A Dependency of the cluster phase on the initial condition

Figure 12: (Color online) Steady-state snapshots (b-c) and (e-f) from two
different initial conditions: (a) random disordered and (d) polar ordered. (b,
e) $\eta=0.1$ and (c, f) $\eta=0.2$. The colorbar represents particle
orientations. Parameters: $q=4$, $L=300$, $v_{0}=0.5$, and $\rho_{0}=6$.

Snapshots in Fig. 12(b-c) and Fig. 12(e-f) illustrate the steady state
behavior of the $q=4$ DVM starting from two different initial conditions: (a)
random disordered and (d) polar ordered. In the top panel, the system exhibits
cluster phase for $q=4$ for both values of $\eta$ when starting from an
unbiased random configuration. The initial coarsening process forms clusters,
but they do not merge to create a single large ordered domain due to the
absence of transverse fluctuations as discussed in the context of Fig. 3. In
the bottom panel, this observation changes for the lowest noise strength when
starting from a polar ordered initial configuration. For $\eta=0.1$, the
steady state remains in an ordered liquid phase (similar to the steady state
behavior of the 4-state APM or ACM at low temperature) signifying that the
fluctuation is weak to alter the broken symmetry phase into a cluster phase.
However, with an increase in the noise ($\eta=0.2$), the steady-state cluster
phase appears again. This also suggests that for DVM with weak anisotropy and
fluctuations, the cluster phase is the non-equilibrium steady-state and the
well-known small $\eta$ or large $\beta$ ordered liquid phase can only be
achieved by taking a strongly polarized ordered initial condition. It is worth
noting that the number of clusters in (c) is higher than in (b) due to more
relaxation via the noise. A further increase of the noise will lead to a
single large ordered domain as shown in the phase diagram of Fig. 15.

##  Appendix B Cluster size analysis for $q=4$

Figure 13: (Color online) Probability distribution $P(S)$ of cluster size $S$
for different noise level. Parameters: $q=4$,
$L_{x}=800,L_{y}=100,\rho_{0}=2$, and $v_{0}=0.5$.

The cluster size analysis for $q=4$ is shown in Fig. 13 for a rectangular
domain of size $800\times 100$. We use a box-counting method to find a cluster
which is described below:

If $L_{x}$ and $L_{y}$ are respectively the linear sizes in a rectangular
domain, we consider $L_{x}\times L_{y}$ as the total number of boxes. For a
box $i$, $n_{i}$ is the number of particles in that box and $c_{i}$ is the
cluster label of the box. If the box does not contain any particles $c_{i}=0$
i.e. it is not part of any cluster. Then we use the Depth-first search (DFS)
algorithm to find the connected boxes which are not void of particles and
label them as a single cluster. For each cluster label $c_{j}$, we calculate
the size of the cluster $S_{c_{i}}$ as following:

$S_{c_{i}}=\sum_{j=1}^{L_{x}L_{y}}\delta_{c_{i},c_{j}}\times n_{j}\,,$ (8)

where $\delta_{c_{i},c_{j}}$ is the Kronecker delta function that equals 1 if
$c_{i}=c_{j}$, and 0 otherwise. Then we calculate the cluster size probability
distribution denoted as $P(S)$:

$P(S)=\frac{\text{Number of clusters with size }S}{\text{Total number of
clusters}}\,.$ (9)

Fig. 13 illustrates that for small noise, the probability of larger cluster
formation is high. This is because reduced fluctuation in the system
facilitates the formation of high-density clusters. However, with noise
($\eta=0.3$), the cluster phase vanishes and the $q=4$ DVM exhibits a
macrophase separation, resulting in a decrease in the probability of obtaining
large clusters. At large noise ($\eta\geqslant 0.4$), the system becomes
disordered, leading to a lesser probability of formation of large clusters. We
would also like to mention that for a fixed noise, $P(S)$ versus $S$ for
various orientations ($\theta=0,\pi/2,\pi,3\pi/2$) shows almost identical
distributions signifying no preference in orientation in the cluster
formation.

##  Appendix C Impact of spatial anisotropy: steady-states for rectangular
domain

Figure 14: (Color online) $\eta-q$ phase diagram of the DVM illustrated by
snapshots on a rectangular domain ($L_{x}=800$, $L_{y}=100$). Parameters:
$\rho_{0}=2$, $v_{0}=0.5$. As a function of $\eta$ and $q$, we observe five
distinct self-organized patterns: cluster ($\eta=0.1$, $q=4$), macrophase
separation ($\eta=0.3$, $q=4$), microphase separation ($\eta=0.4$, $q=7\to
16$), ordered liquid ($\eta=0.2$, $q=7\to 16$), and disordered gas
($\eta=0.5$, $q=4\to 16$).

In Fig. 14, we investigate how spatial anisotropy influences the non-
equilibrium steady-state behavior of the DVM if we switch from a square domain
to a rectangular domain by analyzing late-stage representative snapshots as a
function of noise strength ($\eta$) and discretization parameter ($q$). When
noise is low and $q$ is small, we observe the emergence of a locally ordered
cluster phase similar to our finding for the square domain. However, when $q$
is small and fluctuations are pronounced, those cluster phases relax and
transform into a larger, organized domain (see the snapshot for $q=4$ and
$\eta=0.3$) akin to the Fig. 3. Conversely, when $\eta$ is small and $q$ is
large, weak anisotropy facilitates the merger of the cluster phase into a
larger, well-organized domain (see the snapshot for $q\geqslant 8$ and
$\eta=0.1$) which is observed at $q=10$ in the absence of spatial anisotropy
(see Fig. 3). An increase in the fluctuation for small $q$ might exhibit
multiple bands but those are connected by the periodic boundary conditions and
should be considered a single band (see the snapshot for $q=5$ and
$\eta=0.2,0.3$). However, no cross-sea phase is observed with the rectangular
geometry, which is seen in Ref Kürsten and Ihle (2020) probably due to larger
particle velocity ($v_{0}=1$). An increase in anisotropy $q$ also increases
the no. of bands in the coexistence region at intermediate noise ($\eta=0.4$)
because with $q$, the density fluctuation increases along with the
magnetization fluctuation which prompts the breaking of large domains.

##  Appendix D $(\eta-\rho_{0})$ phase diagram of the DVM

Figure 15: (Color online) Phase diagrams on $\eta-\rho_{0}$ plane computed on
a rectangular domain of $800\times 100$ for various $q$. A macrophase
separation similar to AIM in (a) and (b) transforms to a Vicsek-like phase
separation (microphase) in the coexistence region in (c) and (d) at
intermediate densities. $v_{0}=0.5$.

In Fig. 15, we plot the noise-density ($\eta-\rho_{0}$) phase diagrams
computed on a rectangular domain of $800\times 100$ for various $q$. The
clustering phase is very prominent for small $q$ values ($q=4$, 5) and exists
for high densities. As $q$ is increased, the conventional ordered liquid phase
appears at high densities and for $q\geqslant 8$, the cluster phase disappears
and we notice the emergence of the typical $\eta-\rho_{0}$ diagram observed
for Vicsek-like systems [Fig. 15(d)]. These Vicsek-like phases are
characterized by a more global alignment of the particles, leading to coherent
motion and the absence of distinct clusters or bands. In this regime, the
behavior of the system is predominantly governed by the alignment interactions
between the particles rather than the specific value of $q$. It shows clear
boundaries between different phases based on varying values of $q$ and $\eta$.
At small $q$ and intermediate densities, a macrophase separation similar to
AIM is observed. However, as $q$ exceeds a threshold (e.g., $q\geqslant 8$),
the system transitions to a Vicsek-like phase separation in the coexistence
region.

##  Appendix E Direction of global order in the coexistence region

Figure 16: (Color online) Time evolution of the direction of global order
$\Phi$ in the coexistence region for (a) $L=512$ and various $q$ and (b) $q=9$
and several values of $L$. (c–d) Snapshots showing multiple bands moving
parallelly constituting a microphase-separated coexistence region for $q=9$
and $q=16$, respectively. System size $L=512$. Parameters: $\eta=0.3$,
$\rho_{0}=1.5$, $v_{0}=0.5$.

Here, we analyze the pinned property of the system in the coexistence region
and directly compare the outcome with the steady-state snapshots. In Fig. 16,
the time series of the orientation of global order $\Phi$ is shown for fixed
$L=512$ and various $q$ [Fig. 16(a)] and for fixed $q=9$ and several $L$ [Fig.
16(b)]. When the system size is fixed, we observe pinning to unpinning
transition with $q$ whereas for fixed $q$, the reverse transition happens as
the system size increases. Both observations are similar to the observations
made regarding the time evolution of $\Phi$ in the DVM polar ordered phase
[Fig. 7]. The snapshots shown in Fig. 16(c–d) without any spatial anisotropy
exhibit a coexistence region with multiple bands signifying a microphase-
separated region. Comparing the snapshots with the time evolution of $\Phi$ we
notice that both the unpinning behavior of the DVM for $q=16$ and the pinning
behavior for $q=9$ shows flocking with multiple parallelly moving bands.

##  Appendix F Nature of the DVM ordered phase

Figure 17: (Color online) Zero-activity limit of the order parameter
distributions of the DVM with varying $q$ in the liquid phase. Parameters:
$L=100,\eta=0.3$, and $\rho_{0}=6$. Distinct isolated spots in (a), (b), and
(c) indicate LRO while ring-like distributions in (d), (e), and (f) are
characteristic of the QLRO phase. Figure 18: (Color online) Order parameter
distributions of the DVM at constant activity $v_{0}=0.5$ with varying $q$ in
the liquid phase. Parameters: $L=100,\eta=0.3$, and $\rho_{0}=6$.

In Fig. 17, we show the order parameter (m) distribution on the
$(m_{x},m_{y})$ plane for $v_{0}=0$. At this limit, particles can no longer
move but modify their orientation according to Eq. (3) and the $q$-state DVM
reduces to the two-dimensional $q$-state clock model which shows two distinct
phase transitions, one from disordered to QLRO phase at a higher temperature
and the other from QLRO to LRO phase at a lower temperature for $q\geqslant 5$
Chatterjee  _et al._ (2018). At large $q$, the LRO phase gradually starts to
vanish, and for the two-dimensional XY model ($q\to\infty$), only one phase
transition occurs (Kosterlitz-Thouless phase transition) from the disordered
phase to the QLRO phase. As shown in Fig. 17, at the zero activity limit, the
order parameter distribution shows $q$ distinct isolated spots (signifying
LRO) for small $q$ but as $q$ increases ($q\geqslant 7$), ring-like
distributions characteristic of the QLRO phase appears. These ring-like
distributions become more pronounced as the system size is increased for large
$q$. In contrast, for $v_{0}>0$ (see Fig. 18), the order parameter
distribution exhibits LRO through isolated points of phase ordering as
activity facilitates the broken symmetry phase. For $q>7$, a comprehensible
LRO phase is observed at a large length scale limit as shown in Fig. 9.

Figure 19: (Color online) Order parameter $\langle m\rangle$ versus system
size $L$ for varying $q$ in the ordered liquid phase. (a) $\langle m\rangle$
is independent of $L$ for $v_{0}=0.5$ signifying LRO. (b) $\langle m\rangle$
decays algebraically with $L$ signifying QLRO in the zero-activity limit.
Parameters: $\eta=0.3$, $\rho_{0}=6$.

To understand better the nature of ordering of the DVM liquid phase, we show
the order parameter $\langle m\rangle$ against increasing system size $L$ for
several $q$ in Fig. 19. The data presented are averaged over time and several
initial configurations. We note that, $\langle m\rangle$ remains independent
of the system size $L$ (actually, $m$ scales with $L$ as $m\sim L^{-\lambda}$,
decays much slower than a power law) for all $q$ for $v_{0}=0.5$ [Fig. 19(a)]
signifying LRO. As a result, the liquid phase of the constant-speed DVM is LRO
and the direction of the order parameter exhibits a pinned behavior. For
$v_{0}=0$, shown in Fig. 19(b), however, $\langle m\rangle$ is expected to
algebraically decay to zero for $L\to\infty$ and this effect is more
pronounced for larger $q$ because, for large $q$, the $v_{0}=0$ DVM approaches
the two-dimensional XY model where the low-temperature phase is QLRO.

## References

  * Ramaswamy (2010) S. Ramaswamy, “The mechanics and statistics of active matter,” Annual Review of Condensed Matter Physics 1, 323 (2010).
  * Marchetti  _et al._ (2013) M. C. Marchetti, J. F. Joanny, S. Ramaswamy, T. B. Liverpool, J. Prost, M. Rao, and R. A. Simha, “Hydrodynamics of soft active matter,” Reviews of Modern Physics 85, 1143 (2013).
  * Gompper  _et al._ (2020) G. Gompper, R. G Winkler, T. Speck, A. Solon, C. Nardini, F. Peruani, H. Löwen, R. Golestanian, U. B. Kaupp, L. Alvarez, _et al._ , “The 2020 motile active matter roadmap,” Journal of Physics: Condensed Matter 32, 193001 (2020).
  * Needleman and Dogic (2017) D. Needleman and Z. Dogic, “Active matter at the interface between materials science and cell biology,” Nature Reviews Materials 2, 1 (2017).
  * Vernerey  _et al._ (2019) F. J. Vernerey, E. Benet, L. Blue, A. K. Fajrial, S. L. Sridhar, J. S. Lum, G. Shakya, K. H. Song, A. N. Thomas, and M. A. Borden, “Biological active matter aggregates: Inspiration for smart colloidal materials,” Advances in Colloid and Interface Science 263, 38 (2019).
  * Ghosh  _et al._ (2021) S. Ghosh, A. Somasundar, and A. Sen, “Enzymes as active matter,” Annual Review of Condensed Matter Physics 12, 177 (2021).
  * Bär  _et al._ (2020) M. Bär, R. Großmann, S. Heidenreich, and F. Peruani, “Self-propelled rods: Insights and perspectives for active matter,” Annual Review of Condensed Matter Physics 11, 441 (2020).
  * Shaebani  _et al._ (2020) M. R. Shaebani, A. Wysocki, R. G. Winkler, G. Gompper, and H. Rieger, “Computational models for active matter,” Nature Reviews Physics 2, 181 (2020).
  * Sanchez  _et al._ (2012) T. Sanchez, D. T. N. Chen, S. J. DeCamp, M. Heymann, and Z. Dogic, “Spontaneous motion in hierarchically assembled active matter,” Nature 491, 431 (2012).
  * Peruani  _et al._ (2012) F. Peruani, J. Starruß, V. Jakovljevic, L. Søgaard-Andersen, A. Deutsch, and M. Bär, “Collective motion and nonequilibrium cluster formation in colonies of gliding bacteria,” Physical Review Letters 108, 098102 (2012).
  * Garcimartín  _et al._ (2015) A. Garcimartín, J. M. Pastor, L. M. Ferrer, J. J. Ramos, C. Martín-Gómez, and I. Zuriguel, “Flow and clogging of a sheep herd passing through a bottleneck,” Physical Review E 91, 022808 (2015).
  * Becco  _et al._ (2006) C. Becco, N. Vandewalle, J. Delcourt, and P. Poncin, “Experimental evidences of a structural and dynamical transition in fish school,” Physica A: Statistical Mechanics and its Applications 367, 487 (2006).
  * Calovi  _et al._ (2014) D. S. Calovi, U. Lopez, S. Ngo, C. Sire, H. Chaté, and G. Theraulaz, “Swarming, schooling, milling: phase diagram of a data-driven fish school model,” New Journal of Physics 16, 015026 (2014).
  * Steager  _et al._ (2008) E. B. Steager, C. B. Kim, and M. J. Kim, “Dynamics of pattern formation in bacterial swarms,” Physics of Fluids 20 (2008).
  * Schaller  _et al._ (2010) V. Schaller, C. Weber, C. Semmrich, E. Frey, and A. R. Bausch, “Polar patterns of driven filaments,” Nature 467, 73 (2010).
  * Sumino  _et al._ (2012) Y. Sumino, K. H. Nagai, Y. Shitaka, D. Tanaka, K. Yoshikawa, H. Chaté, and K. Oiwa, “Large-scale vortex lattice emerging from collectively moving microtubules,” Nature 483, 448 (2012).
  * Veigel and Schmidt (2011) C. Veigel and C. F. Schmidt, “Moving into the cell: single-molecule studies of molecular motors in complex environments,” Nature Reviews Molecular Cell Biology 12, 163 (2011).
  * Wong  _et al._ (2016) F. Wong, K. K. Dey, and A. Sen, “Synthetic micro/nanomotors and pumps: Fabrication and applications,” Annual Review of Materials Research 46, 407 (2016).
  * Vicsek  _et al._ (1995) T. Vicsek, A. Czirók, E. Ben-Jacob, I. Cohen, and O. Shochet, “Novel type of phase transition in a system of self-driven particles,” Physical Review Letters 75, 1226 (1995).
  * Toner and Tu (1995) J. Toner and Y. Tu, “Long-range order in a two-dimensional dynamical xy model: how birds fly together,” Physical Review Letters 75, 4326 (1995).
  * Toner and Tu (1998) J. Toner and Y. Tu, “Flocks, herds, and schools: A quantitative theory of flocking,” Physical Review E 58, 4828 (1998).
  * Toner (2012) J. Toner, “Reanalysis of the hydrodynamic theory of fluid, polar-ordered flocks,” Physical Review E 86, 031918 (2012).
  * Solon  _et al._ (2015a) A. P. Solon, H. Chaté, and J. Tailleur, “From phase to microphase separation in flocking models: The essential role of nonequilibrium fluctuations,” Physical Review Letters 114, 068101 (2015a).
  * Chaté  _et al._ (2008) H. Chaté, F. Ginelli, G. Grégoire, and F. Raynaud, “Collective motion of self-propelled particles interacting without cohesion,” Physical Review E 77, 046113 (2008).
  * Solon  _et al._ (2015b) A. P. Solon, J. B. Caussin, D. Bartolo, H. Chaté, and J. Tailleur, “Pattern formation in flocking models: A hydrodynamic description,” Physical Review E 92, 062111 (2015b).
  * Kürsten and Ihle (2020) R. Kürsten and T. Ihle, “Dry active matter exhibits a self-organized cross sea phase,” Physical Review Letters 125, 188003 (2020).
  * Solon  _et al._ (2022) A. Solon, H. Chaté, J. Toner, and J. Tailleur, “Susceptibility of polar flocks to spatial anisotropy,” Physical Review Letters 128, 208004 (2022).
  * Chatterjee  _et al._ (2022) S. Chatterjee, M. Mangeat, and H. Rieger, “Polar flocks with discretized directions: the active clock model approaching the vicsek model,” Europhysics Letters 138, 41001 (2022).
  * Chatterjee  _et al._ (2020) S. Chatterjee, M. Mangeat, R. Paul, and H. Rieger, “Flocking and reorientation transition in the 4-state active potts model,” Europhysics Letters 130, 66001 (2020).
  * Mangeat  _et al._ (2020) M. Mangeat, S. Chatterjee, R. Paul, and H. Rieger, “Flocking with a q-fold discrete symmetry: Band-to-lane transition in the active potts model,” Physical Review E 102, 042601 (2020).
  * Codina  _et al._ (2022) J. Codina, B. Mahault, H. Chaté, J. Dobnikar, I. Pagonabarraga, and X. Shi, “Small obstacle in a large polar flock,” Physical Review Letters 128, 218001 (2022).
  * Benvegnen  _et al._ (2023) B. Benvegnen, O. Granek, S. Ro, R. Yaacoby, H. Chaté, Y. Kafri, D. Mukamel, A. Solon, and J. Tailleur, “Metastability of discrete-symmetry flocks,” Physical Review Letters 131, 218301 (2023).
  * Solon and Tailleur (2015) A. P. Solon and J. Tailleur, “Flocking with discrete symmetry: The two-dimensional active ising model,” Physical Review E 92, 042119 (2015).
  * Benvegnen  _et al._ (2022) B. Benvegnen, H. Chaté, P. L. Krapivsky, J. Tailleur, and A. Solon, “Flocking in one dimension: Asters and reversals,” Physical Review E 106, 054608 (2022).
  * Xue  _et al._ (2023) T. Xue, X. Li, X. Chen, L. Chen, and Z. Han, “Machine learning phases in swarming systems,” Machine Learning: Science and Technology 4, 015028 (2023).
  * Kürsten  _et al._ (2020) R. Kürsten, S. Stroteich, M. Z. Hérnandez, and T. Ihle, “Multiple particle correlation analysis of many-particle systems: formalism and application to active matter,” Physical Review Letters 124, 088002 (2020).
  * Chatterjee  _et al._ (2018) S. Chatterjee, S. Puri, and R. Paul, “Ordering kinetics in the $q$-state clock model: Scaling properties and growth laws,” Physical Review E 98, 032109 (2018).
  * Chatterjee  _et al._ (unpublished) S. Chatterjee, M. Karmakar, M. Mangeat, R. Paul, and H. Rieger, “Metastability of ordered phase in discretized flocking,” (unpublished).