Symmetries and pattern formation in hyperbolic versus parabolic models of self-organised aggregation

The study of self-organised collective animal behaviour, such as swarms of insects or schools of fish, has become over the last decade a very active research area in mathematical biology. Parabolic and hyperbolic models have been used intensively to describe the formation and movement of various aggregative behaviours. While both types of models can exhibit aggregation-type patterns, studies on hyperbolic models suggest that these models can display a larger variety of spatial and spatio-temporal patterns compared to their parabolic counterparts. Here we use stability, symmetry and bifurcation theory to investigate this observation more rigorously, an approach not attempted before to compare and contrast aggregation patterns in models for collective animal behaviors. To this end, we consider a class of nonlocal hyperbolic models for self-organised aggregations that incorporate various inter-individual communication mechanisms, and take the formal parabolic limit to transform them into nonlocal parabolic models. We then discuss the symmetry of these nonlocal hyperbolic and parabolic models, and the types of bifurcations present or lost when taking the parabolic limit. We show that the parabolic limit leads to a homogenisation of the inter-individual communication, and to a loss of bifurcation dynamics (in particular loss of Hopf bifurcations). This explains the less rich patterns exhibited by the nonlocal parabolic models. However, for multiple interacting populations, by breaking the population interchange symmetry of the model, one can preserve the Hopf bifurcations that lead to the formation of complex spatio-temporal patterns that describe moving aggregations.

2004; Gamba et al. 2003). Nevertheless, because of the complexity of these transport models, a common approach to investigate them is to transform them into parabolic models-for which there is an entire analytical toolbox Othmer and Hillen 2000;Degond and Yang 2010;Chalub and Souza 2014;Rousset and Samaey 2013;Bellouquid et al. 2013). The classical approach to transform a 1D hyperbolic equation into a parabolic one is to take a formal parabolic limit using Kac's trick (Kac 1956). This parabolic limit assumes that individuals travel very fast and change direction very quickly (Holmes 1993). In turn, this might lead to a reduced sensitivity to the surroundings (Eftimie 2012), which can affect the exhibited patterns.
Despite these approaches that could be used to investigate various hyperbolic and parabolic models for collective movement, there are not many studies (Hillen 1997;Holmes 1993) that focus on how the patterns exhibited by these models are preserved in the parabolic limits. Also, it is not clear whether all these models can exhibit the same types of patterns. The goal of this article is to address this aspect in a class of 1D nonlocal hyperbolic models for biological aggregations [which were previously introduced in (Eftimie et al. 2007a;Eftimie 2013)] and their parabolic counterparts (obtained via parabolic limits). The original hyperbolic models incorporate one particular communication mechanism, where it is assumed that individuals interact only with neighbours moving towards them. Previous studies have shown that these hyperbolic models, with this inter-individual communication mechanism, can exhibit a large variety of moving and stationary aggregations (Eftimie et al. 2007a;Eftimie 2013). Our aim here is to identify the mathematical and biological mechanisms that can ensure the presence or absence of patterns describing stationary and moving aggregations for the limiting parabolic models. In particular, we use symmetry and bifurcation theory to compare the two types of models in terms of their potential for pattern formation-an approach not attempted before for models of collective movement.
To understand and classify in an unified way the mechanisms behind various spatial and spatio-temporal patterns, one can use the symmetry point of view; see Golubitsky and Stewart (2002). While symmetry breaking behaviours have been described, for example, in swarming models for ants (Altshuler et al. 2005) and bacteria (Leoni and Liverpool 2010;Taylor and Welch 2008), the use of symmetry theory is still not a common approach when investigating pattern formation in models for collective movement. However, the consideration of symmetry is natural in the context of the 1D nonlocal hyperbolic and parabolic models investigated in this article, since the equations are symmetric with respect to translations in the domain because of the invariance of the differential part (constant transport or diffusion) and the translation invariance of the integrals describing nonlocal social interactions. Moreover, since we impose periodic boundary conditions, we obtain a group of spatial translations on the unit circle that sends solutions to solutions of the equation; this symmetry group is denoted by SO(2). In addition, the nature of the nonlinear terms that describe communication-induced social interactions enables the reflection across the midpoint of the interval to also send solutions to solutions. Therefore this reflection along with the SO(2) group generate the group of all symmetries of a circle, which is denoted by O(2). Recently, we showed that this class of O(2) symmetric hyperbolic models can exhibit not only codimension-one bifurcations but also all three types of codimension-two bifurcations, namely Hopf/Hopf, Hopf/steady-state and steadystate/steady-state (Buono and Eftimie 2014a, b). These complex bifurcations are one of the mechanisms that explain the large variety of spatial and spatio-temporal patterns observed in these nonlocal hyperbolic models. For readers unfamiliar with the use of symmetry methods in bifurcation theory, Golubitsky and Stewart (2002) provides all the necessary definitions and results (with particular applications to biology and pattern formation).
In this article, we focus on the nonlocal advection-diffusion equations obtained via parabolic limits from the previously-discussed hyperbolic models, and show that they can exhibit only O(2)-symmetric steady-state bifurcations. Therefore, only stationary aggregations are created via local bifurcations from homogeneous equilibria. In fact, we go further and show that any O(2)-symmetric scalar parabolic equation with nonlocal terms can only exhibit O(2) steady-state bifurcation from a homogeneous equilibrium. This generalizes a result of Britton (1990) who studied a (spatially) nonlocal scalar reaction-diffusion equation with O(2)-symmetry, and showed that Hopf bifurcations cannot occur from the homogeneous equilibrium. In consequence, travelling aggregations (i.e., travelling pulses) cannot emerge from small perturbations of a homogeneous equilibrium because of the absence of these Hopf bifurcations. This does not mean that travelling pulses do not exist in the O(2)-symmetric nonlocal advection-diffusion equation. However, as we explain in Sect. 2.2.4, Matano (1988) gives a negative answer to this question for O(2)-symmetric local scalar parabolic reaction/advection-diffusion equations. For nonlocal reaction/advection diffusion equations, this is still an open question, and addressing it is beyond the scope of this paper.
Since the systems of parabolic equations are usually more rich in terms of pattern formation compared to single parabolic equations, in this study we also focus on nonlocal systems of parabolic equations which are obtained, via the parabolic limit, from hyperbolic systems of multiple interacting populations. We perform a similar analysis for a two-population advection-diffusion model and obtain the occurrence of O(2)-symmetric steady-state bifurcations and Hopf bifurcations. However, the Hopf bifurcations are only possible if the communication mechanisms used by the two populations lead to different social interactions which break the interchange of population symmetry of the models. This again is in contrast with respect to the corresponding hyperbolic model, which has a much wider range of spatio-temporal dynamics. Numerical simulations show some of the emerging patterns for the one-population and two-population parabolic models.
Overall, this analysis not only gives us insight into the types of patterns that can be obtained with nonlocal advective-diffusive models, but also it allows us to compare in terms of symmetry and bifurcations two types of hyperbolic and parabolic models for self-organised biological aggregations (an approach not seen in mathematical models for swarming). Finally, this analysis emphasises that, when we scale our mathematical models to make them more easy to be analysed, even if the symmetry of the initial model is preserved in the scaling, the bifurcation structure might not preserved, leading thus to different types of patterns.
The paper is structured as follows. In Sect. 2 we introduce a one-population hyperbolic model which assumes that individuals interact via different communication mechanisms. We obtain the parabolic limit of this model and then discuss the symmetries and the bifurcation dynamics of the parabolic model, including the preservation of various communication mechanisms and the restrictions on O(2) Hopf bifurcations. In Sect. 3 we introduce a two-population hyperbolic model which assumes that different sub-populations interact using different communication mechanisms. Again, we derive the parabolic limit and discuss the symmetry and bifurcation dynamics of the parabolic two-population system. We also discuss the symmetry-breaking effects of the communication mechanisms necessary to obtain Hopf bifurcations. We conclude in Sect. 4 with a numerical description of some of the patterns exhibited by the parabolic models. Appendix B contains the result that the spectrum of the linear operator at a homogeneous solution for a general O(2)-symmetric scalar parabolic equations with nonlocal terms consists only of real eigenvalues, thus excluding O(2)-symmetric Hopf bifurcation. Finally, Appendix A provides a short summary of some basic concepts and results used in this paper, related to bifurcation theory with symmetry.

One population model
The following hyperbolic model was introduced in Eftimie et al. (2007b) to describe the one-dimensional movement of a population formed of right-moving (u + ) and left-moving (u − ) individuals: with initial conditions u ± (x, 0) = u ± 0 (x). Here γ describes the constant speed, and λ ± describe the turning behaviour in response to the presence of neighbours within the attractive (y ± a ), repulsive (y ± r ) and alignment (y ± al ) interaction ranges (Eftimie et al. 2007a): The turning function f should be positive, increasing and bounded (to make sense biologically). An example of such a function, which is used throughout this article, is f (y) = 0.5 + 0.5 tanh(y − 2). The terms λ 1 + λ 2 f (0) and λ 2 f (y ± ) − f (0) describe the baseline random turning rate and the bias turning rate, respectively. We note that it is biologically realistic to have these different turning rates, since, as mentioned in Lotka (1925), "the type of motion presented by living organisms ... can be regarded as containing both a systematically directed and also a random component". Because f is chosen such that f (0) 1, we can approximate the random turning by λ 1 and Fig. 1 The interaction of reference right-moving individuals u + (x, t) with their neighbours positioned at x + s and x − s, as given by four communication models M2-M5. The interactions of left-moving individuals u − (x, t) can be described in a similar manner. The notation is similar to the one in Eftimie et al. (2007a). Solid black arrows describe the movement direction of individuals, while dotted arrows describe the direction in which information travels. The "X" on the dotted arrows means that information cannot be received from that particular direction. The spatial positioning of neighbours is shown at the bottom of the figure the directed turning-in response to neighbours' density and movement direction-by λ 2 f (y ± ). The opposite signs in front of y ± r and y ± a describe the opposite effects of the repulsive and attractive social interactions. Finally, since individuals in an ecological aggregation can interact-via visual, sound or olfactory signals-with neighbours positioned further away (Larkin and Szafoni 2008b), we assume that the social interaction terms y ± r,a,al are nonlocal. Moreover, we assume that these terms depend not only on the spatial distance between individuals [modelled by spatial kernels K j (s), j ∈ {r,a,al}-given by equation (3)], but also on the communication signals. Figure 1 describes the interaction of a right-moving reference individual u + (x, t) with its neighbours positioned ahead at x + s and behind at x − s, upon the reception of communication signals emitted by these neighbours. The four models shown in Fig. 1, and labeled M2-M5, correspond to four different ways information can be received from neighbours. (Eftimie et al. (2007a) introduced also a fifth toy model denoted by M1. Since that model does not bring anything new to this study, we ignore it here.) Constants q a,al,r describe the magnitudes of the attractive, alignment and repulsive interactions, respectively. Kernels K a,al,r (s), given by Eq. (3), describe the spatial ranges for each of these social interactions Table 1 describes the four communication mechanisms corresponding to models M2-M5. The parameters q r , q al and q a that appear in front of the integrals represent the magnitudes of repulsive, alignment and attractive interactions, respectively. For the nonlocal terms, we choose normalised kernels with s j representing half the length of the interaction ranges and m j = s j /8 representing the width of the interaction kernels [see also Eftimie et al. (2007a)]. We choose m j to ensure that the kernels overlap for less that 2 % of their mass. The integrals in Table 1 are defined on infinite domains. However, for the purpose of numerical simulations, as well as for investigating analytically the problem of pattern formation, we consider a finite domain of length L 0 (i.e., [0, L 0 ]) with wraparound boundary conditions for the nonlocal social interactions. For L 0 large, this can approximate the process of pattern formation on an unbounded domain. The derivation of the model is thus completed by the addition of periodic boundary conditions: Another reason for choosing these types of boundary conditions is to be consistent with experimental studies on collective animal movement on circular domains (Buhl et al. 2006). Such studies have shown "milling" behaviours, described by individuals or groups of individuals [e.g., locusts (Buhl et al. 2006)] moving around the domain. As shown in Sect. 4, these "milling" behaviours are equivalent to our "rotating waves" patterns.
In the following, we take the formal parabolic limit to transform this hyperbolic model into an advection-diffusion model.

Formal parabolic limit
To transform model (1) into a parabolic model, we first differentiate the sum and difference of equations (1a) and (1b) with respect to t and x, respectively (Kac 1956). This leads to two equations for the total density u = u + + u − and for the difference of the densitiesū = u + −u − , which we then use to eliminate γ ∂ 2 u ∂ x∂t [see also Holmes (1993)]. Eventually we obtain the following wave equation for the total population density: Note that here we usedū = −(1/γ ) To transform this wave equation into a parabolic equation, we assume that individuals travel at very large speeds and change direction even faster-for example, during search behaviours. In particular, we fix the space scale (e.g., to meters) and assume that the time scales for the speeding behaviour (τ γ ) and for the turning behaviour (τ λ ) satisfy τ λ τ γ . These assumptions are biologically realistic since, for example, mosquitoes can travel at speeds of 2 km/h (i.e., 20/36 m/s). Hence, on a length scale of 1 m, τ γ ≈ O(1) s. On the other hand, mosquitoes can change their movement direction in just 25-35 ms (Iams 2014) [hence a turning rate of 1/(35 ms)-1/(25 ms)], which means that on a length scale of 1 m, τ λ ≈ O(10 −2 ) s. Therefore, in this case τ λ /τ γ ≈ O(10 −2 ).
With these assumptions, we can rescale the speed and turning parameters as follows in such a way that lim →0 γ 2 2λ 1 = D and lim →0 γ λ 2 2λ 1 = B.
Here, parameter measures the difference in the time scales for the speed and turning behaviours, as we fix the space scale [in the example given before for mosquitoes, fixing the space scale to 1m, we obtained ≈ O(10 −2 )]. Because individuals with very large speeds and turning rates cannot have a perfect perception of the environment (and implicitly of their neighbours), we consider also the scaling f (y) = f * (y), where function f appears in the turning rates λ ± (2). Note that, similarly, we could have chosen to scale the space (x * = x) and time (t * = 2 t) variables, as done by Hillen and Othmer (2000) to model the drift and diffusion time-scales of bacteria such as E. coli.
Substituting the rescaled function f = f * and rescaled parameters γ = γ * / , λ 1,2 = λ * 1,2 / 2 back into Eq. (5), taking the limit → 0 and dropping the "*" (for simplicity), leads to the following scalar parabolic equation: For models M2 and M4, the nonlocal terms are (in the limit → 0) with M2 : For models M3 and M5, the nonlocal terms are (in the limit → 0) Recall that K j , j = r,a,al, are social interaction kernels [defined in (3)], while q j , j = r,a,al, are the magnitudes of the social interactions; see also the caption of Table 1. Also, we observe here that the parabolic equations corresponding to models M2 and M3 do not contain anymore alignment interactions. This is the result of having y M2,M3;± al dependent only onū = −( /γ * ) x 0 ∂u(s,t) ∂t ds [upon substituting u + = (u +ū)/2 and u − = (u −ū)/2 into the y ± al terms in Table 1]. Hence, as → 0, these alignment terms vanish. For models M4 and M5, the terms y M4,M5;± al depend on both u andū, and thus the alignment interactions can persist.
We can conclude from here that the parabolic limit leads to an averaging of the communication mechanisms. The nonlocal advection depends now only on the differences between the density of neighbours to the left/right of the reference individual, and not on the communication mechanisms they use.
To fully define the parabolic model (8), we need to add boundary conditions. To be consistent with the hyperbolic model, we impose again periodic boundary conditions: 2.2 Analysis of the parabolic model In the following we investigate the limiting parabolic model (8) from the point of view of its symmetry (that is, the group of transformations which sends solutions to solutions), as well as the types of local bifurcations it can exhibit. This approach will help us compare the dynamics of this parabolic model and the dynamics of the original hyperbolic model, with the purpose of deciding if both models can exhibit the same types of aggregation patterns.

Symmetries for the one-population hyperbolic and parabolic models
It was shown by Buono and Eftimie (2014b) that the hyperbolic system of equations (1) defined on [0, L 0 ] with periodic boundary conditions is O(2)-symmetric (or O(2)equivariant), where the O(2) action is generated by θ ∈ SO(2), a spatial translation modulo L 0 , and κ, an order two element interchanging left and right moving individuals along with a reflection at mid-domain. Explicitly we have The parabolic limit models (8) obtained in Secti. 2.1 inherit the O(2) symmetry.

Proposition 2.1 Equation (8) is O(2)-equivariant for all communication models, M2, M3, M4 and M5, where the O(2) action is given by
The proof is done in Appendix B, and is a consequence of a result we prove there for general O(2)-equivariant scalar parabolic equations with nonlocal terms.

Remark 2.2
The proof in Appendix B shows that the O(2)-equivariance of (8) holds for all functions f in (2) such that, for y, z ∈ R, f (y) − f (z) is an odd function.

Spatially homogeneous steady states
It was shown in Eftimie et al. (2007b), Eftimie et al. (2007a) that the hyperbolic system (1) can exhibit up to five distinct spatially homogeneous steady states (i.e., states characterised by individuals evenly spread over the whole domain): The different states are: Note that cases (ii) and (iii) can be satisfied by more than one steady state, and hence the five spatially homogeneous states shown in Fig. 2a for models M2, M3 and M4. In Fig. 2b we show an example of stability for the model with communication mechanism M2. We remark here that in certain parameter regions [e.g., q al ∈ (2, 2.6) in Fig. 2b], it is possible that multiple states are unstable, and thus the final pattern depends on the initial condition. For simplicity, throughout this article, we focus only on the spatially homogeneous state (u + , u − ) = (A/2, A/2). For a fixed total population density, the parabolic equation (8) can only have a unique spatially homogeneous state u(x, t) = u * ≥ 0 (as depicted in Fig. 2c). Therefore, the parabolic model cannot display patterns associated with the transitions between different homogeneous states [as it is shown to happen for hyperbolic models (Buono and Eftimie 2014a)].

Spatially heterogeneous steady states
The spatially heterogeneous steady states of the hyperbolic system (1) are solutions of with C some arbitrary constant that connects u + and u − : u + = u − +C. (For simplicity, in (16) we denoted u + = u and u − = u − C.) If we assume that u ± (±∞) = 0 then C = 0. In fact, Buono and Eftimie (2014b) showed a more general result for the case of periodic boundary conditions, namely if (u + (x), u − (x)) is a steady-state solution of the hyperbolic model (1) fixed by the (isotropy) The previous equation can be thus rewritten as: The spatially heterogeneous steady-states of the parabolic equation (8) are given by Here, C 1 is a constant that appears from integration with respect to x, while rates D and B are given by (7). Again, if we assume that u(±∞) = 0 and ∂u ∂ x (±∞) = 0 then C 1 = 0.
Since the spatially heterogeneous steady states for the initial hyperbolic model (1) and the limiting parabolic model (8) coincide, we expect for the parabolic model to support at least stationary pulses, corresponding to resting animal aggregations.

Bifurcation from homogeneous equilibrium: parabolic model
Previous studies of the hyperbolic system (1) with periodic boundary conditions on [0, L 0 ] showed the existence of symmetry-breaking steady-state and Hopf bifurcations, including codimension-two bifurcations (Eftimie et al. 2009;Buono and Eftimie 2014b). These bifurcations are one of the mechanisms leading to the multitude of spatial and spatio-temporal aggregation patterns observed in system (1). Moreover, the bifurcations were shown to occur from the homogeneous equilibrium with O(2) symmetry, namely (u + , u − ) = (u * , u * ), with u * = u + * = u − * (although it is likely that they could occur also from the other four homogeneous equilibria). We now show that the parabolic equation (8) can only support steady-state bifurcations from a homogeneous equilibrium u = u * . We do this analysis for the general O(2) case. In Eq. (8), we let Note that for models M2 and M4, y ± [u] = ±(K + * u) ∓ (K − * u), while for models M3 and M5, y ± [u] = K ± * u. The linearization operator at u = u * is where the gradient is taken with respect to the components of (u, K + * u, K − * u, u x , K + * u x , K − * u x ) considered as variables. In Appendix B, we show a characterization of the spectrum of the linearization of an O(2)-equivariant scalar parabolic equation with nonlocal terms, for which the next result is a special case.
Then, the spectrum of L can only have real eigenvalues. In particular, there cannot be O(2)-symmetric Hopf bifurcations.
In O(2)-symmetric Hopf bifurcation, two branches of periodic solutions emerge simultaneously as a consequence of double purely imaginary eigenvalues forced by symmetry. One is a branch of standing waves and the other is a branch of rotating waves (i.e., travelling pulses on a circle); see Golubitsky and Stewart (2002) for details. (Note that these standing waves and rotating waves correspond to the cases of double mills and single mills obtained with individual-based models on periodic-doughnut shaped-domains Chuang et al. 2007;Carrillo et al. 2009). The non-existence of O(2) Hopf bifurcations strongly restricts the possibility of the models to support patterns of moving animal aggregations. See the end of this section for further discussion.
We now consider the problem of steady-state bifurcations. Using formula (17), we write the linearization for all models where j corresponds to the model chosen with the appropriate K (s), and q = (q r , q a ) or q = (q r , q a , q al ). The linearised equations can exhibit steady-state bifurcations when the following dispersion relations are zero: In addition, we now show that two zero eigenvalues can cross the imaginary axis and thus two modes, corresponding to two different wave numbers k n and k m , can become unstable at the same time. Thus, n : m steady-state/steady-state mode interactions (with n = m) can occur in model (8). Solving for D in Eq. (18), and simplifying the resulting expression yields This condition determines whether n : m mode interactions exist for a given model Mj, with j = 2, 3, 4, 5. Figure 3 shows the dispersion relations (18a) for models M2 and M4, corresponding to steady-state bifurcations (panel a) and steady-state/steady-state bifurcations (panel b). (Similar dispersion relations can be obtained for models M3 and M5-not shown here.) Therefore, the nonlocal parabolic model (8)  mixed-mode solutions appear where the amplitude of the aggregation peaks is modulated along the spatial domain; see Buono and Eftimie (2014b) for examples in the context of model (1).

Existence of rotating waves (travelling pulses) As we show in Appendix B, O(2)equivariant scalar parabolic equations with nonlocal terms do not support O(2)
Hopf bifurcations from a homogeneous equilibrium, and therefore, it is likely that rotating waves (also called "travelling pulses") do not exist for O(2)-equivariant equations. For O(2)-equivariant scalar reaction/advection-diffusion equations with local terms only and satisfying periodic boundary conditions, Matano (1988) showed the nonexistence of rotating waves using a Lyapunov functional approach. Fiedler and Mallet-Paret (1989) showed the same non-existence result as a consequence of a Poincaré-Bendixson theorem.
The question of rotating waves for nonlocal parabolic equations was investigated in Mogilner and Edelstein-Keshet (1999), where the authors studied a family of parabolic equation models (for locusts migration) with nonlocal advection terms defined on R. There, they considered even and odd kernels, and investigated the velocity at the front and back edges of a rectangular-shaped swarm to gain insight into the role of these kernels on the movement of aggregations. In the case of only odd kernels, one can easily show-using the symmetry methods above-that the model in (Mogilner and Edelstein-Keshet 1999) has the Euclidean symmetry E(1) generated by translations x → x + y and the reflection x → −x. With periodic boundary conditions on an interval [−L 0 , L 0 ] their system becomes O(2)-equivariant (in fact, they showed that odd kernels lead to stationary aggregations, which is similar to the situation we describe above). But to induce traveling aggregations, the authors needed to introduce either an even kernel or a local drift term. In both cases, this leads to a forced-symmetry breaking where the reflection symmetry is lost and the equation is then only symmetric with respect to translations on R (if periodic boundary conditions are imposed this means breaking the symmetry from O(2) to SO(2)).
Because systems of parabolic equations usually display a richer dynamics in terms of pattern formation compared to a single parabolic equation, next we focus on models that incorporate two interacting populations (each communicating via a different mechanism).

Two populations model
Eftimie (2013) generalised model (1) to the case when the animal population uses two different communication mechanisms at the same time. This splits the population in two subpopulations, u and v (each communicating via one mechanism). The evolution of densities of left-moving (u − , v − ) and right-moving (u + , v + ) individuals belonging to the two subpopulations, u and v, is described by the following equations: The turning functions are defined in a similar manner as for model (1): with The index, u or v, helps us distinguish the turning interactions (as determined by the communication mechanisms) for population u, from the turning interactions for population v. Note that, to make the model more general, we assume that the two sub-populations have different turning rates λ u,v 1 and λ u,v 2 . In Eftimie (2013), these turning rates incorporated equal strengths of interactions (q a , q r and q al ) and equal interaction ranges (s r , s a and s al ) for both populations u and v. We show below that unequal strengths of interaction or unequal interaction ranges is necessary to obtain periodic dynamics via Hopf bifurcations.
The social interaction terms y u,v,± r (repulsion), y u,v,± a (attraction) and y u,v,± al (alignment), which correspond to the two sub-populations, are described in Table 2. We assume here that the two populations react the same way towards their own members, as well as towards the members of the other sub-population. This assumption is modelled mathematically by the sums (u +v)(x ±s) that appear in y y In Eftimie (2013), this model was investigated only in terms of the spatial and spatiotemporal patterns exhibited. Here, we take the parabolic limit of model (20)-(23), and transform it into a parabolic model for two coupled populations. We note here that the assumptions of unequal turning rates and unequal social interactions that we now incorporate into (20)-(23) break the symmetry of the models, for both the hyperbolic model and its parabolic counterpart. As we see in Sects. 3.3 and 4, this symmetry-breaking aspect has implications on the types of spatio-temporal aggregation patterns obtained.

Formal parabolic limit
As for the one-population model (1), we first consider the formal parabolic limit of model (20) to transform this into a two-population parabolic model. To this end, define ∂t ds. Following the same steps as in Sect. 2.1, we arrive at the following two-dimensional parabolic system: with the diffusion and advection parameters defined as follows: Nondimensionalization of the parabolic system. To simplify the analysis of system (24), we nondimensionalize it. To this end, we introduce the following nondimensional variables and parameters: After dropping the "*", we obtain the following nondimensionalized system: Note that the two parabolic equations are coupled via the advection terms, which depend on the communication mechanisms employed by the two sub-populations. As an example, let us assume that population u communicates with its neighbours via M2 or M4, and population v communicates with its neighbours via M3 or M5. In the limit → 0, with M2, M3 :

Symmetries for the two-populations hyperbolic and parabolic models
The action of O(2) on functions extends diagonally to the second population v. For θ ∈ [0, L 0 ], we have We now have the following result: Table 2 is O(2)-symmetric with respect to the action (29).

then (20) is also symmetric with respect to the order two interchange (u, v) → (v, u). The symmetry group is then
The proof is nearly identical to the one-population case done in Buono and Eftimie (2014b) and we omit it here. As in the one population case, the O(2) symmetry is also inherited by the parabolic limit. The O(2) action is given by Corollary 3.

The same calculation holds for
The invariance of (24) with respect to the θ action is straightforward to show by letting z = x − θ and using the translation invariance of the integrals. We have thus shown the O(2)-equivariance.
If λ ± u = λ ± v , then B u = B v and D u = D v , and the non-dimensional advection and diffusion parameters are B = 1 and D = 1. Therefore u and v can be interchanged, and the symmetry group becomes O(2) × Z 2 .

Local bifurcation analysis for the limiting parabolic model
In Eftimie (2013), it was shown numerically that the hyperbolic model (20) can exhibit rotating waves (i.e., travelling pulses) and stationary pulses, associated with Hopf and steady-state bifurcations. Since the goal of this article is to investigate the bifurcation dynamics of the limiting parabolic models, we focus next on the local bifurcation from the spatially homogeneous steady states corresponding to model (26). In particular we focus on the conditions for the existence of codimension-1 Hopf and steadystate bifurcations, which allow for the existence of rotating waves/standing waves and stationary pulses, respectively. (A detailed linear stability analysis and bifurcation analysis of the hyperbolic model (20) is the subject of an upcoming article.) The parabolic system (26) has a spatially homogeneous equilibrium (u * , v * ). Linearising equations (26) about this equilibrium point leads to By substituting u(x ± s) = u * , v(x ± s) = v * into the integral terms described in Table 2, we obtain f (y u,v, , and thus the last two terms in the linearised equations are cancelled: To investigate the linearized system (32), we use group-theoretic arguments, as in Golubitsky and Stewart (2002). We write (32) as We now consider the spaces C k For a quick reference to the group-theoretic concepts used below, see Appendix A. In Golubitsky and Stewart (2002), it is shown that the subspaces V n = {ae ik n x + c.c | a ∈ C, k n = 2π n/L 0 }, with "c.c" standing for "complex conjugate", are O(2)-irreducible. We consider subspaces X n = {e ik n x a + c.c. | a ∈ C 2 } and a simple calculation shows that if f 1 = (1, 1) T and f 2 = (1, −1) T then where both subspaces on the right are isomorphic to V n . X n is a O(2)-isotypic component of the O(2) action on the space C 2 L 0 . Note that V k , V are non-isomorphic irreducible representations of O(2) if k = . This implies a block diagonalization of L along the X n 's. and we define L n := L| X n . Thus, we compute the linearization operator Here, Q u,v,± depend on the communication mechanisms:  Therefore, The eigenvalues of the 2 × 2 matrix L n (33) are solutions of the characteristic equation Note that b u and b v are purely imaginary terms and we rewrite (34) corresponds to the dispersion relation obtained from (32) by setting u, v ∝ a 1,2 e σ t+ik n x , and setting up the algebraic linear system for a 1 , a 2 .

Equation (34) has a zero eigenvalue if
and this equation can be satisfied forb u andb v having any sign. We now show that this is not the case if one is searching for purely imaginary eigenvalues of (34). Investigating the roots of the dispersion relation (34) for the existence of Hopf bifurcations, reduces to solving the following system: Simple algebra shows that: • for D = 1: if P 1 = 0 then P 2 = −k 4 n < 0, and thus no Hopf bifurcations are possible; Solving (35a) for k n , and then substituting k n into the expression for P 2 leads to This expression is positive (i.e., Hopf bifurcations exist) only if Inequalities (36) and (38) hold true at the same time only ifb u andb v have opposite signs. This cannot be possible if we assume that the two communication mechanisms used by populations u and v lead to the same social interactions, namely same q j and s j , j = r,al,a (see the discussion at the beginning of Sect. 3, and the equations for b u and b v ). Therefore, if the social interactions are the same, system (26) can exhibit only steady-state (ss) bifurcations.
We have the following result: Let us now denote by q u,v j , j = r,al,a, the magnitudes of the social interactions for populations u and v, and s u,v j , j = r,al,a, the interaction ranges for the two populations. Thenb u is defined in terms of q u j , and s u j , whileb v is defined in terms of q v j and s v j , j = r,al,a. To investigate the (D, B) parameter space where Hopf bifurcations can occur, we focus now on model M2&M4 (i.e., population u communicates via M2, while population v communicates via M4). Figure 5a shows the (D, B) parameter space (corresponding to region R 1 in Fig. 4 and P 2 (k 1 ) > 0. In this parameter space, time-periodic solutions can occur. Figure 5b shows the real and imaginary parts of the dispersion relation σ (k n ), for some parameter values where such Hopf bifurcations can occur.
The critical eigenspace corresponding to a purely imaginary eigenvalue iω is obtained by solving where "c.c." stands for complex conjugate. Figure 6 shows contourplots of eigenfunctions corresponding to eigenvalue iω, for the sum of u and v populations. Finally, steady-state/steady-state (ss/ss) bifurcations can occur at the intersection of the neutral stability curves σ (k n 1 ) = 0 and σ (k n 2 ) = 0, for some wavenumebrs k n 1 and k n 2 . Figure 7 shows ss/ss bifurcations that occur for k 3 &k 4 (panels a, c) and for k 1 &k 2 (panels b, d). For model M2&M3 (panel c) we observe that there is an entire line of ss/ss bifurcations (where the curves σ (k 3 ) = 0 and σ (k 4 ) = 0 overlap). On the other hand, model M2&M4 displays a ss/ss bifurcation at (q a , q r )=(1.5, 0.9) (panel d). Note that models M3&M4 and M2&M5 also display a ss/ss bifurcation near the same bifurcation point (q a , q r )=(1.5, 0.9) (not shown here).  Some of the patterns that can be obtained with the nonlocal hyperbolic systems (1) and (20) include stationary pulses, travelling pulses, zigzags, breathers, and we refer the reader to Eftimie et al. (2007a) and Eftimie (2013) for the numerical simulations of these patterns. To discretise the parabolic equations (8) and (24), we use a time-splitting approach (Press et al. 2007). In particular, to solve numerically u t + A[u] + B[u] = 0, with A[u] the advection operator and B[u] the diffusion operator, we first discretise the advection equation (u t + A[u] = 0) using the McCormack method, with periodic boundary conditions. Then, we discretise the diffusion equationū t + B[ū] = 0 with the help of the Crank-Nicholson method with periodic boundary conditions. Here, u denotes the solution of the advection equation solved at the previous step. The integrals are discretised using Simpson's method. The initial conditions are random perturbations of spatially homogeneous steady states u * = v * = c, where c = 1 or c = 2. For the numerical simulations presented in this section, we focus on parameter values that ensure unstable spatially homogeneous steady states. These parameter values (identified via the stability analysis performed in Sects. 2.2.4 and 3.3) lead to the formation of spatial heterogenous patterns (i.e., stationary and travelling pulses, corresponding to resting and travelling aggregations). Figure 8 shows the stationary pulses obtained with the one-population parabolic model (8), in the parameter region where the spatially homogeneous steady states u * = 1 and u * = 2 become unstable via steady-state bifurcations. Figure 9 shows the stationary pulses obtained with the two-population parabolic model (24), in the parameter region where the spatially homogeneous steady state u * = 1 is unstable via steady-state bifurcations. For this case, we do not show simulations obtained via random perturbations of the steady state u * = 2, since the resulting patterns are similar to those in Fig. 8b.
Note that the aggregation patterns shown in Fig. 9, occur at the same position in space for both the u and v populations (i.e, u and v aggregate at x ∈ (1.8, 2.5), x ∈ (4.5, 5.2) and x ∈ (7.2, 8.2)). This is in contrast with the type of stationary aggregations shown in Fig. 10a, a', a", and which are obtained via Hopf bifurcations. Here, populations u and v aggregate at different positions in space: population u can be found near x = 0.9 and x = 9, as well at x ∈ (2.1, 3.1), x ∈ (4.5, 5.5), x ∈ (7, 7.8); population v can be found at x ∈ (1.2, 8.8) with higher densities at x ∈ (1.2, 2.2), x ∈ (3.2, 4.4), x ∈ (5.5, 7) and x ∈ (7.8, 8.8). Figure 10b, b', b" shows rotating waves (travelling pulses) that emerge (through primary bifurcations) from unstable spatially homogeneous steady states. Figure 10c, c', c" shows modulated rotating waves that emerge (through secondary bifurcations) from stationary pulses (which could be unstable). All rotating and modulated rotating waves are unstable, being transient solutions. We stress that in the parameter range where Hopf bifurcations are possible, one could obtain more complex spatial and spatio-temporal patterns (e.g., standing waves characterised by left-moving and right-moving pulses that pass through each other-see Buono and Eftimie (2014a), Buono and Eftimie (2014b); these waves are equivalent to the "double mills" patterns shown by other models for swarming Carrillo et al. 2009). However, the goal of this study is not to investigate all these patterns, but rather to investigate the conditions under which stationary and travelling patterns can occur. For this reason, we also chose not to compare the patterns obtained with the hyperbolic and parabolic models.

Summary and discussion
Pattern formation in self-organised living organisms is a very active research area, with many studies focused on deriving new mathematical models that can reproduce the observed complex spatial and spatio-temporal patterns displayed by these organisms. The difference between these models makes it difficult to compare them in terms of pattern formation. Previous studies on pattern formation in biology have suggested that hyperbolic models could exhibit richer dynamics compared to parabolic models. Here, we investigated this aspect in a class of hyperbolic and parabolic models for selforganised aggregation and movement. To this end, we started with a previously-derived class of nonlocal hyperbolic models that incorporated social interactions mediated by various communication mechanisms. It was previously shown that the use of different communication mechanisms-different ways of receiving/emitting information from/to neighbours-was associated with a large variety of complex spatio-temporal patterns exhibited by these hyperbolic models (Eftimie 2012). To investigate whether these complex patterns could be obtained also with parabolic models, we rescalled the turning rates and the speed and took a formal asymptotic limit to obtain a class of nonlocal parabolic models for self-organised biological aggregations. We performed this parabolic limit for one-population models (where all individuals used the same communication mechanism) and two-population models (where the two populations used different communication mechanisms). While the hyperbolic models incorporated four different communication mechanisms, in the parabolic limit these communication mechanisms were averaged, and individuals interacted only via two possible communication mechanisms. This averaging reduced the types of patterns that one could expect to see in the nonlocal parabolic models.
Furthermore, to investigate these limiting parabolic models, we considered an approach that, to our knowledge, has not been used before to compare and contrast pattern formation in models for collective movement. Namely, we combined linear stability theory, bifurcation theory and symmetry theory. We showed that while both the hyperbolic and parabolic models are O(2) symmetric, the bifurcation dynamics of these models is quite different. In two previous studies (Buono and Eftimie 2014a, b) we showed that the hyperbolic models could exhibit all codimension-1 (Hopf and steady-state) and codimension-2 (Hopf/Hopf, Hopf/steady-state, steady-state/steadystate) bifurcations. Here, we showed that the corresponding parabolic (one-equation and two-equation) models could only exhibit steady-state and steady-state/steadystate bifurcations. Thus, one of the mechanisms for obtaining rotating waves in these parabolic equations, namely Hopf bifurcations with O(2) symmetry, is lost in the parabolic limit. These analytical results were also supported by numerical simulations. To obtain Hopf bifurcations for systems of parabolic equations (but not for single parabolic equations), we needed to impose an extra assumption, namely that different communication mechanisms influenced differently the strength and spatial range of social interactions. Given the heterogeneity of large animal aggregations and the different types of stimuli (visual, chemical, auditory, ...) that individuals use to communicate with each other-stimuli that act on different spatial ranges (Bradbury and Vehrencamp 2011)-we believe that our assumption of different social interactions is biologically realistic. This assumption, which breaks the O(2) × Z 2 symmetry of the models, allows the system of parabolic equations to exhibit rotating waves and modulated rotating waves. However, this assumption increased the complexity of the model, by introducing 6 new parameters (q j , s j , j = r,a,al). Since the purpose of this study was not to describe all possible patterns obtained with nonlocal parabolic equations, but rather to find the conditions under which complex patterns can arise, we did not investigate rigorously all types of spatial and spatio-temporal patterns that could be exhibited by these models.
Even if the parabolic models investigated in this study could not exhibit Hopf bifurcations from spatially homogeneous steady states, theoretically it could be possible for the models to display rotating waves (travelling pulses). The existence of such rotating waves was previously investigated in local parabolic models with O(2) symmetry, and it was shown that such waves could not be exhibited by these local models (Matano 1988). However, for nonlocal parabolic models with O(2) symmetry this is still an open problem. Because of the complexity of the nonlocal social interactions that form the advection term, addressing this problem requires detailed asymptotic investigation of the dynamics of model (8), which is beyond the scope of this article.
Currently, there are not many studies that compare the types of patterns obtained with different models for self-organised biological aggregation and movement. We believe that such studies could help researchers make decisions in regard to what types of models to use when describing particular biological phenomena. Overall, the results in this article support the idea that if we want to model certain biological systems that exhibit complex spatial and spatio-temporal patterns, the use of classical parabolic equations, even with nonlocal terms, does not always lead to the types of patterns that one hopes to see numerically (to match the experimental observations) and this might be a consequence of unaccounted symmetries of the model which precludes bifurcations to those patterns or even their existence. referees and editors for their careful reading of the manuscript and their suggestions which helped improve the scope and readability of the paper.
We now describe the spectrum of the linearization of equation (39). The function u = u * , where u * is a constant, is a homogeneous O(2)-symmetric equilibrium for (39) if h vanishes at the point U * = ( p * 0 , p * 1 , r * 1 , . . . , p * m , r * m , s 0 , s * 1 , t * 1 , . . . , s * n , t * n ), where p * 0 = u * , p * j = p j | u=u * ,r * j = r j | u=u * for j = 1, . . . , m, and s * 0 = s * 1 = t * 1 = · · · = s * n = t * n = 0. The linearization of the operator that appears in the right-hand-side of (39), calculated at u = u * , has the form Theorem 5 Note that K − i (k ) = K + i (k ) and J − i (k ) = J + i (k ). We derive a more explicit form of ∇h(U * ) using the expression in Proposition A. Evaluating the partial derivatives of h and g at U * , we obtain for i = 1, . . . , m and j = 1, . . . , n h p 0 = g p 0 , h p i = g M i , h r i = g M i , h s 0 = 0, h s j = g G j , h t j = −g G j .
Therefore, L (ae ik x ) = (−k 2 ae ik x + (g p 0 , g M 1 , g M 1 , . . . , g M m , g M m , 0, g G 1 , −g G 1 , . . . , g G n , −g G n )V (ae ik x ) Using (41), we see from here that the coefficient of ae ik x must be real. Hence all eigenvalues of L are real.

Appendix B Basic notions of equivariant bifurcation theory
For the convenience of readers, we have collected basic definitions of bifurcation theory of dynamical systems with symmetry (a.k.a equivariant bifurcation theory). The set of translations on the real line given by operators T y (with y ∈ R), defined by T y (x) = x + y, form a group. If one considers a periodic lattice on R, say of length L 0 , translations satisfy T y+L 0 = T y . Therefore, we can only consider the translations T θ for θ ∈ [0, L 0 ], which in fact correspond to the group of translations on the circle. This group is denoted by S 1 , and it has the same algebraic structure as the group SO(2) of 2 × 2 orthogonal matrices with determinant 1; this latter terminology is the one used predominantly. The group consisting of SO(2) along with a reflection symmetry κ satisfying T θ • κ = κ • T −θ , has the same algebraic structure as the group O(2) of 2 × 2 orthogonal matrices. (Note that the symmetry κ depends on the context, as seen in the various models discussed in this paper.) To keep with the tradition, we use the O(2) notation to describe all those groups with identical algebraic structure.
Consider now a partial differential equation for u(x 1 , . . . , x d , t) ∈ R m given by u t = F(x 1 , . . . , x d , u x 1 , . . . , u x d , u x 1 x 1 , . . . , u x 1 ,x Let be a group and u(x 1 , . . . , x d , t) be any solution of (42). Then, (42) is said to be -symmetric or -equivariant if γ.u is also a solution of (42), for all γ ∈ . To each point u ∈ C is attached its symmetry group u := {γ ∈ | γ.u = u}, also called the isotropy subgroup.
One benefit of using symmetry concepts is that it often enables simplifications in the computation of the spectrum of linear operators, via the following properties of group actions and equivariant differential equations.
A subspace V ⊂ C is -invariant if γ.v ∈ V , for all γ ∈ and v ∈ V . The subspace V is a -irreducible representation if it is -invariant and does not contain any proper -invariant subspace. If is a compact Lie group (such as SO(2) and O(2)), then all -irreducible representations are finite dimensional. Two -irreducible representations V 1 , V 2 are called -isomorphic if there exists a nonsingular linear map A : V 1 → V 2 , which commutes with the actions of on V 1 and V 2 respectively. For compact Lie groups, C has a (non-unique) direct sum decomposition into irreducible representations. The direct sum of isomorphic irreducible representations in the decomposition is called a -isotypic component and the decomposition of C into a direct sum of its isotypic components is called the -isotypic decomposition.
A linear mapping L : C → C is -equivariant if L(γ .u) = γ.L(u), for all u ∈ C and γ ∈ . Note that the linearization L of (42) at u * with u * = is always -equivariant. An important result is that L maps -isotypic components V k to themselves; that is, L : V k → V k and so one can "block diagonalize" L along its isotypic decomposition: where L k := L|V k . Therefore, the point spectrum of L can be obtained as the union of the point spectra of L k .