Mean Field Model for Collective Motion Bistability

We consider the Czir\'ok model for collective motion of locusts along a one-dimensional torus. In the model, each agent's velocity locally interacts with other agents' velocities in the system, and there is also exogenous randomness to each agent's velocity. The interaction tends to create the alignment of collective motion. By analyzing the associated nonlinear Fokker-Planck equation, we obtain the condition for the existence of stationary order states and the conditions for their linear stability. These conditions depend on the noise level, which should be strong enough, and on the interaction between the agent's velocities, which should be neither too small, nor too strong. We carry out the fluctuation analysis of the interacting system and describe the large deviation principle to calculate the transition probability from one order state to the other. Numerical simulations confirm our analytical findings.


1.
Introduction. Collective motion has become an emerging topic in biology and statistical physics and has a large potential for the applications in other fields, for example, opinion dynamics in social science [28] and autonomous robotic networks in engineering [32]. Collective motion describes the systemic behavior of a large group of animals, insects or bacteria. In a collective motion model, the system has a large number of individuals, expressed as their locations and velocities, px i , u i q, i " 1, . . . , N . The movement of each individual is affected by the movements of other individuals in the system, mostly attractive interactions. Most biologically realistic and mathematically interesting models suppose local interactions, that is, any individual is affected only by its neighbors within a limited range. In addition, a collective motion model generally adds independent, internal or external randomness to the individuals' movements.
Two approaches are widely used to model collective motion: self-propelled particle (SPP) models and coarse-grained (CG) models. In an SPP model (also known as an individual-based model), the system is viewed as a large number of coupled stochastic processes describing the evolutions of the individuals in the system. In an SPP model we can easily build an elaborated interaction between an individual and its neighbors to match actual biological observations. However, an SPP model consists of numerous nonlinear equations and is difficult to analyze. To obtain the analytical tractability, further simplifications are usually required [41,19]; otherwise the system behavior has to be verified by numerical simulations. On the other hand, a CG model (also known as a macroscopic or mean-field model) assumes infinitely many individuals in the system and considers the evolution of the distribution of the individuals. Most CG models in the literature lead to hydrodynamic equations so that the system behavior can be analyzed. CG models have some limitations. Although some CG models [3,14,23,4,26] are indeed constructed from microscopic dynamics, it is generally difficult to relate the models to actual individuals' properties [1]. Another limitation of CG models is that noise present at the microscopic level can be averaged out by coarse graining. As a result most hydrodynamic equations only describe the deterministic evolution and cannot capture and analyze some important stochastic events like transitions from one order state to another order state observed for instance in [7]. It is possible to incorporate noise in the CG model in a phenomenological way [5] but it is important to model noise at the macroscopic level in a consistent way and to clarify what type of noise in the CG model is generated by noise in the SSP model.
In the recent years, there is a huge amount of literature concerning the modeling of collective motion so we mention only a few papers guiding our own work. Vicsek et al. [39] propose the scalar noise model that is currently the most well-known SPP model. The Czirók model [11] is a one-dimensional model on a torus. It is a simplified version of the Vicsek model and is used to describe the collective motion of locusts for instance (see [40,1] and references therein). Buhl et al. [7] experimentally and numerically study the switching of the alignments of locusts and find that the probability of switching is exponentially decaying with the number of locusts. In [41,19] the authors assume global and uniform interactions in the Czirók model and they derive the nonlinear Fokker-Planck equations for the average velocity and consider overall stochastic behavior. Some interesting features, like the existences of leaders or predators are considered in [10] and [30], respectively. Other various individual-based models related to the Czirók model are discussed in [29,9,31,2], and [36,18,37,3,4,26,25,38,33,8,15,34,35] construct and analyze CG models of collective motion. These papers show that the linear stability analysis of the homogeneous solutions is a simple but powerful technique, and that the formation of "bands" (regions with high density of particles moving in the same direction) can be observed in numerical simulations and related to traveling solutions of CG equations.
One of the most important goals of the collective motion modeling is to mathematically reproduce and analyze the ordering phenomenon where the individuals show a coherence in their movements, such as swarms of locusts, flocks of birds, or schools of fish. Such an ordering phenomenon is called an order state of the system and there are generally multiple order states that represent different directions and rotations. The counterpart of the order states is the disorder state, where the system does not show such a pattern. From the perspective of statistical physics, in this paper we address several interesting questions that are also considered in many of the literature: 1. Existence and number of order states. 2. Quantities that control the phase transition between the disorder state and the order states. 3. Transition probability from one order state to another order state.
Our contributions in this paper are the following. We interpret the Czirók model as an interacting particle system in the McKean-Vlasov framework [12,13]. In this sense, we characterize the state of the system by the empirical probability measure of the locations and velocities of all the individuals, and as the number of individuals goes to infinity, the random, empirical measure converges to a deterministic, probability measure whose density is the solution of a nonlinear Fokker-Planck equation. Therefore, when the number of individuals is large, we can consider the Fokker-Planck equation instead of all the individuals. In this way, we are also able to provide a clear connection between the SPP model and the CG model. This procedure was proved to be very efficient to study the stationary states and their stability properties in the case of opinion dynamics models [21]. It turns out that it is still efficient to address collective motion models, although the framework is different: In the opinion dynamics model addressed in [21], the Fokker-Planck equation is non-degenerate while it is degenerate in the collective motion model that we addresse here, as the external noise affects only the velocities and not the locations. As a result noise plays a very different role which strongly affects the conditions for the existence and stability of the stationary states. For instance, the noise has to be large enough to ensure the existence of the non-trivial stationary states in the opinion dynamics model. In dramatic contrast, we will see that the existence of the order states in the collective motion model does not depend on the noise level and that their stability depends in a nontrivial way of the noise level. Noise has to be large enough to ensure the linear stability of the order states, but large noise also triggers frequent transitions from one order state to another one.
By analyzing the Fokker-Planck equation obtained in the Czirók model, we find that the system may have one disorder state and two order states (clockwise and counterclockwise rotations at constant velocity). We find the necessary and sufficient condition for the existence of the order states. This existence condition quantitatively describes the phase transition between the disorder state and the order states. We then perform the linear stability analysis to examine the stabilities of the order and disorder states, and we find that the order states are stable when the external noises are sufficiently strong. The noise threshold value depends in a nontrivial way of the interaction between individuals, the smallest threshold value is reached when the interaction is neither too weak, nor too strong. Such a phenomenon that the randomness of movement improves the alignment of collective motion coincides with the observation in [41]. When noise is small, the order states are unstable and the complex modulational instability generates moving clusters whose velocities and sizes can be identified. Finally, when the number of individuals is large but finite, the empirical measure is still stochastic and has a small probability to transit from one order state to the other one. This small probability can be described by the large deviation principle. It increases with the noise level. It is exponentially decaying as a function of the number of individuals; this fact is also observed in [7,41]. We confirm our analysis by extensive numerical simulations.
The paper is organized as follows. The interacting particle model and its meanfield limit is presented in section 2. The equilibria of the nonlinear Fokker-Planck equation are analyzed in section 3. In section 4, the linearized stability analysis of the nonlinear Fokker-Planck equation is performed. We provide the fluctuation analysis and the large deviation principle in sections 5 and 6, respectively. In section 7, we verify our analytical findings by numerical simulations. We briefly consider the nonsymmetric Czirók model with normalized influence functions in section 8. We end with a brief summary and conclusions in section 9.
2. Model and Mean Field Limit. We consider the Czirók model, a system of N agents moving along the torus r0, Ls. The position x i and velocity u i of particle i satisfy the following stochastic differential equations: for i " 1, . . . , N , where tW i ptqu N i"1 are independent Brownian motions and u i is a weighted average of the velocities tu j u N j"1 : with the weights depending on the distance }¨} on the torus between the position x i and the positions of the other agents (i.e. }x j´xi } " minp|x j´xi |, L´|x j´xi |q).
Here φpxq is a nonnegative influence function normalized so that and Gpuq is an odd and smooth function. Note that there is no loss of generality in assuming (3) since we can always rescale G to ensure that this condition is satisfied. As we will see in this paper the interesting configuration is when u Þ Ñ u´Gpuq derives from an even double-well potential, with two symmetric global minima and one saddle point at zero. For instance, we can think at Gpuq " 2 tanhpuq (4) and then u´Gpuq derives from the double-well potential u 2 {2´2 logpcoshpuqq, or and then u´Gpuq derives from the double-well potential u 4 {4´u 2 {2. We define the empirical probability measure: Let us assume that µ N p0, dx, duq converges to a deterministic measureρpx, uqdxdu as N Ñ 8. This happens in particular when the positions and velocities of the agents at the inital time t " 0 are independent and identically distributed with the distribution with densityρpx, uq. Then, as N Ñ 8, µ N pt, dx, duq weakly converges to the deterministic measure ρpt, x, uqdxdu whose density is the solution of the nonlinear Fokker-Planck equation starting from ρpt " 0, x, uq "ρpx, uq. Note that the Fokker-Planck equation (7) is degenerate because there is no diffusion in x. The convergence proof is standard and can be found in [12,22].
3. Stationary Analysis for Order and Disorder States. We look for a stationary equilibrium, that is to say, a probability density function stationary in time: Proposition 1. The probability-density-valued solutions ρpx, uq to (8) have the following form: They are uniform in space, Gaussian in velocity, and their mean velocity ξ satisfies the compatibility condition: ξ " Gpξq.
There are, therefore, as many stationary equilibria as there are solutions to the compatibility equation (10). Proof. It is straightforward to check that (9) is a solution of (8). Reciprocally, let ρ be a solution of (8). After integrating (8) from u "´8 to u "`8, we obtaiń which shows that F puq is a Gaussian density function: The condition ξ " L ş uρpx, uqdu requires that To show that ρpx, uq is uniform in x, we first note that ρpx, uq is periodic in x. We expand it as: From (11) and (10), ρ k satisfies

JOSSELIN GARNIER, GEORGE PAPANICOLAOU, AND TZU-WEI YANG
We take a Fourier transform in u: From (16),ρ k satisfies the ordinary differential equatioń´2 πk L`η¯Bρ This equation can be solved. Let k ă 0. For any η P p´8,´2πk{Lq we have which goes to`8 as η Õ´2πk{L ifρ k p0q ‰ 0. However |ρ k pηq| ď ş |ρ k puq|du ď p1{Lq ş ş L 0 ρpx, uqdxdu " 1{L is uniformly bounded, which imposesρ k p0q " 0 and thereforeρ k pηq " 0 for any η P p´8,´2πk{Lq. Similarly, for any η P p´2πk{L, 8q we have πk L˘¯, which goes to`8 as η OE´2πk{L ifρ k p´4πk{Lq ‰ 0. The boundedness ofρ k imposesρ k p´4πk{Lq " 0 and thereforeρ k pηq " 0 for any η P p´2πk{L, 8q. Moreover,ρ k is continuous by the Riemann-Lebesgue lemma and thereforeρ k pηq " 0 for any η P R and k ă 0. Sinceρ´kp´ηq "ρ k pηq, we also find thatρ k pηq " 0 for any η P R and k ą 0. Therefore ρ k puq " 0 for any u P R and k ‰ 0 and the stationary solution is equal to ρpx, uq " ρ 0 puq. l When G is such that u´Gpuq derives from an even double-well potential, such as the two examples (4) and (5), there are three ξ satisfying the compatibility condition (10): 0 and˘ξ e , with ξ e ą 0 the unique positive solution of Gpξq " ξ. In the experiment, the equilibrium ρ 0 px, uq is the disorder state and ρ˘ξ e px, uq are the two order states of the locusts marching on the torus, the clockwise and counterclockwise rotations. The existence of the nonzero solution˘ξ e and hence the stationary states ρ˘ξ e px, uq is independent of the value of the noise strength σ. This observation is unusual and in fact inconsistent with the conclusion in many models in statistical physics and opinion dynamics [12,21]. We will see, however, that the noise strength plays a crucial role in the stability of the stationary states. 4. Linear Stability Analysis. In this section, we assume that the compatibility condition (10) has solutions and we use the linear stability analysis to study the stability of the stationary states. Let ξ be such that Gpξq " ξ and consider for small perturbation ρ p1q . By linearizing the nonlinear Fokker-Planck equation (7) we find that ρ p1q satisfies Our goal is to clarify under which circumstances the linearized system is stable.

4.1.
Expressions of the linearized modes. The perturbation ρ p1q is periodic in x. We expand it as: where the discrete Fourier coefficients are real-valued (because }L´x} " }x}) and bounded by one (because |φ k | ď ş L 0 φp}x}qdx{L " 1). We note that the equations are uncoupled in k and the system is linearly stable if all modes are stable. We then take a Fourier transform in u: From (22),ρ p1q k satisfies the following first-order partial differential equation: starting from an initial condition that is assumed to satisfy By the method of characteristics, we can obtain the following implicit expression forρ p1q k pt, ηq that is useful for the stability analysis: k pt, ηq has the following expression: where D k " 2πk{L, H ξ pηq " ηF ξ pηq and

4.2.
Stability of the 0th-order mode. In this section we find a necessary and sufficient condition for the stability of the mode We use the method of moments. We let m From (28), we get The first moment m p1q 0,1 is stable if and only if G 1 pξq ă 1. More exactly, if G 1 pξq " 1 then the first moment grows linearly in time. If G 1 pξq ă 1 then the first moment is bounded. Once m p1q 0,1 has been shown to bounded, then the solution ρ p1q 0 of (28) can be seen as the solution of a linear parabolic equation which is smooth and bounded as explained in the following proposition.
Proposition 2. The 0th-order mode is stable if and only if G 1 pξq ă 1.
Proof. We know that m with F 0 pηq " G 1 pξqηF ξ pηq. We find that, for any t ě 1, . Therefore, there exists a constant C that depends only on G, ξ and σ such that We can proceed in the same way after differentiating with respect to η the expres- When u´Gpuq derives from an even double-well potential, such as the examples (4) and (5), then, as soon as there exist nonzero solutions˘ξ e to the compatibility equation ξ " Gpξq, we have G 1 p˘ξ e q ă 1 because B u pu´Gpuqq u"˘ξe ą 0. Therefore, from Proposition 1 and Proposition 2, the existence of the order states ρ˘ξ e in (7) is equivalent to the stability of the 0-th mode ρ p1q 0 in (28). In addition, the condition that the equation ξ " Gpξq has the nonzero solutions˘ξ e implies G 1 p0q ą 1 because B u pu´Gpuqq u"0 ă 0, and therefore the disorder state ρ 0 is an unstable equilibrium to the nonlinear Fokker-Planck equation (7) when the order states ρ˘ξ e exist. 4.3. Sufficient condition for the stability of the kth-order modes. In this section we find a sufficient condition for the stability of the kth-order modes for k ‰ 0. We show that all the nonzero modes are stable for sufficiently large σ. From (26) in Lemma 4.1, w k ptq " B ηρ p1q k pt, 0q satisfies where The strategy is to show that ψ k P L 8 p0, 8q and that ş 8 0 |R k ptq|dt ă 1, and thus }w k } L 8 ď p1´}R k } L 1 q´1}ψ k } L 8 by (31). Once it is known that w k is uniformly bounded, it is not difficult to show that the kth mode is stable by the inspection of the formula (26).
1. ψ k defined by (33) belongs to L 8 r0, 8q. 2. Let R k be defined by (34). If |G 1 pξqφ k | ă 1 for all nonzero k, then l Proposition 3. If |G 1 pξqφ k | ă 1 for all nonzero k and σ is large enough, then all the kth-order modes for k ‰ 0 are stable.
Proof. First we writeρ p1q k pt, ηq as (26). By Lemma 4.2, w k P L 8 r0, 8q and thus we can show that }ρ  The linear stability analysis that we have presented reveals that the stability of the k-th mode (k ‰ 0) is determined by the behavior of the function w k solution of (31). The k-th mode is stable if and only if the function w k is bounded. We have shown that the condition ş 8 0 |R k ptq|dt ă 1 is a sufficient condition for the stability of the k-th mode, which is not, however, necessary. It is possible to give a necessary and sufficient condition for the instability of the k-th mode using the Fourier-Laplace transform of (31). This condition is that is non-empty. If this happens, then the Fourier-Laplace transform of w k : blows up when the complex frequency γ goes to an element in C k , so we can conclude that |w k ptq| grows exponentially with the growth rate γ r pkq " suptRepγq, γ P C k u.
If there exists a unique γpkq P C k with Repγpkqq " γ r pkq, then we can predict that the k-th mode grows like exppγpkqtq. The real part γ r pkq is the growth rate of the k-th mode, the imaginary part γ i pkq describes the dynamical behavior of the unstable mode as follows.
The mode with the largest growth rate is the one that drives the instability. If we denote k max " argmax k pγ r pkqq, then the instability has the form (up to a multiplicative constant): expr´ip2πk max x{L´γ i pk max qtqs exprγ r pk max qts.
For a fixed time, its spatial profile has the form of a periodic modulation with the spatial period L{k max . This modulation moves in time at the velocity γ i pk max qL{p2πk max q and grows exponentially in time with the rate γ r pk max q. This shows that the initial uniform distribution of the positions is not stable and that spatial clustering appears.
In Figure 1 we plot the growth rate of the most unstable mode in the situation addressed in the numerical section 7. The double-well potential G defined by (47) is parameterized by h which determines the depths of the wells and one can see that the linear stability is optimal when h is nor too small, neither too large. This will be confirmed by the numerical simulations presented below.  Figure 1. The real (a) and imaginary (b) parts of the complex growth rate of the first mode (the most unstable mode) as a function of σ for h " 5 (solid black), h " 6 (dashed red), h " 8 (magenta dotted), and h " 10 (blue dot-dashed). Here φ is given by (46) and L " 10. G is given by (47) and it derives from the double-well potential plotted in picture c. We address the linear stability of the order state ρ ξe . One can see that the threshold value for the noise level σ to ensure linear stability is 1.8 for h " 5, 0.85 for h " 6, 1.4 for h " 8, and 2.2 for h " 10. The most stable situation is the one corresponding to h " 6.

4.5.
Conclusion. If ξ " Gpξq and G 1 pξq ă 1, then the 0th-order mode of the stationary state ρ ξ is stable. If, additionally, |G 1 pξqφ k | ă 1 for all k ‰ 0, then all the kth-order modes for k ‰ 0 are stable for a sufficiently large σ. Under these conditions the stationary state ρ ξ is stable. When G is such that u´Gpuq derives from an even double-well potential, such as the two examples (4) and (5), this means that the order states are stable equilibria to the nonlinear Fokker-Planck equation (7), while the disorder state is an unstable equilibrium. This shows that the noise strength σ improves the stability of the order states ρ˘ξ e px, uq. Such a phenomenon that the randomness of the movements establishes the alignment of collective motion is also found in [41] both numerically and experimentally.
Note, however, an interesting phenomenon: if G 1 pξq ă 1, then the 0th-order mode is stable, but for some k ‰ 0 the kth-order mode may be unstable. This instability is very different from the instability that results from violating the condition G 1 pξq ă 1, because it manifests itself in an instability in the distribution of the positions of the particles which cannot stay uniform. This situation happens when u´Gpuq derives from a double-well potential and the wells are very deep, as we will see later in the numerical section. The prediction that the order states become unstable when the depths of the wells are too deep seems counterintuitive and indeed this never happens when the interaction is uniform φ " 1, but it happens when the interaction is local.
5. Fluctuations Analysis. In this section, we consider the fluctuation analysis for the interacting system of agents around the stationary state. Recall that µ N is the empirical probability measure of the locations and velocities (6). Here we consider a solution ξ to the compatibility equation ξ " Gpξq. We assume that, at time 0, the initial positions and velocities of the agents are independent and identically distributed with the distribution with density ρ ξ px, uq. Therefore µ N p0, dx, duq converges to the stationary state ρ ξ px, uqdxdu in (9) as N Ñ 8. Moreover, µ N pt, dx, duq converges to ρ ξ px, uqdxdu for all t as N Ñ 8 by (7).
We define, for any test function f px, uq and any measure Xpdx, duq (in the space S 1 of tempered Schwartz distributions on r0, LsˆR): We have We define the fluctuation process around ρ ξ px, uqdxdu: and study its dynamics as N Ñ 8. First, by central limit theorem, at time 0, onverges to a Gaussian random variable with mean zero and variance xρ ξ dxdu, f 2 yx ρ ξ dxdu, f y 2 . This shows that X N pt " 0, dx, duq converges weakly as N Ñ 8 to a Gaussian measure Xpt " 0, dx, duq with mean zero and covariance for any test functions f 1 and f 2 .

JOSSELIN GARNIER, GEORGE PAPANICOLAOU, AND TZU-WEI YANG
The process W N ptq in Cpr0, 8q, S 1 q is such that, for any test function f px, uq, It is a zero-mean martingale whose quadratic variation is, for any test functions f 1 and f 2 , which converges to the deterministic process In other words, as N Ñ 8, W N converges to a Brownian field, and this implies that X N converges to a generalized Ornstein-Uhlenbeck process, as explained in the following proposition.
Proposition 4. Let µ N pt, dx, duq " 1 N ř N i"1 δ pxiptq,uiptqq pdx, duq and ξ be a solution to ξ " Gpξq. We assume that the initial positions and velocities of the agents are independent and identically distributed with the distribution with density ρ ξ px, uq. 1) As N Ñ 8, µ N p0, dx, duq converges to the stationary state ρ ξ px, uqdxdu and X N p0, dx, duq converges to the Gaussian measure-valued process Xp0, dx, duq with mean zero and covariance (39). 2) As N Ñ 8, X N pt, dx, duq converges to the Gaussian measure-valued process Xpt, dx, duq satisfying where W ξ pt, x, uq is a Gaussian process independent of Xp0, dx, duq with mean zero and covariance for any test functions f 1 and f 2 .
The convergences hold in the sense of weak convergence on Cpr0, 8q, S 1 q where S 1 denotes the space of tempered Schwartz distributions on r0, LsˆR. Proof. We only provide the main steps of the derivation of Proposition 4, which follows the one of [12, Theorem 4.1.1]. We first note that the solution X N of (40) is the solution of a martingale problem. To prove the existence of a weak limit, we need to prove that the set tX N u 8 N "1 solution of (40) is weakly compact and thus the sequence tX N , N " 1, 2, 3, . . .u has limits. To prove the uniqueness of the weak limit, we need to prove that any weak limit X is solution of a well-posed martingale problem. Note that the second and third terms in the right-hand side of (40) satisfy when X N weakly converges to X. This shows that X is solution of the martingale problem associated to the Markov diffusion process (41). Then the uniqueness of the limit is ensured by the uniqueness to the solution of the limit martingale problem. l The equation (41) for the fluctuation process X involves the same linearized operator as in (20). So the linear stability analysis carried out in the previous section can be invoked to claim that the system (1) is linearly stable when the initial empirical distribution is initially close to ρ ξ , where ξ is such that Gpξq " ξ, G 1 pξq ă 1, |φ k G 1 pξq| ă 1 for all k ‰ 0, and σ is large enough. When G is such that u´Gpuq derives from an even double-well potential, such as the two examples (4) and (5), this means that the order states ρ˘ξ e are stable equilibria, in the sense that if the initial empirical distribution is close to one of them, then the empirical distribution remains close to it, with small normal fluctuations described by (41) in the limit N Ñ 8, while the disorder state ρ 0 is an unstable equilibrium, in the sense that if the initial empirical distribution is close to it, then the empirical distribution quickly moves away from it. 6. Large Deviations Analysis. Here we assume that two order states ρ˘ξ e exist and are stable and that the initial empirical measure µ N pt " 0, dx, duq converges to ρ ξe . By the mean-field theory, the empirical measure (6) converges to ρ ξe px, uqdxdu as N Ñ 8 and the stability analysis ensures that µ N pt, dx, duq is approximately ρ ξe over the time interval r0, T s. However, when N is large but finite, the empirical measure µ N is still stochastic and there is a very small probability that, eventhough µ N « ρ ξe at time t " 0, µ N « ρ´ξ e at time t " T . Here we use the large deviation principle (LDP) to write down the asymptotic probability of a transition from one stable order state to the other one. The empirical measure µ N satisfies the large deviation principle in the space of continuous probability-measure-valued processes as the following: for a set A of paths of probability measures µ " pµpt, dx, duqq tPr0,T s , we have´i whereÅ andĀ are the interior and closure of A, respectively, with an appropriate topology [13]. The function I is called the rate function: x Bµ Bt pt,¨,¨q´Lμ pt,¨,¨q µpt,¨,¨q, f y 2 xµpt,¨,¨q, p Bf Bu q 2 y dt, where Lν is the differential operator associated to the Fokker-Planck equation: We are interested in the rare event A which corresponds to the transition from one order state ρ ξe px, uqdxdu to the other one ρ´ξ e px, uqdxdu. We therefore assume that the initial conditions correspond to independent agents with the distribution ρ ξe px, uqdxdu and we consider the event A " pµpt, dx, duqq tPr0,T s :~µpT, dx, duq´ρ´ξ e px, uqdxdu~ď δ ( , for some small δ, where~¨~stands for the metric of the space of probability measures. Roughly speaking, for large but finite N , the probability of transitions from one stable state to the other is an exponential function of N whose exponential decay rate is governed by the rate function I: The transition probability (45) of the switching of the alignments of locusts is experimentally and numerically observed in [7].

Remark 1.
In this section we have only provided a formal large deviation description, because the classical large deviation result [13] requires that the noise must be non-degenerate. Therefore, the result of [13] cannot be directly applied to our model and we are not able to well define the space for our case until the rigorous large deviation is proven. However, a rigorous large deviation principle could be possibly constructed by the technique in [6].
7. Numerical Simulations. We use the numerical simulations to illustrate and verify our theoretical analysis. The spatial domain is the torus r0, Ls " r0, 10s, the influence function is φp}x}q " 5ˆ1 r0,1s p}x}q, which is such that ş L 0 φp}x}qdx{L " 1 and φ k " sincpπk{5q, and we let which is such that u´Gpuq derives from the potential V puq " 4´h 10 u 2`h 500 u 4 , that is a double-well potential as soon as h ą h c :" 4. The parameter h allows us to quantify the depths of the two wells of the double-well potential when h ą h c (see Figure 1c).
We discretize the time domain: t n " n∆t, n " 0, 1, 2, 3, . . . and use the Euler method to simulate (1):  Figure 2. The empirical average velocityū n and the square centered L 2 -discrepancy CL 2 2 pnq at each time step t n for h " 2 (a-b) and h " 6 (c-d). The dashed lines in the velocity plots are the solutions of Gpξq " ξ. The dashed lines in the discrepancy plots stand for the value (53) corresponding to a uniform sampling. Here ∆t " 0.1, N " 500, and σ " 2. One can see that the spatial distribution is uniform and the velocity average is 0 when h " 2 and ξ e when h " 6.
7.1. Existence of the order states. From Proposition 1, the number of stationary states (7) is equal to the number of solutions ξ for the compatibility condition (10). When h ď h c " 4, the function u´Gpuq is increasing and 0 is the unique solution of the compatibility equation (10). When h ą h c , the function u´Gpuq derives from a double well-potential and there are three solutions 0,˘ξ e to the compatibility equation (10), with In Figure 2 we test the cases of h " 2 and h " 6 and plot the empirical average velocityū and the square centered L 2 -discrepancy: It measures the discrepancy between the empirical distribution of the positions and the uniform distribution over r0, Ls. By [24] it is given by It is known that, if the positions are independent and identically distributed according to the uniform distribution over r0, Ls, then the square centered L 2 -discrepancy has mean [20]: while its variance is of order N´2. Therefore, as long as the empirical distribution in positions stays uniform, the square centered L 2 -discrepancy stays close to the value (53). In Figure 2 the initial locations and velocities tpx 0 i , u 0 i qu N i"1 are sampled from the distribution ρ 0 px, uqdxdu, the stationary state in (9) with ξ " 0. That is, tx 0 i u N i"1 are uniformly sampled over r0, Ls and tu 0 i u N i"1 are sampled from the Gaussian distribution with mean 0 and variance σ 2 {2.
-When h " 2 ă h c " 4, the only solution of ξ " Gpξq is ξ " 0, which means that ρ 0 is the unique stationary state. Moreover, G 1 p0q " ph`1q{5 ă 1 and the noise level is strong enough to make it linearly stable. Indeed, we can see in Figure 2a that the empirical average velocityū n oscillates around zero and the square centered L 2discrepancy oscillates around the value CL 2 2,u , which indicates that the empirical distribution of locations and velocities stays close to ρ 0 . -When h " 6 ą h c " 4, the order states ρ˘ξ e exist. Moreover, G 1 pξ e q " p13´2hq{5 ă p13´2h c q " 1 and the noise level is strong enough to make the order states ρ˘ξ e linearly stable. In addition, G 1 p0q " ph`1q{5 ą ph c`1 q{5 " 1, therefore, the disorder state ρ 0 is an unstable state. Under these circumstances, the initial empirical distribution that is close to ρ 0 quickly changes as can be seen in Figure 2c-d:ū n oscillates around´ξ e (it could have been ξ e ), while the square centered L 2 -discrepancy oscillates around the value CL 2 2,u . This indicates that the empirical distribution of locations and velocities becomes close to ρ´ξ e . 7.2. Linear stability. Here we assume that h ą h c " 4 so that there are three stationary states ρ 0 , ρ˘ξ e . We have G 1 p0q " ph`1q{5 ą 1, G 1 p˘ξ e q " p13´2hq{5 ă 1, and thus the 0th order mode of ρ 0 is unstable while the 0th-order modes of ρ˘ξ e are stable. An immediate conclusion is that ρ 0 is an unstable state for the nonlinear Fokker-Planck equation (7). The zeroth-order mode of ρ˘ξ e is linearly stable. However, the stability of ρ˘ξ e requires that all the nonzeroth-order modes of ρ˘ξ e are linearly stable. From Proposition 3, we know that for a sufficiently large σ, all the nonzeroth-order modes are stable. Moreover the critical value of σ ensuring the stability can be determined as explained in subsection 4.4.
If h " 6, then the theory in subsection 4.4 predicts that the condition for the stability of the order states ρ˘ξ e is that σ should be larger than 0.85. The simulations shown in Figure 3 confirm these predictions: ρ˘ξ e is stable when σ exceeds the threshold value 0.85. They also allow us to illustrate the instability mechanism exhibited above: when σ is below the threshold value, the average velocity of the a)  Figure 3. The empirical average velocityū n and the square centered L 2 -discrepancy CL 2 2 pnq for σ " 0.5 (a-b), σ " 1 (c-d), σ " 1.5 (e-f). The initial positions tx 0 i u N i"1 are uniformly sampled over r0, Ls and the initial velocities tu 0 i u N i"1 are sampled from the Gaussian distribution with mean ξ e and variance σ 2 {2. Here ∆t " 0.1, N " 2000, and h " 6. One can see that the average velocity is not ξ e and the spatial distribution is not uniform when σ " 0.5 while the average velocity is ξ e and the spatial distribution is uniform when σ " 1 or 1.5. particles is linearly stable, but the uniform distribution of the positions of the particles is not stable and a spatial modulation grows that gives rise to one moving cluster. Note that Figure 3a indeed shows that the uniform distribution of the positions is not stable when σ " 0.5 since the square L 2 -discrepancy takes values much larger than (53). The theory in subsection 4.4 predicts that the most unstable mode is the first one k max " 1, which means that there should be one moving cluster. Figure 1 plots the predicted growth rate γ r p1q of the first mode and the coefficient γ i p1q. For σ " 0.5, we have γ i p1q » 2.15 and the apparent velocity for the cluster is predicted to be approximately equal to γ i p1qˆL{p2πq » 3.4 by (36), which is close to, although slightly larger than, the average velocity of the particles in the stationary distribution ξ e " 2.9. There is no contradiction as the apparent velocity of the cluster can be interpreted as a "phase velocity" while the average velocity of . The initial positions tx 0 i u N i"1 are uniformly sampled over r0, Ls and the initial velocities tu 0 i u N i"1 are sampled from the Gaussian distribution with mean ξ e and variance σ 2 {2. Here ∆t " 0.1, N " 2000, σ " 0.5, and h " 6. One can see that the spatial distribution is not uniform and a cluster is moving with an apparent velocity of 3.6.
the particles can be interpreted as a "group velocity". The growth of the first mode giving rise to the moving cluster can be observed in the simulations (see Figure 4). Moreover, we can see in the right pictures of Figure 4 that the particles in the tail of the moving cluster are transfered to the front of the "next" cluster (which is the same one, by periodicity), which indeed allows the cluster to apparently move faster than the average velocity of the particles: in the simulations, the average velocity of the particles is about 2.8, while the velocity of the cluster is about 3.6, very close to the predicted value 3.4.
The previous results illustrate the fact that the order states ρ˘ξ e are linearly stable when σ exceeds the threshold value. This means that, if the initial conditions are close to ρ˘ξ e , then the empirical distribution remains close to it. This is what we see in Figure 3c-f. When the initial conditions are far from ρ˘ξ e and when the noise level is significantly larger than the threshold value, then the empirical distribution quickly converges to one of the two stationary distributions ρ˘ξ e , as we can see in Figure 5e-f: after a transient period, the empirical distribution in positions becomes uniform and the mean velocity becomes equal to˘ξ e . When the initial conditions are far from the two stable stationary states and when the noise level is larger than but close to the threshold value, then the empirical distribution reaches a state that is not stationary but looks like a moving cluster, as shown in Figure 5a-d.
As revealed by the linear stability analysis carried out in subsection 4.4, the results are not monotoneous as a function of h.  Figure 5. The empirical average velocityū n and the square centered L 2 -discrepancy CL 2 2 pnq for σ " 0.5 (a-b), σ " 1 (c-d), σ " 1.5 (e-f). The initial positions tx 0 i u N i"1 are uniformly sampled over r0, Ls and the initial velocities tu 0 i u N i"1 are sampled from the Gaussian distribution with mean 0 and variance σ 2 {2. Here ∆t " 0.1, N " 2000, and h " 6. One can see that the average velocity is not˘ξ e and the spatial distribution is not uniform when σ " 0.5 or 1 while the average velocity is ξ e and the spatial distribution is uniform when σ " 1.5. stationary states ρ˘ξ e are unstable when σ " 1.0 for h " 5 and for h " 10 while they are stable for h " 6 (see Figure 3c-d), as predicted by the theory in subsection 4.4. This means that the two order states are stable when the depths of the wells of the double-well potential are neither too large nor too small. 7.3. Large deviations. In this subsection we assume that the conditions are fulfilled so that ρ ξe and ρ´ξ e exist and are stable. We consider the transition probability from one stable state to the other:  Figure 6. The empirical average velocityū n (a), the square centered L 2 -discrepancy CL 2 2 pnq (b), and the empirical position distribution smoothed by kernel density estimation (c and d). The initial positions tx 0 i u N i"1 are uniformly sampled over r0, Ls and the initial velocities tu 0 i u N i"1 are sampled from the Gaussian distribution with mean ξ e and variance σ 2 {2. Here ∆t " 0.1, N " 2000, σ " 1, and h " 10. One can see that the average velocity is not ξ e and the spatial distribution is not uniform.
where the rare event A is the set (44) of all possible transitions from ρ ξe at t " 0 to ρ´ξ e at t " T , and the rate function I is defined in (43).
In Figures 8, 9, and 10, we qualitatively examine the effects of N , σ, and h to the transition probability Ppµ N P Aq by checking the frequencies of transitions of the empirical average velocityū n between the two stable order states ρ˘ξ e . In Figure 8, we can observe that the system has less transitions as N becomes larger. In other words, the probability of transition Ppµ N P Aq becomes smaller as N becomes larger; this is consistent with the fact that the large deviation principle predicts that Ppµ N P Aq is an exponential decay function of N . In Figure 9, we can observe that the system experiences more transitions as σ becomes larger. Therefore, Ppµ N P Aq becomes larger as σ becomes larger and this is confirmed by the fact that the rate function is approximately proportional to 1{σ 2 (note, however, that σ is also in the operator L˚). In Figure 10, we find that the probability of transition decreases as h increases. This is qualitatively consistent with the fact that h determines the height of the potential barrier between the two potential wells, in analogy with the classical problem of diffusion exit from a domain [17, Section 5.7].
8. The Nonsymmetric Case. Let us briefly revisit the previous analysis when the average velocity has the form [27] u i "  Figure 7. The empirical average velocityū n (a), the square centered L 2 -discrepancy CL 2 2 pnq (b), and the empirical position distribution smoothed by kernel density estimation (c and d). The initial positions tx 0 i u N i"1 are uniformly sampled over r0, Ls and the initial velocities tu 0 i u N i"1 are sampled from the Gaussian distribution with mean ξ e and variance σ 2 {2. Here ∆t " 0.1, N " 2000, σ " 1, and h " 5. One can see that the average velocity is not ξ e and the spatial distribution is not uniform.
where N i is the weighted number of agents in the neighborhood of the ith agent: If µ N p0, dx, duq converges to a deterministic measureρpx, uqdxdu, then µ N pt, dx, duq converges to the deterministic measure ρpt, x, uqdxdu whose density is the solution of the nonlinear Fokker-Planck equation starting from ρpt " 0, x, uq "ρpx, uq.
is a stationary solution of (56).  Figure 9. The empirical average velocityū n at each time step t n for σ " 4 (a), σ " 4.5 (b), σ " 5 (c), and σ " 5.5 (d). Here ∆t " 0.1, N " 100, and h " 6. The frequencies of the transitions between the two stable order states increases when σ increases. The system has more transitions with a higher σ.  Figure 10. The empirical average velocityū n at each time step t n for h " 5 (a), h " 5.5 (b), h " 6 (c), and h " 6.5 (d). Here ∆t " 0.1, N " 100, and σ " 5. The frequencies of the transitions between the two stable order states decays with h.
Note, however, that we could not prove that any stationary state is of the form (58). It is possible to prove that any stationary state that is uniform in space is of the form (58), but we could not prove that a stationary state must be uniform in space.
The linear stability analysis of the stationary states listed in Proposition 5 follows the same lines as in the symmetric case (2). Let ξ P Ξ and consider ρpt, x, uq " ρ ξ px, uq`ρ p1q pt, x, uq " 1 L F ξ puq`ρ p1q pt, x, uq, for small perturbation ρ p1q . We linearize the nonlinear Fokker-Planck equation (56) Bρ p1q Bt "´u Bρ p1q Bx´B Bu The perturbation ρ p1q is periodic in x and can be expanded as (21). The Fourier coefficientsρ The necessary and sufficient condition for the linear stability of the 0th order mode is G 1 pξq ă 1.
The sufficient condition for the linear stability of the other modes is that |G 1 pξqφ k | ă 1 for all k ‰ 0 and the noise level σ should be larger than a threshold value σ c , which is, however different from the threshold value of the symmetric case (1-2).
9. Conclusion. We have analyzed the Czirók model for the collective motion of locusts. The mean-field theory is used to obtain a nonlinear Fokker-Planck equation as the number of agents tends to infinity. We analyze the phase transition between the disorder and order states by the existence condition for the order states. We then perform the linear stability analysis on the stationary states, and prove that the order states are stable for a sufficiently large noise level and when the interaction between the velocities of the particles is neither too small nor too strong. We provide the fluctuation analysis and the large deviations principle. For a large but finite system we calculate the asymptotic, exponentially small transition probability from one order state to the other. Our analytical findings are verified by the extensive numerically simulations and are found in agreement with the experimental observations.
Therefore a sufficient condition to ensure that ş 8 0 |R k ptq|dt ă 1 is By the fact that |D k | " 2π|k|{L is increasing in |k|, we can have a simplified sufficient condition from (64): If |G 1 pξqφ k | ă 1 for all nonzero k, then this condition is satisfied for all nonzero k if σ is large enough. This completes the proof of the second item of the Lemma.