Martin Parisot
email: martin.parisot@inria.fr

Miroslaw Lachowicz
email: lachowic@mimuw.edu.pl

Miros Law Lachowicz

A Kinetic Model for the formation of Swarms with nonlinear interactions

Keywords: 2010 Mathematics Subject Classification. 47G20, 76N10, 78M35, 82C22, 92C15, 92D25 kinetic formulation, collective behaviour, swarm, self-organisation, orientation interaction, hydrodynamic limit

The present paper deals with the modeling of formation and destruction of swarms using a nonlinear Boltzmann-like equation. We introduce a new model that contains parameters characterizing the attractiveness or repulsiveness of individuals. The model can represent both gregarious and solitarious behaviors. In the latter case we provide a mathematical analysis in the space homogeneous case. Moreover we identify relevant hydrodynamic limits on a formal way. We introduce some preliminary results in the case of gregarious behavior and we indicate open problems for further research. Finally, we provide numerical simulations to illustrate the ability of the model to represent formation or destruction of swarms.

1. Introduction. We are going to model the swarming behavior of an individual population. Let f (t, x, v) be a probability density of individuals at time t ≥ 0 and position x ∈ R d with velocity v ∈ V ⊂ R d . We assume that V, the set of velocities of the individuals, is compact. We model the evolution of populations at the mesoscopic scale by the nonlinear integro-differential Boltzmann-like equation:

∂ t f (t, x, v) + v • ∇ x f (t, x, v) = 1 ε Q [f ] (t, x, v) = 1 ε V (T [f (t, x, .)] (w, v) f (t, x, w) -T [f (t, x, .)] (v, w) f (t, x, v)) dw (1)
with the initial data f (0, x, v) = F (x, v) and the parameter ε that corresponds to the Knudsen number, see [START_REF] Othmer | Models of dispersal in biological systems[END_REF].

The operator Q describes interactions between individuals. We only consider local (in space) interactions that is physically justified at least for large characteristic space scale. Several models of swarming with the long range interactions described by the mean field approach, can be found in the literature, see for example [START_REF] Degond | Continuum limit of self-driven particles with orientation interaction[END_REF][START_REF] Motsch | A new model for self-organized dynamics and its flocking behavior[END_REF][START_REF] Carrillo | A new interaction potential for swarming models[END_REF]. The turning rate T [f ] (v, w) measures the probability for an individual with velocity v to change velocity into w. Linear turning rates with influence of the orientation of individuals were analyzed in details in [START_REF] Geigant | Bifurcation analysis of an orientational aggregation model[END_REF][START_REF] Kang | An integro-differential equation model for alignment and orientational aggregation[END_REF]. In this paper we focus on the following general nonlinear case

T [f (t, x, .)] (v, w) = σ ρ,x β ρ,x (v, w) f (t, x, w) γρ,x (2)
with the macroscopic density of individuals ρ (t, x) = V f (t, x, v) dv. For any ρ and at any x, the interaction rate β ρ,x : V 2 → R + , the attractiveness coefficient γ ρ,x ∈ R + , and σ ρ,x ∈ {-1, 1} characterize the interaction between the individual agents. The interaction rate β ρ,x corresponds to the tendency of individuals to switch to a different velocity. In the following we assume this interaction rate is symmetric, positive, bounded and separated from zero. The attractiveness coefficient γ ρ,x corresponds to the attractiveness or repulsiveness of individuals. Accordingly to observations, the attractiveness depends on the size of the population [START_REF] Simpson | A behavioural analysis of phase change in the desert locust[END_REF], on the local resources [START_REF] Simpson | A comparison of nutritional regulation in solitarious-and gregarious-phase nymphs of the desert locust Schistocerca gregaria[END_REF] or even morphological adaptation [START_REF] Greenwood | Differences numbers of sensilla on the antennae of solitarious and gregarious locusta migratoria l. (orthoptera: Acrididae)[END_REF]. Note that the biologically relevant interactions are such that σ ρ,x = 1. However, in the present paper we are going to study also the negative interactions, i.e. σ ρ,x = -1, because of interesting mathematical properties. A simpler framework with only two possible velocities +1, -1 was considered in [START_REF] Arlotti | On a discrete Boltzmann-type model of swarming[END_REF][START_REF] Banasiak | On a macroscopic limit of a kinetic model of alignment[END_REF]. For simplicity of notation, we do not indicate the density and the position dependence of the collision parameters, as long as this dependence is obvious.

The swarm behavior, for some populations also called herds, flocks, packs, schools, or shoals, is referred to as the self-organization of individual agents. There exists a huge literature related to swarm phenomena. Here we mention few examples and refer to the bibliography therein. Ref. [START_REF] Eftimie | Hyperbolic and kinetic models for self-organized biological aggregations and movement: a brief review[END_REF] reviews hyperbolic and kinetic models for self-organized biological aggregations and traffic-like movement. Book [START_REF] Pareschi | Interacting Multiagent Systems: Kinetic Equations and Monte Carlo Methods[END_REF] provides the mathematical modeling based on a mesoscopic description and the construction of efficient simulation algorithms by Monte Carlo methods for collective phenomena and self-organization in systems composed of large numbers of individuals. Paper [START_REF] Oien | Daphnicle dynamics based on kinetic theory: an analogue-modelling of swarming and behaviour of daphnia[END_REF] deals with plasma kinetic theory to derive the corresponding hydrodynamic equation for the density of Daphnicle. An interesting agent-based stochastic model of vortex swarming in Daphnia has been proposed in [START_REF] Mach | Modeling vortex swarming in daphnia[END_REF]. A cellbased model has been considered in [START_REF] Wu | Social interactions in myxobacterial swarming[END_REF] and the effect of social interactions between cells has been described. The model quantifies the contribution of individual motility engines to swarming. Reference [START_REF] Zohdi | Mechanistic modeling of swarms[END_REF] deals with mechanistic modeling of swarms, the properties of swarm models and the corresponding numerical algorithms. In paper [START_REF] Grünbaum | Non-linear advection-diffusion equations approximate swarming but not schooling populations[END_REF] a numerical scheme has been developed to estimate coefficients in nonlinear advection-diffusion equations from individual based model simulations. The biophysical principles that cause the Proteus mirabilis the swarm phenomena are given in [START_REF] Frénod | An explanatory model to validate the way water activity rules periodic terrace generation in proteus mirabilis swarm[END_REF]. A swarming model on a two-dimensional lattice, where the selfpropelled particles exhibit a tendency to align ferromagnetically, has been studied in [START_REF] Peruani | Traffic jams, gliders, and bands in the quest for collective motion of self-propelled particles[END_REF]. Paper [START_REF] Peruani | Collective motion and nonequilibrium cluster formation in colonies of gliding bacteria[END_REF] has shown that the transition to collective motion in colonies of gliding bacterial cells confined to a monolayer appears through the organization of cells into larger moving clusters. The possible approaches for swarming behavior in locusts have been described and discussed in [START_REF] Edelstein-Keshet | Do travelling band solutions describe cohesive swarm? An investigation for migratory locusts[END_REF]. The paper [START_REF] Erban | From individual to collective behaviour of coupled velocity jump process: a locust example[END_REF] deals with a stochastic individual-based models related to the locusts behavior. In paper [START_REF] Ha | From particle to kinetic and hydrodynamic description of flocking[END_REF] the Cucker-Smale particle model, [START_REF] Cucker | Emergent behavior in flocks[END_REF], for flocking has been discussed and conditions for flocking have been stated. A Vlasov-type kinetic model has been derived and time-asymptotic flocking behavior has been proved. The asymptotic behavior of solutions of the Cucker-Smale model has been studied in [START_REF] Carrillo | Asymptotic flocking dynamics for the kinetic Cucker-Smale model[END_REF]. Reference [START_REF] Carrillo | Double milling in self-propelled swarms from kinetic theory[END_REF] presents the kinetic theory approach for swarming systems of interacting, self-propelled discrete particles. The related macroscopic hydrodynamic equations are derived from the Liouville equation. General solutions include flocks of constant density and fixed velocity and other non-trivial morphologies.

σ = -1 σ = 1 Figure 1.
Properties accordingly to the parameters σ and γ. 1. classification of the evolution system 2. main properties in the space homogeneous case 3. formal hydrodynamic limit

The present paper deals with the mathematical and numerical study of the orientational aggregation. The main novelty with respect to the previous studies [START_REF] Geigant | Bifurcation analysis of an orientational aggregation model[END_REF][START_REF] Kang | An integro-differential equation model for alignment and orientational aggregation[END_REF] is that we consider non-linear turning rate characterized by the coefficient of attractiveness which control the behavior of the solutions. In order to simplify the study, we consider that individuals directly switch to their objective, i.e. using the notation defined in [START_REF] Kang | An integro-differential equation model for alignment and orientational aggregation[END_REF] in the case of deterministic turning G σ = δ 0 (Dirac function) and with optimal reorientation M w (v) = w. The paper is organized as follow. In Section 2, we perform a mathematical analysis of Eq. ( 1) in the space homogeneous case. The main results of this section are summarized in Figure 1. We will see that the behavior of the solutions is mainly controlled by the sign of σ (1 -γ). In the case of σ (1 -γ) > 0, the evolution of the distribution function satisfies an entropy dissipation law and a maximum principle in the velocity domain. From the biological point of view, individuals attempt to disperse and such a behavior is called solitarious, see [START_REF] Simpson | A behavioural analysis of phase change in the desert locust[END_REF][START_REF] Simpson | A comparison of nutritional regulation in solitarious-and gregarious-phase nymphs of the desert locust Schistocerca gregaria[END_REF][START_REF] Greenwood | Differences numbers of sensilla on the antennae of solitarious and gregarious locusta migratoria l. (orthoptera: Acrididae)[END_REF]. In the case of σ (1 -γ) < 0 the evolution of the distribution function is self-organized, i.e. individuals attempt to align themselves along at least one direction. This case corresponds to the so-called gregarious behavior, cf. [START_REF] Simpson | A behavioural analysis of phase change in the desert locust[END_REF][START_REF] Simpson | A comparison of nutritional regulation in solitarious-and gregarious-phase nymphs of the desert locust Schistocerca gregaria[END_REF][START_REF] Greenwood | Differences numbers of sensilla on the antennae of solitarious and gregarious locusta migratoria l. (orthoptera: Acrididae)[END_REF]. The main open questions concern this last case. We are currently not able to predict eventual blow up in finite time for positive self-organized interactions.

In Section 3, we formally describe the spacial hydrodynamic limits in the spacedependent case for the two behaviors. For this purpose, we assume that the main term of distribution, so-called local equilibrium, corresponds to the stable steady state of the homogeneous case. In Section 3.1 we focus on the interactions satisfying the dissipation law. We show that the only stable steady state is the constant function in the whole velocity domain, the so-called diffusive picture. The formal hydrodynamic limit leads to a parabolic equation with a nonlinear coefficient proportional to ρ -γ characteristic to the Carleman-like model, see [START_REF] Salvarani | The diffusive limit for carleman-type kinetic models[END_REF]. We estimate the parameters of the macroscopic model, namely the diffusion parameter for some given velocity domain. The mathematical derivation of an hydrodynamic limit of our model for self-organized interactions is currently an open problem. However, we propose a conjectured solution, namely the aligned picture, based on the asymptotic solution of the space homogeneous equation and in the case of constant interaction rate β.

Finally in Section 4 we present various numerical strategies to solve Eq. ( 1). First, we focus on the interaction operator in Section 4.2. In order to numerically observe the different solutions in the non-Lipschitzian case 0 ≤ γ < 1, we propose several scheme based on linear and nonlinear approximation of the space homogeneous equation. Numerical simulations illustrating the results are presented in Section 5.1. Then, we propose a complete discretization of the PDE (1). The conjectured solution in the gregarious case, namely the aligned picture, is not recover using classical transport scheme due to numerical diffusion. We propose in Section 4.3 an anti-diffusive scheme for Cartesian grid able to recover the aligned picture and we perform numerical simulations in Section 5.2.

2. Mathematical analysis of space homogeneous case. In this section we provide the mathematical analysis of Eq. ( 1) with (2) in the space homogeneous case, i.e. when all functions are x-independent. In addition, we simplify the study by considering restrictive assumptions on the interaction rate Assumption 1. We assume that β is a positive, symmetric and bounded function separated from zero, i.e. there exists 0 < µ ≤ 1 such that for almost any (v, w) ∈ V 2 , we have

µ β L ∞ (V 2 ) ≤ β (v, w) = β (w, v) ≤ β L ∞ (V 2 ) .
The following a priori conservation properties hold Theorem 2.1. The total density of individuals is a priori preserved in time, i.e.

∂ t R d ρ dx = 0 then R d ρ (t, x) dx = R d ρ (0, x) dx . Proposition 1. The operator Q is homogeneous of degree γ, i.e. Q [f ] (t, x, v) = ρ γ Q f ρ (t, x, v) .
Without the lost of generality, we may scale the solution such that the macroscopic density of the scaled solution is 1, i.e. g L 1 (V) = 1. The space homogeneous version of Eq. (1) reads

εσ∂ t g = (β * g) g γ -(β * g γ ) g with g (0, v) = G (v) ( 3
)
where * is the convolution product in the velocity domain V, i.e.

(β * φ) (v) = V β (v, w) φ (w) dw .
Note that we can return to the original variable by setting f (t, x, v) = ρg (ρ γ t, v).

In this section we set ε = 1 for simplicity.

Proposition 2. Let η : R + → R be a bounded convex differentiable function. The solution g of Eq.

(3) a priori satisfies

d dt V η (g) dv ≤ 0, σ (1 -γ) > 0 (entropy dissipation law) ≥ 0, σ (1 -γ) < 0 (self-organization).
Proof. To simplify the notation, we set g v = g (t, v) and g w = g (t, w). By direct estimation we have

d dt V η (g v ) dv = σ V V β (v, w) (g γ v g w -g γ w g v ) η (g v ) dw dv = σ V V β (v, w) g v g w g γ-1 v -g γ-1 w η (g v ) dw dv .
Exchanging v and w and using the symmetric property of β (v, w) we obtain

d dt V η (g v ) dv = σ 2 V V β (v, w) g v g w g γ-1 v -g γ-1 w (η (g v ) -η (g w )) dw dv .
Then the conclusion follows since the application x → x γ-1 is decreasing for 0 ≤ γ ≤ 1 and increasing for 1 ≤ γ.

According to Proposition 2, the solutions of Eq. ( 3) have completely different evolution for different signs of σ (1 -γ). In the following, we are going to analyze Eq. (3) considering separately each domain of parameters. Note that γ = 1 leads to the vanishing interactions and to the trivial solution g (t, v) = G (v). More generally, the non-homogeneous solution is in this case simply shifted with the velocity v, and we have f (t, x, v) = F (x -vt, v). From now on, we assume γ = 1.

We are only interested in positive solutions of Eq. ( 3) for physical reasons. We denote by L p + (V), for 0 < p ≤ ∞, the set of nonnegative functions in the L p -space, i.e. f ∈ L p + (V) iff f ∈ L p (V) and f ≥ 0.

2.1. Entropy dissipative interactions: solitarious behavior.

2.1.1. Case σ = -1 and 1 < γ. In this section we consider Eq. ( 3) with σ = -1 and 1 < γ.

Theorem 2.2. Let σ = -1, 1 < γ and β satisfy Assumption 1. Assume G ∈ L ∞ + (V) be such that there exist nonnegative constants 0 ≤ m ≤ M < +∞ satisfying m ≤ G (v) ≤ M ,
for almost any v ∈ V.

Then there exists an unique global solution g ∈ C 1 (R + ; L ∞ (V)) of Eq. (3) with g (0, v) = G (v). Moreover, the solution g satisfies the maximum principle, i.e. m ≤ g (t, v) ≤ M , for any t > 0 and almost any v ∈ V.

Proof. Consider the auxiliary ODE problem defined by

φ = p -qφ α with φ (0) = φ 0 ≥ 0 , (4)
α > 0 and p and q are given continuous functions on R + . According to the Picard-Lindelöf theorem, there exists an unique local solution and it can be prolonged to R + . In addition φ is nonnegative when both p and the initial data φ 0 are nonnegative.

We approximate the solution of Eq. ( 3) by the sequence starting with g 0 (t, v) = G (v) and g n+1 (t, v) defined as the unique solution of the auxiliary problem (4) with α = γ, p = (β * g γ n ) g n , q = (β * g n ) and the initial data φ 0 = G (v) for any fixed v ∈ V. By the properties of the auxiliary problem (4) we deduce that the sequence (g n ) n is well-defined. In addition, the g 0 satisfies the maximum principle since it is constant in time. Let us now assume that the bounds are verified by g n , i.e.

0 ≤ m ≤ g n (t, v) ≤ M , for any (t, v) ∈ R + × V.
Then we consider the next approximation g n+1 . Using the positive part function 2 (φ) + = |φ| + φ, we write

(g n+1 -M ) + ∂ t (g n+1 -M ) = ((β * g γ n ) g n -(β * g n ) M γ ) (g n+1 -M ) + -(β * g n ) g γ n+1 -M γ (g n+1 -M ) + .
We have the following inequality for 1 < γ

(β * g γ n ) = V β (v, w) g n (w) g γ-1 n (w) dw ≤ V β (v, w) g n (w) dw M γ-1 = (β * g n ) M γ-1 . It follows (β * g γ n ) ≤ (β * g n ) M γ-1 .
Since the function x → x γ is increasing and g n is nonnegative, the RHS is nonpositive. With similar arguments, we treat the lower bound.

The final step of the proof consists in the convergence of the sequence (g n ) n . We compare the time derivative of two consecutive approximations and we multiply by the sign of the difference, i.e.

sgn (g

n+1 -g n ) ∂ t (g n+1 -g n ) = 1 2 sgn (g n+1 -g n ) β * g γ n -g γ n-1 (g n + g n-1 ) + β * g γ n + g γ n-1 (g n -g n-1 ) -sgn (g n+1 -g n ) (β * (g n -g n-1 )) g γ n + (β * g n ) g γ n+1 -g γ n Since 0 ≤ γ the last term is nonpositive: -sgn (g n+1 -g n ) (β * g n ) g γ n+1 -g γ n ≤ 0.
In addition, the function g γ n is γM γ-1 -Lipschitz continuous. Then we integrate over the velocity domain and we obtain

∂ t g n+1 -g n L 1 (V) ≤ (γ + 2) M γ β g n -g n-1 L 1 (V) , with β = β (v, w) L 1 (L ∞ ) = V β (v, w) L ∞ ( dw) dv.
We conclude the convergence by classical contraction arguments. In fact, for any 0 ≤ t ≤ T , we have

g n+1 -g n L 1 (V) (t) ≤ (γ + 2) M γ βT n n! g 1 -g 0 L ∞ (0,T ;L 1 (V)) .
It is clear that the RHS vanishes when n goes to infinity.

PARISOT MARTIN AND LACHOWICZ MIROS LAW 2.1.2. Case σ = 1 and 0 ≤ γ < 1. In this section we consider Eq. ( 3) for σ = 1 and 0 ≤ γ < 1. We show that the result of Theorem 2.2 still holds.

Theorem 2.3. Let σ = 1, 0 ≤ γ < 1 and β satisfy Assumption 1. Assume G ∈ L ∞ + (V) be such that there exist nonnegative constants 0 < m ≤ M < +∞ satisfying m ≤ G (v) ≤ M ,
for almost any v ∈ V.

Then there exists an unique global solution g

∈ C 1 (R + ; L ∞ (V)) of Eq. (3) with g (0, v) = G (v). Moreover, the solution g satisfies m ≤ g (t, v) ≤ M ,
for any t > 0 and almost any v ∈ V.

Proof. We approximate the solution of Eq. ( 3) by the sequence starting from g 0 (t, v) = G (v) and g n+1 (t, v) is the unique solution of the auxiliary problem (4) with α = 1, p = (β * g n ) g γ n , q = (β * g γ n ) and the initial data φ 0 = G (v) for any fixed v ∈ V. We deduce from the auxiliary problem (4) that the sequence (g n ) n is well-defined. Then we show the maximum principle with the same strategy than the proof of Theorem 2.2 with the following inequality for 0 ≤ γ < 1

(β * g n ) = V β (w -v) g γ n (v) g 1-γ n (v) dw ≤ V β (w -v) g γ n (v) (v) dw M 1-γ = (β * g γ n ) M 1-γ .
We conclude following the proof of Theorem 2.2.

Note that we considered initial data separated from zero, i.e. 0 < m. The case not separated from zero is more difficult to treat since the RHS is not Lipschitz continuous for 0 ≤ γ < 1. The results related to the sequence g n (the existence and the maximum principle) could be proved as well. Unfortunately, we cannot prove the convergence of the sequence. In fact, for initial data vanishing in some points, σ = 1 and 0 ≤ γ < 1, the solution of Eq. ( 3) is not unique. Let v such that G (v) = 0, we can easily show that for any time t 1 ≥ 0, there exists a solution g of Eq. (3) satisfying

g (t, v) = 0, if t ≤ t 1 > 0, otherwise.
Remark 1 (Biological context). The solitarious case corresponds to repulsive interactions between individuals. In this context, the individuals try to occupy the whole domain. It leads that the relevant solution is the only solution strictly positive for any time strictly positive.

2.2. Self-organized interactions: gregarious behavior.

2.2.1.

Case σ = -1 and 0 ≤ γ < 1. In this section, we consider Eq. ( 3) for σ = -1 and 0 ≤ γ < 1.

Theorem 2.4. Let σ = -1, 0 ≤ γ < 1, β satisfy Assumption 1 and G ∈ L ∞ + (V).
Then there exists a solution of Eq. (3) in C 1 (R + ; L ∞ (V)) and it is unique among the nonnegative functions. Moreover, the solution vanishes in a finite time on a subdomain of V. In fact, for any v ∈ V and t ∈ R + such that

G (v) < µ G γ L γ (V) 1 1-γ and t ≥ - ln 1 - G γ L γ (V) µ G 1-γ (1 -γ) β L ∞ (V 2 ) G γ L γ (V)
, we have g (t, v) = 0.

Proof. Consider the following auxiliary problem

φ = qφ -pφ γ , with φ (0) = φ 0 ≥ 0, (5)
p and q are given non-negative continuous functions on R + . According to the Picard-Lindelöf theorem, there exits an unique local solution of Eq. ( 5) in the neighborhood of any time t such that φ (t) > 0.

For small enough initial data, the solution is not unique in general. More precisely, for some attractiveness coefficients γ, there exist solutions of the auxiliary problem Eq. ( 5), which are not nonnegative. For example, with β (v, w) = 1 and parameters γ = 2 3 , p and q not depending of time, the following function

φ (t) = p q + φ 1 3 0 - p q e q 3 t
3 is a solution of the auxiliary problem and for φ 0 < p q 3

, is initially positive and become negative for a time large enough. However, the vanishing function is a trivial solution of the auxiliary problem with vanishing initial data. Then, for not vanishing initial data, we can extend the local solution φ on R + by zero from the time that it vanishes, i.e. we define

ψ (t) = 0, if there exists 0 ≤ s ≤ t such that φ (s) = 0 φ (t) , otherwise.
One can check that ψ (t) is a nonnegative solution of Eq. ( 5) in C 1 (R + ). In addition, for the solution nonnegative and small enough, i.e. when 0 ≤ φ < p q 1 1-γ , the solution is non-increasing, i.e. φ ≤ 0. We conclude that the solution ψ is the unique global solution in the set of nonnegative functions of Eq. ( 5) with small initial data, i.e. φ 0 < p q 1 1-γ .

For initial data large enough φ 0 ≥ p q 1 1-γ the solution is increasing. We conclude that the solution is unique since the RHS is Lipschitz-continuous. However, it is not globally defined since it can blow up in a finite time. For parameters p and q that do not depend on time, an explicit solution of Eq. ( 5) can be given by considering the ODE satisfies by φ 1-γ . We have

φ (t) =    p q + φ 1-γ 0 -p q e (1-γ)qt 1 1-γ , for t < τ 0, for t ≥ τ with τ =    - ln(1-q p φ 1-γ 0 ) (1-γ)q , if φ 0 < p q 1 1-γ +∞,
otherwise.

The non-negative solution with parameters p and q not depending on time is globally defined in R + for any nonnegative initial data. We approximate the solution of Eq. ( 3) by the sequence starting from g 0 (t, v) = G (v) and g n+1 (t, v) being the unique non-negative solution of the auxiliary problem (5) with p = (β * g n ), q = (β * g γ n ) and the initial data

φ 0 = G (v) for any fixed v ∈ V. Now assume that g n ∈ L ∞ R + ; L ∞ + (V) which is true for g 0 (t, v) = G (v). We have ∂ t g 1-γ n+1 ≤ (1 -γ) β L ∞ (V 2 ) g n γ L γ (V) g 1-γ n+1 -µ .
Using the Grönwall lemma, we conclude that

g 1-γ n+1 ∈ L ∞ R + ; L ∞ + (V) , then g n+1 ∈ L ∞ R + ; L ∞ + (V) . Since g n (t, v) is nonnegative, we show the convergence of the sequence in L ∞ R + ; L ∞
+ (V) using the classical contraction arguments. Since the function x → x γ is concave and by Proposition (2), g L γ (V) decreases. We have

∂ t g = (β * g γ ) g -(β * g) g γ ≤ β L ∞ (V 2 ) g γ L γ (V) g -µg γ ≤ β L ∞ (V 2 ) G γ L γ (V) g -µg γ .
It follows that the solution of the auxiliary problem ( 5) with the functions

p = µ β L ∞ (V 2 ) and q = β L ∞ (V 2 ) G γ L γ (V)
not depending of time, is an upper bound of the solution of Eq. ( 3). Since the upper-bound is globally defined, we conclude that the solution of Eq. ( 3) exists globally in C 1 R + ; L ∞ + (V) and the estimation of the time after which the solution vanishes in a non-negligible subdomain follows.

2.2.2.

Case σ = 1 and γ > 1. In this section we consider Eq. ( 3) with σ = 1 and γ > 1.

Theorem 2.5. Let σ = 1, γ > 1, β satisfy Assumption 1 and G ∈ L ∞ + (V). Then, there exists an unique solution of Eq. (3) in C 1 ([0, T [ ; L ∞ (V)) with the lower bound of the existence time

T ≥ - ln 1 - µ G γ L γ (V) G γ-1 L ∞ (V) (1 -γ) µ β L ∞ (V 2 ) G γ L γ (V)
. Moreover the solution is nonnegative.

Proof. Consider the following auxiliary problem

φ = pφ γ -qφ, with φ (0) = φ 0 ≥ 0. ( 6
)
The functions p and q are given positive and continuous on R + . By the Picard-Lindelöf theorem there exits an unique local solution of (6). Moreover, φ (t) is nonnegative when φ 0 is nonnegative. Then we approximate the solution of Eq. ( 3) by the sequence starting with g 0 (t, v) = G (v) and g n+1 (t, v) is the unique positive solution of the auxiliary problem (6) with p = (β * g n ), q = (β * g γ n ) and the initial data φ 0 = G (v), for any fixed v ∈ V. We show the convergence of the sequence in the neighborhood of the initial data using the classical contraction arguments.

We denote by g the solution of the auxiliary problem ( 6) with the functions

p = β L ∞ (V 2 ) and q = µ β L ∞ (V 2 ) G γ L γ (V)
not depending of time and the initial data φ 0 = G (v), for any fixed v ∈ V. Since the function x → x γ is convex, by Proposition (2), g L γ (V) is increasing. We have

∂ t g = (β * g) g γ -(β * g γ ) g ≤ β L ∞ (V 2 ) g γ -µ g γ L γ (V) g ≤ β L ∞ (V 2 ) g γ -µ G γ L γ (V) g .
It follows that g is an upper bound of the solution of Eq. ( 3). We estimate the solution of Eq. ( 6) for p and q that do not depend on time by setting ψ = φ 1-γ . We obtain

φ (t) =   q pφ γ-1 0 + q -pφ γ-1 0 e (γ-1)qt   1 γ-1 φ 0 ,
and the time of blowup for any initial condition such that φ 0 > q p 1 γ-1 given by

T = - ln 1 -q pφ γ-1 0 (γ -1) q ,
which ends the proof.

Note that assuming the initial data close to an equilibrium is not enough to conclude the global existence. In fact, even in the case β = 1, the lower bound of the existence time tends to infinity iff the initial data satisfies

G γ L γ (V) = G γ-1 L ∞ (V)
with G L 1 (V) = ρ = 1. We will see that this condition corresponds to the steady state of Eq. ( 3), see Proposition 3. The next step consists in identifying the existence time of the solution of Eq. ( 3) with the parameters σ = 1 and γ > 1. We indicate that the mathematical result on the global existence or blowup in a finite time is an open problem. To give a meaning of the model, the solution has to be at least in L γ (V). One may prove the global existence of a similar problem with an additional diffusive operator, modeling the individual reflection, see [START_REF] Ben-Naim | Opinion dynamics: rise and fall of political parties[END_REF][START_REF] Toscani | Kinetic models of opinion formation[END_REF][START_REF] Boudin | A kinetic approach to the study of opinion formation[END_REF].

3. Formal hydrodynamic limits. In this section we are going to discuss a formal derivation of hydrodynamic limits of Eq. ( 1) in the case of a constant attractiveness coefficient, i.e. γ (x, ρ) = γ. This is a preliminary step and the rigorous results are still open problems, in particular in the case of self-organized interactions σ (1 -γ) < 0. The motivation of the hydrodynamic limit is to identify particular asymptotic solutions in the perspective to be compared to numerical solutions, see Section 4. Since the interactions conserve the mass accordingly to Theorem 2.1,

∂ t ρ + ∇ x • (ρu) = 0 with ρ (t, x) u (t, x) = V vf (t, x, v) dv . (7)
In order to estimate the macroscopic velocity u at least for small ε, we expand the solution of Eq. (1) using the Hilbert expansion

f (t, x, v) = ∞ i=0 ε i f i (t, x, v)
and we identify the terms of the same order of ε. The first step of the derivation of the hydrodynamic limit consists in determining the distribution f 0 , the so-called local equilibrium, such that Q (f 0 ) = 0. Equation (1) corresponds to an infinite set of local equilibria, which may lead to different hydrodynamic limits. The local equilibria of Eq. (1) are constant functions on their support, i.e. a nonnegative function

f 0 ∈ L ∞ + (V) satisfies Q (f 0 ) = 0 if and only if there exists a measurable subset W (t, x) ⊆ V such that f 0 (t, x, v) = ρ |W| , if v ∈ W, 0, otherwise where |W| = W 1 dw. Proof. Let I a,b (g) be I a,b (g) = W W β (v, w) g a w g b v dw dv (8)
with W = Supp (g). Using the Cauchy-Schwarz inequality on the product measure yields

I 1,1 (g) = W W β g γ /2 w g 1-γ /2 v β g 1-γ /2 w g γ /2 v dw dv ≤   W W βg γ w g 2-γ v dw dv   1 2   W W βg 2-γ w g γ v dw dv   1 2 ≤ W W βg γ w g 2-γ v dw dv = I γ,2-γ (g) , for any nonnegative function g ∈ L ∞ + (V). Since f 1-γ 0 Q (f 0 ) = 0 on the support of f 0 , we have I 1,1 (f 0 ) = I γ,2-γ (f 0 ).
It follows that the local equilibrium satisfies the equality in the Cauchy-Schwarz inequality, i.e. there exists a constant C such that for almost any (v, w)

∈ V 2 we have β (v, w) f γ 0 (v) f 2-γ 0 (w) = Cβ (v, w) f 2-γ 0 (v) f γ 0 (w)
. We conclude that the condition for an equilibrium solution of Eq. ( 1) to be constant on its support is required. Finally, we easy check that the condition is sufficient.

From now on, we define the relevant local equilibrium of Eq. ( 1) by the steady state of the homogeneous case Eq. (3) which is Lyapunov stable. Such a choice of a local equilibrium is motivated by considering Eq. ( 1) along the trajectory, see [START_REF] Kaniel | The Boltzmann equation[END_REF]

∂ t f # = 1 ε Q # [f ] , where f # (t, x, v) = f (t, x -vt, v) . (9)
Then for the Hilbert expansion we obtain

Q (f 0 + εf 1 ) = O (ε) .
If the equilibrium is not Lyapunov stable, the solutions can leave the vicinity of the equilibrium and the term f 1 is not controllable with respect to ε. In such a case, the Hilbert expansion cannot be justified.

Macroscopic evolution for entropy dissipative interactions.

Proposition 4. Let β satisfy Assumption 1. In both entropy dissipation cases σ (1 -γ) > 0 (i.e. σ = -1, γ > 1 as well as σ = 1, 0 ≤ γ < 1), the only steady state of Eq. (3) that is Lyapunov stable in L ∞ (V) is the constant function in the domain V.

Proof. By Proposition 3, the steady states of Eq. ( 3) are constant functions on their supports. We will show that if the support W of the steady state g 0 is not the whole domain V, then the steady state is not stable. We introduce the perturbation g ∈ L ∞ (V) such that V g dv = 0; the perturbed state g ζ = g 0 + ζ g, with ζ > 0, is nonnegative for ζ small enough and its support is strictly larger than W, i.e. W = Supp (g ζ ) \ W has a nonzero measure. Such a perturbation could always be obtained by setting

g (v) = |W| -|V|, if v ∈ W, |W|, otherwise .
Then for any w ∈ W, we have g 0 (w) = 0, g (w) > 0 and

∂ t g ζ (w) = σ (β * g ζ ) g γ ζ -β * g γ ζ g ζ = σ (β * g ζ ) g γ ζ γ -σ β * g γ ζ gζ + o (ζ) .
When 0 ≤ γ < 1, the first term of the RHS dominates for ζ small. Similarly, when 1 < γ, the second term of the RHS dominates for ζ small. It follows ∂ t g ζ (w) > 0 and then we conclude that the steady state cannot be stable. On the other hand, the maximum principle ensures that the constant function in the whole domain V is stable.

Then, the relevant local equilibrium is given by

f d (t, x, v) = ρ(t,x)
|V| , the so-called diffusive picture -see [START_REF] Banasiak | On a macroscopic limit of a kinetic model of alignment[END_REF]. Introducing f d in Eq. ( 7) we get u = 1

|V| V v dv + O (ε).
The average velocity u is characterized by the set V of possible velocities. For physical reasons, assume that the domain V is symmetric, i.e. for any v ∈ V, we have -v ∈ V. It leads u = O (ε) and then ∂ t ρ = O (ε). From the Hilbert expansion we obtain

D f d Q (f 1 ) = σ (γ -1) f γ d   V β (v, w) dwf 1 -(β * f 1 )   = v • ∇ x f d . ( 10
)
where

D f d Q (f )
is the Fréchet derivative with respect to the equilibrium. Note that V f 1 dv = 0 according to the Hilbert expansion. We denote by

L ∞ 0 (V) the set of function g ∈ L ∞ (V) such that V g dv = 0 and L : L ∞ 0 (V) → L ∞ 0 (V) g 1 → V β (v, w) (g 1 (v) -g 1 (w)) dw. Proposition 5. Let β satisfy Assumption 1 and G 1 ∈ L ∞ 0 (V). (i) There exists a unique solution g 1 ∈ L ∞ 0 (V) of Lg 1 = G 1 . (ii) V g 1 Lg 1 dv is positive for any nonzero g 1 .
Proof. We prove (ii) using the following Cauchy-Schwarz inequality

I 1,1 (g) = W W β (v, w) g 1 (v) β (v, w) g 1 (w) dw dv ≤ W W β (v, w) (g 1 (v))
2 dw dv = I 2,0 (g 1 ) .

PARISOT MARTIN AND LACHOWICZ MIROS LAW

with I a,b defined by Eq. ( 8). We conclude that V g 1 Lg 1 dv is nonnegative. In addition, it vanishes only with the equality case of the Cauchy-Schwarz inequality, corresponding to the constant function almost everywhere. Since the average of the solution vanishes, it follows that the solution of V g 1 Lg 1 dv = 0 vanishes almost everywhere. Then, we prove (i) using the Fredholm theory. Since β ∈ L ∞ V 2 , we only have to check that the homogeneous equation Lg 1 = 0 does not admit any nontrivial solution. Multiplying the homogeneous equation by g 1 and integrating yields V g 1 Lg 1 dv = 0. Using (ii) we conclude that the only solution of the homogeneous equation is the zero solution.

Accordingly to Eq. ( 10), the macroscopic equation for the entropy dissipation interactions (i.e. σ (1 -γ) > 0) at the second order is given by the following nonlinear parabolic equation

∂ t ρ -ε∇ x • κ ρ γ ∇ x ρ = O ε 2 with κ = V v • L -1 v dv |1 -γ| |V| 1-γ . ( 11
)
By Proposition 5 the diffusion coefficient κ is positive since it can be written as

κ = V χ • Lχ dv |1 -γ| |V| 1-γ with χ = L -1 v.
In addition, the nonlinearity ρ -γ of the diffusion coefficient is classical for the generalized Carleman-type kinetic models, see [START_REF] Salvarani | The diffusive limit for carleman-type kinetic models[END_REF]. The limit case when ρ vanishes, the so-called regime of very fast diffusion is well-known and described in [START_REF] Daskalopoulos | On nonlinear parabolic equations of very fast diffusion[END_REF][START_REF] Vásquez | Smoothing and decay estimates for nonlinear diffusion equations[END_REF][START_REF] Salvarani | The diffusive limit of carleman-type models in the range of very fast diffusion equations[END_REF].

Note that the nonlinearity of the diffusion coefficient can be modified since the interaction rate β can be a function of ρ.

The diffusion coefficient κ, for the classical simple domains V with β (v, w) = 1, is given by

V κ S d V = v ∈ R d : |v| = V V 3-γ-(1-γ)d |1 -γ| S 1-γ d B d V = v ∈ R d : |v| ≤ V d 2-γ V 2-(1-γ)d |1 -γ| (d + 2) S 1-γ d , (12)
where V > 0 is a (given) maximal speed, S d is the surface area of the unit sphere

in dimension d, i.e. S d-1 = 2π d 2 Γ( d 2 )
and Γ is the Euler Gamma function. In some particular cases, for example V = B 3 V , β (v, w) = 1 and γ = 1 3 , the diffusion coefficients is not a function of the maximal speed V . Then the hydrodynamic limit is defined for unbounded domains V.

3.2.

Macroscopic evolution for self-organized interactions. Proposition 6. Let β (v, w) = 1. In both self-organized cases σ (1 -γ) < 0 (i.e. σ = -1, 0 ≤ γ < 1 as well as σ = 1 and γ > 1), there is no steady states of Eq. (3) that are stable in L ∞ (V).

Proof. In both self-organized cases it is obvious that for any point v ∈ V such that

G (v) < G γ γ-1
L γ (V) , g (t, v) tends to zero when t tends to infinity. By Proposition 3 the steady states of Eq. ( 3) are constant functions g 0 on their supports W, i.e.

g 0 (v) = g 0 γ γ-1 L γ (V)
. We consider g ζ the solution of Eq. ( 3) with the initial data G ζ defined as a perturbation of the steady state

G ζ = g 0 1 + ζ 1 W+ -1 W-with the two subsets W -⊂ W and W + ⊂ W such that |W -| = |W + | (and ζ small). For any v ∈ W -we have G ζ < G ζ γ γ-1
L γ (V) then it follows that g ζ (t, v) tends to zero for any v ∈ W -. We conclude that the steady state is not stable.

By Proposition 6 the local equilibrium have to be search among distributions or measures. We may expect the local equilibrium in the form of a Dirac function

f a (t, x, v) = ρ (t, x) δ (v -u (t, x))
, where u ∈ V is the mean velocity, see [START_REF] Degond | Continuum limit of self-driven particles with orientation interaction[END_REF][START_REF] Motsch | A new model for self-organized dynamics and its flocking behavior[END_REF][START_REF] Carrillo | A new interaction potential for swarming models[END_REF] -the so-called aligned picture. This kind of limit can be established considering a discrete set of possible velocities, i.e. for V a finite subset of Z, see [START_REF] Banasiak | On a macroscopic limit of a kinetic model of alignment[END_REF]. Note that in this framework, the problem of global existence disappear since every norm are equivalent, and the l 1 -norm is preserved. Considering a continuous set of possible velocities it is clear that Eq. ( 1) has no meaning in the framework of distributions or measure. Therefore we cannot derive macroscopic limit using classical tools. The derivation of a macroscopic limit for Eq. ( 1), even in the formal way, is an open problem.

In the following, we formulate a conjecture for the macroscopic limit which seems to correspond to the limit of Eq. ( 1) in the case of self-organized interactions when ε tends to zero. We discuss the relevance of the conjecture in Section 5.2.2.

By Proposition 6, there is no local equilibrium in L ∞ (V) but still the L 1 (V)norm is finite. We assume that the local equilibrium is mono kinetic (see [START_REF] Carrillo | A new interaction potential for swarming models[END_REF]), i.e. f (t, x, v) = ρ (t, x) δ (v -u(t, x)). The main problem is the definition of the mean velocity u (t, x).

Unfortunately, we are not able to solve the problem in the general case. However, in the case where β (v, w) is constant, the mean velocity is the velocity v such that the initial datum F (x, v) is maximal, if it is unique. Let us explain the strategy in one dimensional framework for simplicity. It is clear that the asymptotic solution of the homogeneous cases is the mono kinetic function with the mean velocity given by the velocity v such that the initial condition F (x, v) is maximal, if it is unique. More precisely, for the initial data F , we define U (x) such that

F (x, U (x)) = max v∈V F (x, v) .
Then the solution ρ (t, x) can be defined using the characteristic method before a shock, i.e. as long as the trajectories do not intersect. In particular, it is possible to create vacuum, i.e. ρ = 0, if the velocity at left is smaller than the velocity at right. At each point where N trajectories intersect, the solution at the intersection point is given by the solution of the homogeneous equation ( 3) with an initial data composed by the N Dirac functions, i.e.

G (v) = N i=1 ρ i δ (v -u i ) .
Note that it is exactly the case where the interaction operator Q does not have a sense since the function G is not in L γ (V). We may however assume that the asymptotic solution is given by the velocity characterized by the largest population, when this is unique.

More precisely the conjectured solution after a shock is given by

g (v) = j∈M N i=1 ρ i Card (M) δ (v -u j ) with M = 1 ≤ i ≤ N | ρ i = max 1≤j≤N (ρ j ) .
To illustrated it, we construct the conjectured solution in the case of initial data that are piecewise constant in Section (5.2.2). In the general case, i.e. β is not constant, we are not able to define the macroscopic limit since we are not able the define the initial mean velocity or the velocity after a shock.

4. Numerical solutions. The present section is devoted to the numerical solutions of Eq. ( 1). We are going to solve independently the interactions between individuals and the transport in space, based on the formulation along the trajectory [START_REF] Daskalopoulos | On nonlinear parabolic equations of very fast diffusion[END_REF]. Such a splitting strategy is classical for the numerical solutions of the kinetic equations. We present the numerical results in the 1-dimensional case. We propose several schemes for the interaction operator (Section 4.2). The objective is to design a scheme that is able to recover each solution of Eq. ( 3) and not to restrictive in the limit ε goes to zero. All these schemes can be easily extended to the multi-dimensional case.

In Section 4.3 we propose a solution of the transport operator in the 1-dimensional case. For the Cartesian grid in space and in velocity, the strategy does not introduce the numerical diffusion, which is required to describe the formation of swarms. The approximation of the transport operator cannot be extended to the multidimensional case or to the general grids. However, several more sophisticated approaches to design an anti-dissipative schemes can be found in the literature -see e.g. [START_REF] Després | An Antidissipative Transport Scheme on Unstructured Meshes for Multicomponent Flows[END_REF].

4.1. Numerical integration in velocity. We start by introducing a numerical integration with respect to variable v. We introduce a Cartesian grid of the velocity space that is symmetric V = [-V, V ], i.e. V Nv = {v i } -Nv≤i≤Nv such that v i = i dv with dv = V Nv . We propose to use Newton-Cotes formulas, with positive weigh ω j > 0, denoted • , to approximate the nonlocal operators. The Newton-Cotes formula of degree 1 reads

φ = dv Nv j=-Nv ω j φ j , with ω j = 1 2 , if j = -N v or j = N v 1, otherwise , with φ i = φ (v i ).
In addition, we introduce the notation of the discrete interaction rate

β ij = β (v i , v j ) > 0. The nonlocal operator V β (v, w) φ (w) dw for the velocity v = v i is naturally approximated by φ i = β i φ = dv
[G] i = ρ |I| if i ∈ I 0 otherwise. with |I| = dv i∈I w i .
Proof. The proof is similar to the continuous case Proposition 3.

Numerical approximations of the interactions.

In the present section we analyze several numerical schemes for solution of the space homogeneous case [START_REF] Ben-Naim | Opinion dynamics: rise and fall of political parties[END_REF]. We indicate the dependence with respect to the Knudsen number ε. Our main objective is to provide a numerical scheme stable for long time step in order to use it at the macroscopic regime, i.e. when ε tends to 0.

4.2.1. Euler explicit scheme. The simplest way to solve Eq. ( 3) is probably the Euler explicit time scheme, i.e.

[g e ]

n+1 i = [g e ] n i + σ dt ε [g e ] n i ([g e ] n i ) γ -([g e ] n ) γ i [g e ] n i , (13)
with the time step dt, [g e ] n i is an approximation of the solution g (t n , v i ) and the discrete initial data [g e ]

0 i = [G] i = G (v i ).
Proposition 8 (Properties of the Euler explicit scheme ( 13)). Let the initial data be nonnegative and integrable, i.e.

C max e = dv max -Nv≤i≤Nv Nv j=-Nv [ge] n j =[ge] n i ω j β ij ([g e ] n i ) γ [g e ] n j -[g e ] n j γ [g e ] n i [g e ] n i -[g e ] n j
, the approximation [g e ] satisfies the maximum principle

min -Nv≤i≤Nv [G] j ≤ [g e ] n i ≤ max -Nv≤i≤Nv [G] j .
Proof. The properties i), ii) and iii) are obvious. In order to prove iv) we assume that the approximation in iteration n is nonnegative which is true at the initial data [G] i ≥ 0. The scheme [g e ] (13) can be written in the form

[g e ] n+1 i =     1 + σ dt dv ε Nv j=-Nv [ge] n j =[ge] n i ω j β ij ([g e ] n i ) γ-1 -[g e ] n j γ-1 [g e ] n i -[g e ] n j [g e ] n i [g e ] n j     [g e ] n i -σ dt dv ε Nv j=-Nv [ge] n j =[ge] n i   ωjβij ([g e ] n i ) γ-1 -[g e ] n j γ-1 [g e ] n i -[g e ] n j [g e ] n i [g e ] n j    [g e ] n j .
For the solitarious cases, i.e. σ (1 -γ) > 0, for any (a, b) ∈ R + , the function

-σ a γ-1 -b γ-1 a-b
ab is nonnegative. Then, under the CFL condition iv), the approximation in time iteration n + 1 is a convex function of the approximation at time n. We conclude the maximum principle of the discrete approximation. 4.2.2. Euler semi-implicit scheme. The Euler explicit scheme ( 13) is stable under the CFL condition, see Proposition 8. The CFL conditions could become restrictive, in particular considering macroscopic regime, i.e. when ε tends to zero. To avoid this difficulty, it is classical to use Euler implicit scheme. The nonlinear term requires using of iterative process based on the linearization. For simplicity, we look for the numerical strategy that does not require the solution of matrix systems. However, it is not possible to find a semi-implicit scheme satisfying all the properties of the continuous solution, in particular the discrete version of the mass conservation Theorem 2.1. This conservation is the key point of stability of the numerical scheme since for the self-organized models, i.e. σ (1 -γ) < 0, it is the only norm decreasing in time. In order to fix the mass conservation at the discrete level, we propose a correction step using an operator P ρ : R 2Nv+1 + \ {0} → R 2Nv+1 . Lemma 4.1. Assume that the operator P ρ satisfies the following properties:

i) The image of the operator

P ρ is {φ j } ∈ R 2Nv+1 + | φ = ρ .
ii) For the gregarious cases, i.e. σ (1 -γ) < 0, the operator P ρ is nonnegative.

Moreover, for any -N v ≤ i ≤ N v such that φ i = 0, we have P ρ i = 0. iii) For the solitarious cases, i.e. σ (1 -γ) > 0, any interval of R + is stable by the operator P ρ , i.e. for any 0 ≤ a ≤ b and {φ i } ∈ [a, b] 2Nv+1 , we have

P ρ i {φ j } j i ∈ [a, b] 2Nv+1 .
iv) For any {φ i } i ∈ R 2Nv+1 \ {0} such that φ = ρ and the perturbation {ψ i } i ∈ R 2Nv+1 such that ψ i = O (ν) and {φ i + ψ i } i ∈ R 2Nv+1 \ {0}, we have

P ρ i {φ j + ψ j } j = φ i + O (ν) .
In the present paper, we set P ρ i {φ j } j = 0 if φ i = 0 and

P ρ i {φ j } j = ρ 1 φ =0 + (φ-ρ 2V ) -1 φ =0 (φ-ρ 2V ) + -(φ-ρ 2V ) + 1 φ =0 (φ-ρ 2V ) - max (φ-ρ 2V ) -1 φ =0 , (φ-ρ 2V ) + 1 φ =0 (14)
otherwise. The proofs of the properties of Lemma 4.1 using the formula ( 14) are obvious. Note that other definitions of the operator P ρ are possible.

We design the semi-implicit scheme

[g i ] n,r+1 i = P [G] [g i ] n,r j j
, where

[g i ] n,r i = ε [g i ] n,0 i + dt (σ) + [g i ] n,r i ([g i ] n,r i ) γ + (σ) -([g i ] n,r ) γ i [g i ] n,r i ε + dt (σ) + ([g i ] n,r ) γ i + (σ) -[g i ] n,r i ([g i ] n,r i ) γ-1 , (15)
where

[g i ] n+1,0 i = lim r→∞ [g i ] n,r i and [g i ] 0,0 i = [G] i = G (v i ). The discrete unknown [g i ] n,0 i
is an approximation of the solution of Eq. ( 3) at the point (t n , v i ). If the approximation at the previous time step vanishes, i.e. [g i ] n,r i = 0, we have [g i ] n,r i = 0 except in the case of negative gregarious interactions, i.e. σ = -1 and 0 ≤ γ < 1, where [g i ] n,r i is not defined. By the continuous analysis, see Theorem 2.4, we known that the only nonnegative solution in this case is the vanishing solution. Then we can extent the numerical scheme setting

[g i ] n,r i = 0 if [g i ]
n,r i = 0 for any σ and γ.

Proposition 9 (Properties of the Euler semi-implicit scheme [START_REF] Frénod | An explanatory model to validate the way water activity rules periodic terrace generation in proteus mirabilis swarm[END_REF]). Let the initial data be integrable, i.e. [G] < ∞ and nonnegative, i.e. [G] i ≥ 0. The scheme [g i ] (15) satisfies the following properties:

i) The sequence [g i ] is an approximation of the solution of Eq. (3) in order 1 in time.

ii) The approximation [g i ] satisfies the mass conservation, i.e.

[g i ] n,r = [G] . iii) The approximation [g i ] is nonnegative, i.e. [g i ] n,r i ≥ 0. Moreover, for any -N v ≤ i ≤ N v such that the initial data [G] i are positive, the approximation [g i ]
n,r i is positive. iv) For the solitarious case, i.e. σ (1 -γ) > 0, the approximation [g i ] satisfies the maximum principle, i.e.

min -Nv≤i≤Nv [G] j ≤ [g i ] n,r i ≤ max -Nv≤i≤Nv [G] j . Proof. It is clear that [g i ] n,r i
( dt) is an approximation of order 1 in time of the solution of Eq. (3). Then using Lemma 4.1.iv), we conclude the property i). The property ii) is a corollary of the Lemma 4.1.i). Using the recursive argument, the non-negativity (positivity) iii) is obvious since all the term are nonnegative (positive).

Assume that the previous approximation [g i ] n,r satisfies the maximum principle min -Nv≤j≤Nv

[g i ] n,0 j = m n ≤ [g i ] n,r i ≤ M n = max -Nv≤j≤Nv [g i ]
n,0 j .

Note that the initial approximation [g i ]

n,0 satisfies this assumption. We show that [g i ] n,r satisfies the property in the case of positive interaction, i.e. σ = 1. For any X ∈ R we have

[g i ] n,r i -X = ε [g i ] n,0 i -X + dt [g i ] n,r i ([g i ] n,r i ) γ -([g i ] n,r ) γ i X ε + dt ([g i ] n,r ) γ i
.

The first term of the nominator is clearly nonpositive with X = M n and nonnegative with X = m n . Since 0 ≤ γ ≤ 1 the functions x γ-1 is decreasing and x γ is increasing then

[g i ] n,r i ([g i ] n,r i ) γ -([g i ] n,r ) γ i M n ≤ [g i ] n,r i ([g i ] n,r i ) γ -M γ n ≤ 0 .
The similar result holds for the lower bound. We conclude m n ≤ [g i ] n,r i ≤ M n . We proceed similarly for the negative interaction case. More precisely, we write

[g i ] n,r i -X = ε [g i ] n,0 i -X + dt ([g i ] n,r ) γ i [g i ] n,r i -[g i ] n,r i ([g i ] n,r i ) γ-1 X ε + dt [g i ] n,r i ([g i ] n,r i ) γ-1
.

Since 1 ≤ γ, the function x γ-1 is increasing and

([g i ] n,r ) γ i [g i ] n,r i -[g i ] n,r i ([g i ] n,r i ) γ-1 M n ≤ [g i ] n,r i M γ-1 n [g i ] n,r i -([g i ] n,r i ) γ-1 M n ≤ 0
The similar result holds for the lower bound. Finally we conclude the property iv) using Lemma 4.1.iii).

Note that the semi-implicit scheme ( 15) can be used in the limit ε = 0. Thanks to the iterative process, when r tends to infinity, the solution with ε = 0 tends to the asymptotic solution of Eq. ( 3) which is given by the steady state Proposition 7.

In practice, we perform several tests (positive or negative interactions, solitarious or gregarious behavior, convex or concave interaction rate β) and the iterative process converges in any cases using the semi-implicit scheme [START_REF] Frénod | An explanatory model to validate the way water activity rules periodic terrace generation in proteus mirabilis swarm[END_REF]. 4.2.3. Non-linear scheme. The explicit and implicit Euler schemes are based on a linearization of the interaction operator. Then the previous schemes ( 13), ( 15) conserve the vanishing points, i.e. for any points

v i such that G (v i ) = 0, we have [g e ] n i = [g i ] n,0 i = 0.
Even if this solution is relevant from mathematical point of view, it is not satisfactory at the biological level, see Remark 1. In addition, none of the previous schemes are able to vanish in a finite time in the case of negative gregarious interaction, see Theorem 2.4. The objective of this section is to design a numerical scheme able to recover the solution in the case of no Lipschitz interaction coefficient 0 ≤ γ < 1, i.e. the solution ables to leave zero in case of positive solitarious interaction and the solution ables to vanish in a finite time in case of negative gregarious interaction.

As we mentioned, this solution required the non-linear structure of the equation. We propose to base the numerical scheme on the following nonlinear ODE problem

εσφ = (pφ γ -qφ) with φ (0) = φ 0 . (16)
Assuming γ = 1, we use the substitution of unknown ψ = φ 1-γ to get the linear ODE problem

εσψ = (1 -γ) (p -qψ) with ψ (0) = ψ 0 . ( 17
)
Then it is possible to use the explicit or implicit Euler schemes to get an approximation of solution of Eq. ( 17). However, this schemes lead to restrictive conditions on the time step in the case of gregarious interaction σ (1 -γ) ≤ 0. In fact, the time step ensuring the positivity of the numerical solution vanishes when the solution vanishes in some points. To overcome this drawback, we use the exact solution of the ODE [START_REF] Greenwood | Differences numbers of sensilla on the antennae of solitarious and gregarious locusta migratoria l. (orthoptera: Acrididae)[END_REF] assuming that the parameters p and q are two given constants. From now on, we denote by ψ t, ψ 0 , p, q the solution of the linear ODE (17) and we have the following explicit formula

ψ t, ψ 0 , p, q = p q + ψ 0 - p q exp -σ (1 -γ) q t ε . (18)
For the nonnegative initial data ψ 0 and positive parameters p and q, the solution Eq. ( 18) is nonnegative in the solitarious case, i.e. σ (1 -γ) > 0. For the negative gregarious case, i.e. σ = -1 and 0 ≤ γ < 1, the only nonnegative solution of the nonlinear ODE ( 16) is given by the extension of ( 18) by zero, see Section (2.2.1).

For the positive gregarious case, i.e. σ = 1 and 1 < γ, the solution of the nonlinear ODE vanishes only if φ 0 = 0. In addition, the solution of the nonlinear ODE ( 16) blows up, see (2.2.2).

Finally, we design the nonlinear scheme [g n ] n,r+1 i

= P [G] [g n ] n,r j j
, where

[g n ] n,r i = max 0, ψ dt, [g n ] n,0 i 1-γ , [g n ] n,r i , ([g n ] n,r ) γ i 1 1-γ , (19) where
[g n ] n+1,0 i = lim r→∞ [g n ] n,r i and [g n ] 0,0 i = [G] i = G (v i ).
For the positive gregarious case, i.e. σ = 1 and 1 < γ, the numerical scheme ( 19) is defined only if, for any

-N v ≤ i ≤ N v such that [g n ] n,0 i > 0, the time step satisfies dt ≤ - ε ln 1 - ([gn] n,r ) γ i [gn] n,r i ([gn] n,0 i ) γ-1 (γ -1) ([g n ] n,r ) γ i . (20)
This time step corresponds to the time of blow up of nonlinear ODE [START_REF] Geigant | Bifurcation analysis of an orientational aggregation model[END_REF].

Proposition 10 (Properties of the nonlinear scheme [START_REF] Ha | From particle to kinetic and hydrodynamic description of flocking[END_REF]). Let the initial data be integrable and nonnegative:

[G] < ∞ and [G] i ≥ 0. The scheme [g n ] (19)
, under the CFL condition [START_REF] Kang | An integro-differential equation model for alignment and orientational aggregation[END_REF], satisfies the following properties: i) The sequence [g n ] is an approximation of the solution of Eq. (3) in order 1 in time.

ii) The approximation [g n ] satisfies the mass conservation

[g n ] n,r = [G] . iii) The approximation [g n ] is nonnegative [g n ] n,r i ≥ 0. iv) For the solitarious case, i.e. σ (1 -γ) > 0, the approximation [g n ] n,0 i satisfies the maximum principle, i.e. min -Nv≤i≤Nv [G] j ≤ [g n ] n,r i ≤ max -Nv≤i≤Nv [G] j . Proof. It is clear that [g n ]
n,r i is an approximation in order 1 in time of the solution of Eq. (3). Then using Lemma 4.1.iv) we obtain i). The property ii) is a corollary of the Lemma 4.1.i). The property iii) is obvious by the using the formula ( 19) and Lemma 4.1.ii).

For the solitarious cases, i.e. σ (1 -γ) > 0, for a fixed value of

{[g n ] n,r i } i , the new approximation [g n ] n,r i
[g n ] n,0 j = m n ≤ [g n ] n,r i ≤ M n = max -Nv≤j≤Nv [g n ]
n,0 j .

Note that the initial approximation [g n ]

n,0 satisfies this assumption. Then for 0 ≤ γ < 1, we have the following inequalities

([g n ] n,r ) γ i m 1-γ n ≤ [g n ] n,r i = ([g n ] n,r ) γ ([g n ] n,r ) 1-γ i ≤ ([g n ] n,r ) γ i M 1-γ n .
Similarly, for 1 < γ, we have the following inequalities

[g n ] n,r i m γ-1 n ≤ ([g n ] n,r ) γ i = [g n ] n,r ([g n ] n,r ) γ-1 i ≤ [g n ] n,r i M γ-1 n .
Using the recursive argument and Lemma 4.1.iii), we obtain iv).

Spacial evolution.

In the following, we propose a numerical strategy for the space-depending PDE (1) in the whole space of 1-dimension. The boundary conditions are not considered in this work and the numerical simulation are performed using the periodic domain. We propose the Eulerian discretization on the Cartesian grid x k = k dx, with dx the space step. The choice of the Eulerian discretization time step accordingly to the velocity. In addition, to treat the interaction between the direct neighbors, we have to considerate the interaction during the transport step.

Let us decompose the anti-diffusive scheme in five steps. The steps 1, 2, and 4 are similar to those in the strategy proposed in the previous Section 4.3.1. The transport step 3 is similar as in the case of a transport with the CFL parameter set to 1. In addition, a special treatment of the part of the population which does not move, i.e. i = 0, is required. Finally, the step 5 is an actualization of the solution when the different time step are synchronized.

Note that this synchronization is realizable by the fact that the Cartesian grid is used. The generalization to a general grid required more sophisticated treatment, see [START_REF] Després | An Antidissipative Transport Scheme on Unstructured Meshes for Multicomponent Flows[END_REF].

1. We construct an approximation of the distribution at the interfaces x k + dx 2 using the distribution at the left for the positive velocity and at the right for the negative velocity, i.e. we set

[f a z ] n,l k+ 1 /2,i =      [f a z ] n,l k+1,i if i < 0 [f a z ] n,l k,i +[f a z ] n,l k+1,i 2 if i = 0 [f a z ] n,l k,i if 0 < i .
2. We estimate the distribution at the interface after interaction, i.e.

[f a z ]

n,l

k+ 1 /2,i = Q z,i dx 2 dv , [f a z ]
n,l k+ 1 /2,j j .

3. We transport the distribution at the interface considering that only the fast enough velocity, i.e. |v i | > l dv, passes through the interface, i.e. we set

[f a z ] n,l+1 k,i =                    [f a z ] n,l k+ 1 /2,i if i ≤ l -N v [f a z ] n,l k-1 /2,i if l -N v < i < 0 1 2 [f a z ] n,l k,0 [f a z ] n,l k-1 / 2 ,0 [f a z ] n,l k-1 / 2 ,0 + [f a z ] n,l k,0 [f a z ] n,l k+ 1 / 2 ,0 [f a z ] n,l k+ 1 / 2 ,0 if i = 0 [f a z ] n,l k+ 1 /2,i if 0 < i < N v -l [f a z ] n,l k-1 /2,i if N v -l ≤ i . (22)
Because of the reconstruction at the interface, the distribution is split and leads to numerical diffusion (only for i = 0) without this special treatment.

To overcome this difficulty we reconstruct the approximation at the grid point x k considering the ratio of the distribution at the interface coming from the point x k . 4. We estimate the distribution at the grid point after interaction, i.e.

[f a z ] n,l+1 k,i = Q z,i dx 2 dv , [f a z ]
n,l+1 k,j j .

5. We iterate the steps 1 to 4 until all the velocities pass trough the interface, i.e. for 1 ≤ l ≤ N v and we set the approximation at the new time step

[f a z ] n+1,0 k,i = [f a z ] n,Nv
k,i . The numerical scheme can be interpreted as a time splitting with N v loops with the sub time step dx V , with V the larger speed. In each loop, we have split the transport (step 1 and 3) and the interaction (step 2 and 4) into 2 steps. It follows that the new approximation [f a z ] n+1,0 k,i is an approximation of the solution of Eq.

(1) at the point t n + dx dv , x k , v i . In addition, the numerical unknown [f a z ] n,l k,i can be considered as an approximation in the point t n + l dx V , x k , v i .

5. Numerical simulations and validation.

5.1. Numerical simulations of the space homogeneous case. In the present section we compare the numerical schemes of the interaction operator described in Section 4.2. We consider in the space homogeneous case and with V = [-V, V ], V = 0.9, dv = 4.5 • 10 -3 . The initial data are defined by

G (v) = max 0, 0.75 sin 4 0.1 + v 2 π e -(v-0.75) 2 2 , ( 23
)
the interaction rate is β (v, w) = e -(v-w) 2 and the Knudsen number is ε = 1. In Figure 2 we represent the numerical solutions in the case of negative solitarious interactions, i.e. σ = -1, γ = 2. The results of the proposed schemes ( 13), ( 15) and [START_REF] Ha | From particle to kinetic and hydrodynamic description of flocking[END_REF] are almost the same. The vanishing areas are preserved along the simulation and the steady state is the constant solution on the support of the initial data. Note that in this case, the solution is unique.

In Figure 3 we present the numerical solutions in the case of positive solitarious interactions, i.e. σ = 1, γ = 0.5. The results of the linear schemes [START_REF] Eftimie | Hyperbolic and kinetic models for self-organized biological aggregations and movement: a brief review[END_REF] and [START_REF] Frénod | An explanatory model to validate the way water activity rules periodic terrace generation in proteus mirabilis swarm[END_REF] are almost the same. The vanishing areas are preserved along the simulation and the steady state is the constant solution on the support of the initial data. The numerical solutions using the nonlinear scheme [START_REF] Ha | From particle to kinetic and hydrodynamic description of flocking[END_REF] lead to the constant steady state on the whole domain V.

Both solutions, constant on the initial data support as well as constant on the whole domain, are steady states of Eq. ( 1), but only the latter is Lyapunov stablesee Proposition 4. It follows that only solution satisfying Remark 1 is the solution approached by the scheme [START_REF] Ha | From particle to kinetic and hydrodynamic description of flocking[END_REF]. We will see that this property leads to important consequences at the macroscopic level see Section 5.2.1.

In Figure 4, we represent the numerical solutions in the case of negative gregarious interactions, i.e. σ = -1, γ = 0.5.

The CFL condition of the explicit scheme (13) (Proposition 8) is too restrictive in the velocity point grid such that the solution is small. In fact the time step satisfying the CFL condition tends to vanish with the solution.

The results with the two other schemes ( 15) and ( 19) are almost the same. Note that the asymptotic solution is a Dirac function. However, the support of the Dirac function (v = 0.2) is not the velocity for which the initial datum is maximal (v = 7.2). The correspondence between the support of the asymptotic solution and the localization of the maximum of the initial datum is only valid when the interaction rate is constant, i.e. β = 1 which is not the case in the present simulation.

In Figure 5 we represent the numerical solutions in the case of positive gregarious interactions, i.e. σ = 1, γ = 2. The results of the two linear schemes [START_REF] Eftimie | Hyperbolic and kinetic models for self-organized biological aggregations and movement: a brief review[END_REF] and [START_REF] Frénod | An explanatory model to validate the way water activity rules periodic terrace generation in proteus mirabilis swarm[END_REF] are almost the same. We can clearly identify a time (t ≈ 5.2) after which the solution is concentrated in a single velocity grid point. In addition, the iterative process of the nonlinear scheme [START_REF] Ha | From particle to kinetic and hydrodynamic description of flocking[END_REF] does note converge at this time. From now on we refer to this time as the time of numerical blow up. Using a smaller velocity step, the time of numerical blow up seams to be unchanged. This observations are even more clear since the numerical diffusion is limited. In the case of positive solitarious interactions we discredit this solution considering biological arguments, see Remark 1. However, in the case of negative solitarious interactions, the solution is unique, see Theorem (2.2), and the scheme of the interaction operator leads to the unstable equilibrium. The hydrodynamic limit is thus not defined in the negative solitarious interactions case, except when the initial datum is separated from zero.

In the fourth and fifth lines of Figure 6 we present the numerical solutions obtained with the nonlinear scheme for the interactions [START_REF] Ha | From particle to kinetic and hydrodynamic description of flocking[END_REF]. More precisely, the fourth line of Figure 6 is obtained using the up-wind scheme in space [START_REF] Kaniel | The Boltzmann equation[END_REF] [f u n ] and the fifth line of Figure 6 is obtained using the anti-diffusive in space [START_REF] Mach | Modeling vortex swarming in daphnia[END_REF] [f a n ]. The results are qualitatively similar to the macroscopic solution, i.e. first line of Figure 6. However, the solution is significantly more diffusive than the macroscopic solution in particular using the up-wind scheme [f u n ]. In addition, the diffusion in the areas with large density is clearly faster, which characterize that the nonlinearity of the diffusion is not well-represented. In the case of the negative solitarious interactions, i.e. σ = -1 and γ = 2, the nonlinear scheme for interaction [START_REF] Ha | From particle to kinetic and hydrodynamic description of flocking[END_REF] preserves the points in the velocity grid such that the distribution vanishes, as well as the semi-implicit scheme [START_REF] Frénod | An explanatory model to validate the way water activity rules periodic terrace generation in proteus mirabilis swarm[END_REF]. The solution obtained in this case is qualitatively the same as presented in the second line of Figure 6 for the up-wind scheme in space [START_REF] Kaniel | The Boltzmann equation[END_REF] or the third line of Figure 6 for the anti-diffusive scheme in space [START_REF] Mach | Modeling vortex swarming in daphnia[END_REF]. 5.2.2. The aligned picture. In the following section we consider the negative gregarious case, i.e. σ = -1 and γ = 0.5, with β = 1 in the macroscopic limit ε = 0, for the initial datum [START_REF] Oien | Daphnicle dynamics based on kinetic theory: an analogue-modelling of swarming and behaviour of daphnia[END_REF]. The following results are similar to the positive gregarious interactions case, i.e. σ = 1 and γ = 2, in the macroscopic limit ε = 0, except that the non-linear fixed point for interactions converges slowly. First we estimate the conjectured solution following the rules design in Section 3.2 for an initial condition that is piecewise constant. At the initial time, two distributions meet at point x 4 . Since u 2 > u 3 , the trajectory starting from x 4 is a shock. Since ρ 2 < ρ 3 , the shock is moving with the velocity u 3 and the density from the left part gathers on this trajectory: we have a space Dirac function S 4 δ (x -(x 4 + u 3 t)) with the mass S 4 = (u 2 -u 3 ) tρ 2 for a time t small enough to do not intersect another trajectory. Similarly, at the initial time, two distributions meet at point x 7 . However, since u 4 < u 5 , there is no interactions and the two distributions have their velocities, creating an empty area between them. We proceed similarly until all the populations tend to the same direction. In Figure 7 we present the conjectured solution for the initial datum [START_REF] Oien | Daphnicle dynamics based on kinetic theory: an analogue-modelling of swarming and behaviour of daphnia[END_REF] with x in abscissa and t in ordinate. i ] with a space step dx = 10 -3 (second line), and antidiffusive space scheme [START_REF] Mach | Modeling vortex swarming in daphnia[END_REF] [f u i ] with a space step dx = 1 999 (third line).

In Figure 8 we present the numerical solutions obtained with the different schemes proposed above and the initial datum [START_REF] Oien | Daphnicle dynamics based on kinetic theory: an analogue-modelling of swarming and behaviour of daphnia[END_REF]. The space x is plotted on the abscissa and the time t is plotted on the ordinate. We use the semi-implicit scheme [START_REF] Frénod | An explanatory model to validate the way water activity rules periodic terrace generation in proteus mirabilis swarm[END_REF] to treat the interaction operator. The nonlinear interaction scheme (19) cannot be used in the asymptotic limit ε = 0 in the case of gregarious interaction σ (1 -γ) < 0. However, similar results were obtained for the nonlinear interaction scheme [START_REF] Ha | From particle to kinetic and hydrodynamic description of flocking[END_REF] with the small Knudsen numbers, i.e. ε < 10 -3 . The velocity step set to dv = 0.1. The results with the other velocity steps are exactly the same as long as the initial velocities are represented in the velocity grid.

In the first line of Figure 8, we use the up-wind scheme in space [START_REF] Kaniel | The Boltzmann equation[END_REF] [f u i ] with a space step dx = 10 -3 . The numerical diffusion makes the numerical approximation useless because the part of the distribution able to pass through an interface is proportional to the time step. Since the time step is small for the stability reasons, the distribution, initially in a mesh point, is generally larger than the part incoming from the neighbors. Thus after application of the interaction operator, the population tends to the direction of the initial data in each point, except in the part where positive gregarious interactions (σ = 1 and 1 < γ) the existence and uniqueness of the solutions are proved only locally. In addition, the macroscopic limit of the gregarious cases, even if the solutions globally exist, cannot be realized using the classical strategies. In fact, the local equilibrium of the solution seems to be a Dirac function, which is not in the domain of definition of the interaction operator. Nevertheless, we propose a conjecture of the macroscopic limit in a simplified case. This conjectured solution cannot be recover using classical transport scheme and we propose a numerical anti-diffusive scheme in time to recover it.

The present work leads to several important open questions. The most important of them is the global existence or the blow up of the solutions in the case of positive gregarious interactions (σ = 1 and 1 < γ). In addition, the derivation of the macroscopic limit in the case of the gregarious interactions is an open question, even in a simplified case β = 1.

At the numerical level the present work highlights the main difficulties in the numerical approximation of the non-linear interaction operator. The main property after the mass conservation of a numerical scheme is the entropy monotony of the solution. It is clear that the proposed numerical strategies do not satisfy the entropy monotony in both gregarious and solitarious cases. The design of a numerical scheme to deal with the two behaviors is also an open question. The space discretization realized in the present work cannot be extended to the multi-dimensional case. The adaptation of anti-diffusive schemes in multi-dimensional case describe in the literature is in progress.

Proposition 3 .

 3 Let γ = 1 and β satisfy Assumption 1.

Nv j=-Nv ω j β ij φ j . Proposition 7 (

 7 Discrete local equilibrium). The local equilibria of the discretization of Q using the Newton-Cotes formulas are the sets of values ([G] i ) -Nv≤i≤Nv constant on their supports, i.e. there exists a set I ⊂ {-N v , . . . , N v } and ρ ≥ 0 such that

1 ,

 1 [G] i ≥ 0 and [G] < ∞. The scheme [g e ] (13) satisfies the following properties: i) The sequence [g e ] is an approximation of the solution of Eq. (3) in order 1 in time. ii) The approximation [g e ] satisfies the mass conservation, i.e. [g e ] n,r = [G] . iii) Under the CFL condition C pos e the approximation [g e ] is non-negative, i.e. [g e ] n i ≥ 0. iv) For the solitarious cases, i.e. σ (1 -γ) > 0, and under the CFL condition C max e dt ≤ ε with

  is a monotonous function of the time step ∆t. More precisely, the new approximation [g n ] when ∆t tends to infinity. In order to prove the maximum principle, we have to show that the two bounds are larger than the minimum m n = min -Nv≤j≤Nv [g n ] n,0 j and smaller than the maximum M n = max -Nv≤j≤Nv [g n ] n,0 j . Assume that the previous approximation [g n ] n,r satisfies the maximum principle min -Nv≤j≤Nv

Figure 2 .

 2 Figure 2. Numerical result of the negative solitarious interaction case, i.e. σ = -1, γ = 2. Explicit scheme (13) [g e ] (first line), semiimplicit scheme (15) [g i ] (second line), and non-linear scheme (19) [g n ] (third line).

Figure 3 .

 3 Figure 3. Numerical result of the positive solitarious interaction case, i.e. σ = 1, γ = 0.5. Explicit scheme (13) [g e ] (first line), semi-implicit scheme (15) [g i ] (second line), and nonlinear scheme (19) [g n ] (third line).

Figure 4 .Figure 5 .

 45 Figure 4. Numerical result of the negative gregarious interaction case, i.e. σ = -1, γ = 0.5. Semi-implicit scheme (15) [g i ] (first line), and non-linear scheme (19) [g n ] (second line).

Figure 6 .

 6 Figure 6. Numerical result of the diffusive picture, i.e. σ = -1, γ = 0.5. Nonlinear parabolic equation (11) (first line), antidiffusive space scheme (22) with nonlinear scheme for the interactions (19) [f a n ] (second line), upwind space scheme (22) nonlinear scheme for the interactions (19) [f u n ] (third line), anti-diffusive space scheme (22) with semi-implicit scheme for the interactions (15)[f a i ] (fourth line), upwind space scheme (22) with semi-implicit scheme for the interactions (15)[f u i ] (fifth line).

Figure 7 .

 7 Figure 7. Conjecture of the hydrodynamic gregarious solution with initial data that is piecewise constant given by[START_REF] Oien | Daphnicle dynamics based on kinetic theory: an analogue-modelling of swarming and behaviour of daphnia[END_REF].

Figure 8 .

 8 Figure 8. Numerical result of the aligned picture, i.e. σ = -1, γ = 0.5 and using the semi-implicit interaction scheme[START_REF] Ha | From particle to kinetic and hydrodynamic description of flocking[END_REF]. Upwind space scheme[START_REF] Kaniel | The Boltzmann equation[END_REF] [f u i ] (first line), anti-diffusive space scheme[START_REF] Mach | Modeling vortex swarming in daphnia[END_REF] [f u i ] with a space step dx = 10 -3 (second line), and antidiffusive space scheme[START_REF] Mach | Modeling vortex swarming in daphnia[END_REF] [f u i ] with a space step dx = 1 999 (third line).

Acknowledgements. M.L. acknowledges a support from the University of KwaZulu-Natal Research Found (South Africa). The work of M.P. was carried out during the tenure of an ERCIM "Alain Bensoussan" Fellowship Programme. The research leading to these results has received funding from the European Union Seventh Framework Programme (FP7/2007-2013) under grant agreement n o 246016.

M.L. acknowledges a support from the University of KwaZulu-Natal Research Found (South Africa). The work of M.P. was carried out during the tenure of an ERCIM "Alain Bensoussan" Fellowship Programme. The research leading to these results has received funding from the European Union Seventh Framework Programme (FP7/2007-2013) under grant agreement n o 246016. xxxx 20xx; revised xxxx 20xx.

is motivated by the diffusion picture in the case of solitarious behavior, see Section 3.1. The diffusion equation [START_REF] Després | An Antidissipative Transport Scheme on Unstructured Meshes for Multicomponent Flows[END_REF] are classically solved using Eulerian discretization. In addition, we set ∆t the time step for the transport scheme. The question of the numerical solution in the multi-dimension case and with the general grid will be addressed in further works. [START_REF] Boudin | A kinetic approach to the study of opinion formation[END_REF].3.1. Up-wind scheme. Once the velocity domain is discretized on a grids, it is classical to use an explicit up-wind scheme to treat the transport with constant velocity. We define by [f u z ] n k,i the approximation of solution f of Eq. ( 1) at point (t n , x k , v i ) using the up-wind scheme and one of the interaction schemes defined in Section 4.2, i.e. z ∈ {e, i, n}.

We decompose the up-wind scheme in the four steps:

1. We construct an approximation of the distribution at the interfaces x k + dx 2 using the distribution at the left for the positive velocity and the velocity at the right for the negative velocity. More precisely we set

2. We estimate the distribution at the interface after interaction, i.e.

.

with the interaction step Q z,i T, {φ j } j that is the solution of the scheme z ∈ {e, i, n} proposed in Section 4.2 at the final time T using the initial data [G] i = φ i . 3. We transport the distribution at the interface using a classical up-wind scheme for the constant velocity. More precisely, we set

It is well known that the up-wind scheme is stable under the CFL condition |v i | ∆t dx ≤ 1. 4. We estimate the distribution at the grid point after interaction, i.e.

Unfortunately we will see from the simulation in Section 5.2.2 that the numerical diffusion is too large to represent the swarms. In fact the numerical diffusion is proportional to the CFL parameter |v i | ∆t dx . In the following Section 4.3.2, we propose a numerical scheme based on a time splitting accordingly to the velocity, such that the CFL parameter links to the considered velocity is set to 1 and the numerical diffusion vanishes. 4.3.2. Anti-diffusive scheme. We define by [f a z ] n k,i the approximation of f solution of Eq. (1) at point (t n , x k , v i ) using the anti-diffusive scheme and one of the interaction schemes defined in Section 4.2, i.e. z ∈ {e, i, n}. The strategy is based on the fact that when |v i | ∆t dx = 1, the trajectories intersect the grid points and the numerical diffusion vanishes. To set the CFL parameter to 1 for any velocity, we use different (non-mathematical) arguments showing a probable blow-up of the solution in the positive gregarious case. 5.2. Numerical simulations of the hydrodynamic regimes. The following section is devoted to the numerical validation of the asymptotic regimes identified in Section 3. We compare the different schemes proposed in Section 4.2 and Section 4.3 and discuss the relevance of the solution. Due to the CFL condition function of the Knudsen number ε, we do not consider the explicit scheme for interaction [START_REF] Eftimie | Hyperbolic and kinetic models for self-organized biological aggregations and movement: a brief review[END_REF]. We consider the mono kinetic initial data

in a [0, 1]-periodic domain.

5.2.1.

The diffusive picture. In the following section, we consider the positive solitarious case, i.e. σ = 1 and γ = 0.5 with β = 1 of the initial data [START_REF] Oien | Daphnicle dynamics based on kinetic theory: an analogue-modelling of swarming and behaviour of daphnia[END_REF]. In the present section, we assume that the agent can move until a maximum velocity set to V = 1.1, i.e. V = B 1 1.1 . To get a relevant numerical solution of the interaction operator, the discrete space step dx should be smaller than the mean fee path . Thus the Knudsen number has to be larger than the dimensionless space step, i.e. ε = L ≥ dx L . Since the solution is 1-periodic, the characteristic length scale is L = 1. All the numerical results given in this section are obtained using the space step dx = 10 -3 and the (fixed) Knudsen number ε = dx L = 10 -3 . In Figure 6, we present the numerical solutions obtained with the different schemes proposed above and the initial datum [START_REF] Oien | Daphnicle dynamics based on kinetic theory: an analogue-modelling of swarming and behaviour of daphnia[END_REF]. The space x is plotted on the abscissa and the time t is plotted on the ordinate. Accordingly to Section 4 the macroscopic regime is given by the nonlinear parabolic equation [START_REF] Després | An Antidissipative Transport Scheme on Unstructured Meshes for Multicomponent Flows[END_REF]. In the first line of Figure 6 we present the numerical solution obtained with the nonlinear parabolic equation [START_REF] Després | An Antidissipative Transport Scheme on Unstructured Meshes for Multicomponent Flows[END_REF] neglecting the term in ε 2 and using the classical centered implicit scheme. Using the estimation realized in [START_REF] Edelstein-Keshet | Do travelling band solutions describe cohesive swarm? An investigation for migratory locusts[END_REF], the diffusion coefficient of the asymptotic nonlinear parabolic equation is set to κ = V 3 / 2 3π . In the second line of Figure 6, we present the numerical solution obtained with the up-wind scheme in space ( 21) and the semi-implicit scheme for the interactions [START_REF] Frénod | An explanatory model to validate the way water activity rules periodic terrace generation in proteus mirabilis swarm[END_REF]

The approximated solution is very far from the macroscopic solution, i.e. first line of Figure 6. As we have seen in Section 5.1 the local equilibrium of the scheme ( 15) is not the constant solution in the velocity domain. It does not converge to the Lyapunov stable equilibrium so-called diffusive picture used to estimate the nonlinear parabolic equation [START_REF] Després | An Antidissipative Transport Scheme on Unstructured Meshes for Multicomponent Flows[END_REF]. More precisely, the set of velocity of the agent is stable using the semi-implicit scheme [START_REF] Frénod | An explanatory model to validate the way water activity rules periodic terrace generation in proteus mirabilis swarm[END_REF], thus the solution vanishes in each velocity grid points on which the initial condition vanishes in the whole space domain and at any time, i.e.

This drawback clearly links to the use of the semi-implicit scheme [START_REF] Frénod | An explanatory model to validate the way water activity rules periodic terrace generation in proteus mirabilis swarm[END_REF] and we obtain similar result using the semi-implicit scheme [START_REF] Frénod | An explanatory model to validate the way water activity rules periodic terrace generation in proteus mirabilis swarm[END_REF] with the anti-diffusive scheme in space [START_REF] Mach | Modeling vortex swarming in daphnia[END_REF] in the third line of Figure 6. The velocity of the agent is the initial data take value 0. In addition the numerical diffusion is not the same for all the velocities : the higher is the velocity, the less is the numerical diffusion.

In the second line of Figure 8 we present the numerical solutions obtained with the anti-diffusive scheme in space [START_REF] Mach | Modeling vortex swarming in daphnia[END_REF] [f u a ] with a space step dx = 10 -3 . The conjectured solution is almost recover. The only difference is in the representation of the swarms in space. Instead of the Dirac function the population is distributed in a volume of size dx. This drawback can have very important impact since the population size that is not in the swarm can be comparable or even larger than the population size in the swarm distributed in a volume of size dx. For example, in the considered test case, the interactions between the trajectory, denoted S 3 , delivered from the point x 3 = 0.4 with the velocity u 2 = 0.2 intersect the trajectory, denoted S 6 , delivered from the point x 6 = 0.7 with the velocity u 4 = -0.3 at time t 36 = x6-x3 u2-u4 . At the intersection time the swarm S 3 gathers the whole population between x 1 and x 2 , i.e. S 3 (t 36 ) = (x 2 -x 1 ) ρ 1 = 0.25 whereas the swarm S 6 does the whole population between x 3 and x 5 , i.e. S 6 (t 36 ) = (x 4 -x 3 ) ρ 2 + (x 5 -x 4 ) ρ 3 = 0.23. Since S 3 (t 36 ) > S 6 (t 36 ), the swarm after interaction goes with the velocity u 2 . However, the populations of the discrete scheme are S 3 (t 36 ) = (x 2 -x 1 ) ρ 1 + dxρ 2 and S 6 (t 36 ) = (x 4 -x 3 -dx) ρ 2 + (x 5 -x 4 ) ρ 3 + dxρ 4 . It follows that for the discrete schemes ( 15) and ( 19), the swarm after interaction has the velocity u 2 only if

In the third line of Figure 8, we present the numerical solutions obtained with the anti-diffusive scheme in space [START_REF] Mach | Modeling vortex swarming in daphnia[END_REF] [f u a ] with a space step dx = 1 999 . Note that the swarm after the interaction at time t 36 goes with the wrong velocity. Then the approximation after this time is completely different. It is very hard to predict the space step required in general to obtain the good approximation because it required the estimation of all the population of the swarm at any time. The discretization which may overcome this drawback will be the subject of the further works. A possible improvement can be a discretization of the swarm at the interface, which is a space Dirac function. However, it is not possible to distinguish the swarm a priori out of the strictly asymptotic limit ε = 0. Then it is not possible to estimate the part of the population that should be discretized at the interface, leading to the swarm, and the part discretized at the control volume.

Conclusion.

The present paper is devoted to the modeling of the formation and destruction of swarms using nonlinear Boltzmann-like equation. A new model is proposed using mainly two parameters to characterized the interaction between agents, namely the sign of interactions σ and the interaction coefficient γ. We highlight that the solutions can be self-organized of satisfy an entropy dissipation law as a function of these two parameters. These mathematical properties are directly linked to the biological behavior, namely gregarious or solitarious.

First, a preliminary mathematical analysis is realized in the space homogeneous case. For the solitarious interactions, the existence of the global solutions is shown. The solution in the case of interaction coefficient 0 ≤ γ < 1 is not unique. However, we identify the solution which seems relevant for biological point of view and we provide a numerical strategy to recover it. In addition, we identify the macroscopic limit of the model, namely in the case of solitarious interactions the diffusive picture.

In the case of negative gregarious interactions (σ = -1 and 0 ≤ γ < 1) we show the global existence and uniqueness of the non-negative solutions. In the case of