Rémi Sentis

Didier Paillard

Céline Baranger

Patricia Seytor

Modelling and numerical simulation of plasma flows with two-fluid mixing

Keywords: Mixing model, Multi-material two-temperature Euler system, Concentration diffusion equation, Mixing kinetic energy equation

For the modelling of plasma flows at very high temperature such the ones produced by laser beams, one must account for a bi-temperature compressible Euler system coupled to electron thermal conduction and radiative conduction. Moreover, mixing of two different fluids can occur, the two fluids occupying the same volume. For modelling such a phenomenon, instead of dealing with the conservation of mass, momentum and energy for each fluid, we propose here a simplified model which will be easier to implement in a multi-physics Lagrangian 2D code. The principle is to use a closure for expressing the relative velocity between the two fluids with the help of the gradient of the concentration. So, besides the classical system, the final model consists in a non-linear diffusion equation for the concentration and an equation for the mixing kinetic energy (analogous to the one used in turbulence models). We present also first numerical 2D simulations using this model.

Introduction. In the multicomponent plasma flows at very high temperature (for instance in the plasma produced by laser beams), mixing of two different fluids can occur. In the physical interesting cases, the mixing phenomena occur in narrow regions of the simulation domain, where the relative velocity between the two fluids is large if compared to the thermal one. A classical model for such phenomena consists in a system of six equations which corresponds to the conservation of mass, momentum and energy for each fluid. The two fluids are assumed to occupy the same volume (the global pressure is the sum of the pressures of the two fluids). The aim of this work is to perform numerical simulations dealing with this mixing problem ; this has been addressed by numerous earlier works, see [START_REF] Scannapieco | A multifluid interpenetration mix model[END_REF], [START_REF] Stewart | Two-phase Flows : Models and methods[END_REF], [START_REF] Ishii | Thermo-Fluid dynamic theory of two phase flows[END_REF], [START_REF] Resnyansky | A two-phase Shock-wave model[END_REF].

Here, our aim is to perform global two-dimension simulations taking into account the mixing phenomena and the coupling with an electron temperature equation (with thermal conduction). The related theoretical model has been written for instance in [START_REF] Decoster | Modeling of collisions[END_REF]: it consists of the previous system of six equations with an evolution equation for the electron temperature. Recall that in the single species case, the simplest model is the so-called two-temperature Euler system which reads as follows (see for instance [START_REF] Bowers | Numerical Modeling in Applied Physics[END_REF], [START_REF] Desprès | Lagrangian system of conservation laws[END_REF]) 1 corresponding author, remi.sentis@cea.fr 2 permanent address : CEA, DAM, CESTA, 33114 Le Barp, France where we denote ρ the density, U the velocity, E e the electron energy and E 0 the total ion energy (kinetic plus internal energy). Moreover, Υ 0 denotes a relaxation term between the ion and electron temperatures. Notice that other physical phenomena are to be taken into account such as the electron thermal conduction, a radiative transfer model and a laser energy deposition model. All these phenomena correspond to added terms in the evolution equation for E e which are coupled to other equations. This leads to heavy complex simulations, in particular in a multi-dimension framework, which are generally performed using Lagrangian codes, see for instance FCI2 [START_REF] Buresi | Laser Program Development at CELV[END_REF], [START_REF] Dautray | La Fusion par Confinement Inertiel[END_REF]. For all these reasons, reduced models in which the mixing process is roughly simulated are very usefull. This paper is aiming at proposing such a reduced model which may be easily implemented in an existing Lagrangian code.

In the first section we recall the classical model for mixing phenomena which consists in a system of three equations for each fluid (conservation of mass, momentum and energy), besides an equation for the electron energy, as it is written in the reference paper [START_REF] Decoster | Modeling of collisions[END_REF]. In the second section, we emphasize the closures for getting a five-equation model which consists in the conservation equations of mass, momentum and energy for the averaged fluid coupled with an equation of concentration and an equation for the mixing kinetic energy. These closures are made under the assumption that the friction coefficient becomes large. The relative velocity is given by a closure using the gradient of the concentration ; it allows to state a diffusion approximation model for the concentration (it looks like a diffusion model for combustion or turbulence problem, see [START_REF] Youngs | Numerical simulations of turbulent mixing[END_REF]). In the second part of this section we introduce the supplementary terms due to the coupling with the electron energy balance equation. In a certain way, our model is in the same family of one of the mixing models introduced in [START_REF] Scannapieco | A multifluid interpenetration mix model[END_REF] or [START_REF] Wilson | Degradation of radiatively driven ICF capsula[END_REF].

In the last section, we give some enlightenments on the numerical methods for solving this mixing model. Some preliminary numerical results are also presented firstly in a one-dimension toy problem and secondly in a two-dimension problem where we account for electron temperature coupling as well as a simple laser energy deposition model.

We thank sincerly Gilles Carré for his work on the implemention of the numerical method and his help for tuning the numerical code.

Initial model for mixing of plasmas

We recall here the detailed model given in [START_REF] Decoster | Modeling of collisions[END_REF] (see particularly § §1.4 ,1.5). In the sequel, we do not account for the ion thermal conduction flux and the thermo-electric terms, for the sake of simplicity. We first write this model without the coupling between the electron and ion temperatures for focusing on its main features. The coupling between the electron and ion will be introduced in the second subsection.

Model without the coupling with the electron temperature

∂ ∂t ρ q + ∇. (ρ q u q ) = 0, (1)
∂ ∂t (ρ q u q ) + ∇. (ρ q u q u q + P q ) = Ξ q ,

∂ ∂t (ρ q E q ) + ∇. (ρ q u q E q + u q P q ) = Ω q (3)

knowing that the friction or drag term may read as

Ξ a = -Ξ b = νρ a ρ b (u b -u a ).
(uu denotes as usual a tensor). In the sequel, we assume that each fluid obeys a perfect gas law with a polytropic coefficient γ q ; that is to say we have P q = (γ q -1)ρε q .

In [START_REF] Decoster | Modeling of collisions[END_REF], the author gives the expression of the relaxation coefficient ν which depends on the temperatures of the two species and on the physical characteristics of the particles (and on the relative velocity when it is not small compared with sound speed). On the other hand, we have

Ω a = -Ω b = Ω # + νρ a ρ b (u b -u a ) m a u a + m b u b m a + m b , Ω # = ν # ρ a ρ b (T b -T a )
where m a , m b are the mass of the two ion species and coefficient ν # corresponds to the relaxation between the temperatures of the two species. As we will see in the sequel, the order of magnitude of the coefficient ν is crucial in order to determine whether the mixing of the two plasmas can occur or not.

Model taking into account the electron temperature

We account here for the coupling of the previous model with the electron population. Let us denote the electron density, temperature, energy per unit of volume and pressure by N e , T e , E e and P e . Firstly, if η a and η b are the ionization levels divided by the atomic mass, according to the electric neutrality of the plasma, we get

N e = η a ρ a + η b ρ b
Notice also that the electron velocity U e is not equal to the averaged ion velocity but is characterized by

N e U e = η a ρ a u a + η b ρ b u b
(that is to say the electric current is zero). Secondly, we know that

E e = 3 2 N e T e , P e = N e T e .
Of course the continuity equation (1) still holds.

The ion momentum and energy balance equations have to be modified to take into account the electron pressure, that is to say a force proportional to 1 Ne ∇P e (corresponding to the electrostatic field); so we get

∂ ∂t (ρ q u q ) + ∇. (ρ q u q u q + P q ) + η q ρ q N e ∇P e = Ξ q , (4)
∂ ∂t (ρ q E q ) + ∇. (ρ q u q E q + u q P q ) + η q ρ q N e u q .∇P e = Υ q + Ω q .

(5)

here Υ q denotes a relaxation term between the ion and electron temperatures and reads as

Υ a = ν E a ρ a N e (T e -ε a /C v,a ), Υ b = ν E b ρ b N e (T e -ε b /C v,b ).
where C v,q is the heat capacity and ν E q is a positive coefficient. Let us state now the electron energy balance equation accounting for the electron thermal conduction; so we set q th the Spitzer thermal conduction flux equal to -κ(T e )∇T e (where κ is roughly speaking proportional to T 5/2 e ) and we have

∂ ∂t E e + ∇. (E e U e ) + P e ∇.U e + ∇.q th = -(Υ a + Υ b ). (6)
Notice that the friction phenomenon between ions and electrons has not been accounted here since it is generally less important than the friction between the two ion species (notice that this is not always true , particularly at the very beginning of the collision process, but it is sufficient here to explain the main ideas of our method).

Of course, one can check easily that the total energy E e +ρ a E a +ρ b E b satisfies a balance equation

∂ ∂t (ρ a E a + ρ b E b + E e ) + ∇.q th = -∇. (ρ a E a u a + ρ b E b u b + u a P a + u b P b + U e (P e + E e ))
2 Reduced model

For the sake of simplicity, we withdraw in the first sub-section the coupling terms with the electron energy equation. The correction terms accounting this coupling are given in the second sub-section.

Model without the coupling with the electron temperature

Let us define the averaged physical quantities:

ρ = ρ a + ρ b , U = (ρ a u a + ρ b u b )/ρ (7)
are the mass density and the mean ion velocity V = u a -u b , the relative velocity c = ρ a /ρ, the mass concentration of material a (so we have ρ(1 -c) = ρ b ) P 0 = P a + P b , the global ion pressure.

The global specific internal energy ε is

ε = (ρ a ε a + ρ b ε b )/ρ = cε a + (1 -c)ε b (8)
Let us also denote

K = 1 2 |V| 2 c(1 -c)
which is the specific mixing kinetic energy. It is easy to check that the total energy E 0 satisfies

E 0 = ρ a E a + ρ b E b = ρε + 1 2 ρ|U| 2 + ρK.
One can show after standard calculus that :

∂ ∂t ρ + ∇. (ρU) = 0, (9)
∂ ∂t (ρU) + ∇. (ρUU) + ∇P 0 = -∇.(ρVVc(1 -c)), (10)
∂ ∂t (ρc) + ∇. (ρcU) + ∇. (ρc(1 -c)V) = 0, (11)
Moreover, if one defines σ = ρ ν, one may state an equation for the relative velocity

∂ ∂t V + ∇. (UV) -∇. |V| 2 2c -1 2 + B + A = -σV, (12)
B = -∇. VV 2c -1 2 + ∇. |V| 2 2c -1 2 + 2c -1 2 [(∇.V)V -(V.∇)V] + (V.∇)U -(∇.U)V. with A = (γ a -1) 1 ρ a ∇(ρ a ε a ) -(γ b -1) 1 ρ b ∇(ρ b ε b )
In the monodimensional case, B is equal to zero, so in the sequel, we always withdraw B.

Afterwards, the ion total energy balance equation reads as

∂ ∂t E 0 + ∇. (E 0 U) + ∇. (P 0 U + ρc(1 -c)VG) + ∇. (2ρKU -(2c -1)ρKV) = 0. with G = 1 ρ a (P a + ρ a ε a ) - 1 ρ b (P b + ρ b ε b ) = γ a ε a -γ b ε b .

The closures i)

We first assume that the mixing pressure tensor ρVVc(1 -c) reduces to a scalar one; that is to say in the r.h.s. of equation [START_REF] Resnyansky | A two-phase Shock-wave model[END_REF] we set -∇(2ρK).

ii) Secondly, we assume that there exists an increasing function Ψ of the variable c satisfying

Ψ(1) = γ a , Ψ(0) = -γ b (13)
with the following approximation for the two terms G and

A G =γ a ε a -γ b ε b ≃ εΨ(c) A = 1 ρ 1 c ∇P a - 1 (1 -c) ∇P b ≃ ε∇Ψ(c)
iii) We now make the crucial assumption that the friction coefficient σ is large enough if compared to the inverse of the characteristic evolution time. Thus, in equation ( 12), we withdraw all the terms but the two last ones, so it leads to the following approximation

V ≃ - ε σ ∇Ψ(c) = -D∇c, with D = ε σ Ψ ′ .
Remark 1. The closure for A may be justified by the following way. Near c = 1, if the flow is assumed to be isentropic, we have P a = P r (cρ) γa /ρ γa r and P r ,ρ r constant; so we may state

1 cρ ∇P a = γ a Pr ρ γa r ρ γa-1 c γa-2 ∇c + γ a Pr ρ γa r
c γa-1 ρ γa-2 ∇ρ; then after withdrawing the second term and

setting ρ γa-1 ≃ ρ γa-1 r , we get A ≃ 1 cρ ∇P a ≃ ε∇Ψ(c) with Ψ ′ (1) = γ a (γ a -1).
In the case where γ a = γ b = γ, the simplest choice for a odd function Ψ satisfying [START_REF] Richtmyer | Difference methods for initial value problems[END_REF] and this last condition is

Ψ(c) = γ 7 -γ 2 (c - 1 2 ) + 2(γ -3)(c - 1 2 ) 3 .
Otherwise, we may choose

Ψ(c) = γ 7-γ 2 (c -1 2 ) + 2(γ -3)(c -1 2 ) 3 + (γ a -γ b )/2, where γ = (γ a + γ b )/2.
Remark 2. If the flows for both species are assumed to be isentropic, then the closure for A comes from the closure for G, after withdrawing the gradient of ε.

According to the approximation of G and the relation [START_REF] Ishii | Thermo-Fluid dynamic theory of two phase flows[END_REF], we get

ε a = ε γ b + (1 -c)Ψ cγ b + (1 -c)γ a , ε b = ε γ a -cΨ cγ b + (1 -c)γ a ( 14
)
If the polytropic coefficients γ a and γ b are equal for the two ion species, we have P 0 = ρ(γ -1)ε. Otherwise, we may give an equation of state for the averaged fluid. Indeed, according to [START_REF] Youngs | Numerical simulations of turbulent mixing[END_REF], we get the following relation

P 0 =ρ(Γ(c) -1)ε, with Γ(c) = (γ a γ b + c(1 -c)(γ a -γ b )Ψ(c))/(cγ b + (1 -c)γ a );
it is a good equation of state for the averaged fluid, since Γ(c) -1 is strictly positive. Now, using the previous closures for A and V, the evolution equation for concentration and for the total energy reads as

∂ ∂t (ρc) + ∇. (ρcU) -∇. (ρc(1 -c)D∇c) = 0 ∂ ∂t E 0 + ∇. (E 0 U) + ∇. ((P 0 + 2ρK)U) -∇. [(c(1 -c)Ψρε -(2c -1)K)ρD∇c] = 0 ( 15
)
and the velocity equation reads as

∂ ∂t (ρU) + ∇. (ρUU) + ∇(P 0 + 2ρK) = 0
Now, to get a full model, it suffices to determine the mixing kinetic energy K. Of course, one could claim that K has to be equal to expression c(1 -c)D 2 |∇c| 2 /2, but this leads to complicate balance equation for the internal energy. Thus, we prefer to use a simple form for the internal energy balance and to deduce an evolution equation to determine K given by E 0 ρ -1 -ε -1 2 |U| 2 . Notice that according to initial equations ( 2 [which is justified when σK is large compared to ∇(c(1 -c)εD∇c). ] Then, from (15) one gets a coherent evolution equation for K which is

∂ ∂t (ρK) + ∇. (ρKU) + 2ρK∇.U + ∇ ((2c -1)ρKD∇c) = ∇. (ρc(1 -c)ΨεD∇c) -2σρK.
But since εΨ ′ = σD we may use the following relation

∇. (ρc(1 -c)ΨεD∇c) = σρc(1 -c)D 2 |∇c| 2 + Ψ∇. (ρc(1 -c)εD∇c)
Therefore, using the Lagrangian derivative D t = ∂ ∂t +U.∇, we may summarize the system which we address besides the continuity equation ( 9)

ρD t U + ∇(P 0 + 2ρK) = 0 ( 16
)
ρD t ε + P 0 ∇.U = 2σρK (17) ρD t c -∇. (ρc(1 -c)D∇c) = 0 ( 18
)
ρD t K + 2ρK∇.U + ∇ ((2c -1)ρKD∇c) -Ψ∇. (ρc(1 -c)εD∇c) = σρ c(1 -c)D 2 |∇c| 2 -2K . ( 19
) Of course, there is an energy balance relation associated to this system which is [START_REF] Wilson | Degradation of radiatively driven ICF capsula[END_REF]. Notice that when coefficient σ goes to infinity, coefficient D goes to 0 then this system reduces to the classical Euler equations for two fluids without mixing ; indeed if D would be set to zero and if the initial value of K would be zero, then K(t, .) is always equal to 0.

In equation ( 19), one sees that when coefficient 1/σ is small if compared to the characteristic time of evolution of the plasma, then K will be close to the quanity c(1 -c)D 2 |∇c| 2 /2. Remark 3. If the initial value c in of c is a Heavyside function, then the function c in remains a solution of equation ( 18), but it is unstable solution. On the other hand, it may be proved easily that for any bounded regular function D of x, the solution of (18) remains positive.

Remark 4 This diffusion concentration equation ( 18) has a form analogous to the one written in [START_REF] Landau | Fluids Dynamics[END_REF], (chap. VI, but in this last one, there are a pure diffusion term and a thermo-diffusion term). This kind of equation is also introduced in turbulence models, in motor combustion or in classical model for plasmas (see [START_REF] Hirshman | Transport of a multiple-ion species plasma[END_REF] , [START_REF] Youngs | Numerical simulations of turbulent mixing[END_REF] for instance). It may be also compared to the one introduced [START_REF] Decoster | Modeling of collisions[END_REF] ( §9.5). In this last case (besides a thermo-diffusion term), one can find an equation which reads with our notations

ρD t c -∇. (ρD 0 ∇c) = 0, (20)
with D 0 depending on the temperature and on the plasma components but independent of c.

Notice that if one wants to follow the expansion of a mixing domain with a model like (20), the relative velocity becomes

V = - 1 c(1 -c) D 0 ∇c = -D 0 ∇ log c 1 -c ,
which is infinite near the front of the mixing domain (when c tends to 0 or 1). Our model overcomes this difficulty and may be interpreted as a limitation of the mass flux at the boundary of the mixing model. In a certain way, we have changed the unbounded function log c 1-c by a bounded one ψ(c).

Remark 5. For modeling the mixing, we use the formula given by [START_REF] Decoster | Modeling of collisions[END_REF] (p. 119) and we assume that γ a = γ b = 5/3. Then we get

σ = ρβ 0 9π 2α 2 0 1/3 2 3 (ε a + ε b ) + V 2 r -3/2 , with β 0 = 4πη 2 a η 2 b (m a + m b )q 4
e log Λ ab where q e is the electron charge, log Λ ab a Coulomb logarithm (which is in the order of some units), α 2 0 is a dimensionless coefficient in the order of 1 depending on the plasma and V r is a corrector term in the order of the relative velocity. Then, by using the expressions of ε a and ε b (see ( 14)), D has the following form

D(c) ≃ ε β 0 Ψ ′ (c) ρ 9π 2α 2 0 1/3 4 3 + (2 -4c) Ψ(c) 5 ε + V 2 r 3/2 . ( 21
)
One has to evaluate the corrector term V r 2 , but it is possible to check that this added term is not negligible only when the ion temperature is not large, that is to say at the beginning of mixing phenomenon, so it may be seen that a very crude evaluation of this term is sufficient. As a matter of fact, we set

V r 2 = ε 5 /(Lρβ 0 ) 2
where L is a constant which is a characteristic value of the width of the mixing zone (it depends on the actual cases and in our hot plasma cases it is set to some micrometers).

According to the fact that |Ψ(c)| ≤ 5/3, we can see that D(c) ≥ C ε 5/2 ρβ 0 for some constant C.

Model accounting for the electron temperature

This technical section is devoted to the statement of the extension of the the previous mixing model, when we take into account the coupling between the ion and electron temperatures (as described in section 1.2). For the sake of simplicity, we assume here that γ g = γ l = γ. Recall that

N e = ρ(cη a + (1 -c)η b ).
Let us first evaluate the quantity Υ 0 = Υ a + Υ b , we get

Υ 0 = ρ 2 (cη a + (1 -c)η b ) (cν E a + (1 -c)ν E b )T e -cν E a ε a C v,a -(1 -c)ν E b ε b C v,b .
Recall that according to the closure we have made, the relations ( 14) express ε a and ε b with the help of ε ; therefore, it leads to

Υ 0 = ρ 2 (cη a + (1 -c)η b )(cν E a + (1 -c)ν E b ) T e - ε C v,0 ,
where C v,0 is given by

cν E a + (1 -c)ν E b C v,0 = ν E a c(γ + (1 -c)Ψ) γC v,a + ν E b (1 -c)(γ -cΨ) γC v,b . (22)
Thus, quantity Υ 0 has a classical form Ω(T e -ε/C v,0 ), with Ω positive. Now, defining ζ = (η a -η b ) /(cη a + (1 -c)η b ), we get a modified expression for the ion energy balance

ρD t E 0 ρ + ∇. (P 0 U + 2ρKU) + U.∇P e -∇. ((c(1 -c)Ψρε -(2c -1)K)ρD∇c) = Υ 0 -ζc(1 -c)D∇c.∇P e (23)
Moreover, for the electron energy, equation ( 6) reads by a classical way and the full model becomes

ρD t (ρ -1 ) -∇.U = 0, ( 24
)
ρD t U + ∇(P 0 + P e + 2ρK) = 0, ( 25
)
ρD t ε + P 0 ∇.U = 2σρK + Υ 0 , (26)
ρD t c -∇. (ρc(1 -c)D∇c) = 0, (27)
ρD t K + 2ρK∇.U + ∇ ((2c -1)ρKD∇c) = Ψ∇. (c(1 -c)ρεD∇c) +σρ[c(1 -c)D 2 |∇c| 2 -2K] + ζc(1 -c)D∇c.∇P e ( 28
)
ρD t (E e /ρ) + P e ∇.U -∇.(E e c(1 -c)ζD∇c) -P e ∇.(c(1 -c)ζD∇c) + ∇.q th = -Υ 0 . ( 29
)
Remark The structure of these three energy equations is satisfying because we can observe the coupling terms 2ρσK, 2ρ 2 ν 0 e K and Υ 0 appear with opposite signs in the right hand side of these equations. On the other hand, one can check that the global energy E 0 + E e satisfies a conservative balance equation.

Numerical issues

The aim of this work was the implementation of the previous model in a 2D multi-material hydrodynamics code for plasma flow. Without going into the details of the numerical methods for the simulations without the mixing model, we only recall briefly that for the hydrodynamics part of the model, one uses a Wilkins scheme which is based on a staggered grid method; that is to say, the density, the energy, the pressures and the concentration are evaluated in the center of each cell and the velocity at the cell vertices (in a one-dimension framework, it reduces to the classical Von Neuman-Richtmyer scheme, cf. [START_REF] Richtmyer | Difference methods for initial value problems[END_REF]).

Some remarks on the numerical methods

At each time step, there are two stages. At the first stage of the leapfrog method, one updates the velocity and one moves each node with its velocity. At the second stage, one updates the density and one solves the internal energy equation.

The modifications due to the mixing model are the following. The mixing kinetic energy K and the concentration c are evaluated at the center of each cell. First one has to deal with the concentration equation (18): we make a time implicit discretization of the diffusion term

∂ ∂x c(1 -c)D ∂ ∂x c
the key point is the evaluation of the mean value of [c(1 -c)D] i+1/2 on each interface i + 1/2 between the cells i and i + 1. We state in the one-dimension framework

c(1 -c)D ∂c ∂x i+1/2 = 2 c i+1 -c i (δx i+1 + δx i ) (c i+1 (1 -c i+1 )D i+1 + c i (1 -c i )D i ) ,
where the c i are the unknown values; we have the analogous for the two-dimension framework. So at each time step, one has to solve a non linear matrix system by an iterative method. Then, using the value of the normal relative velocity -Dn.∇c at the interfaces between two neighbouring cells (n is the unit normal to the interface), one solves the mixing kinetic energy equation (28). It is done with a finite volume scheme; the only technical point is the advection term (2c -1)ρKD∇c which is handled with an upwind technique which is explicit with respect to the time (in each cell, it depends on the sign of (2c -1)n.∇c in the neighbouring cells). Of course, the time step has to satisfy a classical stabilty criterium for this upwind technique which is of the type |(2c -1)Dn.∇c|δt/δx ≤ 1, but this criterium is generally less strong than the standard one related to the Wilkins scheme v ther δt/δx ≤ 1, where v ther is the sound speed.

Moreover, for the ion internal energy equation ( 26) , one has to take into account extra terms in the right hand side : no difficulty for the term 2σρK and the coupling term Υ 0 which reads as Ω(T e -ε/C v,0 ) is handled classically by implicit evaluation of T e and ε.

Numerical results

a. One-dimension toy problem.

In a first stage, we consider here a toy one-dimension problem and we wish to compare the numerical results obtained with the proposed model and the ones obtained with a multifluid code which handles a model of the type of [START_REF] Decoster | Modeling of collisions[END_REF], (cf. MULTIF code [START_REF] Chenais-Popovics | Kinetic to thermal energy transfer and interpenetration in the collision of laser-produced plasmas[END_REF]). We consider here a simple case, that is to say without laser and without radiative loss: there are two slabs occupied initially by a pure gold plasma which are separated by vacuum. Both fluids are characterized by the same polytropic coefficient γ (equal to 5 3 ) and the same heat capacity coefficient; their initial densities are equal to 0.0029 g/cm 3 and their initial temperatures equal to 1 keV. The initial velocities are 9.5 10 7 cm/s for the left fluid and -9.5 10 7 cm/s for the right one. For the evaluation of D defined by (21), we take α 0 = 0.84 (since the two fluids are identical [START_REF] Decoster | Modeling of collisions[END_REF]) and β 0 = 2.22 10 36 (in CGS units) which corresponds to log Λ ab = 5 (a mean value of the Coulomb logarithm in the mixing region during the time of simulation).

With these values, the slowing down distance of the ions on the other ion population in the mixing region is in the order of 130 µm and the corresponding slowing down time can reach 160 ps.

The spatial grid is 1000 µm long and composed of about 1000 meshes. As the slowing down distance represents several meshes, a strong collision with mixing between the two fluids occurs.The ionization level fixed to approximately 58, that is to say the electron equation of state is For the closure function Ψ, we have chosen as above. Moreover, in the above mentioned Spitzer formula for the electron thermal flux, we take κ(T e ) = k(x)T 5/2 e , where k(x) is a function depending weakly on the space variable (due to the so-called electron Coulomb logarithm).

E e = 3 2 58ρT e .
We show on figure 1, the density profiles at 300 ps, in the standard case without accounting for the mixing model and the same density profiles in the case where the mixing model is taken into account by using the above value of σ 0 . Notice that 300 ps corresponds to a few relaxation times. We can observe that the level of the maximum of the density is about 0.015 g/cm 3 without the mixing model but only 0.012 g/cm 3 with the mixing model, it is closer to the reference value which is 0.011 g/cm 3 given by the MULTIF code.

Notice that in the simulations with the mixing model, we have to give an ad hoc value for the length L, in the formula defining the corrector term V 2 r ; but we may check that if we change this value by a factor 2, the above values of density do not vary by more than 5%. Moreover, choosing other functions for Ψ (satisfying [START_REF] Richtmyer | Difference methods for initial value problems[END_REF]) leads also to small variations for these values of density.

Two dimension simulations

In this simulation, we have two gold diskes which are separated by vacuum. Each disk is heated by some laser beams; it is expanded and a collision between the two diskes occurs. We have performed numerical simulations with the two-dimension axi-symmetric code FCI2 [START_REF] Buresi | Laser Program Development at CELV[END_REF] with an arbitrary Lagrange-Euler framework accounting for laser heating (by a ray-tracing method) and electron thermal conduction. We have compared the results of a simulation made without mixing model and another one accounting for mixing model from the time 1.9 ns (which is the time of collision between the two diskes).

Thus, on figure 2, one may compare the density without and with mixing model at time t = 2 ns. Moreover, on figure 3, one compares firstly the ion temperature without and with mixing model at time t = 2 ns and secondly, on figure 4, the map of electron temperature at the same time.

It may be noticed that the ion temperature quickly increases up to 150 10 6 K in the collision region with or without mixing. This is due to the fact that almost all the initial kinetic energy of the two fluids is very quickly converted into internal energy.

Conclusion

We have proposed a model for mixing of two plasmas; it has been extended easily to the framework of the two-temperature Euler system. The underlying assumption is that the relative velocity between the two fluids is small compared to the sound speed of the flow. This model consists in coupling the classical Euler equations with a concentration equation (which is a non linear diffusion equation) and an evolution equation for the mixing kinetic energy K. The implementation of this model has been performed in a two-dimension Arbitrary-Lagange-Euler code where electron thermal conduction, radiative conduction and laser enegy deposition are taken into account.

Without the mixing model, when a strong collision between the two fluids occurs, a bump of the density ρ may be observed in the collision region. Accounting the mixing model, the numerical results show a decrease of the level of this bump; this leads to a notable decrease of the radiative emission which is, roughly speaking, proportional to the square of ρ. Notice that our model may be also generalized to mixing accounting for three fluids, using two concentration diffusion equations.

  + ∇. (ρUU) + ∇(P 0 + P e ) = 0, ∂ ∂t E 0 + ∇.(UE 0 ) + ∇.(UP 0 ) + U.∇P e = Υ 0 , ∂ ∂t E e + ∇.(UE e ) + P e ∇.U = -Υ 0 .

For

  each fluid (a and b ) denote by ρ a , ρ b the densities, u a , u b the velocities and by E a , E b the total specific energies. So ε a = E a -1 2 |u a | 2 and ε b = E b -1 2 |u b | 2 are the specific internal energies. The initial model is the following (with q = a or b)

  )(3) and the previous closures, the balance equation for the internal energy ρ a ε a + ρ b ε b would read ∂ ∂t (ρε) + ∇. (ρεU) + P 0 ∇.U = -Ψ∇.(ρc(1 -c)εD∇c) + 2σρK. But, we prefer to set simply ∂ ∂t ρε + ∇. (ρεU) + P 0 ∇.U = 2σρK.

Figure 1 :

 1 Figure 1: Profile of the density at time t = 300 ps, without the mixing model (solid line), with the mixing model (doted line), compared to the result given the MULTIF code (the lowest dashed line).

Figure 2 :

 2 Figure 2: Map of the density at 2. ns, without (on the left) and with the mixing model (on the right).

Figure 3 :

 3 Figure 3: Map of the ion temperature at 2. ns, without (on the left) and with the mixing model (on the right).

Figure 4 :

 4 Figure 4: Map of the electron temperature at 2. ns, without (on the left) and with the mixing model (on the right).


Katia Chirkova Crlao

On principles and practices of language classification

This article is an overview of issues in language classification, in particular in connection with three subgroups of the Sino-Tibetan language family: Tibetic, Sinitic and Qiangic. First, I discuss the practical application of currently prevalent classifications of Tibetic, Sinitic and Qiangic languages. Then, with reference to insights from classification practices in biology, I review alternatives for a practical classification in linguistics.

Introduction

It is a great privilege and pleasure to contribute to the Festschrift in honor of Alain Peyraube. I gladly and admiringly join this celebration of his scientific achievements, the diversity of his interests, his renowned erudition, and his foresight on issues of great impact to the field. In this article, I reflect on one of Alain Peyraube's most recent interests: linguistic classification from a general and more specific, East Asian, perspective. His contributions to the issue of classification [START_REF] Peyraube | L'homme aux six mille langues[END_REF](Peyraube , 2007) ) touch upon major trends in the field of East Asian linguistics, namely:

(1) unabating interest in classification issues, from periodic revisions of the precise composition of traditional groupings of languages in China and East Asia (such as Sino-Tibetan, Austronesian or Austro-Asiatic) to the construction of new macrofamilies (such as Austric, Austro-Tai, Sino-Tibetan-Austronesian or Proto-East-Asian)

(2) strong emphasis on interdisciplinary studies, whereby linguistics is coupled with other disciplines (such as genetics, evolutionary biology, archaeology or anthropology) to test relevant phylogenetic hypotheses

(3) fruitful exchange between biology and linguistics, inspired by a deeply entrenched view of parallelism between the development of living organisms and the development of languages These major trends guide and shape daily linguistic practice. Since classification is fundamental to all aspects of linguistic research, my main concern in this article is the practical application of existing classifications and the impact of the aforementioned trends on ongoing synchronic investigations. Given my area of expertise, the ensuing discussion will be limited to the field of Sino-Tibetan languages.

The theory of evolution and genetic classification in linguistics

In a metaphor that is nearly as old as linguistics itself, language is often viewed as an organism. Ranging from a direct equation1 to a more conventional symbolic portrayal, the conceptualisation of language as an organism has long accompanied the development of the discipline. 2 Indeed, parallels between languages and organisms abound (e.g. [START_REF] Janda | On language, change, and language change -or, of history, linguistics, and historical linguistics. The Handbook of Historical Linguistics[END_REF]Joseph 2003:59-81, Peyraube 2007), allowing for profitable exchange between linguistics and biology.

Fundamentally, both biology and linguistics are grounded in the theory of evolution, which both disciplines hold to be the most adequate and unitary hypothesis to which a great variety of both biological and linguistic phenomena can be related. As a result, the most widespread and influential type of classification in linguistics is genetic, based on the assumption of common ancestry of languages and using basic vocabulary, sound correspondences and, whenever available, grammatical (essentially morphological) evidence as classification criteria. This type of classification is so basic to linguistics that no explicit justification of its validity and legitimacy is normally deemed necessary. Characteristically, the recent edition of the Encyclopaedia of Languages and Linguistics (2006) names genetic classification "the most satisfying way to group languages" for "professional linguists and general readers" (Blake 2006:446), but does not explain what is so satisfying about it. One of the few scholars to explicitly address this issue, Michael Noonan (2008:4), points out that genetic classification is not always useful in providing information about the structure of a language, especially on higher taxonomic levels. Instead, he notes, information about where in the world a language is spoken provides more useful information about grammatical structure. On the other hand, Noonan argues, genetic classification has proven a boon to historical linguistics and theories of language change. Ultimately, Noonan concludes, the most satisfying aspect of genetic classification is its sheer familiarity, since "What family does it belong to?" is typically "the first thing a linguist will ask on being told of an unfamiliar language." In other words, genetic classification is popular due to its familiarity, but its practical value with regard to concrete linguistic data is questionable. Nonetheless, genetic classification is considered superior to other types of linguistic classification, for example, those based on typological features (see §3 for examples). This assumed superiority arguably results from the fundamental assumption that a classification based on evolution, which is in turn responsible for the existence and structure of natural languages, is imbued with explanatory power. Due to the deeply entrenched connection between explanation and theory in science, explanatory (viz. theoretical) analyses and frameworks are routinely valued higher than descriptive ones (cf. Dryer 2006:212-214). 3The dominance of genetic classification is the result of the sweeping success of Indo-European historical linguistics, which, in turn, is due to the clear evidence of a common origin for groups of languages in Europe (e.g. Romance languages) and the extensive amount of written records. The main metaphor that is used by Indo-European linguistics to explain the historical relationships, the genealogical family tree, is based on the assumption that the phylogenetic relationship between groups of languages is analogous to the genealogical relationship between individual languages, i.e. descent along a single genetic line (inheritance or vertical transmission). 4 The straightforwardness of the model and its appeal in the Indo-European context stimulated its export to other linguistic areas.

Uncritical application of the family tree model to languages outside the Indo-European family has repeatedly engendered criticism, most recently and forcefully in [START_REF] Dixon | The Rise and Fall of Languages[END_REF] or Aikhenvald and Dixon (2001). The conclusion of these two studies is that the family tree model is appropriate and can be proved through reconstruction in many cases (as in the case of Indo-European languages), but is problematic or even unsuitable in those cases where historical evolution is obscured by contact-related phenomena (as in the case of Australian languages). In other words, linguistic convergence and areal diffusion obscure historical development, so that originally genealogically unrelated languages may develop surface resemblances, forming groups which are synchronically similar, but genetically polyphyletic (that is, with multiple ancestors).

Borrowing and contact-induced change (i.e. horizontal transmission) as factors affecting the linear evolution of languages (i.e. vertical transmission) are now commonly taken into account in mainstream historical linguistics. Nonetheless, the impact of contact on linear descent remains subject of heated debate in linguistics, as in biology (cf. Dagan and Martin 2006:118.1). At present, awareness of the importance to enrich the dominant family tree model by taking into account "all the processes that affect language formation and development" (Chappell 2001:354) is on the increase. Much research is being undertaken to find objective criteria for assessing the balance between the effects of vertical and horizontal transmission on language change. Some of the most acute challenges include distinguishing loans from cognates, and common innovations from retentions.5 Attempts to "undo" the effects of contact include a methodology for separating innovations resulting from genetic inheritance from those diffused through contact, and a fully articulated theory of the relative diffusibility of features in a contact situation (e.g. [START_REF] Curnow | What languages features can be 'borrowed'? Areal Diffusion and Genetic Inheritance: Problems in Comparative Linguistics[END_REF][START_REF] Dench | Descent and diffusion: The complexity of the Pilbara Situation. Areal Diffusion and Genetic Inheritance: Problems in Comparative Linguistics[END_REF]. To date, all these attempts have been proven futile [START_REF] Curnow | What languages features can be 'borrowed'? Areal Diffusion and Genetic Inheritance: Problems in Comparative Linguistics[END_REF][START_REF] Dench | Descent and diffusion: The complexity of the Pilbara Situation. Areal Diffusion and Genetic Inheritance: Problems in Comparative Linguistics[END_REF] and are even argued to be doomed to failure [START_REF] Thomason | On the unpredictability of contact effects[END_REF].

In those cases where the linear development is known to be obscured by contact-related phenomena, as in the case of many East Asian languages (e.g. [START_REF] Benedict | Sino-Tibetan: A Conspectus[END_REF][START_REF] Matisoff | Genetic versus contact relationship: prosodic diffusibility in South-East Asian languages. Areal Diffusion and Genetic Inheritance: Problems in Comparative Linguistics[END_REF], and where the horizontal dimension of language change is supposedly as important for a coherent assessment of data as the vertical dimension, the emphasis is invariably on the latter (viz. genetic descent). The large number of relevant publications characteristically reflects the dominance of the genetic approach and the associated genetic classification. Also symptomatic of this dominance is the fact that general studies in language classification are limited to genetic classifications (e.g. Ruhlen 1994, Campbell and[START_REF] Campbell | Language Classification: History and Method[END_REF]. Such studies mostly focus on higher taxonomic levels: families and macro-families. While major debates revolve around distant genetic relationships of considerable time depth (such as Nostratic or Amerind), classification issues at more empirically bound lower taxonomic levels have all along attracted less attention. This is remarkable, because it is at these lower taxonomic levels, where data are most voluminous and versatile, that the need for meaningful classifications to organize and store empirical knowledge is at its greatest. I have been confronted with the practical application of these low level classifications of languages (Tibetan and Chinese dialects (or Tibetic and Sinitic languages) and Qiangic languages) and will focus on this issue in the main part of this article.

Classifications of Tibetan and Chinese dialects and of Qiangic languages

Classification is a highly practical activity. It provides a practical means by which objects can be identified and compared, and the knowledge about them organized. In addition, a meaningful classification enables the researcher to make inductive generalizations concerning the classified data and to discover new knowledge of the variation and distribution of the examined attributes. In the words of Ernst Mayr (1976:427), a good classification, like a good scientific theory, "has a high predictive power with respect to the assignment of newly discovered species and the pattern of variation of previously unused characters. That classification is the best which is at least affected by such new discoveries".

How do currently prevalent classifications of Tibetic, Sinitic, and Qiangic languages square with this desideratum for classification? I deal with this question as a practitioner who would like to use existing classifications as a descriptive tool and as a basis for inductive generalizations, while profiting from their predictive power with respect to adding new dialects or languages to the existing classification schemes.

Importantly for the ensuing discussion, the three groups considered (Tibetic, Sinitic, Qiangic) can be divided into two groups, based on whether the genetic relationship among the languages within each group is presupposed or established by comparative linguistics.

The first group, comprising Tibetan and Chinese dialects, consists of languages whose respective genetic relatedness is a matter of oral and written memory, hence the designation "dialects", even though both Tibetan and Chinese dialects are, respectively, synchronically as diverse as the Romance or Germanic languages within the Indo-European language family. The knowledge of their respective genetic relatedness and even the accepted division into distinct groups by and large derive from non-linguistic sources: self-awareness, culture, history, and geographical distribution of the groups in question. In this context, existing linguistic classifications have mostly sought to translate this established division into linguistic terms. Chinese and Tibetan dialectal studies commonly hold that most substantial differences between dialects lie in phonology and lexicon, whereas syntactic differences are relatively negligible. Hence, the favoured criteria for classification are for the most part phonological. Several classifications have been advanced for both Tibetan and Chinese dialects (the mostly widely accepted classifications are considered in some detail below).

The second group, Qiangic, brings together languages whose genetic relatedness is not known, but hypothesized.

In terms of a genetic argument, all considered cases (Tibetic, Sinitic, Qiangic) are subgroups. That is to say, they are known or hypothesized monophyletic groups. A subgroup admits as evidence only shared innovations as distinguished from retentions or parallel developments (e.g. Harrison 2003:232-239).

The most widespread classification of Tibetan dialects spoken in China divides all dialects into three groups: Dbus-gtsang (or Central), Khams, and Amdo, which correspond to the three namesake historical provinces of Tibet. 6 In linguistic terms, this classification is based, among others, on the presence of tone, on the simplification of consonant clusters, on the presence of voiced obstruent initials, on the development of Old Tibetan clusters with medials -w-, -y-and -r-, on the number of consonant codas, and on the presence and the number of diphthongs. Viewed from the principles of classification outlined above, this classification yields little practical value. Overall, the underlying classification criteria are neither mutually exclusive nor representative for each respective group. As a result, they fail to single out homogeneous groups and to do justice to the actual diversity of Tibetan dialects. Consequently, the resulting classification does not allow for generalizations about the formed groups and it has no predictive power with respect to the assignment of newly described varieties. Take, for instance, the Báimǎ language (also known as the Báimǎ Tibetan dialect) of Northern Sìchuān and Southern Gānsù provinces in the People's Republic of China (PRC). The linguistic affiliation of Báimǎ has long been disputed, due to the popular belief among the Báimǎ people that they are descendents of the Dī people, an ancient Tibeto-Burman ethnic group who inhabited roughly the same area until their gradual assimilation into the Hàn and the Tibetans during the Táng Dynasty (Chirkova 2008c). If judged solely on linguistic grounds, the language of the Báimǎ people is predominantly Tibetan-like in all its linguistic sub-systems, even though it exhibits some non-Tibetan features in its lexicon, morphology, and syntax (Chirkova 2008b). Accordingly, Báimǎ is argued to be a Tibetan dialect (Huáng and Zhāng 1995:104;Zhāng 1997:134-135, 140). The problem, however, is that this new Tibetan dialect cannot be straightforwardly accommodated by the dominant tripartite dialect classification of Tibetan dialects spoken in China. 8 Based on the exhibited characteristics, Báimǎ is at the same time Khams-like (given the drastic simplification of the syllable structure, the elimination of Old Tibetan syllable codas, and the presence of tones) and Amdo-like (given, for instance, the development of Old Tibetan velar stop clusters with medial -r-into alveopalatal affricates and of Old Tibetan dental and labial stop clusters with the same medial into retroflex affricates). 9 In sum, the dominant classification of Tibetan dialects of China fails to accommodate Báimǎ. In fact, the same holds true for any of the newly discovered dialects in its neighbourhood, such as Zhongu (Sun 2003a), Chos-rje (Sun 2003b), or the Tibetan dialects of Western Sichuan (Suzuki 2008).

In his discussion of the popular tripartite classification of the Tibetan dialects of China, Jackson Sun (2003a:794-797) points out that this classification mixes possibly convergent, but independent developments (such as the presence of tone) with shared retentions (such as the preservation of Old Tibetan voiced obstruent onsets). This, in Sun's analysis, explains its inadequacy. The envisaged solution is to upgrade this routinely accepted (but overall unrevealing) classification into a skye rgu, Nángqiàn 囊謙 nang chen Counties of Yùshù Prefecture in Qīnghǎi Province, and Bāqīng 巴 青 sbra chen, Dīngqīng 丁青 steng chen and Suǒ 索 Counties in Tibet). 8 Presumably for this reason, Báimǎ is not included in any of the comprehensive outlines of Tibetan dialects of China (e.g. Zhāng 1996). 9 Chirkova 2008c, 2008d. meaningful (i.e. genetic) classification, which, in turn, should be based on shared innovations and phonological isoglosses to delimit precise dialect boundaries. Sun does not discuss the feasibility of establishing a classification along these lines and, to my knowledge, such a classification has not yet been attempted. However, a close investigation of a possibility of such a stricto sensu genetic classification of Chinese dialects by Laurent [START_REF] Sagart | On distinguishing Hakka and non-Hakka dialects[END_REF] suggests that such a classification is not viable.

The currently widely accepted classification of Chinese dialects is based on one single criterion: the development of the Middle Chinese voiced stops into the modern dialects [START_REF] Fang-Kuei ; | Chinese languages and dialects[END_REF], Yuán et al. 1961).10 While this criterion works fairly well in the majority of cases, as Norman (1988:181) points out, it is clearly valid only for a relatively recent historical period, given that the voiced initials probably persisted in the North down to at least the tenth or eleventh centuries. This traditional classification of Chinese dialects suffers from the same flaws as its Tibetan counterpart. It fails to do justice to the actual diversity of dialects and to form a reliable basis for generalizations regarding the formed groups. This problem is particularly acute in relation to dialect syntax (e.g. Chirkova 2008a for the Mandarin group).

The insufficiency of one criterion for a meaningful classification of Chinese dialects in general, and Hakka dialects in particular, has been pointed out by [START_REF] Sagart | On distinguishing Hakka and non-Hakka dialects[END_REF], who further notes that a choice of some other criterion would result in an entirely different classification.

Similar to Sun's assessment of the current classifications of Tibetic languages, Sagart views the non-conformity of the traditional classification of Chinese dialects to the standards of genetic classification as one reason for its shortcomings. He notes that the traditional classification of Chinese dialects defines dialect zones whose historical statuses are very different: some are defined by innovation, others by retention, whereas it is generally admitted that only innovations are useful in linguistic classification. However, in stark contrast with Sun's conclusion, Sagart finds a stricto sensu genetic classification of Chinese dialects impossible and even unnecessary. It is impossible, Sagart argues (1998:299), because innovations may be obliterated or reversed through contact. Hence, dialect groups that are currently recognized in Chinese dialectology are fuzzy entities that owe much of their make-up to contact (i.e. horizontal transmission) as opposed to vertical inheritance.11 This is not to say that a genetic classification of Chinese dialects is altogether unattainable. The linguistic history of the group can be reconstructed to the extent possible from written sources. However, as Sagart argues, a genetic dialect classification is unnecessary, because such a classification has at best a minor role to play in the furthering of knowledge on Chinese dialects. Instead, Sagart points out, geographical approaches to Chinese dialectology hold more promise.

The use of a small number of criteria and the use of common innovations as a valid criterion are disfavoured as basis for a meaningful classification in the case of well-documented dialects or languages known to be related (Sinitic). Curiously, these criteria gain in explanatory value in the case of little-known languages that are not known to be related. Take the Qiangic languages of the Chinese South-West as an example.

Qiangic is a putative group of eleven Sino-Tibetan languages of Sìchuān and Yúnnán provinces in the PRC. After Sūn (1983,2001), these eleven languages are often grouped in a northern branch and a southern branch. 12 The Qiangic subgroup in its current understanding was proposed in the 1970s by Chinese linguists. 13 It brings together geographically adjacent, non-written, and under-researched languages, which are not known to be related. This subgroup relies primarily on shared typology: the widespread use of measure words, directional prefixes and similar case marking, to which some, equally general phonological features such as the presence of uvulars and tones have also been added (Sūn 1983). Seen from mainstream historical linguistics, these are type-identifying features (Nichols 1996:48), in other words, features that are found in enough unrelated language families and that are consequently low in identificational value as to the genetic relationship. Overall, the percentage of common vocabulary among the languages of the Qiangic subgroup (the feature which, in contrast to the typological characteristics above, could provide more reliable support for the hypothesis of a common origin of the languages) is relatively small. Altogether it makes up less than 20% between any two Qiangic languages (Huáng 1991:355). Notably, the type-identifying features in Qiangic are not found in their totality in languages from the neighboring subgroups. At the same time, languages of these neighboring (and phylogenetically better-understood) subgroups spoken in the Qiangic area, appear to share with Qiangic languages many features that are held specific of the Qiangic subgroup. An example of such a language is Kami, a dialect of Tibetan spoken in the southern end of the Qiangic area. Kami has typically Qiangic and at the same time strikingly non-Tibetan features, such as directional prefixes, measure words, or uvular phonemes (Chirkova forthcoming). The presence of these 'Qiangic' features in a Tibetic language suggests that Qiangic is rather a Sprachbund: a geographically delimited area that includes languages that are genetically unrelated, but share a fair number of highly distinctive traits.

Even though the grounds for its classification are demonstrably weak, the Qiangic subgroup has found its way in all present classifications of the Sino-Tibetan language family. This is perhaps because even despite their diversity, Qiangic languages are possibly more closely related to each other than to anything else in Sino-Tibetan (Bradley 1997:36).

Formulated on such loose grounds, the Qiangic subgroup fails the test of a meaningful classification: it cannot serve as basis for generalizations or single out features specific for the group. In fact, the precise composition of this subgroup has been debated since its very proposal. In the 1980s, the aforementioned Báimǎ language was argued to be Qiangic, because it exhibits some features that are held to be typical of this subgroup, such as the presence of directional prefixes (Sūn 1980). Conversely, in the 1990s, one of the Qiangic languages proper, Nàmùyī, was claimed to be genetically related to the Ngwi languages instead [START_REF] Cíwò | Nàmùyīyǔ zhīshǔ yánjiū 納木依語支屬研究 [On the genetic affiliation of the Nàmùyī language[END_REF]Huáng 1997:13-15), based on a large amount of basic vocabulary shared between Ngwi and Nàmùyī.

Recent years have witnessed an upsurge of interest in Qiangic languages and linguistics. Many Qiangic languages are currently being investigated. Attempts are also undertaken to place this methodologically problematic grouping on a scientific footing: by focusing on common innovations, lexical comparisons and, eventually, 12 Northern Qiangic languages include Qiāng, rGyalrong (Jiāróng), Ěrgōng, Prinmi (Púmǐ) and Minyak (Mùyǎ). Southern Qiangic languages comprise Zhābà, Choyo (Quèyù), Ěrsū, Nàmùyī, Shǐxīng and Guìqióng. Qiangic is a recent subgroup, located at the intersection of Tibetic, Ngwi-Burmese, and Na languages. 13 For the history of the Qiangic subgroup, see [START_REF] Sun | Review of Zangmianyu Yuyin He Cihui / Tibeto-Burman Phonology and Lexicon[END_REF]. phonological reconstructions. 14 Viewed in the context of the above discussion, such a classification is bound to be based on few criteria, which are furthermore spurious, as no objective criteria to distinguish common features (in lexicon or elsewhere) from accidental similarity in form, borrowing, and/or genetic drift have yet been found (cf. Thomason and Kaufman 1988:6). 15 In fact, there is not necessarily a correlation between the number of underlying criteria and the quality of the resulting classification. Even one character may be sufficient, if it is known to be correlated with many others of varied nature. If the correlation has been established, the use of the one feature will amount to the use of all correlated features together. Unfortunately, to my knowledge, lexical sharings are not correlated with any other linguistic feature. The use of this one single criterion is hence bound to result in a biased classification, whose validity is furthermore contingent on the acceptance of lexical sharings as probative of genetic relatedness.

What would be the gains of the genetic approach in the case of Qiangic languages, if pursued rigorously? A genetic argument in the subgroup context can ultimately lead to either of these two possibilities:

First, the comparative method could show that the Qiangic subgroup is polyphyletic. The fact that Qiangic languages have grown to share a number of salient features is then to be explained through effects of convergence and hybridisation. Notably, the same conclusion can be directly arrived at by a closer investigation of empirical data, without resorting to the formal apparatus of lexical comparisons, loanword stratification, search for cognate etyma, establishment of sound correspondences, and phonological reconstructions.

Second, the comparative method could show that the Qiangic group is monophyletic. In that case, the extreme synchronic diversity of the languages forming the group, and the relatively meagre amount of shared cognate vocabulary, would have to be explained. The comparative method cannot provide this explanation and does not even posit such among its explicit goals. Consequently, a closer investigation of empirical data will again be required if any sense is to be made of the linguistic situation and the historical development of the relevant languages.

The following three conclusions can be made on the basis of the three cases reviewed in this section.

First, the dominant status of genetic classification and the expectations associated therewith are at odds with the actual practical merits of genetic classification. After all, genetic classification reflects known or hypothesized lines of descent, and in so doing, it restricts the available data in a way that excludes much information.

Second, in those cases where the line of descent is unknown (e.g. Qiangic), to proceed with an explanatory approach (genetic classification) before assembling a representative set of empirical data, runs the risk of applying potentially subjective criteria and of involving circular reasoning. 14 The establishment of a genetic relationship in the Sino-Tibetan context must rely-faute de mieuxon shared derivational morphology and lexical comparisons, in spite of the fact that both can be borrowed, given a sufficient degree of contact. Because few inflectional paradigms can be identified across the family, shared and especially irregular inflectional paradigms cannot be relied on to prove the existence of the family (Handel 2008: 435). 15 Subgrouping is known in mainstream comparative linguistics to be methodologically and factually problematic. However, its limitations are argued away on the basis of the practical necessity of subgrouping assumptions for the application of the comparative method (Harrison 2003 : 239). Third, as suggested by studies of well-documented languages that are known to be related (Sinitic, Tibetic), it is impossible to combine information on genetic descent with maximum information on synchronic similarities between the classified varieties into one classification scheme. Therefore, at least two approaches (and classifications) are warranted to organize these two types of information, for example, one genetic and one geographical (cf. Sagart for Sinitic).

Overall, given the little practical value of existing classifications of Sino-Tibetan languages, genetic or otherwise, and given the overall pressure to reflect on phylogenies, a linguist working with Sino-Tibetan languages cannot but get involved in classification work. Indeed, the lion's share of production in Sino-Tibetan linguistics deals with genetic classification issues: from individual languages to subgroups, subgroups to families, and families to macro-families (cf. [START_REF] Handel | What is Sino-Tibetan? Snapshot of a field and a language family in flux[END_REF]. While existing linguistic classifications are far from satisfactory and no immediate linguistic solutions are in sight, it may be worthwhile, given the success of the "language as an organism" metaphor, to turn to biology for inspiration and possible solutions.

From classifications in biology to classifications in linguistics

In its simplest form, classification in biology is the ordering of entities into groups on the basis of their relationships, that is, of their associations by contiguity, similarity, or both (cf. Simpson 1961:9, quoted from Sokal andSneath 1967:3). 16 Notably, the term "relationship" may imply relationship by ancestry or, alternatively, it may indicate the overall similarity (or affinity) as judged by the characters of the organisms without any implication as to their relationship by ancestry. It has been argued that no classificatory scheme is able simultaneously to yield information on the degree of synchronic similarity and genetic descent [START_REF] Michener | A quantitative approach to a problem in classification[END_REF][START_REF] Michener | A quantitative approach to a problem in classification[END_REF]Sneath 1967:55-57). Consequently, two major types of classification exist: empirical (synchronic) and genetic (diachronic). The empirical classification is often referred to as "natural", i.e. denoting a system that occurs naturally in the empirical world. The two types of classification, natural and genetic, differ in scope, in purpose, and in the number of underlying criteria.

A natural classification is one that endeavors to utilize many attributes of the classified entities and to accurately reflect most of the various natures of the objects. While it is not explanatory as to the biological phenomena that it describes, a natural classification forms solid basis for generalizations concerning the described data and it is useful for a wide range of purposes (e.g. Gilmour 1940:472). An example of an empirical classification in biology, that groups animals and kinds according to their shared physical characteristics, is the familiar Linnaean taxonomy widely used in the biological sciences.17

A natural system of classification in biology is contrasted with more or less 'artificial' systems, which have no natural occurrence prior to their construction by the researcher. Such an 'artificial' classification is based on a limited range of attributes, selected with an eye to a specific classification purpose. Consequently, it is useful for limited purposes (Blackwelder 1967:186).

Given that evolution is responsible for the existence and structure of the natural system, classification in biology (as in linguistics) is possible because evolution has produced diverse kinds (or languages). From this it is possible to generalize that classification is based on evolution, but it must be recognized that this means that it is based on the results of evolution: empirical facts, synchronically observable concrete living organisms. Evolution can naturally be seen as the central cause in biology. Nonetheless, it cannot be used as an explanatory concept or an independent basis for classification. The main reason for that is that common ancestry is never known as a fact, but merely hypothesized. Furthermore, the phylogeny is based, on the one hand, on previous classifications and, on the other hand, on known empirical facts, so that there are no examples of classifications that are actually based on what were taken to be phylogenetic facts as distinct from the features of living organisms (Blackwelder 1967:191-193, Sokal andSneath 1967:56-57). The genetic type of classification, which focuses on providing a summary of existing knowledge of phylogeny, is hence a subsidiary classification (e.g. Gilmour 1940:473).

To return to the field of linguistics, we can observe that linguistics does not have in its arsenal any analogy to natural classification in biology. Instead, genetic classification has come to be used as an explanatory concept and a general linguistic classification. It is, however, a secondary classification, because it selects few criteria (essentially basic vocabulary in the domain of Sino-Tibetan languages) and it pursues one specific domain of inquiry (the common ancestry of languages). As such, it falls short of accounting for the versatility of empirical data, and it is unable to accommodate new data, as discussed above on the basis of Tibetic, Sinitic, and Qiangic languages. While a significant and intriguing facet of linguistic inquiry, indispensable for the purposes of historical linguistics, it is unsuitable as general linguistic classification.

Successful in those cases where written records (linguistic "fossil record") are available and linguistic similarities self-evident (through history or mutual intelligibility), genetic classification reaches its limits in those cases where relationships between languages are less evident and written records are scarce or absent. Since no reliable method to ascertain characters probative of historical relatedness has been developed, much caution is required when using genetic arguments for classification purposes. Of course, there is no harm in speculating about possible relationships between languages in the latter type of situation, but in linguistics, as in biology (cf. Sokal and Sneath 1967:8), it appears more informative to concentrate on collection and analysis of the maximum of empirical facts, because facts inform and advance progress in linguistic theory and historical linguistics.

Discussion

To return to my departing point as a practitioner in need of a practical linguistic classification, I cannot but conclude that a genetic classification might not be "the most satisfying way to group languages" for "professional linguists and general readers", as it is praised to be.

There is need for another general linguistic classification, which will be satisfying in the sense that it can give justice to the versatility of the subject matter of linguistics, while providing a practical means to organize our knowledge of various facets of the continuously growing body of data. To serve this purpose, this classification is bound to be empirical, i.e. based on a variety and quantity of features.

On a more abstract level, in classification as elsewhere in linguistics, there is no fundamental contrast between description and explanation, as there is naturally need for both descriptive and explanatory approaches [START_REF] Dryer | Descriptive theories, explanatory theories, and basic linguistic theory[END_REF]. While there is no intrinsic incompatibility in describing a language and then explaining the things described, it appears desirable that both types of approaches enjoy equal footing, so that explanatory approaches do not channel descriptive ones. The primary goal in linguistics, as in any area of scientific inquiry, remains first and foremost to describe a set of facts without any particular theoretical implication. This set of facts is further open to synthesis, hypothesis generation and interpretation, including, for instance, grading of observed similarities in terms of their respective historical status (retentions, innovations, parallel developments), and assessing the balance between horizontal and vertical lines of transmission in language change. All in all, my claim is that we need both empirical and genetic classifications, as well as many other types of classification (for instance, typological) for various specific purposes.

In terms of close parallels between linguistics and biology, the two disciplines continue to go hand in hand. Both owe their success to the theory of evolution with the tree model as its central concept. This theory and the tree concept have all along exerted strong influence on the development of the two disciplines. For over a century, linguistics and biology were mainly concerned with fleshing out a universal tree in their respective disciplines. The ongoing revision of this concept in biology, where evidence begins to mount up that a substantial component of the evolutionary process is non-tree-like [START_REF] Doolittle | Phylogenetic classification and the universal tree[END_REF][START_REF] Dagan | The tree of one percent[END_REF][START_REF] Lawton | Why Darwin was wrong about the tree of life[END_REF], is bound to have a major impact on linguistics. In fact, linguistics (where the need to incorporate horizontal transmission into ideas about linguistic evolution has long been evident) is already undergoing a similar paradigm shift.

In sum, further progress is possible only if a broad range of varied data is taken into consideration, and this broad range of data requires well-grounded theories and adequate classifications to be well-understood.

  Figure 1. Classification of Tibetan dialects(adapted from Hú 1991: 177)

Recently, language as a parasitic species[START_REF] Kortlandt | A parasitological view of non-constructive sets[END_REF]

[START_REF] Harrison | On the limits of the comparative method[END_REF][START_REF] Mufwene | The Ecology of Language Evolution. Cambridge approaches to language contact[END_REF] or language as a symbiont [[START_REF] Driem | Studies on Sino-Tibetan Languages: Papers in Honor of Professor Hwang-cherng Gong on this Seventieth Birthday[END_REF]van Driem , 2008]. 2 See van Driem 2008 for an overview.

In the words of[START_REF] Popper | Conjectures and Refutations: The Growth of Scientific Knowledge[END_REF][START_REF] Sokal | Principles of Numerical Taxonomy[END_REF]: 103): "The scientist aims at finding a true theory or description of the world[…], which shall also be an explanation of the observable facts."

The drawbacks of the group to an individual analogy have been pointed out, among others, byGilmour (1940: 469-470) for systematics and, recently,Mufwene (2001: 16-21) for linguistics.

A relevant example from the field of Sino-Tibetan linguistics is the Bái language, whose genetic affiliation has long been the subject of scholarly speculation. The ongoing controversy is due to the so far unresolved issue of whether numerous lexical similarities between Chinese and Bái reflect longterm contact between the two languages (loans) or inherited features (cognates).

Cf. Gésāng 1964[START_REF] Ǎitáng | Zàngyǔ fāngyán de yánjiū fāngfǎ 藏語方言的研究方法 [Methods in the study of the Tibetan dialects[END_REF][START_REF] Péng | Zàngyú Jiǎnzhì 《 藏 語 簡 志 》 [Outline of the Tibetan language[END_REF]:114-145, Hú 1991:173-177, Zhāng 1993[START_REF] Nichols | The comparative method as heuristic[END_REF]

Alternative classifications of Tibetan dialects spoken in China also exist. Targeting the presence of tone as the sole criterion, Hú[START_REF] Tǎn | Zàngyǔ 《藏語》 [The Tibetan language]. Hàn-Zàng Yǔ Gàilùn 《漢藏語概論》 [A general introduction to Sino-Tibetan languages[END_REF] arrives at a division of all dialects into just two groups: Khams-gtsang and Amdo (Figure1). Increasing the number of criteria to five (presence of tone, voiced obstruent initials, consonant clusters, consonant codas and diphthongs),[START_REF] Bùfán | Yùshù Zàngyǔ de yǔyīn tèdiǎn hán lìshǐ yǎnbiàn guīlǜ 玉樹藏語的語音特點 含歷史演變規律 [Characteristics in Yushu Tibetan phonology and its rules of historical change[END_REF] arrive at four groups: Dbus-gtsang, Khams, Amdo and Yùshù 玉樹 (spoken in Záduō 雜多 rdza rdo, Jiégǔ 結古

Other classifications of Chinese dialects exist, based on different sets of criteria, e.g.[START_REF] Norman | Chinese[END_REF] and[START_REF] Wurm | Language Atlas of China[END_REF]. These classifications by and large correspond to the currently widely accepted division of Chinese dialects into groups in the tradition classification by[START_REF] Fang-Kuei ; | Chinese languages and dialects[END_REF].

This is a received view in Chinese dialectology, e.g.[START_REF] Pulleyblank | Chinese dialect studies. Languages and Dialects of China[END_REF]Pulleyblank ( , 1998)).

The compilation in this section is based on Huxley (1940),[START_REF] Sokal | Principles of Numerical Taxonomy[END_REF],[START_REF] Blackwelder | Taxonomy: A Text and Reference Book[END_REF],[START_REF] Jardine | Mathematical Taxonomy[END_REF],[START_REF] Salthe | Evolutionary Biology[END_REF],[START_REF] Mayr | Evolution and the Diversity of Life: Selected Essays[END_REF]Mayr ( , 1988)),[START_REF] Bailey | Typologies and Taxonomies: An Introduction to Classification Techniques[END_REF], and[START_REF] Pigliucci | Making Sense of Evolution: The Conceptual Foundations of Evolutionary Biology[END_REF].

A phenetic view of classification in biology, according to which the sole aim of classifications is to reflect as accurately as possible the relative similarities or dissimilarities of populations, is frequently said to be pioneered by Michel Adanson in his Familles des Plantes (1763) (cf.Jardine and Sibson 1971: 136-137).


Aimé Lachapelle

Human Crowds and Groups Interactions: A Mean Field Games Approach

This paper is devoted to some Mean Field Games (hereafter MFG) modeling of human crowds behavior. More precisely, we study 2population dynamics (each of whom consisting of a continuum of individuals) with opposit interests. We focus on the crowd aversion case inside the group and also onto the other group (xenophobia in some cases). We write a macro-Nash problem between the two populations, then we give an existence and uniqueness result and characterize optimal points as MFG solutions. Finally we provide a simple gradient descent method to approximate the solutions and show some simulations.

I. INTRODUCTION MFG have been recently introduced by Lasry &

Lions ( [START_REF] Lasry | Jeux à champ moyen. I. Le cas stationnaire. (French) [Mean field games. I. The stationary case[END_REF], [START_REF] Lasry | Jeux à champ moyen. II. Horizon fini et contrôle optimal. (French) [Mean field games. II. Finite horizon and optimal control[END_REF], [START_REF] Lasry | Mean field games[END_REF]) and seem to be very useful to model big groups interactions with "intelligent" individuals (control aspects). To the best of our knowledge there exist very few papers on the topic (MFG and human crowds), e.g. the works of the author (see [START_REF] Lachapelle | Jeux à champ moyen, simulation numérique et application à un modèle de choix technologiques[END_REF]) or Guéant (see [START_REF] Guéant | A reference case of Mean field Games[END_REF]). The present paper is based on [START_REF] Lachapelle | [END_REF] and is mainly motivated by the macroscopic crowd motion models (see the recent work of Buttazzo, Jimenez & Oudet ( [START_REF] Buttazzo | An optimization problem for mass transportation with congested dynamics[END_REF]), Carlier & Salomon ([3]), Hughes ([5]), or Maury, Roudneff-Chupin & Santambrogio ( [START_REF] Maury | A macroscopic crowd motion model of gradient flow type[END_REF]) for a gradient flow setting). Nevertheless, in our case, the MFG use corresponds to discrete-continuous fundaments. We adopt, as in [START_REF] Lachapelle | Computation of Mean Field Equilibria in Economics[END_REF], the optimal control point of view of MFG and take fully advantage of it to prove existence and to approximate solutions, in particular using the transformation originally due to Benamou & Brenier (see [START_REF] Benamou | A computational fluid mechanics solution to the Monge-Kantorovich mass transfer problem[END_REF]). But the main point of this work is that we particularly focus on the modeling of interactions between two groups in some xenophobia situation. The paper is organized as follows. In section II we briefly motivate the use of MFG for human crowds and present the problem we consider. Section III is devoted to the study of optimality conditions (that is in fact the link between MFG and optimal This work was supported by the french A.N.R, "Programme blanc OTARIE" A. Lachapelle, CEREMADE, UMR CNRS 7534, Université Paris IX Dauphine, Pl. de Lattre de Tassigny, 75775 Paris Cedex 16, France. lachapelle@ceremade.dauphine.fr control) and existence. In the next part (section IV) we provide a numerical strategy to approximate the solutions, based on a gradient descent method. Finally, in section V we close with some simulations.

II. TWO-GROUP DYNAMICS: THE PROBLEM

A. MFG & human crowds

As mentioned in the introduction, it seems to be very natural to think of MFG to model human crowds. The MFG setting is well adapted for several reasons. First it is a micro-macro approach, being an approximation of N -player differential games when N tends to infinity (huge systems of individuals). Secondly it is useful to model interactions between "intelligent" individuals (or agents), that is: they are rational players with rational expectations. The next reason is that MFG describe non-cooperative equilibrium configurations (Nash point asymptotics). Note also that since MFG deal with atomized and anonymous agents (this last assumption is quite natural in human crowds), the mathematical framework (optimal control of Fokker-planck, MFG system) enables numerical simulations. To fix ideas, we now turn to recall briefly what are MFG in our crowd motion setting.

We work on the d-dimensional Torus Ω := T d to avoid boundary difficulties (note that they could be treated, as in [START_REF] Lachapelle | Computation of Mean Field Equilibria in Economics[END_REF] for instance). In the present work, we focus on the finite horizon framework. Define then the time-space domain, Q = [0, T ] × Ω. To simplify, we consider in this paragraph a unique group (or population) composed with a continuum of individuals. It is fully characterized during the time period by the evolutive measure (m t ) 0≤t≤T (we will abusively refer to it as a density) with initial situation m 0 given. Then, the agent starting at x ∈ Ω at the beginning of the period evolves controlling the drift of the following stochastic process

dX x t = α t dt + σdW t , X x 0 = x,
where α t is the control parameter, W t a standard Brownian motion and σ the given noise.

Her individual problem is then to minimize (over a certain class of controls α) the quantity

E T 0 L(X x t , α t ) + V [m t ](X x t )dt + Ψ(X x T ) ,
where the Lagrangian L stands for a control and position cost. In what follows we take the quadratic cost L(x, α) := |α| 2 2 (note that it does not depend on the position). Function Ψ is a final cost (incentive to reach a certain area) and V is a state cost (that is a criteria depending on the mean field created by the others, the density of agents). This last cost is the key point of the MFG modeling: individuals anticipate the crowd evolution (m t ) 0≤t≤T and take it into account in their optimization problem. It can be found in [START_REF] Lasry | Jeux à champ moyen. II. Horizon fini et contrôle optimal. (French) [Mean field games. II. Finite horizon and optimal control[END_REF], [START_REF] Lasry | Mean field games[END_REF], that the solution of the continuum of individual problems satisfies the so called MFG system

∂ t m -σ 2 2 ∆m + div(m∂ p H(x, ∇v)) = 0 ∂ t v + σ 2 2 ∆v + H(x, ∇v) = V [m]
, with in the one hand initial and transversality conditions: m| t=0 = m 0 v| t=T = Ψ, in the other hand α := -∇v and H is the Legendre transform of L. It is well known (see [START_REF] Lachapelle | [END_REF], [START_REF] Lachapelle | Computation of Mean Field Equilibria in Economics[END_REF], [START_REF] Lasry | Jeux à champ moyen. II. Horizon fini et contrôle optimal. (French) [Mean field games. II. Finite horizon and optimal control[END_REF]) that when V is the derivative (e.g. the Gâteaux derivative) of a potential Φ on bounded measures on Ω, i.e. V = Φ ′ , then the critical points of the optimal control problem of Fokker-Planck

inf α Q |α| 2 2 m + T 0 Φ(m t )dt + Ω Ψm(T ) ∂ t m -σ 2 2 ∆m = -div(αm), m(0) = m 0 , (1)
are solutions of the MFG system. Note also that it is a sufficient condition as soon as Φ is convex, which is the case in this paper. The convexity corresponds in fact to the crowd aversion (contrary to attraction situations, see for instance [START_REF] Lachapelle | [END_REF], [START_REF] Lachapelle | Computation of Mean Field Equilibria in Economics[END_REF]).

In what follows, we call problem (1) a global optimization problem of a MFG.

B. Writing the Nash problem

Let us focus on the case where two populations interact inside Ω. We want to study equilibria between the two groups (typically Nash points as suggested by Lasry and Lions in [START_REF] Lasry | Jeux à champ moyen. II. Horizon fini et contrôle optimal. (French) [Mean field games. II. Finite horizon and optimal control[END_REF]). Before giving the problem of group i, i = 1, 2, let us recall a classical notation. For any point x = (x 1 , x 2 ) ∈ R 2 , and for a fixed coordinate i, we denote by

x -i the element of R, x -i := {x 1 , x 2 } -{x i }.
Formally, the global optimization problem (linked to a continuum of individual problems) of group i, given the control and the mass evolution of the other group (i.e. (α -i t , m -i t )), reads as:

inf α i J i λ (α)
where

J i λ (α) := Q |α i | 2 2 m i + T 0 Φ i λ (m 1 t , m 2 t ) + Ω Ψ i m i (T ). (2)
From now on, m i depends on α i , more precisely it is viewed as a bounded nonegative measure (i.e. belonging to the set M b (Q, R + )) which is a weak solution of the Fokker-Planck equation:

∂ t m i - σ 2 2 ∆m i = -div(α i m i ), m i (0) = m i 0 .
(3) We distinguish the populations by considering different initial densities m i 0 (.) and different final incentive costs Ψ i . However we study the simple case where the Brownian motion and the noise are similar for both groups. In the definition of the criterion (2), the coupling cost we consider is

Φ i λ (m t ) := Ω (m i t ) 2 + λm 1 t m 2 t ,
for a nonnegative real constant λ. This models of a situation involving individual's aversion to the members of her own group and onto the other group, so that λ stands for a "xenophobia" parameter. To see that, note that the individual mean field criterion is:

V i [m t ](x) = 2m i t (x) + λm -i t (x).
In fact we consider that xenophobia is significant when the aversion λ to the other group is higher than the one to the own group (equal to 2). The Nash problem between the groups is then:

Find α = (α 1 , α 2 ) such that: (N ) J i λ (α) = inf α i ∈M b (Q,R d ) J i λ (α i , α -i ), for i = 1, 2.

III. EXISTENCE & OPTIMALITY

A. Optimality conditions

In this part we present a characterization of the Nash equilibria. To do so, let us introduce the MFG system for two groups: for i = 1, 2,

∂ t m i - σ 2 2 ∆m i + div(m i ∇v i ) = 0, m i (0) = m i 0 , (4)
∂ t v i + σ 2 2 ∆v i + |∇v i | 2 2 = Φ i λ (m) ′ , v i (T ) = Ψ i , (5)
and the joint minimization problem (for the two groups)

(Q) inf α=(α 1 ,α 2 ) J λ (α) := J 1 λ/2 (α) + J 2 λ/2 (α),
under the constraints: m i is a solution of (3), for i = 1, 2. Note that J λ is convex if the xenophobia parameter is not too large, that is λ ≤ 2.

We are now in the position to give the optimality conditions. Proposition 1: If λ ≤ 2 then the following assertions are equivalent:

1) α ∈ M b (Q, R d ) is a solution of (N ) and m satisfies (3) for α = α, 2) α ∈ M b (Q, R d ) is a solution of (Q) and m satisfies (3) for α = α, 3) (m, v
) is a solution of the MFG system (4)- [START_REF] Hughes | A continuum theory for the flow of pedestrian[END_REF], with α = α = ∇v.

If λ > 2 then it is only necessary i.e. 2. ⇒ 1.,3.

A proof of this statement is based on classical arguments of differential calculus and can be found in [START_REF] Lachapelle | [END_REF].

B. Existence

If there exists a solution of problem (Q), then proposition 1 ensures the existence of a Nash point between the two groups. Before giving an existence result, we may reformulate the problem in a more rigorous way, following [START_REF] Benamou | A computational fluid mechanics solution to the Monge-Kantorovich mass transfer problem[END_REF], [START_REF] Buttazzo | An optimization problem for mass transportation with congested dynamics[END_REF], [START_REF] Lachapelle | Computation of Mean Field Equilibria in Economics[END_REF]. We adopt a vectorial point of view, and use the following notations:

• M b (Q, R d ) is the set of bounded d-vectorial measures on Ω. • m = (m 1 , m 2 ) ∈ M b (Q, R 2
+ ), and for all

x = (x 1 , x 2 ) ∈ R 2 + * , 1 x := ( 1 x 1 , 1 x 2 ), • q = (q 1 , q 2 ) ∈ M b (Q, R 2d ), and for all y = (y 1 , y 2 ) ∈ R d × R d , |y| 2 := (|y 1 | 2 , |y 2 | 2 ), • A := {(q, m) ∈ M b (Q, R 2d ) × M b (Q, R 2 + ) : Q (∂ t u + σ 2 2 ∆u)dm i + Q ∇u.dq i = Ω Ψ i m i (T ) -u 0 m i 0 , ∀u ∈ C ∞ (Q) , pour i = 1, 2}.
Let us introduce the (q = αm, m) formulation, i.e. following [START_REF] Buttazzo | An optimization problem for mass transportation with congested dynamics[END_REF], the functions

ϕ 1 (a, b) := |a| 2 2b if (a, b) ∈ R 2d × R 2 + * +∞ else, ϕ 2 (b) := |b| 2 + λb 1 b 2 if b ∈ R 2 + *
+∞ else, and

K(q, m) := K 1 (q, m) + K 2 (m), K 1 (q, m) := Q ϕ 1 ( dq dL d+1 , dm dL d+1 )dL d+1 , K 2 (m) := Q ϕ 2 ( dm dL d+1 )dL d+1 + Ω Ψ.dm| t=T ,
where L d+1 denotes the Lebesgue measure in R d+1 . Then we can rewrite K in a simpler form:

K(q, m) = J λ (α), if q << m and q = αm +∞ else.
In this setting, a rigorous formulation of

(Q) is inf (q,m)∈A K(q, m), (6)
and we can give the result as announced before. Proposition 2: If λ ≤ 2 and m 1 0 , m 2 0 ∈ L 2 (Q), then problem (Q) possesses a solution (which is unique as soon as λ < 2). Moreover there exists a Nash point i.e. a solution of (N ).

Remark 1: Note that existence does not fail when adding a constraint of the type m ≤ constant as in [START_REF] Maury | A macroscopic crowd motion model of gradient flow type[END_REF].

Given the reformulation (6), the proof is a simple adaptation of the one obtained by Buttazzo, Jimenez and Oudet in [START_REF] Buttazzo | An optimization problem for mass transportation with congested dynamics[END_REF]. A complete proof of this proposition is provided in [START_REF] Lachapelle | [END_REF]. In the next section we deal with defining a numerical procedure to approximate the solution.

IV. NUMERICAL SETTING

In this part we introduce the discretization and a gradient descent method in order to approximate the solution(s) of problem (N ). More precisely, we distinguish the cases when the joint problem (Q) is convex from when it is not. In the convex setting (i.e. when λ ≤ 2) we describe the gradient descent that we apply to the joint functional. The non-convex case λ > 2 (in which the xenophobia is significant) is more involved but interesting (we expect non-uniqueness). We then provide an alternating directions method taking advantage of the convexity of group i's problem, given group (-i)'s evolution.

Gradient First of all, let us write the gradient formula of the functional. We look at the reformulated problem given by [START_REF] Lachapelle | Jeux à champ moyen, simulation numérique et application à un modèle de choix technologiques[END_REF]. We slightly modify the point of view considering that the density m is an affine function of the momentum q. To fix ideas, the joint problem reads as:

inf q F (q) := i=1,2 Q |q i t | 2 2m i t +Φ i λ/2 (m t )+ Ω Ψ i m i T , (7)
where m i , i = 1, 2 solves :

∂ t m i - σ 2 2 ∆m i = -div(q i ) , m i (0, .) = m i 0 (.). (8
) Thanks to a classical differential calculus (see for instance [START_REF] Lachapelle | [END_REF]), it is easy to get the explicit formula of the gradient

∀(q, m) ∈ A, ∀w = (w 1 , w 2 ) ∈ M b (R 2d ), ∇F (q).w = Q q i m i + ∇θ i .dw i i=1,2 , (9)
where θ i satisfies, for i = 1, 2,

-∂ t θ i - σ 2 2 ∆θ i = - |q i | 2 2(m i ) 2 + (2m i + λm -i ), θ i | t=T = Ψ i . ( 10
)
Algorithm for the convex case Since problem ( 7) is convex when λ ≤ 2, we decide to apply a gradient descent method. We focus on the 2Dcase (d = 2) and take Ω = [0, 1] 2 with periodic boundary conditions. Let M and N be two positive integers, we define the time and space steps by dt = 1

N and dx = 1 M . For (i, j, k) ∈ A := {0, ..., N } × {0, ..., M } 2 , for a given function f defined on Q, f i j,k denotes the numerical approximation of f (idt, jdx, kdy). Equations ( 8) and [START_REF] Lasry | Jeux à champ moyen. II. Horizon fini et contrôle optimal. (French) [Mean field games. II. Finite horizon and optimal control[END_REF] are iteratively solved by using finite differences, after initializations (m 0 j,k = m 0 (jdx, kdy) and θ N j,k = Ψ(jdx, kdy)), using the following approximations

∂ t f (idt, jdx, kdy) = f i+1 j,k -f i j,k dt , ∆f (idt, jdx, kdy) = f i j+1,k -2f i j,k + f i j-1,k (dx) 2 + f i j,k+1 -2f i j,k + f i j,k-1 (dy) 2 . At step n, let f (n) := f i,(n) j,k (i,j,k)∈A
. Then the gradient descent method (hereafter GDM) is the following: 1) Initialization: Choose q (0) then compute m (0) by solving [START_REF] Lachapelle | Computation of Mean Field Equilibria in Economics[END_REF] with the finite difference scheme.

2)

Step n:

• Compute θ (n) by solving numerically [START_REF] Lasry | Jeux à champ moyen. II. Horizon fini et contrôle optimal. (French) [Mean field games. II. Finite horizon and optimal control[END_REF] with q (n-1) and m (n-1) , then compute the discretized gradient ∇F (q (n-1) ) (formula ( 9)), using θ (n) . • Compute the descent: q (n) = q (n-1)ρ n ∇F (q (n-1) ). • If ||q (n)q (n-1) || <Tol1, then stop the algorithm (Tol1 is a tolerance threshold defined by the user).

Else, n = n + 1. Note that ρ n above is the descent step size, it is chosen optimal, i.e. minimizing the following: ρ ∈ [0, 1] → F (q (n-1) -ρ∇F (q (n-1) )).

Alternating directions method for the nonconvex problem The case where aversion to the other group is significant (λ > 2), for which we have less theoretical results, also seems interesting. One of the main goal of the present work is to obtain numerical simulations in this situation. Consequently, it is convenient to describe the way we try to approximate the Nash points between the groups when λ > 2. To do so we choose an alternating directions method, provided the convexity of both group i's problem, given group (-i)'s evolution:

inf q i F i (q) := T 0 Ω |q i t | 2 m i t + Φ i λ (m t )dt + Ω Ψ i m i T ,
where m i solves (8) for q i , i = 1, 2. One can easily get the formula of the gradient of F i looking at the joint case ( 9)- [START_REF] Lasry | Jeux à champ moyen. II. Horizon fini et contrôle optimal. (French) [Mean field games. II. Finite horizon and optimal control[END_REF].

In what follows, by writing that we compute q we also mean that we compute the corresponding m solution of the discretized versions of the Fokker-Planck equations [START_REF] Lachapelle | Computation of Mean Field Equilibria in Economics[END_REF], for i = 1, 2.

The strategy to approximate the Nash points is to apply GDM successively to each group. Note that above, the upper index refers to the group number and the lower one to the iteration. 1) Initialization: Choose q 1 0 then compute q 2 1 with GDM and q 1 0 .

2)

Step k 1:

We know q 2 k . • Compute q 1 k then q 2 k+1 by using GDM (with, respectively, q 2 k and q 1 k ).

• If ||q i k -q i k-1 || <Tol2 for i = 1, 2, then stop the procedure. Else, k = k + 1.

V. SIMULATIONS

The GDM shows good convergence results when the initials densities of individuals are significantly positive (i.e. m i 0 > constant > 0). This last section is devoted to the exhibition of some very first tests. In the next simulations we take T = 1 and σ 2 2 = 0.01.

A. Test 1: crowd aversion in a single group

In the first example we focus on a case involving only one population (m 2 0 = 0), i.e. a similar framework as the one studied by Buttazzo, Jimenez and Oudet in [START_REF] Buttazzo | An optimization problem for mass transportation with congested dynamics[END_REF]. Fig. 1 shows the initial density of agents (centralized around the point (0.1, 0.1)) and the final cost, modeling a strong incentive for individuals to be in some neighborhood of (0.5, 0.8) and (0.8, 0.5) at instant T .

B. Test 2: groups interactions

Let us now look at the more interesting case with two populations. Recall that we look for Nash equilibria between two groups whose global optimization problem is [START_REF] Lachapelle | Jeux à champ moyen, simulation numérique et application à un modèle de choix technologiques[END_REF]. We use the procedure detailed before (starting with group 2). In order to emphasize the xenophobia behavior we choose λ = 20 in the definition Φ i λ (m 1 , m 2 ) = Ω (m i ) 2 + λm 1 m 2 . We consider a symmetric configuration and represent the graphs of the initial densities (m i 0 , i = 1, 2) and final costs (Ψ i , i = 1, 2) in Fig. 5. Group 1 is initially centralized around (0.35, 0.5), group 2 around (0.5, 0.35). Concerning the final costs the situation is still symmetric since they model incentives to reach positions in the neighborhood of (respectively for group 1 and 2) (0.65, 0.5) and (0.65, 0.5). With such a situation we are interested in crossing phenomenon. The graphs of both densities are depicted in Fig. 6. We can notice the same spreading over we observed in Test 1. However, the most interesting evolution period is described in Fig. 7. Indeed, we can see that group 1 gives the priority to group 2 to go to its attractive area passing through the center of the domain (the shortest road for the euclidian metric). Some of the individuals of group 1 wait, some others go through the border (periodic conditions), and the lasts go through the center (the most congested area). Anyway we note that group 2 reaches quicker than group 2 its goal. Looking at Fig. 8, we can check that both group are finally centralized around the points (0.65, 0.5) and (0.5, 0.65). The last remark is that we obtain the opposite (or symmetric) situation when starting to optimize on group 1. Then, the symmetry break seems to be a consequence of this starting choice.

VI. CONCLUSIONS AND FUTURE WORKS

A. Conclusions

In this work we provide a macroscopic model for human crowds and groups interactions for intelligent individuals. It is a MFG model also inspired by [START_REF] Buttazzo | An optimization problem for mass transportation with congested dynamics[END_REF]. We take advantage of the optimal control point of view of MFG in order to prove existence and to develop a numerical approach which consists of a simple gradient descent method. We also test the algorithm in a case involving xenophobia between two groups.

B. Future works

To the best of our knowledge this paper is one of the first works exploring both theoretically and numerically a MFG approach to model groups interactions and human crowds, so that many things remain to be done. For instance we think of: proving the gradient method convergence, testing its robustness when σ → 0, or taking into account congestion effects (i.e. considering moving cost of the type L = |α| 2 m b with b > 0).

5 Fig. 2 . 1 Fig. 3 .

 5213 Fig. 1. Data

Fig. 4 .

 4 Fig. 4. Value of F for each iteration

  (a) Initial densities (b) Final costs Ψ i

Fig. 5 . 4 Fig. 6 .

 546 Fig. 5. Data

( a ) 75 Fig. 7 .Fig. 8 .

 a7578 Fig. 7. Group 2 go through the middle of the domain as some agents of group 1, nevertheless most of group 2 individuals transit by the border (periodic conditions)

VII. ACKNOWLEDGMENTS

The author wishes to express his thanks to G. Carlier, J. Salomon & F. Santambrogio for helpful discussions and G. Buttazzo for his active interest on the topic.