Table Of ContentPhase transitions in systems of self-propelled agents and related network models
M. Aldana,1,2,∗ V. Dossetti,2 C. Huepe,3 V. M. Kenkre,2 and H. Larralde1
1 Instituto de Ciencias F´ısicas, Universidad Nacional Auto´noma
de M´exico. Apartado Postal 48-3, Cuernavaca, Mor.62251, M´exico
7 2 Consortium of the Americas for Interdisciplinary Science,
0 University of New Mexico. 800 Yale Boulevard NE, Albuquerque, NM 87131, USA
0 3 614 N. Paulina St., Chicago, IL 60622-6062, USA
2
An important characteristic of flocks of birds, school of fish, and many similar assemblies of self-
n propelled particles is the emergence of states of collective order in which the particles move in the
a
same direction. When noise is added into the system, the onset of such collective order occurs
J
throughadynamicalphasetransition controlled bythenoise intensity. Whileoriginally thoughtto
0 be continuous, the phase transition has been claimed to be discontinuous on the basis of recently
3
reportednumericalevidence. Weaddressthisissuebyanalyzingtworepresentativenetworkmodels
closely related to systems of self-propelled particles. We present analytical as well as numerical
]
h resultsshowingthatthenatureofthephasetransition dependscrucially ontheway inwhich noise
c is introduced into thesystem.
e
m PACSnumbers: 05.70.Fh,87.17.Jj,64.60-i
-
t
a
The collective motion of a group of autonomous par- from the one presented in [4] to the one described in [5].
t
s ticles is a subject of intense research that has potential The first network model, which we will refer to as
.
t applications in biology, physics and engineering [1, 2, 3]. the vectorial network model, consists of a network of
a
m One of the most remarkable characteristics of systems N 2D-vectors (represented as complex numbers), {σ1 =
- such as a flock of birds, a school of fish or a swarm of eiθ1,σ2 = eiθ2,...,σN = eiθN}, all of the same length
d locusts, is the emergence of ordered states in which the |σ | = v and whose angles {θ (t),θ (t),...,θ (t)} can
n 1 2 N
n particles move in the same direction, in spite of the fact change in time. Each vector σ interacts with a fixed
n
o
that the interactions between the particles are (presum- set of K other vectors, {σ ,...,σ }, randomly chosen
c n1 nK
ably)ofshortrange. Giventhatthesesystemsaregener- from anywhere in the system. We will call this set of K
[
ally out of equilibrium, the emergence of ordered states vectors the inputs of σn. Once each vector σn has been
1
cannotbeaccountedforbythestandardtheoremsinsta- providedwithafixedsetofK inputconnections,thedy-
v
tistical mechanics that explain the existence of ordered namics of the network are then given by one of the two
3
3 states in equilibrium systems typified by ferromagnets. following interaction rules:
7
A particularly simple model to describe the collective
1 K
0 motionofagroupofself-propelledparticleswasproposed θ (t+1) = Angle 1 σ (t)+ηξ(t), (1)
7 by Vicsek et al. [4]. In this model each particle tends to n vK X nj
j=1
0 move in the average direction of motion of its neighbors
/ while being simultaneously subjected to noise. As the K
mat amplitude of the noise increases the system undergoes a θn(t+1) = Anglev1K Xσnj(t)+ηeiξ(t), (2)
phase transition from an orderedstate in which the par- j=1
- ticles move collectively in the same direction, to a disor-
d
n deredstate in whichthe particles moveindependently in where for any vector ~v = |v|eiφ we define the function
o random directions. This phase transition was originally Angle(~v) = φ, and ξ(t) is a random variable uniformly
c thought to be of secondorder. However,due to a lackof distributed in the interval [−π,π]. The the dynamics of
:
v ageneralformalismtoanalyzethecollectivedynamicsof the network is fully deterministic for η =0 and becomes
i the Vicsekmodel, the natureofthe phasetransition(i.e. more random as the parameter η increases. In what fol-
X
K
r whether itissecondorfirstorder)hasbeenbroughtinto lows, we will refer to the quantity (1/vK)Pj=1σnj as
a question [5]. the average contribution of the inputs of σn.
To quantify the amount of order in the system we de-
In this letter we show that the nature of the phase
fine the instantaneous order parameter ψ(t) as
transition can depend strongly on the way in which the
noise is introducedinto these systems. We illustrate this
N
by presenting analytical results on two different network ψ(t)= lim 1 (cid:12) σ (t)(cid:12). (3)
systemsthatarecloselyrelatedtothe self-propelledpar- N→∞vN (cid:12)(cid:12)X n (cid:12)(cid:12)
(cid:12)n=1 (cid:12)
ticle models. We show that in these two network models (cid:12) (cid:12)
the phase transition switches from second to first order In the limit t → ∞, the instantaneous order parameter
when the way in which the noise is introduced changes ψ(t) reaches a stationary value ψ [4, 5, 6, 7]. Thus, in
2
the stationary state all the vectors are aligned if ψ ∼ 1, 1 Vicsek model v = 0.1
Vicsek model v = 0.5
whereas if ψ ∼0 the vectors point in random directions. 0.8 Vicsek model v = 1.0
Vicsek model v = 5.0
The interaction rules given in Eqs. (1) and (2) were Network model K = 5
proposed by Vicsek et al. in Ref. [4], and by Gr´egoire ψ0.6
Network Model
and Chat´e in Ref. [5], respectively. The difference be-
0.4
tween these two interaction rules consists in the way in Vicsek model
which the noise is introduced: In Eq. (1) the noise is 0.2 with different speeds
added outside the Angle function, i.e. after the Angle
0
function has been applied to the averagecontribution of 0 0.2 0η.4 0.6 0.8
the inputs. On the other hand, in Eq. (2) the noise is
added inside the Angle function, i.e. it is added directly
FIG. 1: (Color online.) Phase diagram of the Vicsek model
to the average contribution of the inputs. In Ref. [5],
andthevectorialnetworkmodelforthecaseinwhichthenoise
Gr´egoire and Chat´e posed the question as to whether
isaddedasinEq.(1). Whenthespeedvoftheparticlesinthe
thesetworulesleadtothesametypeofphasetransition.
Vicsekmodelincreases,thephasetransitionconvergestothat
In this letter we show that the interaction rules in of the vectorial network model. The numerical simulations
Eq.(1)andEq.(2)producedifferenttypesofphasetran- werecarriedoutforsystemswithN =20000particlesandan
sitionsinthenetworksystemsunderconsideration,which average numberof interactions perparticle K =5.
suggests that a similar effect is being observed in [5] for
the self-propelled systems.
served when the noise is introduced as in Eq. (1), albeit
Obviously, in the self-propelled particle models the el-
the finite size effects observed near the critical point.
ements do notinteractthrougha network. Instead, they
moveina2Dspace,eachparticleinteractinglocallywith The probability distribution function (PDF) of the
the particles that fall within a certain radius. This mo- sum v1K PKj=1σnj(t)+ηeiξ(t) that appears in Eq. (2) is
tion allows particles that are initially far apart to meet, computed as for a random walk assuming that all the
interact, and separate again, giving rise to effective long terms are statistically independent. By projecting this
range interactions. On the other hand, in our vectorial PDF ontothe unit circle we canestablisha recursionre-
network model the particles are fixed to the nodes of a lationfortheorderparameter,whichforK ≫1becomes
network. The long-range correlations produced by the ψ(t+1)=Mη(ψ(t)), where
motion of the particles in the self-propelled models are
proxiedinournetworkmodelthroughrandomlychoosing ψ(t) F 1,1,2,[ψ(t)]2 if ψ(t)<η
2η 2 1(cid:16)2 2 η2 (cid:17)
theinputsofeachelementfromanywhereinthenetwork.
An underlying assumption of our work is that the exis- Mη(ψ(t))≡
F 1,−1,1, η2 if η <ψ(t)
toenncteheanodccnuartruerneceofotfhseupchhalsoentgr-arnasnigtieonindteerpaecntdiosnms,oastnldy 2 1(cid:16)2 2 [ψ(t)]2(cid:17) (4)
less crucially on whether they are produced by the mo- and F (a,b;c,x) is the Gauss hypergeometric function.
2 1
tion of the particles or by the network topology [6, 7]. ThefunctionMη(ψ)isshowninFig.2afordifferentval-
While the exact relation between these two ways of es- ues of η (solid curves). This figure also displays with
tablishinglong-rangeinteractionsisnotyetknown,ithas symbolsthe numericaldynamicalmappingcomputedfor
been shown that a strong parallelcan be established be- the self-propelled model with the interaction rule given
tween them [7, 9]. Further, below we show that there in Eq. (2), N = 20000 particles, and an average num-
are atleast two limits in which they arefully equivalent: ber of interactions per particle K = 100. Clearly, the
forlargeparticlespeedsandforhighdensities(seeFig.1 numerical mapping coincides with the theoretical result
and Fig. 2). for Mη(ψ), showing that the network and self-propelled
InRef.[7]itwasproventhat,asthenoiseamplitudeη systemsare alsoequivalentin the highdensity limit case
increases, the vectorial network model with the interac- considered here. The numerical mappings for the self-
tionrulegivenasinEq.(1)undergoesacontinuousphase propelled system were obtained by placing the particles
transitionfromorderedstateswhereψ >0,todisordered in various random initial conditions constrained to pro-
states where ψ = 0. Fig. 1 shows this phase transition duce everyorder parametervalue ψ(t) in the x-axis,and
obtained numerically for N = 20000 and K = 5. It also thencomputing one time stepusing Eq.(2)to obtainthe
displays the phase transition in the Vicsek model for a corresponding value ψ(t+1) in the y-axis.
systemwiththe sameN,adensitysuchthattheaverage The fixedpoints ofthe dynamicalmapping ψ(t+1)=
number of interactions per particle is also K = 5, and M (ψ(t)) give the stationary values of the order param-
η
increasing particle speeds. As can be seen from Fig. 1, eter. From Eq. (4) it is clear that ψ = 0 is always a
the Vicsek model curves approach continuously the net- fixed point. However, the stability of this fixed point
work model curve as v → ∞. This supports the idea changes depending upon the value of η. By numerically
that in both cases a second order phase transition is ob- solving Eq.(4) to obtain the fixed point, we find that
3
1.0 1
(a)
0.8
0.8
0.6 esis
0.4 ηη == 00..46 ψ0.6 Theory yster
ψ) 0.2 ηη == 00..6872 0.4 SSiimmuullaattiioonn of h
( 0.0 n
Mη 1.0 0.2 egio
(b) R
0.8 0
0.6 0 0.2 0.4 η 0.6 0.8 1
η = 0.2
0.4 η = 0.4
0.2 η = 0.61 FIG. 3: (Color online.) Phase diagram of the vectorial net-
η = 0.8
0.0 work model for the case in which the noise is added as in
0 0.2 0.4 0.6 0.8 1
ψ Eq.(2). Thesolidlinecorrespondstothepredictionobtained
from Eq. (4). The dashed and dotted-dashed curves are the
results of the numerical simulation starting out the dynam-
FIG. 2: (Color online.) Graph of the dynamical mapping
Mη(ψ). (a) The interaction rule is as in Eq. (2). The solid ics from initial conditions for which ψ(0) ≈ 1 and ψ(0) ≈ 0,
respectively. The phase transition in this case is clearly dis-
curvescorrespond to theanalytical solution given in Eq. (4),
continuous.
andthesymbolstothenumericalsimulationcarriedoutfora
systemwithN =20000andanaveragenumberofinteractions
perparticleK =100. (b)TheinteractionruleisasinEq.(1).
limits of metastability and decay at values of η slightly
ThecurveswerecomputednumericallyforasystemwithN =
20000 particles and K = 5. In (a) the non-zero stable fixed above 1/2 and below 0.672, as observed in the graph.
point appears discontinuously as η decreases, whereas in (b) Additionally, Eq. (4) was obtained in the limit of large
it appears continuously. K, however, already for K =20 the agreement is good.
The secondmodel thatwe consideris a majority voter
model inwhichthenetworkelementsσ canacquireonly
n
for 0.672 < η the only stable fixed point is ψ = 0. As two values, +1 or −1. We can think of this system as
η →0.672fromabove,the graphof Mη(ψ) moves closer a society in which every individual σn has to make a
to the identity and eventually another non zero stable decision about an issue with two possible alternatives,
fixed point ψ′ appears discontinuously when η ≈ 0.672 either +1 or −1. Again, each element σn receives inputs
(see the point indicated with an arrow in Fig. 2a). For from a set of K other elements randomly chosen from
1/2 < η ≤ 0.672 there are actually two stable fixed anywhere in the system. Let us first consider the case
points. In this region of bi-stability the system shows in which the interaction between σn and its K inputs,
hysteresis. Finally, when η < 1/2 the fixed point ψ = 0 {σn1,σn2,...,σnK}, is given by
′
becomes unstable and only the non zero fixed point ψ
K
remains stable. Contrary to this, when the noise is in- 1 ξ(t)
σ (t+1)=Sign Sign σ (t) + (5)
troduced as in Eq. (1) the non-zero stable fixed point n K X nj 1−η
appears continuously, as can be seen from Fig. 2b. j=1
Thevalidityoftheseresultsiscorroboratedbynumer- where Sign[x] = 1 if x > 0, Sign[x] = −1 if x < 0, η is
ical simulations carried out for networks with N = 105 a parameter that takes a constant value in the interval
andK =20. Fig.3showsthefixedpointψofEq.(4)asa [0,1],andξ(t)isarandomvariableuniformlydistributed
functionofthenoiseintensityη(solidline). Thedisconti- between [−1,1]. For the Sign function to be well defined
nuity ofthe orderparameterψ atη =0.672andη =1/2 we choose K as an odd integer. Eq. (5) is similar to
is apparent. The dashed and dotted-dashed curves are Eq.(1) in thatthe noise is addedto the signof the aver-
the plotsofthe resultsfromthe numericalsimulationfor age contribution of the inputs. Since in this case σ is a
n
thecasesinwhichallthe vectorswereinitiallyalignedin discrete variablethattakesonly the twovalues +1 or-1,
thesamedirection(ψ(0)=1),andwhenthevectorswere the Sign function has to be applied again. This interac-
initiallyorientedinrandomdirections(ψ(0)≃0),respec- tion rule reflects the fact that an individual in a society
tively. In the region of bi-stability 1/2 < η < 0.672, the usually tends to be of the same opinion as the majority
system reaches one or the other of the two stable fixed ofhis“friends”(inputs), though,withprobabilityη/2he
points depending upon the initial condition. can have the opposite opinion.
The theoretical curves presented in Fig. 3 show the Theinstantaneousorderparameterψ(t)isdefinedasin
“limits of metastability” for the system; i.e. the maxi- Eqs.(3),butnowtheverticalbarsrepresenttheabsolute
mal and minimal values of η for which the system has value instead of the norm of a vector.
bistable behavior (hysteresis). Clearly, specific realiza- In Ref. [8] it has been shown that the majority voter
tions of the system cannot be driven all the way to the model with the interaction rule given by Eq. (5) under-
4
1 In summary, we have analyzed numerically and an-
alytically the phase transition from ordered to disor-
0.8
dered states in two network models that capture some
ψ0.6 of the main aspects of the interactions in systems of
self-propelled particles. In particular, the self-propelled
0.4
modelbecomesequivalenttothevectorialnetworkmodel
0.2 in the limit of large speeds or high densities. We have
0 shown that for the two network models, the phase tran-
0 0.1 0.2 0.3 η0.4 0.5 0.6 0.7 sitionchangesfromsecondordertofirstorderdepending
on the way in which the noise is introduced into the sys-
tem. This change is consistent with the results reported
FIG. 4: (Color online.) Phase transitions in the majority
by Vicsek et al. in [4], and by Gr´egoire and Chat´e in
voter model when the noise is added as in Eq. (5) (solid line
[5]. This consistency suggests that a similar effect is be-
with circles), and as in Eq.(6) (dashed line with diamonds).
The numerical simulations were carried out for systems with ing observed in the self-propelled model and motivates
N =105 and K =3. a deeper analysis in order to determine the nature of its
phase transition. Clearly, the two ways of introducing
noise correspondto different physical situations. On the
goes a continuous phase transition as the value of η is
one hand, with the Vicsek type of noise the uncertainty
varied. In Fig. 4 we reproduce this phase transition for
falls onthe decisionmechanism. Onthe other,introduc-
networks with N =105 and K =3 (solid curve with cir-
ing the noise a` la Gr´egoire-Chat´e, the decision function
cles). It is clear from this figure that, when the noise is
is perfectly determined and the uncertainty falls on the
added as in Eq. (5), the phase transition in the majority
the argumentsofthis function. Thereisno reasonto ex-
voter model is indeed continuous.
pect,apriori,similarbehaviorsunderthesetwodifferent
Let us now consider the case in which the interaction
physical situations.
rule between σ and its inputs is given by
n
Noteaddedinproof: While this manuscriptwasbeing
K reviewed, Ref. [10] appeared, in which it is shown that
1
σn(t+1)=Sign K Xσnj(t)+2ηξ(t) (6) thefirst-orderphasetransitionfoundinRef.[5]bymeans
j=1 of the Binder cumulant is a numerical artifact.
This work was partly supported by NSF Grant No.
where ξ(t) is a random variable uniformly distributed in
INT-0336343and CONACyT Grant No. P47836-F.The
the interval [−1,1] and η ∈ [0,1]. This rule is similar
work of C.H. was supported by the National Science
to that given in Eq. (2) in that the noise is added to the
Foundation under Grant No. DMS-0507745.
averagecontributionoftheinputsandthentheSignfunc-
tion is evaluated. Again, the PDF of the sum appearing
in Eq. (6) can be computed as for a random walk as-
suming that all the terms are statistically independent.
Integrating the PDF over all positive values of the sum ∗
Electronic address: max@fis.unam.mx
we obtain a recursion relation for the order parameter [1] C.M. Topaz and A.L. Bertozzi, SIAM J. Appl. Math.
(see [8]), which for K =3 becomes 65, 152 (2005). M.R. D’Orsogna, Y.-L. Chuang, A.L.
Bertozzi, and L. Chayes, Phys. Rev. Lett. 96, 104302
3ψ(t)− 1[ψ(t)]3, for 0<η < 1 (2006).
2 2 6 [2] J.K. Parrish and L. Edelstein-Keshet, Science 248, 99
ψ(t+1)= 1+8η6η ψ(t)− 1−8η2η [ψ(t)]3, for 61 <η < 21 (B1i9o9l.9)B.uJll..K2.0P2a,r2r9is6h,(2S0.0V2.).Viscido, and D. Gru¨nbaum,
[3] I.D.Couzin,J.Krause,N.R.Franks,andS.A.Levin,Na-
21η ψ(t), for 12 <η ture433,513(2005).J.Buhl,D.J.Sumpter,I.D.Couzin,
J. Hale, E. Despland, E. Miller, and S.J. Simpson, Sci-
Althoughψ =0 is alwaysa fixedpoint ofthe previous
ence 312, 1402 (2006).
equation,astabilityanalysisrevealsthatfor0<η <1/2
[4] T. Vicsek, A. Czir´ok, E. Ben-Jacob, I. Cohen, and O.
the solution ψ = 0 is unstable and the only stable fixed Shochet, Phys.Rev. Lett.75, 1226 (1995).
point is ψ = 1. In the region 1/2 < η the fixed point [5] G. Gr´egoire and H. Chat´e, Phys. Rev. Lett. 92, 025702
ψ = 1 disappears altogether and the only fixed point (2004).
that remains is ψ = 0, and it is stable. Therefore, at [6] J. Toner and Y. Tu,Phys. Rev.Lett. 75, 4326 (1995).
[7] M.AldanaandC.Huepe,J.Stat.Phys.112,135(2003).
the critical value η = 1/2 the system undergoes a dis-
[8] C.HuepeandM.Aldana,J.Stat.Phys.108,527(2002).
continuous phase transition from a totally ordered state
M. Aldana and H. Larralde, Phys. Rev. E 70, 066130
(ψ = 1) to a fully random state (ψ = 0). Fig. (4) shows
(2004).
the results ofthe numericalsimulationfor a systemwith [9] J. D. Skufcaand E. M. Bollt, Math. Boisci. and Eng. 1,
N =105 and K =3 (dashed curve with diamonds). 347 (2004).
5
[10] M. Nagy, I. Daruka and T. Vicsek, Physica A 373, 445
(2007).
