Flocking and turning: a new model for self-organized collective motion  Andrea Cavagna a,b,c , Lorenzo Del Castello a,b , Irene Giardina a,b,c , Tomas Grigera d,e , Asja Jelic a,b , Stefania Melillo a,b , Thierry Mora f   , Leonardo Parisi a,b , Edmondo Silvestri a,b,g   , Massimiliano Viale a,b , and Aleksandra M. Walczak h  a Istituto Sistemi Complessi (ISC–CNR), Via dei Taurini 19, 00185 Roma, Italy  b Dipartimento di Fisica, “Sapienza” Universit ́ a di Roma, P.le Aldo Moro 2, 00185 Roma, Italy  c Initiative for the Theoretical Sciences, The Graduate Center, 365 Fifth Avenue, New York, NY 10016 USA  d   Instituto de Investigaciones Fisicoqu ́ ımicas Te ́ oricas y Aplicadas (INIFTA) and Departamento de F ́ ısica, Facultad de Ciencias Exactas, Universidad Nacional de La Plata, c.c. 16, suc. 4, 1900 La Plata, Argentina  e CONICET La Plata, Consejo Nacional de Investigaciones Cient ́ ıficas y T ́ ecnicas, Argentina  f   Laboratoire de Physique Statistique de l’  ́ Ecole Normale Sup ́ erieure, CNRS and Universites Paris VI and Paris VII, 24 rue Lhomond, 75231 Paris Cedex 05, France  g   Dipartimento di Matematica e Fisica - Universit ́ a Roma Tre, Via della Vasca Navale, 84, 00146 Roma, Italy, and  h Laboratoire de Physique Th ́ eorique de l’  ́ Ecole Normale Sup ́ erieure, CNRS and University Paris VI, 24 rue Lhomond, 75231 Paris Cedex 05, France  (Dated: January 22, 2015) Birds in a flock move in a correlated way, resulting in large polarization of velocities.   A good understanding of this collective behavior exists for linear motion of the flock. Yet observing actual birds, the center of mass of the group often turns giving rise to more complicated dynamics, still keeping strong polarization of the flock. Here we propose novel dynamical equations for the collec- tive motion of polarized animal groups that account for correlated turning including solely social forces. We exploit rotational symmetries and conservation laws of the problem to formulate a the- ory in terms of generalized coordinates of motion for the velocity directions akin to a Hamiltonian formulation for rotations.   We explicitly derive the correspondence between this formulation and the dynamics of the individual velocities, thus obtaining a new model of collective motion. In the appropriate overdamped limit we recover the well-known Vicsek model, which dissipates rotational information and does not allow for polarized turns.   Although the new model has its most vivid success in describing turning groups, its dynamics is intrinsically different from previous ones in a wide dynamical regime, while reducing to the hydrodynamic description of Toner and Tu at very large length-scales.   The derived framework is therefore general and it may describe the collective motion of any strongly polarized active matter system.  Flocks of birds represent paradigmatic cases of emer- gent collective behavior, where short range interactions between individuals lead to nontrivial movement patterns on the large scale[1–5]. Several models have been devel- oped by biologists [6–10], physicists [11–17] and control theorists [23, 24] to describe flocking behavior.   Some of these models greatly improved our understanding of collective motion, offering a powerful description of the large scale properties of active systems [18–21].   How- ever, recent experiments on real flocks [25] revealed new surprising features in the way collective turns are per- formed that are still unaccounted for, urging for a novel theoretical explanation. The most conspicuous feature of flocking is the pres- ence of collective order:   flocks are strongly polarized groups, where individuals all move in approximately the same direction.   This kind of orientational order can be naturally explained in terms of alignment interactions between the birds velocities and basically all models of flocking include a ‘social’ force term, describing the ten- dency of individuals to adjust their flight direction to those of neighbors. Recently, using inference techniques on real data of flocking birds [26], it was in fact shown that simple models based on pairwise alignment topolog- ical interactions are able to explain the long-range corre- lations between flight directions and the large degree of coherence exhibited by natural flocks [27]. The role of alignment is emphasized in the simplest model of flocking, the Vicsek model [11]. In the Vicsek model flocking individuals are described as self-propelled particles of constant speed, where the velocity (i.e.   the flight direction) of each particle is updated from one time step to the next by computing the average direction of motion of neighbors. Several variants of the model have been analyzed in the literature, where additional cohesive terms are added [9, 13, 14, 18], and interactions are cho- sen topologically rather than metrically (as in the orig- inal version) [17, 28]. In all these cases the structure of the dynamical equations always remains the same, where alignment forces due to neighbors directly act on the in- dividual velocity.   In other terms, the velocity evolves according to a first order Langevin equation, its time increments being determined by the social force and ran- dom noise. This kind of dynamics generates long-range order at low noise and scale-free correlations in the polar- ized phase, both experimentally observed features of nat- ural flocking. Additionally, at very large scales all these models are well described by the hydrodynamic theory of flocking, whose intriguing predictions of giant density fluctuations, propagating sound modes and anomalous  arXiv:1403.1202v2 [cond-mat.stat-mech] 21 Jan 2015
2 diffusion have been observed in a variety of active inani- mate and animate fluids [18, 20–22]. However, the Vicsek model does not correctly describe the way collective turns occur in natural flocks. Recent experiments on starling flocks [25] show that collective turns involve linear propagation of purely directional in- formation (the direction of motion of the birds) and cur- vature (of the individual trajectories), without any local change in density. Turns start locally in space (one sin- gle individual starts the turn) and occur very quickly, the turning front propagating with a speed that does not depend on the density of the flock.   These ‘turning modes’ travel extremely fast and undamped through the group, with a linear dispersion law analogous to super- fluid transport [25]. This is different from what happens in the Vicsek model: the hydrodynamic theory of flock- ing, that describes the long-wavelength behavior of the Vicsek model, predicts propagating sound modes that couple orientational and density disturbances [18, 22]. These are not, however, the kind of modes observed in experiments, where density does not play any role and curvature propagates together with direction.   In fact, turns typically occur in finite groups and on short scales, outside the hydrodynamic regime described by [18, 22]. On such scales the Vicsek model does not account for directional propagation:   if a directional disturbance is created in the system, such as a few individuals willing to turn and change direction of motion, the turning in- formation is transmitted diffusively and attenuated, i.e. ‘turning’ modes are non-propagating overdamped modes (see Fig. 1 and sec. V). In Ref.   [25] some of us showed that there are two crucial and connected ingredients missing in the micro- scopic dynamics of the Vicsek model, that must be taken into account to theoretically reproduce the dispersion law observed in turning flocks: 1) the presence of a conser- vation law associated to the rotational symmetry of the system; 2) the existence of a behavioral inertia mediating the effect of the social force. In this paper we introduce a new microscopic model for flocking where these two ingredients are appropriately incorporated.   This model provides a unified general set of equations to describe sys- tems of interacting self-propelled individuals, where all the salient traits of real flocking - not only order and cor- relations, but also turns and information propagation - are correctly reproduced. Interestingly, the Vicsek model can be retrieved as the overdamped limit of this more general model. Additionally, in the long wavelength limit the model is well described by the hydrodynamic theory of flocking of Toner and Tu [18] and is therefore consis- tent with its predictions on the large scale behavior. Because the derivation of our new model will be some- what lengthy, we start by stating the final equations. As we shall see, the crucial ingredient of the new dynamics is the conservation of the internal angular momentum, namely of the spin and the introduction of a spin inertia. For this reason we will refer to these new equations as start of turn  start of turn  Vicsek model real !ock  FIG. 1. Schematic representation of turn propagation in real and Vicsek flocks.   a . In real flocks when a bird starts a turn the turning information is propagated unattenuated through the whole group and the flock as a whole performs a collec- tive turn.   b.   In the Vicsek model the turning information is damped and does not propagate. The representation concerns scales (group sizes and time) where collective turns occur. In the hydrodynamic regime the Vicsek model exhibits propa- gating sound modes - see Sec.IV A.  the   inertial spin model :  d~ v i  dt   =   1  χ~ s i   ×   ~ v i  d~ s i  dt   =   ~ v i ( t )   ×      J v 2 0  ∑  j  n ij ~ v j   −   η v 2 0  d~ v i  dt   +   ~ ξ i  v 0      (1)  d~ r i  dt   =   ~ v i ( t ) with noise correlator  〈 ~ ξ i ( t )   ·   ~ ξ j   ( t ′ ) 〉   = (2 d )   η T δ ij   δ ( t   −   t ′ )   .   (2) In these equations   ~ r i   represents the position of particle  i ,   ~ v i   its velocity (with fixed speed   | ~ v i |   =   v 0 ), and   ~ s i   is a new variable associated to each particle, representing a generalized momentum (the spin), connected to the in- stantaneous curvature of the particle’s trajectory.   The parameter   χ   is a generalized moment of inertia,   η   a fric- tion coefficient, and   J   the strength of the alignment force to neighbors. The connectivity matrix   n ij   describes who is the neighbour of whom.   Finally   T   is a generalized temperature. These equations might look unfamiliar and difficult to interpret, but we will show that they are quite natural. Most of the paper is devoted to explaining how the equa- tions can be derived from general symmetry considera- tions, what is their meaning, what are the interesting limiting cases (the Vicsek model and the deterministic limit), and what are the predictions of the model for in- formation propagation in the ordered phase. We hope the reader will go through this whole exercise, but for those less interested in the details, the last section of the paper presents in a self-contained way numerical simulations of model (1) and illustrates its behavior.
3  I.   SETTING UP THE STAGE  Our aim in this work is to introduce a model of flocking able to describe collective turning and information prop- agation in flocks. Since these phenomena occur on finite scales of space and time, we need to go beyond hydrody- namics.   As a consequence, the universality of the long length scale behavior is lost, and a degree of arbitrariness arises in the way we chose the terms of the microscopic dynamical equations.   One way to address this issue is to propose a model and check a posteriori its validity. The reader can look at Eqs.   (1) as our guess, and then directly skip to the last section, where we show with nu- merical simulations that model (1) correctly describes the experimental results of [25]. In the following sections we will follow a longer route, and show that it is possible to construct a microscopic model of collective motion starting from effective vari- ables, symmetries, conservation laws and Hamiltonian dynamics.   This may seem a bit formal when we think about systems as complex and highly dissipative as real flocks. However, our point will be that all the nontrivial dissipative nature of animals’ active motion can be cap- tured by an effective constraint on their velocity.   Once we find the right set of canonical variables enforcing this constraint, we are able to give a Hamiltonian description of flocks. In a second part of the paper, we then gener- alize this description to include non-Hamiltonian contri- butions and noise. In Sec.   I A we want to make this formal description easier to understand by going back to a very simple (in fact, quite trivial) example, showing that the use of effec- tive variables and symmetries in presence of a constraint is quite standard in physics. For the sake of simplicity we will assume in this section that bird velocities lie on a plane. Although this is true in most experimental cases [25], it is not an essential mathematical condition, and in the next sections we will analyze the fully three-dimensional case.   We will also disregard the reciprocal movement of birds (which is very low during real turns [25]) when deriving analytically an expression for the dispersion law.   Later in Sec.   V we then perform numerical simulations of model (1) where the complete dynamics of positions, velocities and spins is implemented, and network rearrangements are fully taken into account.  A.   Circular motion in the external space of position  Consider a point on a plane, with coordinate   ~ r   and linear momentum   ~ p , subject to the constraint,  | ~ r   |   =   r 0   .   (3) The constraint indicates that the point is confined on a circle of radius   r 0 . How is the circular constraint imple- mented in practice?   We do not know, but probably in a complicated way.   We know there must be forces act- ing on the particle that are hard to describe explicitly. We also know that the standard canonical variables ( ~ r, ~ p ) are not very convenient to describe the system. The stan- dard thing to do is to introduce a new coordinate that automatically enforces constraint (3). This coordinate is the polar angle   θ , whose relation to the coordinate of the particle is simply  ~ r   =   r 0   exp ( iθ )   .   (4) The variable   θ   parametrizes the rotation in the (external) space of coordinates. The momentum conjugated to   θ   is the generator of this rotation [29, 30], namely the orbital angular momentum,   l . Hence, ( θ, l ) are a pair of canon- ical conjugated variables and, in absence of forces other than the ones enforcing the constraint, the Hamiltonian is given by  H   =   l 2  2 I ,   (5) where   I   is the moment of inertia – the resistance of the point to change its angular momentum.   The Hamilton equations of motion are  dθ dt   =   { θ, H }   =   ∂H ∂l   =   l/I ,   (6)  dl dt   =   { l, H }   =   −   ∂H ∂θ   = 0   ,   (7) where the Poisson bracket is defined as  { A, B }   =   ∂A ∂θ ∂B ∂l   −   ∂A ∂l ∂B ∂θ   .   (8) The fact that the angular momentum   l   is the generator of the rotation parametrized by   θ   is expressed by the relation,  dA dθ   =   { A, l }   .   (9) The   rotational   symmetry   of   the   problem   ensures (through Noether’s theorem) that the generator of the rotation,   l , is conserved, so that the particle performs rotational uniform motion with constant angular veloc- ity,    ̇ θ .   However, if a dissipative term proportional to    ̇ θ  were present in (7) the particle would eventually stop ro- tating, unless there is some external injection of angular momentum. The message of this repetition of classical mechanics is that, regardless of the complexity by which the constraint is enforced in practice (3), once we introduce the right canonical variables we can completely forget about these complications. Additionally, the symmetry of the Hamil- tonian with respect to the transformation parametrized by the effective coordinate ensures the conservation of the generator of the transformation itself, that is the momen- tum.
4  B.   Circular motion in the internal space of velocity  Flocks are out-of-equilibrium systems, where energy is continuously injected and dissipated at the individ- ual level. For this reason, linear momentum is not con- served and it is not possible to use a Hamiltonian de- scription in the standard canonical variables (positions and velocities).   However, the only crucial consequence of this very complicated aerodynamics and energetics of bird motion (including the interaction with the surround- ing medium) is that individuals move at approximately constant speed. Experimental findings show that fluctu- ations in the speed are indeed very small, both during straight flight [27] and during turns [25].   The simplest way to describe the active nature of the system is there- fore to model individuals as self-propelled particles mov- ing at constant speed, as in the original Vicsek model,  | ~ v   |   =   v 0   .   (10) This constraint on the velocity of a single bird is analo- gous to the constraint on the position of standard circular motion which we have seen above, equation (3). Hence, as in that case, we must find a generalized coordinate which automatically enforces this constraint.   This co- ordinate is the phase,   φ , namely the orientation of the velocity on the plane,  ~ v   =   v 0   exp ( iφ )   .   (11) Note that the phase   φ   is the parameter of the rotation in the   internal   space of velocity, while   θ   was the param- eter of the rotation in the   external   space of positions. This difference is crucial: the generator of the rotation parametrized by   φ   is   not   the orbital angular momentum  l , but rather the internal angular momentum, or spin,   s . The variables ( φ, s ) are the canonically conjugated pair enforcing constraint (10), exactly as ( θ, l ) are those en- forcing constraint (3). For this reason, the single-particle Hamiltonian in absence of forces other than those enforc- ing the constraint is given by  H   =   s 2  2 χ ,   (12) where   χ   is a generalized moment of inertia quantifying the resistance of a bird to a change of its spin.   The equations of motions are  dφ dt   =   { φ, H }   =   ∂H ∂s   =   s/χ , ds dt   =   { s, H }   =   −   ∂H ∂φ   = 0   .   (13) The symmetry of the Hamiltonian under rotations in the internal space ensures that the spin is conserved. Hence, this motion is simply one where the velocity rotates at constant angular velocity,    ̇ φ   =   s/χ . Here we are consid- ering the idealized case where noise and dissipation are absent in the system, so that the spin is fully conserved. Later on we will introduce dissipation, which causes even- tually the spin, and thus the angular velocity, to decay to zero in absence of external injection of spin.  C.   From one to many birds  To describe flocks we must generalize Hamiltonian (13) to the case of many interacting particles, each one with its phase and spin, ( φ i , s i ).   The obvious generalization is the following [25]:  H   =   V   ( { φ i } ) +   ∑  i  s 2  i  2 χ .   (14) In this expression,   V   ( { φ i } ) is a potential describing the interaction between particles.   The Hamilton equations of motions are  dφ i  dt   =   ∂H ∂s i  =   s i  χ ds i  dt   =   −   ∂H ∂φ i  =   −   ∂V ∂φ i  .   (15) The only thing we need to know about the interaction po- tential   V   is that it is rotationally symmetric: because all directions of motion are   a priori   equally likely at the indi- vidual and collective level, it is reasonable to expect that interactions between individuals will respect this symme- try. This means that   V   (and therefore the Hamiltonian) is invariant under a global rotation of all the velocity vec- tors   ~ v i , that is under the transformation   φ i   →   φ i   +   δφ . Thanks to Noether’s theorem, this global symmetry im- plies a global conservation law. Although the individual spin   s i   is now   not   constant, due to the force   ∂V /∂φ i , the  total   spin,   S   =   ∑  i   s i , is conserved.   This is clear from Hamilton’s second equation,  dS dt   =   −   ∑  i  ∂V ∂φ i  = 0   ,   (16) where the r.h.s. vanishes as a consequence of the symme- try,   V   ( { φ i } ) =   V   ( { φ i   + δφ } ). We will see later on that in the continuum space limit, equation (16) becomes a con- tinuity equation for the spin field. We will also see that this conservation law is the crucial mathematical ingre- dient that allows for the propagation of undamped turns in natural flocks [25]. It is very important to understand that in the case of many particles, a uniform rotation of the polar angle   θ  and a uniform rotation of the phase   φ   give two completely different results. The   θ -rotation, generated by the orbital angular momentum   l , acts on the positions of the points and therefore it gives rise to   parallel paths   trajectories: these all have the same origin as a centre of rotation, but different radii of curvature.   On the other hand, the   φ - rotation, generated by the spin   s , acts on the velocities of the points, giving rise to   equal radius   trajectories: these have different centres of rotation, but the same radius
5 of curvature. Turning flocks of birds are known to move along equal radius trajectories [25, 31], as this is the only way to keep cohesion at constant speed. This experimen- tal observation confirms that phase and spin are indeed the correct canonical variables. Let us summarize this introductory Section.   We de- scribe birds’ motion through a set of effective conjugated variables, the phase   φ i   and the spin   s i , that automati- cally enforce the constraint of constant speed,   | ~ v i   |   =   v 0 . This allows us to disregard the highly complex and dissi- pative mechanism enforcing this constraint and to switch to a Hamiltonian description. Once the speed is constant, the only transformation that the velocity can undergo is a rotation, which is parametrized by the phase. The gen- erator of this rotation is the spin (internal angular mo- mentum), so that in presence of rotationally symmetric forces the total spin of the system is conserved. The way real birds turn is of course complicated, involving the de- tails of flight aerodynamics.   However, as long as birds are able to control their motion keeping their speed fixed and performing equal radius turning - as they do - they can be described in a simple way within this framework.  II.   DETERMINISTIC EQUATIONS FOR THE VELOCITY AND THE SPIN  We have written the dynamical equations for the phase and the spin, which are the canonical variables, but not for the velocity. However, in order to actually run a sim- ulation of a system of self-propelled particles, we need to know how the velocities evolve in time.   Because the velocity is   not   a canonical variable, the update equations for   ~ v i   will be nontrivial. Moreover, we need to study the fully three-dimensional case, in order to have a model as general as possible. In this Section we will write the dy- namical equations for the velocity, in absence of noise and dissipation, which amounts to describing a deterministic flock. Noise and dissipation will be introduced later on.  A.   From the phase to the velocity  In the context of Hamiltonian dynamics, there is a standard method to retrieve the dynamical evolution of any observable once the Hamiltonian is given in terms of canonical variables (here phase and spin). This method is to calculate the Poisson brackets with the Hamilto- nian, which is the time-evolution operator. In the simple planar case discussed in the previous section, we have  v i   = exp( iφ i ), hence  dv i  dt   =   { v i , H }   =   ∂v i  ∂φ i  ∂H ∂s i  −   ∂v i  ∂s i  ∂H ∂φ i  = i   v i  s i  χ .   (17) In the planar case the velocity is a two dimensional vector on the ( x, y ) plane, whereas the spin (the generator of the rotations of   v   in the plane) can be seen as a vector along the   z   direction orthogonal to the plane. In this vectorial description, the term i   v i s i   at the r.h.s of (17) is a vector orthogonal to both   v i   and   s i , hence it is simply the cross (vector) product between   v i   and   s i .   We therefore can write (with a slight abuse of notation)  d~ v i  dt   = 1  χ~ s i   ×   ~ v i   .   (18) The derivation of equation (18) can be generalized to the case of a fully three-dimensional velocity.   In this case we must introduce three phases,   { φ a i   }   ( a   =   x, y, z ), each one parametrizing the rotation around a different carte- sian axis.   Accordingly, the generator of the symmetry, i.e.   the spin, becomes a full three-dimensional vector,  ~ s i   =   { s a i   } . With this notation, for example,   φ z   is the an- gle parametrizing rotations around the   z   axis (i.e. acting in the   xy   plane) and   s z   is the corresponding spin (ori- ented along   z ). Rotations around the   z   axis leave the   z  component of   ~ v   unaltered while effectively rotating the projection of   ~ v   in the plane orthogonal to   z .   Thus   φ z  plays the role of a polar coordinate in the   xy   plane, and it is equivalent to   φ   in the planar case discussed in the previous section. We have,  v y   =  √  v 2 0   −   v 2  z   sin( φ z   )   (19)  v x   =  √  v 2 0   −   v 2  z   cos( φ z   )   .   (20) Analogous relations hold connecting   φ x   and   φ y   to the components of   ~ v .   We can now explicitly compute the equations for the velocities by using Poisson’s brackets,  dv a i  dt   =   { v a i   , H }   =   ∑  c  ∂v a i  ∂φ c i  ∂H ∂s c i  −   ∂v a i  ∂s c i  ∂H ∂φ c i  =    abc v b i  s c i  χ ds a i  dt   =   { s a i   , H }   =   ∑  c  ∂s a i  ∂φ c i  ∂H ∂s c i  −   ∂s a i  ∂s c i  ∂H ∂φ c i  =    abc v c i  ∂H ∂v b i  (21) where    abc   is the fully antisymmetric Levi-Civita sym- bol, and we exploited the fact that   ∂v a i   /∂φ b i   =   −  abc v c i   . We can rewrite Eqs. (21) in a more compact form using vectorial notation, and finally get the set of dynamical equations of the inertial spin model (deterministic case),  d~ v i  dt   =   − ~ v i   ×   δH δ~ s i  = 1  χ~ s i   ×   ~ v i   (22)  d~ s i  dt   =   − ~ v i   ×   δH δ~ v i  =   ~ v i   ×   ~ F i   (23)  d~ r i  dt   =   ~ v i   ,   (24) where we have also added the kinematic relationship be- tween the velocity   ~ v i   of a particle and its position in space  ~ r i . We note that these equations are formally analogous to what one would get for rotations and circular motion in real space, with the external variables   ~ r   and   ~ l   playing the role of   ~ v   and   ~ s .   In this last case, though, equation (24) would of course not be present.
6 s i  v i   v j  s j  FIG. 2.   Schematic representation of velocity, spin, and tra- jectory for two particles.  In (22-24)   H   is the same Hamiltonian as in (14) (but written in terms of velocities rather than phases) and  ~ F i   =   − δH/δ~ v i   is therefore the force acting on particle   i  due to interactions with other particles.   For alignment interactions we have  H ( ~ v, ~ s ) =   −   J  2 v 2 0  ∑  ij  n ij ~ v i   ·   ~ v j   +   1 2 χ  ∑  i  s 2  i   ,   (25) where   n ij   is the connectivity matrix (being 1 if   j   is a particle interacting with   i   and 0 otherwise), and   J   is the alignment strength. Therefore   ~ F i   = ( J/v 2 0   )   ∑  j   n ij ~ v j   . Note that, contrary to Eqs. (15), Eqs.(22-24) do not have a Hamiltonian form, precisely because ( ~ v i , ~ s i ) are not canonical variables ( ~ s i   is the conjugated moment to the phase   ~ φ i , not to the velocity vector   ~ v i ).   Still, they retain a pseudo-Hamiltonian structure where the deriva- tives of the Hamiltonian are combined with a vectorial product. Thanks to this structure, the speed of the indi- vidual particles is automatically conserved by the dynam- ics. Besides, the equations also conserve   ~ s i   ·   ~ v i , which we fix equal to zero (the only solution in absence of forces). Finally, the Hamiltonian   H   itself is conserved by these equations. We can now explore how Eqs.(22-24) determine the motion of particles in space (see Fig. 2). When the spin  ~ s i   is equal to zero, the direction of motion of particle   i  always remains the same and the particle therefore fol- lows a straight path in real space.   When   ~ s i   is different from zero but constant, the flight direction obeys uniform circular motion and the particle performs a turn in real space with a constant radius of curvature   R   ∼   v 0 χ/ | ~ s i | . Hence, the spin   ~ s i   has a clear kinematic meaning being related to the instantaneous curvature   κ   = 1 /R   of the trajectory.   When there are forces acting on the parti- cle the local spin/curvature   ~ s i   changes in time, inducing a variation in the rate of direction changes.   The paths followed by the particles in real space depend on the in- stantaneous realization of the forces.   However, the im- portant point is that - whatever these forces are - in our new model they act on the spin   ~ s i   and not directly on the velocity,   ~ v i . In other terms, forces cannot change the di- rection of motion abruptly, but there is an inertial effect mediated by   ~ s i . In this respect, we see that the general- ized inertia   χ   measures the resistance of the particle to change the instantaneous radius of curvature of its tra- jectory. In the case of birds,   χ   can be explicitly related to the resistance of the bird to change its banking angle (see [25]).  B.   Spin conservation and information propagation  Another crucial consequence of the presence of a global continuous symmetry is that the generator of the symme- try is a constant of motion of the dynamics. In our case, this implies that the global spin   ~ S   =   ∑  i   ~ s i   is conserved. Interestingly, for a polarized flock this conservation law also has local implications (i.e.   at the level of individ- ual particles): if a strong misalignment occurs between a particle   i   and its neighbors, this causes an excess of curvature and an excitation of the local spin   ~ s i ; since the global spin is a conserved quantity, the local excitation cannot be re-absorbed or dissipated, but must be carried away.   For this reason we expect propagating modes in this system. To better understand how this occurs, let us consider a polarized flock of particles all moving approximately in the same direction. We stick to the simpler planar case for convenience, but the argument is general.   For an ordered flock, the rotational symmetry is spontaneously broken since a specific global direction of motion is se- lected by the system. Still, a residual symmetry remains in the space orthogonal to this direction: we can choose the phase   φ i   of each particle as the angle with respect to the global flock velocity, and use Eqs.(15) to describe the evolution of such phases. Since the flock is ordered, phases are small and the potential can be written as   V   =   J/ 4   ∑  ij   n ij   [ φ i   −   φ j   ] 2 . If phase variations are smooth and we look at distances larger than the typical distance   a   between interacting neighbors, we can consider the continuum limit in space  φ i ( t )   →   φ ( ~ r, t ),   s i ( t )   →   s ( ~ r, t ).   The phase variation between a particle and its neighbors   ∑  j   n ij   [ φ i   −   φ j   ] 2  can be conveniently approximated by   Jn c a 2 [ ∇ φ ] 2 , where  n c   = (1 /N   )   ∑  i  ∑  j   n ij   is the average number of interact- ing neighbors. Then, Eqs.(15) become  ∂φ ∂t   =   s χ   (26)  ∂s ∂t   =   n c a 2 J ∇ 2 φ .   (27) Interestingly, the second of these equations can be writ- ten as a continuity equation for the field   s ( ~ r, t ),  ∂s ∂t   − ∇ ·   ~ j   = 0   ,   (28) with current   ~ j   =   n c a 2 J ~ ∇ φ .   This conservation law de- scribing spin flow is precisely the local counterpart of global spin conservation that we mentioned above. If we now take a second derivative of   φ   with respect to time, we immediately get d’Alembert’s equation,  χ ∂ 2 φ ∂t 2   =   Jn c a 2 ∇ 2 φ .   (29)
7 This equation can be easily solved in Fourier space to get the dispersion relation, namely the law describing how directional information travels through the system.   We find  ω   =   c s k ,   (30) with   c s   =   √ Jn c a 2 /χ .   A real value for the frequency  ω   corresponds to propagating modes:   what we find is therefore that a deterministic flock exhibits undamped propagating modes of the phase. Besides, the dispersion law is linear, meaning that propagating modes travel at a well defined speed   c s , which is a function of the alignment strength and the generalized inertia. Linear propagation is reminiscent of sound propagation in a medium. Here, however, the modes that we are looking at are not related to density fluctuations, but to phase fluctuations. In fact the same equations for the phase we are looking at, and the same dispersion law, would also hold for a fixed net- work of particles, e.g. a regular lattice. Phase propagat- ing modes mean that if in a flock a particle starts turn- ing, this change will affect through the alignment term the spin/curvature of nearby particles, which will start turning themselves.   The information will travel with a speed   c s   which is larger the more ordered the flock, and the whole system will perform a collective turn.   As we shall discuss, this is similar to what happens in natural flocks of birds, indicating that far from being an idealized abstract description of collective motion, the determinis- tic model introduced so far has some crucial ingredients occurring in real systems. The field equations we have written for   φ ( ~ r, t ) and  s ( ~ r, t ) in this Section make one important approxima- tion, namely that the interaction matrix   n ij   is fixed in time.   Without this approximation we cannot use the continuous Laplacian, nor obtain d’Alembert’s equation. This condition is, in general, violated, and this is what marks the difference between a lattice model and an off- lattice model of moving particles, where   n ij   =   n ij   ( t ) de- pends on time. However, if the time scale of propagation in d’Alembert’s equation (29) is shorter than the typi- cal time scale of reshuffling of the network   n ij   , then the equations are justified. This is indeed the case in natural flocks of starlings [25], where phenomena like collective turns happen on a time scale indeed significantly shorter than the time needed for an individual to change its inter- acting neighbours. Nevertheless, at the end of this work we will run numerical simulations of the full model with varying   n ij   ( t ) and show that when   c s   is large enough, the propagation phenomena are quantitatively in agreement with the field equations we have written above.  C.   Dynamical universality class  The relationship between symmetries of the system, equations of motions and conservation laws is a well known feature of Hamiltonian dynamics and a partic- ular manifestation of Noether’s theorem.   The presence of a continuous symmetry implies the existence of a con- served quantity, and the concomitant linear propagation of local disturbances.   What is novel in our analysis is that we have studied these properties in the context of the dynamics of flight directions of flocking particles. Interestingly,   the   same   structure   of   equations   and dispersion law that we have introduced for flocks also hold in other kinds of systems, where a non-conserved order   parameter   (the   flight   direction   in   the   case   of flocks) co-evolve with a conserved quantity (e.g.   the spin/curvature).   In fact, Eqs.   (27)(28) are formally identical and belong to the same dynamical universal- ity class as the planar ferromagnet [36] and superfluid Helium [34, 37] (model F in the Halperin-Hohenberg clas- sification [33]). In liquid Helium, for example, the quan- tum phase plays the role of the flight direction, the (con- served) superfluid component plays the role of the spin,  ~ j   is the superfluid current, and the propagating wave is the so-called ‘second sound’ [34, 35, 37]. In the same way, the full three-dimensional model Eqs (22-24) is equiva- lent to the isotropic antiferromagnet (model G dynamical universality class), with exactly the same equations, the staggered magnetization playing the role of the flight di- rection, and the (conserved) total magnetization playing the role of the spin. In both these cases (model F and model G of [33]), propagating modes and a linear dispersion law are the consequence of the rotational symmetry and of the con- servation of the spin deriving from it.   This is not a generic feature of spin models. For example, in the case of the isotropic ferromagnet, the conserved field is the or- der parameter itself and this conservation does   not   derive from any symmetry of the Hamiltonian. In this case the linear propagation of the phase is lost, and a dissipative, diffusive dispersion law holds [33]. Therefore, the essen- tial ingredient of linear spin propagation is not simply the presence of a conserved quantity, but rather the fact that this conserved quantity must be the generator of a spon- taneously broken continuous symmetry.   This condition is the true hallmark of superfluid transport, namely of the undamped and linear propagation of the spin. In the case of Helium, the spin is mapped onto the superfluid component, which propagates through second sound [33], whereas in the isotropic antiferromagnet one has actual spin superfluid propagation [33].   In flocks, spin carries the turning information.  III.   A NEW MODEL OF COLLECTIVE MOTION  The deterministic inertial spin model we have intro- duced in the previous Section conserves spin (through the continuity equation), but it does so a bit too effec- tively, as in absence of forces the angular velocity remains forever constant. Here we want to cure this anomaly. An example may help: although linear momentum conserva- tion and inertia is essential to propagate standard sound
8 in a crystal, we certainly do not expect sound to propa- gate forever and for infinite distance, as some dissipation will eventually dampen the signal.   The same must be true for second sound, which is carried by spin conser- vation. Therefore, we must introduce a dissipative term coupled to   ~ s . Moreover, the dissipative term will be as- sociated to a noise term, making the model stochastic. It is important to understand that the dissipation we need to introduce is   not   the standard linear one, propor- tional to  ̇ ~ x . Standard linear dissipation would drive the modulus of the velocity to zero, while this cannot happen here: the constraint   | ~ v i |   =   v 0   must be automatically en- forced by our equations.   Instead, spin dissipation must account for the natural and rather intuitive fact that in absence of spin injection and of interaction with neigh- bours, a bird is unlikely to turn forever with constant angular velocity. As we shall see, regulating the strength of this dissipative term will make the difference between propagating and nonpropagating information.  A.   The Inertial Spin Model  As we would do for a physical system, we start from the deterministic structure of the equations (22-24) and add noise and dissipation to the ‘conservative’ forces   ~ F i   =  − δH/δ~ v i . In this way we obtain the final equations of the inertial spin model of collective motion,  d~ v i  dt   =   1  χ~ s i   ×   ~ v i   (31)  d~ s i  dt   =   ~ v i ( t )   ×      J v 2 0  ∑  j  n ij ~ v j   −   η v 2 0  d~ v i  dt   +   ~ ξ i  v 0      (32)  d~ r i  dt   =   ~ v i ( t )   ,   (33) with   ~ v i   · ~ s i   = 0. Here   η   is a generalized viscous coefficient and   ~ ξ i   is an i.i.d vectorial noise with variance  〈 ~ ξ i ( t )   ·   ~ ξ j   ( t ′ ) 〉   = (2 d )   η T δ ij   δ ( t   −   t ′ )   ,   (34) where we have introduced the generalized temperature   T   , in analogy to physical systems.   Note that, consistently with our entire description, the constraint   | ~ v i |   =   v 0   is satisfied by these equations. The inertial spin model (31-33) represents a full de- scription of the system’s dynamics in terms of the vari- ables ( ~ v i , ~ s i ).   Being a first order system it is relatively easy to implement numerically, and we shall discuss re- sults of simulations later in the paper. For the time being, we note that there are four important parameters enter- ing Eqs.(31-33): the alignment strength   J , the moment of inertia   χ   of the spin, the viscous coefficient   η , and the temperature   T   . The spin inertia,   χ , is the new ingredient in the dynamics. To better understand the effect of dissipation, let us consider the case where there are no social forces acting on the birds and the noise is zero (i.e.   F i   = 0,   T   = 0).   We remind that under these conditions, if   η   = 0 a bird conserves its local spin, i.e.   it keeps turning with uniform angular velocity.   When   η   6   = 0, on the other hand, Eqs. (31-32) give,  d~ v i  dt   =   1  χ~ s i   ×   ~ v i   (35)  d~ s i  dt   =   −   η χ~ s i   .   (36) The spin exponentially decays to zero, so that the ve- locity vector stops rotating.   A bird that was originally turning ( ~ s i   6   = 0) reduces its spin/curvature until reach- ing a state of straight motion ( ~ s i   = 0).   Thus   η   acts as a damping on the spin and plays the role of a rotational dissipation.  B.   A closed equation for the velocity  By taking a further time derivative of Eq. (31) and exploiting Eq. (32), one can get a closed second order equation in the velocities,  χ d 2 ~ v i  dt 2   +   χ ~ v i  v 2 0  (   d~ v i  dt  ) 2  +   η d~ v i  dt   = =   J v 2 0    ~ v i   ×   ∑  j  n ij ~ v j      ×   ~ v i   +   v 0   ~ ξ ⊥  i   ,   (37) with  ~ ξ ⊥  i   =   ~ ξ i   −   ( ~ ξ i   ·   ~ v i  v 0  )   ~ v i  v 0  ,   (38) and  〈 ~ ξ ⊥  i   ( t )   ·   ~ ξ ⊥  j   ( t ′ ) 〉   = 2( d   −   1)   η T δ ij   δ ( t   −   t ′ )   .   (39) Equation (37) is a generalized Langevin equation with re- spect to the velocity vectors. The first term in the l.h.s is an inertial second order term, where   χ   is   not   the conven- tional inertia (i.e. the mass of the particle) but measures the resistance of a particle to change the curvature of its trajectory.   The second term at the l.h.s.   is a con- tribution arising from the constraint on the speed, the analogue of what would be the centripetal acceleration for circular motion in   ~ r   space. Its presence is a signature of the underlying Hamiltonian structure of the equations and is therefore essential to recover the appropriate de- terministic limit when   η   →   0. The last term in the l.h.s. is the dissipative term. In the r.h.s. we find the effect of the social force and the noise. As compared to other second order equations that have been considered in the literature [38–42], in Eq.(37) the connection between inertial terms and symmetries is au- tomatically implemented. We stress that this is in fact a crucial point, since - as we showed in the previous sections - the way information propagates is strongly influenced by this connection.
9  C.   Overdamped limit: the Vicsek model  Interestingly, the original Vicsek model happens to be the overdamped limit of the inertial spin model, that is the limit in which the spin inertia becomes negligible compared to the spin dissipation.   Equation (37) is the full Langevin equation for the velocity, including both inertial terms (proportional to   χ ) and a dissipative term (proportional to   η ).   As usual in the Langevin equation [43, 44], the overdamped limit is obtained by taking the limit   χ/η 2   →   0, which in our case gives  η d~ v i  dt   =   J v 2 0    ~ v i   ×   ∑  j  n ij ~ v j      ×   ~ v i   +   v 0   ~ ξ ⊥  i   .   (40) Note that a rescaling of time,   t   →   t ′   =   t/η , would allow us to get rid of   η   both in the equation for   ~ v   and in the noise correlator (this is the reason why the overdamped limit is   χ/η 2   →   0, rather than   χ/η   →   0).   This time rescaling is usually taken for granted in the overdamped limit and for this reason   η   is normally set equal to 1 in this limit. Equation (40) is identical to the Vicsek model (in the continuous time limit version). The double cross product in the r.h.s.   has a clear interpretation.   Calling   ~ F i   = ( J/v 2 0   )   ∑  j   n ij ~ v j   the social force acting on particle   i , we can exploit the properties of the cross product and write ( ~ v i   ×   ~ F i )   ×   ~ v i   =   v 2 0   [   ~ F i   −   (   ~ F i   ·   ~ v i ) ~ v i /v 2 0   ] =   v 2 0   F   ⊥  i   =   J ( ∑  j  n ij ~ v j   ) ⊥   .   (41) Only the perpendicular component of the social force contributes to changing the velocity vector, as it should be since the speed is constant due to the constraint. Thus, the double cross product is a formal and clean way to ensure that the norm of   ~ v i   is conserved during the dynamics.  D.   Steady state distribution of velocity and spin  We have seen that the dynamical equations of the in- ertial spin model have their deep roots in the existence of a symmetry, and its associated conservation law. In fact, the rotational symmetry and its spontaneous breaking also have implications when looking at the steady state distribution of velocity and spin. In general, one might wonder whether a steady state distribution exists at all in flocks, since the system is completely out of equilib- rium.   In fact many of the most intriguing features of active systems precisely come from their off-equilibrium nature [21].   Still, if we look at a flock on timescales smaller than the typical swapping time between an indi- vidual an its neighbors, then we can consider the veloci- ties and spins as dynamical variables evolving on a fixed interaction network.   In this case we can show that the joint distribution probability for velocities and general- ized momenta obeys a Fokker-Planck equation [43, 44]. The steady state solution of such an equation reads  P   ( ~ v, ~ s )   ∼   exp  [  −   1  T H ( ~ v, ~ s )  ]  =   (42) exp    −   1  T    −   J  2 v 2 0  ∑  ij  n ij ~ v i   ·   ~ v j   +   ∑  i  s 2  i  2 χ        .  As usual in statistical physics, this probability distribu- tion factorizes between coordinates and momenta, that is velocities and spins. This result tells us that velocities are in fact amenable to a statistical description (on scales where network rearrangements can be disregarded):   a measure exists, which can be used to compute expec- tation values. This distribution is a Boltzmann distribu- tion, where the ratio between noise and dissipation - what we called temperature - in fact acts as a real tempera- ture in the statistical sense.   Additionally, inertia does not enter in the velocity distribution.   For this reason, all the static properties of our model (including the ve- locity correlation functions and the steady-state ordering properties), are the same as the Vicsek model, which is described by the very same configurational part. If we now look at the explicit form of the Hamilto- nian in (42), we see that the velocities are described by a Heisenberg model on a random euclidean lattice (the one realized by the particles positions). In the ordered phase, where the rotational symmetry is spontaneously broken, this distribution generates soft Goldstone modes and long range correlations corresponding to fluctuations of the in- dividual velocities perpendicular to the global velocity of the flock [45]. Long-range correlations in velocity fluctua- tions have been experimentally observed in natural flocks of birds [27]. Besides, in [26] using inference techniques on empirical data, we showed that the equal time prop- erties of the velocity fluctuations can be fully described by a Boltzmann distribution with local pairwise align- ment interactions. The model we have discussed so far is therefore fully consistent with all experimentally known dynamic   and   static features of real flocks of birds.  IV.   INFORMATION PROPAGATION IN THE ORDERED PHASE  Equations (31-33) give rise at low noise to polarized flocks, much as in the Vicsek model. In Section II we have seen that for a deterministic flock the rotational sym- metry and the concomitant conservation law have strong implications on the properties of the ordered phase. Now we would like to understand whether this scenario holds in presence of noise and dissipation. Let us consider the system in the strongly polarized phase.   In this case, the individual velocities are very aligned to each other and the flock as a whole has a non- zero collective direction of motion   ~ V   . To fix ideas, let us
10 assume that the group velocity is along the   x   direction,  ~ V   =   V ~ n x , with   ~ n z   = (1 ,   0 ,   0). The individual velocities can be conveniently written in terms of a longitudinal component along the direction of motion,   ~ n x , and a per- pendicular component, i.e.   a two dimensional vector   ~ π  in the two-dimensional ( y, z ) plane,  ~ v i   =   v L i   ~ n x   +   ~ π i   .   (43) Because of the large polarization all individual velocities will be very close to   ~ V   , so that   π 2  i      1 and  v L i   =  √  v 2 0   −   π 2  i   ∼   v 0 [1   −   π 2  i   / (2 v 2 0   )]   .   (44) Using the phases defined in the previous sections, we can write  π y   =   v 0   sin( φ z   )   ∼   v 0   φ z   ,   (45)  π z   =   v 0   sin( φ y   )   ∼   v 0   φ y   .   (46) However confusing these relations may seem, they are correct: a nonzero   y   component of   ~ π   is obtained by ro- tating   ~ v   around the   z   axis, a rotation parametrized by  φ z   .   Similarly, a nonzero   π z   is generated by a rotation parametrized by   φ y   . The crucial point is that, once we have spontaneous symmetry breaking, only   two   relevant phases are left. By using (44-46) into equation (37) and expanding at the first order (spin wave expansion), we obtain for each one of the two phases,   φ z   and   φ y   , the same equation,  χ d 2 φ i  dt 2   =   J   ∑  j  Λ ij   φ j   −   η dφ i  dt   +   ξ ⊥  i   ,   (47) where Λ ij   =   − n ij   +   δ ij  ∑  k   n ik   is the discrete Laplacian. In general, solving equation (47) is not an easy task, because the dynamics of the   ~ φ ’s is coupled in a non-trivial way with the movement of the particles in space, i.e. the network   n ij   changes with time. Hydrodynamic the- ories of flocking address this problem by coarse-graining the microscopic dynamics and looking at the large-scale long-time behavior.   Here, however, we are interested in addressing the dynamics of the system on shorter timescales.   We can therefore assume that   n ij   does not change significantly, as indeed happens in flocks during collective turns [25]. Later in Sec. V we perform numer- ical simulations of model (1) where network rearrange- ments on the scale of turns are fully taken into account. For a fixed interaction network, the equation for   φ   can be solved analytically to obtain the dispersion relation. To do this, we must choose a diagonal representation both in space and time, i.e.   we must diagonalize the Laplacian Λ, write the equations in terms of its eigen- modes, and Fourier transform with respect to time. This can be done exactly numerically. Yet a simple analytical approximation gives us some intuition about the nature of the dispersion relation. To do this we can proceed as in the previous section: if we look at spatial scales larger than the nearest-neighbor distances, we can approximate the discrete Laplacian with its continuous counterpart, i.e.   J   ∑  j   Λ ij   →   Jn c a 2 ∇ 2   (where, as before,   a   is the typical distance of the interacting neighbors and   n c   their number). In this case we can write  χ ∂ 2 φ ∂t 2   =   Jn c a 2 ∇ 2 φ   −   η ∂φ ∂t   +   ξ ⊥   .   (48) The propagator (Green function) of this differential equa- tion can be easily computed in Fourier space giving the following dispersion law  χω 2   −   iηω   −   Jn c a 2 k 2   = 0   .   (49) This dispersion law tells us how local disturbances in the flight direction of a particle affect the rest of the flock. If a bird/particle changes its flight direction, will this perturbation propagate through the entire group or be damped?   The answer to this question depends on the values of the parameters entering the equations, and in particular on the balance between dissipation and inertia. In the deterministic limit   η   →   0 we get the linear dis- persion law   ω   =   c s k   discussed in the previous sections with  c s   =  √  Jn c a 2  χ   .   (50) This is the perfectly undamped case of propagating ori- entational modes. At the opposite end of the spectrum there is the over-damped limit,   χ/η 2   →   0, which, as we have seen, corresponds to the Vicsek model. In this case  ω   is purely imaginary with a quadratic diffusive disper- sion law,   ω   =   i ( Jn c a 2 /η ) k 2 .   We conclude that in the Vicsek case there is no propagation, but pure exponen- tial damping. We will see through numerical simulations that this is indeed what happens: in the Vicsek model it is not possible to locally initiate a turn that propagates to the entire flock. In the general case, where both dissipation and inertia are different from zero, we obtain  ω   =   i/τ   ±   ω 0  √  1   −   k 2 0   /k 2   ,   (51) with  ω 0   ≡   c s k   ,   k 0   ≡   η  2 √ Jn c a 2 χ   ,   τ   ≡   2 χ/η .   (52) Here   ω 0   is the zero dissipation frequency, while   k 0   and   τ  are the two relevant scales, respectively in wave number and time, related to the effect of dissipation. With zero dissipation, we get   k 0   = 0,   τ   =   ∞   and   ω   =   ω 0 : the time scale needed for dissipation to have an effect is infinite and linear propagation occurs on all spatial scales.   For  η   6   = 0, on the other hand, we have two regimes, accord- ing to the value of the friction coefficient and of the wave number   k .   For   k   ≥   k 0   we have   attenuated   propagating waves, as the frequency has both a real and an imaginary
11 part.   For   k < k 0   we have   evanescent   waves with expo- nential decay and we recover the overdamped Vicsek-like behavior. Consistently with the above analysis we see that the large scale   k   →   0 behavior of a system of self-propelled particles is well accounted for by the Vicsek overdamped case and is therefore well described by the hydrodynamic theories of flocking.   There is however a crucial ingredi- ent to be considered when dealing with natural or even artificial groups: flocks are finite. The smallest value of  k   in the system is   k min   ∼   1 /L , where   L   is the linear size of the flock. Hence, even if dissipation is present, if it is sufficiently small it does not affect the scales relevant for information to travel through the group. More precisely, if,  η <  √ Jn c χa L   :   underdamped limit   ,   (53) then   there   is   linear   propagation   of   the   information throughout the whole flock.   Indeed, for these values of   η   the time scale of the exponential decay is   τ   > √ χ/ ( Jn c a 2 )   L   =   L/c s . Therefore, small dissipation im- plies that the damping time constant is larger than the time the information takes to travel across the flock. In other words, the signal is effectively very weakly damped across the length scale of interest. We conclude that even when a small dissipation is present, propagation of infor- mation is qualitatively the same as that described by the zero dissipation theory. Natural flocks of birds appear to be precisely in this underdamped, quasi-deterministic regime [25].   Exper- iments on large flocks of starlings performing collective turns show that, once the turn (i.e. a direction change) is initiated by some individuals at the edge of the flock, the turning front propagates through the group with a linear dispersion law with very limited attenuation.   The rela- tionship (50) between the experimentally measured speed of information propagation   c s   and the alignment strength  J , predicted by the model in the deterministic limit, is also remarkably verified [25]. Thus, far from being an ide- alized construction, the deterministic model introduced previously truly captures some important mechanisms occurring in real biological groups.  A.   Transition to the hydrodynamic regime  The natural flocks investigated in [25] live in the quasi- deterministic regime   k > k 0 . Other systems, however, (or even very large flocks) might well explore larger length- scales. This is certainly the case for a variety of fluids or large assemblies of active individuals.   As mentioned in the previous section, hydrodynamic theories of flocking offer a sophisticated and detailed description of the   k   →   0 regime, which has been investigated in numerical simula- tions and experiments on granular active matter. While we do not want to give a complete account of this descrip- tion (we refer the reader to [18]) we briefly comment on what happens as these large length (and time) scales are approached. On the one hand, as we noted in the previ- ous section, when   k   →   0 the frequency Eq. (51) becomes purely imaginary indicating that there are not any more purely directional propagating modes. The physical rea- son is that for times scales   t      τ   = 2 χ/η   the global spin relaxes due to dissipation and does not act as a quasi- conserved variable, thus the conservation law that deter- mined the propagating nature of the modes is lost.   In other terms, the second order terms of the microscopic dynamics become irrelevant, as well as the coupling be- tween direction and curvature.   On the other hand, on such large time scales the movement of individuals be- comes relevant, and crucial:   orientational fluctuations (which would be evanescent on a fixed network) couple with fluctuations in density. As a consequence, there are propagating hybridized longitudinal sound modes trans- porting both directional and density disturbances.   The speed of propagation of such modes usually depends on the average density (as in standard sound), contrary to the purely directional propagating modes of the under- damped regime, where the propagation speed   c s   does not depend of the flock density and is determined by the global alignment strength   J   (i.e.   the degree of po- lar order). Finally, if the interaction between the moving particles and the surrounding medium becomes impor- tant the hydrodynamic theory also predicts shear-bend modes [18, 46] that propagate directional changes. These modes are transverse (i.e. they do not exist in   d   = 2) and have limited speed in a co-moving reference frame, and are therefore well distinct from the turning modes ob- served in flocks, which are longitudinal and have very fast propagation speeds [25]. Summarizing, in the quasi-deterministic regime   k > k 0 , where the spin is a quasi-conserved variable, there are propagating modes thanks to the coupling between velocity and spin and the consequent conservation law; these turning modes transport changes of direction   and  curvature and describe collective turns in flocks. In the hydrodynamic regime, due to the movement of the net- work, changes in direction couple with changes in local density generating propagating sound modes; this regime is universal, i.e.   independent of the details of the mi- croscopic dynamics, and does not describe transport of curvature.  V.   NUMERICAL SIMULATIONS  In the previous sections we have used analytical argu- ments to derive the predictions of the inertial spin model on how orientational information propagates in the or- dered phase. We did that, however, working under two assumptions:   we disregarded effects due to the reshuf- fling of the interaction network, and we considered the continuous space limit. While these assumptions are rea- sonable in certain cases, we want to test the behavior of the model in full generality.   We performed numerical
12 simulations for a wide range of the relevant parameters  χ ,   η , and   J   and explored how the values of these param- eters determine the collective dynamics of the system. We implemented a time discretized version of the in- ertial spin model,  ~ v i ( t   +   dt ) =   ~ v i ( t ) + 1  χ~ s i ( t )   ×   ~ v i ( t )   dt   (54)  ~ s i ( t   +   dt ) =   ~ v i ( t ) +   ~ v i ( t )   ×   J v 2 0  ∑  j  n ij ~ v j   ( t )   dt   (55)  −   η χ s i ( t )   dt   +   ~ v i ( t )   ×   ~ ξ i ( t )  v 0  √ dt ~ r i ( t   +   dt ) =   ~ r i ( t ) +   ~ v i ( t )   dt   (56) where we exploited the equation for   ~ v i   in order to replace  d~ v i /dt   with   ~ s i , and where we recall   ~ v i   ·   ~ s i   = 0. The noise variance is  〈 ~ ξ i ( t )   ·   ~ ξ j   ( t ′ ) 〉   = (2 d )   η T δ ij   δ t,t ′   .   (57) Our objective is twofold: on the one hand we want to check the model’s behavior in the ordered phase in the various possible regimes; on the other hand we want to compare it to results obtained in natural flocks of birds [25]. To this end, we chose the connectivity network be- tween particles according to a topological rule, as in real flocks [26, 27]:   each particle interacts with the first   n c  neighbors so that   n ij   = 1 if   j   is among the first   n c   neigh- bors of   i , and zero otherwise. We focus on finite groups of particles with open boundary conditions, but quali- tatively similar results also hold for periodic boundary conditions. We are interested in understanding how, for different values of   η ,   χ   and   J , directional information propagates in the system, and whether the simulated flock is able to sustain coherent collective turns as natural flocks do. In real flocks collective turns start locally:   an individ- ual bird, the initiator, starts to turn and then the in- formation about turning propagates through the system with a linear dispersion law [25].   To mimic this situ- ation, we consider the system at low enough values of noise as to generate an ordered flock with large polar- ization (Φ   >   0 . 95 as in natural flocks), and we initialize all the   ~ s i   = 0 so that the flock as a whole performs a straight motion.   Then, we chose at random a particle inside the flock and make it artificially turn changing its flight direction with a ramp-like time dependence (red line in Fig. 3d,e,f).   The change in direction of the ini- tiator affects - through Eqs. (56) - the dynamics of the spins of nearby particles, which might or might not cause a collective movement. We then look at the behavior of the whole flock in time to assess whether the signal has propagated, and what kind of dispersion law is obeyed. According to the analytical arguments given in the pre- vious section, we expect two different regimes, according to the values of   η   and   χ : 1)   Overdamped regime:   η 2 /χ > n c J ( a/L ) 2 .   In this regime, given a system of size   L , some attenu- ated propagation occurs up to certain spatial scales ( k   ≥   k 0 , with   k 0   = 1 / (2 a √ n c Jχ )   >   1 /L ).   On larger scales ( k < k 0 ) however, dissipation takes over leading to an exponential decay of the signal. The extreme case occurs for   η 2 /χ   → ∞   (or   χ/η 2   →  0), corresponding to the Vicsek model, when propa- gation of orientational perturbations does not occur (unless we reach hydrodynamic lenght-scales where density and orientational modes couple, see [18, 22] and Sec. IV A). 2)   Underdamped regime:   η 2 /χ      n c J ( a/L ) 2 .   In this regime there is linear propagation of the signal throughout the whole system with negligible atten- uation.   The speed of information propagation is determined solely by the ratio of alignment to in- ertia,   c s   =   a √ n c J/χ . We performed simulations for values of the parameters corresponding to these two different regimes. As can be seen from Fig. 3, what we observe is fully consistent with what is described above.   The top panels in this figure display the trajectories of the individual particles in the flock for three decreasing values of   η 2 /χ , the trajectory of the initiator being represented by the thick black line. The intermediate panels show, for the same events, the cosine of the velocity of each particle with the original flight direction of the flock.   Here the initiator is repre- sented as a red line, and one can see the ramp-like depen- dence of the perturbation.   The lower panels reproduce the individual acceleration profiles, and are particularly useful to determine if and how particles turn, and how the turning signal is attenuated during time. The first two cases in Fig. 3 (panels a,b) belong to the overdamped regime. For values of   η 2 /χ   much larger than  n c J ( a/L ) 2   (which is of order 1 for our simulations) the turn of the initiator does not trigger any collective turn, since the signal does not propagate and is quickly damped (panel a).   We can see from the cosine curves (panel d) that some of the neighbors feel the disturbance caused by the initiator, but this effect is very small: no parti- cle follows the initiator, no collective turn is triggered and the flock ultimately keeps its original flight direction (see also attached Video 1). For smaller values of   η 2 /χ  (panel b) there is some propagation, but it is progres- sively attenuated as the signal travels through the group and eventually decays (panels e,h).   Some particles are reached by the directional perturbation before the signal is damped (corresponding to scales   k > k 0 ), and follow the initiator changing their flight direction and perform- ing a turn.   Most of the particles, however, are not sig- nificantly affected by the perturbation, consistently with a strong damping of the signal on large enough scales. This effect is clearly seen by looking at the acceleration profiles of the particles (Fig. 3h). Some particles display a profile similar to the one of the initiator with a clear peak marking the occurrence of a turn.   The height of the peak, however, decreases as the turning information is transmitted (i.e. particles start turning one after the other) and at some point it abruptly drops to zero (see
13  FIG. 3. Information propagation in different regimes.   a.   Strongly overdamped regime,   η 2 /χ   = 28 . 8   ×   10 2 . The signal is quickly dissipated and there is no collective turn. All the particles keep following the initial direction of motion without following the initiator. The trajectory of the initiator is displayed as a thick black line.   η   = 60,   χ   = 1 . 25,   J   = 0 . 8.   b.   Overdamped regime,  η 2 /χ   = 1 . 8   ×   10 2 . Some propagation occurs, but the signal is strongly attenuated before reaching throughout the whole group. As a consequence, the flock looses cohesion and coherence while turning.   η   = 15,   χ   = 1 . 25,   J   = 0 . 8.   c.   Underdamped regime,  η 2 /χ   = 7 . 2   ×   10 − 2 . The signal is propagated unattenuated through the whole flock, which performs a neat turn retaining shape and cohesion.   η   = 0 . 3,   χ   = 1 . 25,   J   = 0 . 8.   d ,   e ,   f   Cosine of the individual velocities with the original flight direction of the flock for the three cases displayed in panels a,b,c. The cosine curve of the initiator is displayed as a red curve: at time   t   =   t 0   the flight direction of the initiator is tuned from its original direction (coherent with the flock motion) to a final direction following a ramp-like pattern. The inset in panel   f   shows the threshold used to compute the ranking of the particles (see text).   g ,   h ,   i  Individual acceleration profiles for the three cases displayed in panels a,b,c. The acceleration curves have been smoothed with a low-pass filter to cure noise and high frequency oscillations. The insets display the intensity of the peak as a function of time. Particles are ranked according to the order of reaction to the initiator (see text). The other parameters of the simulations are  N   = 512,   T   = 8   ×   10 − 5 ,   n c   = 6,   v 0   = 0 . 1. The integration time is chosen as   dt   = 0 . 1 √ J/χ   to ensure proper simulation time for all values of   J   and   χ .  inset). In this case, the flock breaks into two subgroups: the largest component keeps the original flight direction, while a smaller fragment is dragged by the initiator (see attached Video 2). In this respect we note that, due to active nature of the particles, the effect of overdamping is even stronger than discussed in sec. IV: not only the turning information fails to percolate through the group, but the flock breaks up and looses cohesion. In the underdamped regime (panel c), on the con- trary, there is robust propagation of the signal through the whole group:   all the particles follow the initiator changing their flight direction with negligible attenua- tion. Both the cosine curves and the acceleration profiles show that the particles turn in a highly coordinated way: the curves are very similar one to the other with a shift in time corresponding to the propagation of the turn from the initiator to the last particle (Fig. 3f,i). The flock per- forms a neat collective turn retaining shape and cohesion during motion (Fig. 3c, and Video 3). In this regime, we can quantify the dispersion law using the method introduced in [25]. We rank all the particles according to their order of turning: the initiator has rank  r   = 0, then the particle that first starts turning after the initiator has rank   r   = 1 and so on.   To determine the rank of a particle we follow two different procedures. We look at the cosine curve of the particle and determine its absolute turning delay from the initiator (and therefore its rank) as the time when the cosine reaches a threshold value (we choose 0.9 - see inset in Fig. 3f). Alternatively, we proceed as in [25]:   from the acceleration curves we
14  FIG. 4. Propagation curve and dispersion law.   a.   Propagation curve. Distance   x   traveled by the turning front vs time, for three different values of the parameters   η ,   χ ,   J   in the small dissipation regime. The distance   x   traveled in a time   t   is proportional to the radius of the sphere containing the first   r ( t ) birds in the rank, namely   x ( t ) = [ r ( t ) /ρ ] 1 / 3 .   t   = 0 corresponds to the time when the initiator starts the turn. The three curves correspond respectively to   η   = 0 . 3,   χ   = 0 . 83,   J   = 1 . 2 (red circles);   η   = 0 . 3,  χ   = 1 . 25,   J   = 0 . 8 (green squares);   η   = 0 . 3,   χ   = 2 . 50,   J   = 0 . 4 (blue diamonds). The rank is computed using the cosine curves (see text). The speed of propagation,   c s , is the slope of the propagation curve in linear regime. The colored straight lines show the linear fits for the three different curves.   b.   Speed of propagation   c s , normalized by the typical distance   a   of the interacting neighbors, as a function of   √ J/χ .   For each value of   J/χ   we run several simulations and estimate   c s   from the slope of the propagation curve, the point corresponds to the average value and bars to standard error. The speed of propagation depends on the ratio   J/χ   as predicted by the analytic arguments in the previous section:   more ordered flocks transmit the turning information quicker.   Inset:   same curve as in main panel, but using accelerations curves to compute ranks and propagation curves.   b.   Speed of propagation   c s , normalized by   a , as a function of 1 / √ χ (1   −   Φ). Inset: same curve as in main panel, but using accelerations curves to compute ranks and propagation curves.  compute the relative time delays between pairs of parti- cles as the shift necessary to superimpose their curves. Once the relative delays are found for all particles we use a standard algorithm to rank them all.   This method is more robust when it is not known a priori who is the ini- tiator (as in experiments [25]). The two procedures give equivalent results. If we now look at the rank as a function of delay time we get the ranking curve, describing how the turn is ex- ecuted through the group.   Since the turn starts locally and the flock moves in three dimensions,   r 1 / 3 ( t ) is a mea- sure of the distance   x   travelled by the turning wave in time   t . As can be seen from Fig. 4a, the distance vs. time curve displays a clear linear dependence corresponding to a linear dispersion law, as predicted by the analytic argument in previous section. We can also check quanti- tatively the prediction for the propagation speed  c s   =   a  √  n c J χ   (58) To do so, we run numerical simulations for several values of the parameters in the underdamped regime.   We fix  η   = 0 . 3 and vary   J   and   χ   so as to obtain propagation curves with different slopes, but keep   Jχ   constant in or- der to have the same value of   k 0   <   1 /L . Then, for each simulation, we compute the distance vs. time curve as in Fig.   4a and we estimate   c s   from a fit of the linear part of the curve. As can be seen from Fig. 4b, the dependence of   c s   on   J/χ   is very well reproduced by numerical data. An alternative way to check the speed dependance - which is particularly useful when comparing to real data - is to plot it as a function of the polarization. We know from Section III D that the statistical properties of the velocities are described through the Boltzmann measure Eq. (42). If we use this distribution, we can compute the expected value of the polarization in the ordered phase. This gives Φ = 1   −   T J − 1 TrΛ − 1   (see [26] for the details of the computation).   We can therefore eliminate   J   in favour of Φ, and obtain,  c s   =   a  √  n c   T   (TrΛ − 1 ) (1   −   Φ) χ   .   (59) The parameters   T   and   n c   are kept constant in our simu- lations; moreover, the trace of the inverse discrete Lapla- cian Λ ij   fluctuates little from sample to sample. There- fore, to check equation (59) we plot the speed propa- gation   c s /a   as a function of 1 / √ (1   −   Φ) χ   (see Fig. 4c). Also in this case the numerical data agree very well with the analytical prediction. The behavior displayed in Fig. 4 is analogous to what is observed in natural flocks of birds.   The propagation curves look indeed very similar to the ones computed for real flocks and shown in [25].   Even more importantly, Fig. 4c is remarkably similar to what found in real data, and the model correctly reproduces the experimentally found relationship between propagation speed and de- gree of order in the system.   We can therefore conclude that Eqs.   (1) in the underdamped limit fully describe the correlated collective turning exhibited by real flocks. We remind that in this limit inertia and deterministic effects dominate over dissipation: it is the Hamiltonian structure of the dynamical equations, and the connected conservation law, that cause the linear dispersion law and the highly non-trivial relationship between speed of prop- agation and degree of order in the system.
15  VI.   CONCLUSIONS  The main point of our work is that momentum mat- ters.   More precisely, the powerful mathematical entan- glement between symmetry, momentum and conserva- tion, is the crucial ingredient giving rise to a dispersion relation, namely to a propagation law, identical to the one observed in real flocks of birds during global changes of direction. If momentum is disregarded, the dispersion relation changes dramatically and linear propagation col- lapses: the practical result is that no collective change of direction can be achieved by the group. But what momentum is that? We hope we have con- vinced the reader that the momentum essential to our new model is   not   linear momentum,   nor   orbital angular momentum.   Accordingly, the generalized inertia in our equation is not the mass, nor the moment of inertia. Lin- ear momentum conservation is clearly violated in active matter systems: this is because the particles are, indeed, active.   Momentum (and energy) is injected and dissi- pated continuously at the individual level. For this rea- son describing active matter systems in terms of Hamil- tonian dynamics in the ( x, p ) variables is impossible. But once we wrap all the complicated active mechanism into an effective constraint for the velocity, we can define an effective angular variable, the phase, which parametrizes the rotation of the velocity and satisfies the constraint. The relevant momentum then is the generator of this in- ternal rotation, that is the spin.   The resistance of the particle to changes in its spin is the turning inertia. This set of variables, ( φ, s, χ ), has nothing to do with the more familiar set, ( θ, l, I ), of polar angle, orbital angular mo- mentum and moment of inertia. This unusual definition of momentum is the precise reason why the propagation phenomenon described by the inertial spin model is so akin to superfluidity.   One could object that linear propagation of information is in fact the most mundane thing of all, sound waves be- ing its most obvious example.   The point is that sound waves transport density fluctuations, the symmetry in- volved is translation, and the relevant momentum is   p . On the other hand, second sound waves, those associated to superfluidity, transport phase fluctuations, the sym- metry is internal rotation, and momentum is   s , the spin. This transport of phase and spin   is   typical of superfluid systems that transport quantum phase.   However weird this connection between flocks and quantum systems may seem, it is mathematically exact. What is surprising is that the new phenomenology of linear propagation is obtained at the rather low price of introducing just one more parameter to the standard Vic- sek model, that is the spin inertia   χ . Moreover, because of the standard Hamiltonian decoupling between coordi- nates and momenta, the static (short time) properties of the inertial spin model are exactly the same as the Vicsek model, so that we do not have to worry about how to re- cover the rich physical description of static correlations previously achieved.   In the overdamped limit momen- tum and inertia become irrelevant and ones recovers the original Vicsek dynamics. Finally, the large length-scale behavior is described by the standard hydrodynamics of flocking [18]. Hence, the inertial spin model introduced here achieves quite a lot in terms of new and experi- mentally accurate phenomenoogy, at a very little cost in terms of modeling complexity. At the more concrete numerical level, we have clearly shown that the inertial spin model describes turns accu- rately. A change of direction can successfully be initiated by just   one   bird and the turn propagates very efficiently to the rest of the flock. The speed of propagation is re- lated to the strength of the alignment interaction and to the turning inertia exactly in the way predicted by our mathematical equations and observed in real exper- iments [25]. On the contrary, if spin conservation is ne- glected (overdamped or Vicsek case) no collective turn is achieved. Interestingly, due to the active nature of flock- ing, the effect of overdamping is even more dramatic than predicted in the fixed-lattice case: instead of a diffusive and damped propagation of the signal, one observes no propagation at all, as the failure to transport the turn- ing information ultimately produces the break up of the flock. An interesting question is what happens when consid- ering very large systems and time scales.   This is the realm addressed by the hydrodynamic theories of flock- ing [12, 18, 20]. On these scales the motion of the network cannot be disregarded and a non-trivial hybridization oc- curs between density modes (caused by the relative move- ment of the particles) and orientational modes (due to alignment between velocities).   As we already discussed in the paper, if some rotational dissipation is present in the system (i.e.   η   6   = 0), the asymptotic long wavelength behavior of our model is the same as the overdamped Vicsek limit, and we therefore expect the predictions of current hydrodynamic theories to describe it appropri- ately. One can however wonder whether a hydrodynamic description exists that also accounts for the underdamped behavior occurring at smaller scales, where many relevant phenomena like collective turns in flocks take place. To do so, one would need to modify the hydrodynamic equa- tions to include a weakly damped spin field. This opens new interesting perspectives at the theoretical level. What we have done here and in Reference [25] is to push one step further the ambitious program to tame the vast richness of biological phenomena using the pow- erful conceptual framework of theoretical physics.   This program is of course not new. Kinetic approaches [49, 50] and hydrodynamic theories of active matter [18, 20, 21] used tools as kinetic theory,   Navier-Stokes equations and dynamical renormalization group to determine the long distance and large time scaling properties of flocks and to explain why the Mermin-Wagner theorem is vi- olated.   At the static level, the maximum entropy ap- proach [26, 51, 52] endeavors to define a thermodynamic description of biological collective systems based exclu- sively on the experimental data. In this approach, theo-
16 retical physics tools, as the Goldstone theorem and scal- ing relations, have powerful consequences [27]. Here, however, we do something different.   We bring into the arena the abstract power of symmetries and con- servation laws in linking static properties to dynamical predictions at the level of the microscopic description of the system.   Symmetry is not necessarily a dynamical concept: it is a very general and profound condition that a system may or may not enjoy. However, the conserva- tion laws predicted by symmetry through Noether’s the- orem and Hamilton’s description have deep dynamical implications, which, as we have seen, can be experimen- tally verified in natural system and in numerical simula- tions alike. To do this the right set of   effective   variables must be found - a recurrent and quite natural theme in physics. We believe that this kind of approach may have some generality in biological systems other than flocks.  Acknowledgments.   We thank William Bialek, Ser- ena Bradde, Paul Chaikin and Dov Levine for discus- sions.   Work in Rome was supported by grants IIT– Seed Artswarm, ERC–StG n.257126 and US-AFOSR - FA95501010250 (through the University of Maryland). Work in Paris was supported by grant ERC–StG n. 306312.  [1] Camazine, S., Deneubourg, J.-L., Franks, N. R., Sneyd, J., Theraulaz, G., and Bonabeau, E.   Self-Organization in Biological Systems   (Princeton University Press, Prince- ton, 2001) [2] Self-organization and collective behavior in vertebrates Couzin, I. D. Krause, J.   Adv Study Behavi   32 :   1–75 (2003.) [3] Giardina I. Collective behavior in animal groups:   the- oretical models and empirical studies   HFSP Journal   2 : 205–219 (2008) [4] Sumpter DJT. Collective Animal Behavior (Princeton University Press, Princeton 2010) [5] Cavagna A, Giardina I, Bird Flocks as Condensed Mat- ter,   Ann. Rev. Cond. Matt. Phys.   DOI: 10.1146/annurev- conmatphys-031113-133834 (2014). [6] Aoki I., A simulation study on the schooling mechanism in fish.   Bull. Jpn. Soc. Sci. Fish.   48 : 1081-1088 (1982) [7] Reynolds CW. 1987.   Computer Graphics   21 : 25-33. [8] Huth, A., & Wissel, C. The Simulation of the Movement of Fish Schools.   J. Theor. Biol.   156 , 365–385 (1992). [9] Couzin, I. D., Krause, J., James, R., Ruxton, G. D., and Franks, N. R. Collective memory and spatial sorting in animal groups   J. Theor. Biol.   218 , 1–11 (2002) [10] Hildenbrandt, H., Carere, C., and Hemelrijk, C. Self- organized aerial displays of thousands of starlings:   a model   Behav. Ecol.   21 : 1349–1359. (2010) [11] Vicsek, T., Czir ́ ok, A., Ben-Jacob, E., Cohen, I., and Shochet, O. Novel type of phase transition in a system of self-driven particles   Phys. Rev. Lett.   75 , 1226–1229 (1995). [12] Toner,   J.   and   Tu,   Y.   Long-range   order   in   a   two- dimensional dynamical XY model: how birds fly together  Phys. Rev. Lett.   75 : 4326–4329 (1995) [13] G Gr ́ egoire, H Chat ́ e & Y Tu, Moving and staying to- gether without a leader Physica D   181 : 157-170 (2003) [14] Gr ́ egoire, G. and Chat ́ e, H. Onset of collective and cohe- sive motion   Phys. Rev. Lett.   92 , 025702 (2004) [15] DOrsogna MR, Chuang YL, Bertozzi AL, and Chayes LS. Self-propelled particles with soft-core interactions: patterns, stability, and collapse Phys Rev Lett 96 , 104302. (2006) [16] Aldana A, Dossetti V., Huepe C, Kenkre VM, Larralde H. . Phase transitions in systems of self-propelled agents and related network models   Phys. Rev. Lett.   98 : 095702 (2007) [17] Ginelli, F. & Chat ́ e, H. Relevance of metric-free inter- actions in flocking phenomena. Phys. Rev. Lett.   105 : 168103 (2010) [18] Toner, J. and Tu, Y. Flocks, herds, and schools: A quan- titative theory of flocking   Phys. Rev. E   58 , 4828–4858 (1998) [19] Vicsek T, Zafeiris A. Collective motion   Physics Reports  517 : 71140 (2012) [20] Ramaswamy S. The Mechanics and Statistics of Active Matter   Ann. Rev. Cond. Matt. Phys.   1, 301 (2010). [21] Marchetti MC, Joanny JF, Ramaswamy S, Liverpool TP, Prost J et al. Hydrodynamics of soft active matter   Rev. Mod. Phys.   85   1143 (2013). [22] Tu Y., Toner J., Ulm M. Sound waves and the absence of galilean invariance in flocks.   Phys. Rev. Lett.   80   4819 (1998). Toner,   J.   and   Tu,   Y.   Long-range   order   in   a   two- dimensional dynamical XY model: how birds fly together  Phys. Rev. Lett.   75 : 4326–4329 (1995) [23] Justh, E.W., & Krishnaprasad, P.S.   Equilibria and steer- ing laws for planar formations , Systems & Controls Let- ters   52 , 25–38 (2004) [24] Tanner, H.G.; Jadbabaie, A.; Pappas, G.J., Flocking in Fixed and Switching Networks,   Automatic Control, IEEE Transactions   52 , 863–868 (2007) [25] Attanasi A, Cavagna A, Del Castello L, Giardina I, Jelic A, Melillo S, Parisi L, Shen E, Viale M. Information transfer and behavioural inertia in starling flocks   Nat. Phys.   10 , 691–696 (2014) [26] Bialek, W., Cavagna, A., Giardina, I., Mora, T., Silvestri, E., Viale, M., and Walczak, A. M. Statistical mechanics for natural flocks of birds   Proc. Natl. Acad. Sci. USA  109 , 4786–4791 (2012) [27] Cavagna, A., Cimarelli, A., Giardina, I., Parisi, G., San- tagati, R., Stefanini, F., and Viale, M. Scale-free cor- relations in starling flocks.   Proc. Natl. Acad. Sci. USA  107 :11865–11870 (2010) [28] Ballerini, M., Cabibbo, N., Candelier, R., Cavagna, A., Cisbani,   E.,   Giardina,   I.,   Lecomte,   V.,   Orlandi,   A., Parisi, G., Procaccini, A., et al. Interaction ruling an- imal collective behavior depends on topological rather than metric distance: Evidence from a field study   Proc. Natl. Acad. Sci. USA   105 , 1232–1237 (2008). [29] Goldstein H., Classical Mechanics (Addison-Wesley Pub- lishing Company, Reading, MA, 1980)
17  [30] Fetter AL, Walecka, JD. Theoretical Mechanics of Parti- cles and Continua ( Courier Dover Publications, 2012) [31] Pomeroy, H, Heppner F. Structure of turning in airborne rock dove (Columba livia) flocks   Auk   109 : 256267 (1992) [32] Cavagna A., Duarte Queirs SM, Giardina I, Stefanini F and Viale M. Diffusion of individual birds in starling flocks   Proc. R. Soc. B   280 , 20122484 (2013) [33] Hohenberg, P. C. and Halperin, B. I., Theory of dynamic critical phenomena,   Rev. Mod. Phys.   49 , 435–479 (1977) [34] Halperin, B.I. & Hohenberg, P.C. Hydrodynamic Theory of Spin Waves,   Phys. Rev.   188 , 898–918 (1969) [35] Sonin, E.B. Spin currents and spin superfluidity,   Ad- vances in Physics   59 , 181–255 (2010) [36] Matsubara, T. & Matsuda, H., A lattice model of Liquid Helium, I,   Prog. Theor. Phys.   16 , 569–582 (1956) [37] Lane, C.T., Fairbank, H.A. & Fairbank, W.M. Second Sound in Liquid Helium II,   Phys. Rev.   71 , 600–605 (1947) [38] Justh, E.W. & Krishnaprasad, P.S. Equilibria and steer- ing laws for planar formations,   Systems & Controls Let- ters   52 , 25–38 (2004). [39] Szabo, P., Nagy, M., & Vicsek, T. Transitions in a self- propelled-particles model with coupling of accelerations  Phys. Rev. E   79 , 021908 (2009) [40] Hemelrijk C.K. & Hildenbrandt H. Some causes of the variable shape of flocks of birds   PLoS ONE   6 , e22479 (2011) [41] Gautrais, J., Ginelli, F., Fournier, R., Blanco, S., Soria, M., Chate ́ e, H. & Theraulaz, G. Deciphering interactions in moving animal groups,   Plos Comp. Biol.   8   e1002678 (2012) [42] Sumino   Y.,   Nagai   K.H.,   Shitaka   Y.,   Tanaka   D., Yoshikawa K., Chat ́ e H., Oiwa K. Large-scale vortex lat- tice emerging from collectively moving microtubules.   Na- ture   483 , 448–452 (2012) [43] Zwanzig, R. Nonequilibrium statistical mechanics (Ox- ford University Press, Oxford, 2001). [44] Gardiner, C. W. Handbook of stochastic methods. Vol. 3.(Springer, Berlin, 1985). [45] Goldstone, J. Field theories with Superconductor solu- tions.   Il Nuovo Cimento   19 , 154–164 (1961) [46] Ramaswamy S., Simha R.A., Hydrodynamics Fluctua- tions and Instabilities in Ordered Suspensions of Self- Propelled Particles,   Phys. Rev. Lett. 89 , 058101 (2002) [47] Camperi M, Cavagna A, Giardina I, Parisi G, Silvestri E, Spatially balanced topological interaction grants op- timal cohesion in flocking models   Interface Focus 2 , 715– 725 (2012) [48] FJ   Dyson,   General   theory   of   spin–wave   interactions.  Phys Rev   102 , 1217–1230 (1956). [49] Bertin E, Droz M, and Gregoire G.   Phys. Rev. E   74 , 022101 (2006) [50] Ihle T, Kinetic theory of flocking: Derivation of hydro- dynamic equations.   Phys. Rev. E 83   030901 (2011) [51] Bialek, W., Cavagna, A., Giardina, I., Mora, T., Pohl O., Silvestri, E., Viale, M., and Walczak, A. M. Social interactions dominate speed control in poising natural flocks near criticality.   Proc. Natl. Acad. Sci. USA   111 , 7212-7217 (2014) [52] Cavagna, A., Giardina, I., Ginelli I, Mora, T., Piovani D., Tavarone R. and Walczak, A. M. Dynamical maximum entropy approach to flocking.   Phys. Rev. E   89   042707 (2014)