Luigi Cantini
email: luigi.cantini@cyu.fr

Ali Zahra
email: ali.zahra@cyu.fr

Hydrodynamic behavior of the 2-TASEP

. We find that the auxiliary variables, introduced previously in the literature to express the density dependence of particle currents, turn out to be the Riemann variables of the conservation equations. This allows us to work out quite explicitly the rarefaction and shock solutions and to completely solve the associated Riemann problem. Our theoretical results are confirmed by Monte Carlo simulations.

Introduction

The asymmetric simple exclusion process (ASEP) is a minimal model of transport in (quasi) one-dimensional systems. It consists of particles which occupy the sites of a one dimensional lattice with only one particle allowed on each lattice site. These particles hop under the effect of an external driving force which breaks detailed balance and creates a stationary current. This model was introduced in the late 60s in biology to model translation in protein synthesis [START_REF] Macdonald | Kinetics of biopolymerization on nucleic acid templates[END_REF] and independently in probability [START_REF] Spitzer | Interaction of Markov processes[END_REF] and afterwards it has found a wide spectrum of applications, ranging from theoretical and experimental studies of biophysical transport [START_REF] Chou | Non-equilibrium statistical mechanics: from a paradigmatic model to biological transport[END_REF] to modeling traffic flow [START_REF] Chowdhury | Statistical physics of vehicular traffic and some related systems[END_REF][START_REF] Evans | Bose-einstein condensation in disordered exclusion models and relation to traffic flow[END_REF]. As soon one considers models which are more suited for physical/biological systems, one will encounter variants of ASEP containing localized or mobile defects and several species of particles, which have different behaviors. As a result, typically these models are not exactly solvable and even for some of the most basic questions, like the study of large scale behavior of the system (which in the case of ASEP is known to be described by the Burgers equation [START_REF] Rost | Non-equilibrium behaviour of a many particle process: Density profile and local equilibria[END_REF][START_REF] Benassi | Hydrodynamical limit for the asymmetric simple exclusion process[END_REF][START_REF] Rezakhanlou | Hydrodynamic limit for attractive particle systems on Z d[END_REF][START_REF] Kipnis | Scaling limits of interacting particle systems[END_REF]), approximations schemes like mean-field are necessary.

In this paper we address the question of the large scale or hydrodynamic behavior of an exactly solvable multispecies generalization of ASEP, consisting of two kinds of particles, •-particles and •-particles, moving in opposite directions. One can think of them as opposite charged particles moving under the influence of an external electric field or as cars moving on two opposite lanes. Each site of a one-dimensional lattice is either empty or occupied by one of the two kinds of particles. For convenience, empty sites can be treated as a third species of particles, the * -particles. In continuous time, a •-particle jumps forward on empty sites with rate β, while a white •-particle jumps backward on empty sites with rate α. On top of this, an adjacent pair •• swaps to •• with rate 1.

• * → * • rate β * • → • * rate α • • → • • rate 1 (1)
This model has appeared in the literature under different names. It has been first considered in [START_REF] Derrida | Statphys-19: 19th IUPAP Int. In Conf. on Statistical Physics[END_REF][START_REF] Mallick | Shocks in the asymmetry exclusion model with an impurity[END_REF], where the stationary measure on a finite periodic lattice was written in a matrix product ansatz form [START_REF] Derrida | Exact solution of a 1-D asymmetric exclusion model using a matrix formulation[END_REF][START_REF] Blythe | Nonequilibrium steady states of matrix-product form: a solver's guide[END_REF]. It is also a particular case (q = 0) of the so called AHR model [START_REF] Arndt | Spontaneous breaking of translational invariance in one-dimensional stationary states on a ring[END_REF][START_REF] Arndt | Spontaneous breaking of translational invariance and spatial condensation in stationary states on a ring. I. the neutral system[END_REF][START_REF] Rajewsky | Spatial particle condensation for an exclusion process on a ring[END_REF], in which the swap •• → •• is allowed with rate q. Being a natural 2-species generalization of TASEP we shall call this model 2-TASEP.

It turns out that the 2-TASEP is Yang-Baxter integrable [START_REF] Cantini | Algebraic Bethe ansatz for the two species ASEP with different hopping rates[END_REF] for arbitrary values of the parameters α and β. It belongs indeed to a larger family of integrable multispecies exclusion processes introduced in [START_REF] Cantini | Inhomogenous Multispecies TASEP on a ring with spectral parameters[END_REF]. Bethe ansatz techniques can be used to solve exactly for the long time limit behavior of the generating function of the currents [START_REF] Cantini | Algebraic Bethe ansatz for the two species ASEP with different hopping rates[END_REF][START_REF] Derrida | Bethe ansatz solution for a defect particle in the asymmetric exclusion process[END_REF]. More recently, in the case α + β = 1, the transition probabilities as well as the joint current distribution for some specific initial distribution of a finite number of • and •-particles have been obtained [START_REF] Chen | Exact confirmation of 1d nonlinear fluctuating hydrodynamics for a two-species exclusion process[END_REF][START_REF] Chen | Limiting current distribution for a two species asymmetric exclusion process[END_REF], and an asymptotic analysis of these results has allowed to prove that the joint current distribution is given by a product of a Gaussian and a GUE Tracy-Widom distribution in the long time limit, as predicted by non-linear fluctuating hydrodynamics [START_REF] Van Beijeren | Exact results for anomalous transport in one-dimensional hamiltonian systems[END_REF][START_REF] Spohn | Nonlinear fluctuating hydrodynamics for anharmonic chains[END_REF][START_REF] Ferrari | Coupled kardar-parisi-zhang equations in one dimension[END_REF].

When α + β = 1 the stationary measure factorizes and the currents have a simple expression as function of the densities. In [START_REF] Fritz | Derivation of the leroux system as the hydrodynamic limit of a twocomponent lattice gas[END_REF][START_REF] Tóth | Perturbation of singular equilibria of hyperbolic two-component systems: a universal hydrodynamic limit[END_REF] the hydrodynamic limit of the 2-TASEP for α = β = 1 2 has been studied and proven to converge to the classical Leroux system of conservation laws [START_REF] Leroux | Analyse et approximation de problèmes hyperboliques non linéaires[END_REF][START_REF] Serre | Existence globale de solutions faibles sous une hypothèse unilaterale pour un système hyperbolique non linéaire[END_REF]. The Leroux system is a notable example of a Temple class system i.e. a 2-components conservation law whose shock and rarefaction curves coincide [START_REF] Temple | Systems of conservation laws with invariant submanifolds[END_REF]. The theory of Temple class systems shares several common features with the theory of single component conservation laws [START_REF] Serre | Systems of Conservation Laws 2: Geometric Structures, Oscillations, and Initial-Boundary Value Problems[END_REF], in particular well-posedness results for Temple systems are available for a much larger class of initial data compared to general systems of conservation laws.

For arbitrary α and β only numerical results based on mean field approximation are available [START_REF] Arndt | Spontaneous breaking of translational invariance and spatial condensation in stationary states on a ring. ii. the charged system and the two-component burgers equations[END_REF]. In the present paper we study the exact hydrodynamic equations of the 2-TASEP and show that they are Temple class for arbitrary α and β. This allows us to compute their rarefaction and shock solution, as well as to solve completely the Riemann problem, which consists in determining the density profile starting from a domain-wall initial data.

The paper is organized as follows. In Section 2 we review and expand on results about the 2-TASEP currents obtained in [START_REF] Cantini | Algebraic Bethe ansatz for the two species ASEP with different hopping rates[END_REF]. The core of the paper is Section 3 where the conservation laws are studied. We derive the rarefaction waves as well as the shock solutions and finally we solve the full Riemann problem. In Section 4 we compare the prediction of the hydrodynamic equations with Monte Carlo simulations. Conclusions and some outlooks for further works are discussed in Section 5.

Currents

In this section we reproduce and expand the results of the analysis in [START_REF] Cantini | Algebraic Bethe ansatz for the two species ASEP with different hopping rates[END_REF] in a convenient way, which makes manifest the symmetries of the model. In order to compute the particle currents as functions of the local densities, we consider our model on a periodic ring with a fixed number of particles of each species. Call M i the number of particles of species i, they are related to N , the length of the ring, by N = M • + M • + M * . Let the system evolve starting at time t = 0 from an arbitrary fixed configuration and call n i,j (t) the number of swaps of consecutive ordered pairs of particles of type i, j up to time t. This number increases by +1 each time two consecutive ordered particles of species i, j exchange their position i, j → j, i. The average rate of swaps i, j → j, i in the steady state is just given by lim t→+∞ 1 t E [n i,j (t)], irrespectively of the initial state. The particle currents in the steady state are hence given by

J i = lim t→+∞ 1 N t E j =i n i,j (t) -n j,i (t) (2)
In our case, it is convenient to introduce the following quantity

Φ(ν •,• , ν •, * , ν * ,• ) = lim t→+∞ 1 t E [ν •,• n •,• (t) + ν •, * n •, * (t) + ν * ,• n * ,• (t)] . (3)
The currents are obtained as specialization of Φ(ν

•,• , ν •, * , ν * ,• ) J • = Φ(1, 1, 0) N , J • = Φ(-1, 0, -1) N , J * = Φ(0, -1, 1) N . (4)
In [START_REF] Cantini | Algebraic Bethe ansatz for the two species ASEP with different hopping rates[END_REF] the function Φ(ν •,• , ν •, * , ν * ,• ) was shown to be given by the solution of the following equation det

G(Φ(ν •,• , ν •, * , ν * ,• ), ν •,• , ν •, * , ν * ,• ) = 0. ( 5
)
where the matrix

G(Φ, ν •,• , ν •, * , ν * ,• ) is given by G(Φ, ν •,• , ν •, * , ν * ,• ) =    Φ F α [M • , M • , M * ] F β [M • , M • , M * ] ν •,• M • + ν * ,• M * F α [M • + 1, M • , M * ] -F β [M • , M • + 1, M * ] ν •,• M • + ν •, * M * -F α [M • , M • + 1, M * ] F β [M • + 1, M • , M * ]    (6) with F γ [a, b, c] := 0 dz 2πi 1 z a (z -1) b (z -γ) c . (7)
When one of the particle species is strictly absent (i.e. when one among M • , M • , M * vanishes) the model reduces to a single species TASEP and it is not difficult to see that one of the currents vanishes, while the others boil down to the usual TASEP current. On the other hand in the following we shall assume that at least one particle per species is present (M i = 0) and we shall be mainly interested in the thermodynamic limit of these quantities as N → ∞, with lim N →∞ M i N = ρ i . We shall see, as already found in [START_REF] Derrida | Statphys-19: 19th IUPAP Int. In Conf. on Statistical Physics[END_REF][START_REF] Mallick | Shocks in the asymmetry exclusion model with an impurity[END_REF], that the presence of even a single particle of a given species (i.e. an infinitesimally vanishing but not strictly zero density) can affect the macroscopic behavior of the system. With this in mind, we consider the limit a -→ ∞, with b/a and c/a fixed, of the function

F γ [a, b, c], that behaves like ‡ F γ [a, b, c] ∼ 1 z a γ (z γ -1) b (z γ -γ) c
, where z γ is the zero of the saddle point equation a z + b z-1 + c z-γ = 0, belonging to the interval [0, min [1, γ]]. Applying this expression in eq.( 5) we get in the thermodynamic limit lim

N →∞ Φ(ν •,• , ν •, * , ν * ,• ) N = (ν •,• ρ • + ν * ,• ρ * )z α (1 -z β ) + (ν •,• ρ • + ν •, * ρ * )z β (1 -z α ) (8)
where with z α ∈ [0, min(1, α)] and z β ∈ [0, min(1, β)] are the solution of the saddle point equations The result for the currents then reads

ρ • z α + ρ • z α -1 + 1 -ρ • -ρ • z α -α = 0 (9) ρ • z β + ρ • z β -1 + 1 -ρ • -ρ • z β -β = 0. ( 10
J • = z α (z β -1) + ρ • (z α -z β ) (11)
J • = z β (1 -z α ) + ρ • (z α -z β ) (12)
J * = ρ * (z α -z β ) (13)
Notice that eqs. [START_REF] Benassi | Hydrodynamical limit for the asymmetric simple exclusion process[END_REF][START_REF] Rezakhanlou | Hydrodynamic limit for attractive particle systems on Z d[END_REF] are invariant under exchange ρ • ↔ ρ • , α ↔ β and z α ↔ z β . This implies as expected, that under exchange ρ • ↔ ρ • and α ↔ β we have J • ↔ -J • . Let us finish this section by showing how some known results fit in the analysis above.

♦ β = 1. In this case •-particles do not distinguish •-particles from * -particles and so they behave just as particles in a single species TASEP. This is reflected in eq.( 10), where ρ • disappears and one finds z β = ρ • , which replaced in eq.( 12) gives

J • = ρ • (1-ρ • ). The case α = 1 is completely analogous: z α = ρ • and J • = ρ • (ρ • -1). ♦ α + β = 1.
In this case it is known that the stationary measure takes a factorized form [START_REF] Rajewsky | Spatial particle condensation for an exclusion process on a ring[END_REF]. At the level of the currents, we have indeed

J • = -ρ • (ρ • + αρ * ) and J • = ρ • (ρ • + βρ * ).

The z variables

In our analysis the variables z = (z α , z β ) will play a prominent role, it is therefore important to work out their domain of definition D z (α, β) corresponding to the physical domain D in the variables ρ = (ρ

• , ρ • ), ρ • , ρ • ≥ 0, ρ • + ρ • ≤ 1.
First of all we have already seen above that z has to satisfy z α ∈ [0, min(1, α)] and z β ∈ [0, min(1, β)]. At fixed z, the system of equations (9,10) is just the crossing of two lines in the ρ plane: α coming from eq.( 9) and β coming from eq. [START_REF] Rezakhanlou | Hydrodynamic limit for attractive particle systems on Z d[END_REF]. In Fig. 1 we show by a simple geometrical argument that these lines cross inside D whenever z α + z β ≤ 1. So in conclusion the domain D z (α, β) is given by z α ∈ [0, min(1, α)], z β ∈ [0, min(1, β)] and z α + z β ≤ 1. The geometrical reasoning explained in Fig. 1 allows also to conclude that at fixed z α , ρ • is an increasing function of z β , while at fixed z β , ρ • is an increasing function of z α .

In Figs. [START_REF] Mallick | Shocks in the asymmetry exclusion model with an impurity[END_REF] we have reported on the left the physical domain D in the densities plane and on the right the corresponding domain D z in the z variables plane for the cases α, β > 1 (2a) , α, β < 1 and α + β > 1 (2b) and α + β < 1 (2c). Notice that in Fig. 2a the thick red segment on the ρ • axis is mapped to the point (z α = 1, z β = 0) and the thick blue segment on the ρ • axis is mapped to the point (z α = 0, z β = 1). In Fig. 2c: the overlap of the red and blue segments on the boundary ρ * = 0 is mapped to the point (z α = α, z β = β). This indicates that the mapping ρ → z can be singular on the boundary of the physical domain D, where at least one of the densities vanishes. Let's analyze the different possibilities and work out the portion of the boundary where the mapping is singular. ρ • → 0 In this case from eq.( 9) we deduce that z α = 0 while the two solutions of eq.( 10) are z β = 1, z β = ρ • β and we have to retain the smallest one. If β ≤ 1 then on the ρ

• = 0 1 1 ρ • ρ • β z β β (1 -z β , z β ) α zα α (zα, 1 -zα)
C α of coordinates ρ • = z α , ρ • = 1 -z α . Similarly given z β ∈ [0, min(1, β)],
the line β (blue) intersects the boundary of D at the left side and at the diagonal side at the point

C β of coordinates ρ • = 1 -z β , ρ • = z β .
For the two lines to cross we need the point C α to lie above the point C β , which is the case iff

z α + z β ≤ 1.
axis the map is 1-to-1 and there are no singularities, on the other hand, if β > 1 then all the points ρ • ≥ β -1 are mapped to the same point (z α = 0, z β = 1).

ρ • → 0 This case is treated similarly to the previous one: we have z β = 0 and the two solutions of eq.( 9) are z α = 1, z α = ρ • α. If α ≤ 1 then on the ρ • = 0 axis the map is 1-to-1 and there are no singularities, on the other hand, if α > 1 then all the points ρ • ≥ α -1 are mapped to the same point (z α = 1, z β = 0).

ρ * → 0 Eqs.(9, 10) have solutions, z α = α, z α = ρ • , z β = β, z β = 1 -ρ • . Whenever α + β ≤ 1, all the points on the line ρ • + ρ • = 1 such that α ≤ ρ • ≤ 1 -β are mapped to the single point z α = α, z β = β.

Behavior at the boundary of the physical domain

The singularities of the mapping ρ → z reflect some important features of the model. Let's consider the currents of non zero density particles at the boundary of D ρ • → 0 We have for the current

J • (ρ • ) = βρ • (1 -ρ • ) 0 ≤ ρ • ≤ β -1 (1 -ρ • ) β -1 ≤ ρ • ≤ 1. ( 14
)
ρ • → 0 We have for the current ρ * → 0 We have for the current

J • (ρ • ) = -αρ • (1 -ρ • ) 0 ≤ ρ • ≤ α -1 -(1 -ρ • ) α -1 ≤ ρ • ≤ 1. ( 15
) ρ• ρ• 1 β 1 α 1 1 zα z β 1 1 (a) α, β > 1 ρ• ρ• (1 -β, β) (α, 1 -α) 1 1 zα z β α β 1 1 (b) α, β < 1, α + β > 1 ρ• ρ• 1 1 zα z β α β 1 1 (c) α + β < 1
J • (ρ • ) = ρ • (1 -ρ • ) 0 ≤ ρ • ≤ β, 1 -α ≤ ρ • ≤ 1 β(1 -α) + (α -β)ρ • β ≤ ρ • ≤ 1 -α. ( 16
)
These results have to be compared to the situation in which we have strict absence of a species of particles (and not just vanishing density). Consider for example a system without •-particles. Such a system is effectively a single species TASEP with jump rates equal to β and hence with current just

J • (ρ • ) = βρ • (1 -ρ • ).
Comparing this with eq.( 14) we see that for β > 1 this behavior holds only for 0 ρ • > β -1 the presence of even a single •-particle affects the macroscopic behavior of the system, giving rise to a modified current. At the boundary of D, the average speed of the zero density particles displays also an interesting behavior.

≤ ρ • ≤ β -1 , while for ρ• J•(ρ•) β 1 β -1 (a) ρ• - J•(ρ•) α 1 α -1 (b) ρ• J•(ρ•) 1 β 1 -α (c)
ρ • → 0 Speed of •-particles v • (ρ • ) = -α+β(1-α)ρ•(1-ρ•) 1+(α-1)ρ• -βρ • 0 ≤ ρ • ≤ β -1 -1 β -1 ≤ ρ • ≤ 1, (17)
ρ • → 0 Speed of •-particles v • (ρ • ) = β+α(1-β)ρ•(1-ρ•) 1+(β-1)ρ• -αρ • 0 ≤ ρ • ≤ α -1 1 α -1 ≤ ρ • ≤ 1, ( 18
)
ρ * → 0 Speed of * -particles v * (ρ • ) = z α -z β (19)
with

z α = α ρ • ≤ 1 -α 1 -ρ • ρ • ≥ 1 -α , z β = β ρ • ≥ β ρ • ρ • ≤ β (20)
The results in eqs. [START_REF] Cantini | Algebraic Bethe ansatz for the two species ASEP with different hopping rates[END_REF][START_REF] Cantini | Inhomogenous Multispecies TASEP on a ring with spectral parameters[END_REF][START_REF] Derrida | Bethe ansatz solution for a defect particle in the asymmetric exclusion process[END_REF] have been obtained in the literature by considering systems with a single particle of either species: the cases ρ * → 0, eqs. [START_REF] Rajewsky | Spatial particle condensation for an exclusion process on a ring[END_REF][START_REF] Derrida | Bethe ansatz solution for a defect particle in the asymmetric exclusion process[END_REF] first appeared in [START_REF] Mallick | Shocks in the asymmetry exclusion model with an impurity[END_REF][START_REF] Derrida | Bethe ansatz solution for a defect particle in the asymmetric exclusion process[END_REF], the cases ρ • → 0 or ρ • → 0, eqs. [START_REF] Arndt | Spontaneous breaking of translational invariance in one-dimensional stationary states on a ring[END_REF][START_REF] Cantini | Algebraic Bethe ansatz for the two species ASEP with different hopping rates[END_REF] and eqs. [START_REF] Arndt | Spontaneous breaking of translational invariance and spatial condensation in stationary states on a ring. I. the neutral system[END_REF][START_REF] Cantini | Inhomogenous Multispecies TASEP on a ring with spectral parameters[END_REF] first appeared in [START_REF] Lee | Two-way traffic flow: Exactly solvable model of traffic jam[END_REF].

Conservations laws

Under Euler-scaling (where site position and time scale as -1 n, -1 t for → 0) the density profiles are expected to evolve deterministically as solutions of a system of conservation laws. Consider initial data ρ (0) i (x) and to such data associate a family of initial conditions of the 2-TASEP of product Bernoulli form, with local probability at site n given by

E[χ i (n, t = 0)] = ρ (0) i ( n),
where χ i (n, t) is the i-th species indicator function at time t and site n. We expect that the random variable χ i ( -1 x , -1 t, ) converges for → 0 to a deterministic density profile. More precisely we expect that lim -→0 n:a≤ n≤b

χ i (n, -1 t) = b a ρ i (x, t)dx, a.s. ( 21
)
where ρ = (ρ • , ρ • ) is the solutions of a system of conservation laws

∂ t ρ + ∂ x J = 0. ( 22
)
with initial condition ρ(t = 0) = ρ (0) = (ρ

(0) • , ρ (0)
• ). By making the usual hypothesis of local stationarity we identify the local currents with the stationary currents at density ρ, given by eqs. [START_REF] Kipnis | Scaling limits of interacting particle systems[END_REF][START_REF] Derrida | Exact solution of a 1-D asymmetric exclusion model using a matrix formulation[END_REF]. A more precise statement and proof of this result for the case α = β = 1 2 can be found in [START_REF] Fritz | Derivation of the leroux system as the hydrodynamic limit of a twocomponent lattice gas[END_REF]. While the approach developed in [START_REF] Fritz | Derivation of the leroux system as the hydrodynamic limit of a twocomponent lattice gas[END_REF] can be extended to the full α+β = 1 line (for which the stationary measure is product), it is not expected to work for arbitrary values of α and β [33]. In the present paper we take eqs. [START_REF] Chen | Limiting current distribution for a two species asymmetric exclusion process[END_REF][START_REF] Van Beijeren | Exact results for anomalous transport in one-dimensional hamiltonian systems[END_REF] as a working hypothesis. Eqs. [START_REF] Van Beijeren | Exact results for anomalous transport in one-dimensional hamiltonian systems[END_REF] form a system of coupled conservation laws, whose non-linearity is known to be at the origin of characteristic phenomena such as shocks formation in finite time, and rarefaction waves, i.e. self-similar solutions, which present regions expanding in time at constant speed where the densities interpolate between two boundary values. In the following we shall analyze in detail eqs. [START_REF] Van Beijeren | Exact results for anomalous transport in one-dimensional hamiltonian systems[END_REF]. We shall show that the variables z are Riemann variables for this system. On general grounds we know that the rarefaction fans can be expressed in an implicit form involving the Riemann variables. What is more surprising is that for our system also the shock solutions are explicitly written in terms of the Riemann variables: they correspond to a discontinuity of only one of the two Riemann variables (the other one being continuous). Putting together fans and shock we can explicitly solve the Riemann problem. In the Section 4 these theoretical results are compared to Monte Carlo simulations of the 2-TASEP.

3.1. The cases α = β = 1 and α + β = 1 Before discussing the system of equation [START_REF] Van Beijeren | Exact results for anomalous transport in one-dimensional hamiltonian systems[END_REF] in full generality let us start with some comments on two particular cases: α = β = 1 and α + β = 1

• α = β = 1.
When β = 1, the •-particles don't distinguish •-particles from * -particles. This means that the •-particles evolve as in a single species TASEP. At the level of currents, for β = 1 we have indeed J • = ρ • (1-ρ • ). In this case the conservation law for ρ • completely decouples from that of ρ • and takes the usual form of the non-viscous Burgers equation

∂ t ρ • + (1 -2ρ • )∂ x ρ • = 0.
Analogously, for α = 1, the conservation law for ρ • completely decouples from that of ρ • and takes the form

∂ t ρ • -(1 -2ρ • )∂ x ρ • = 0.
So for α = β = 1, system [START_REF] Van Beijeren | Exact results for anomalous transport in one-dimensional hamiltonian systems[END_REF] just decouples completely into two Burgers equations.

• α + β = 1.
In this case, thanks to the factorization of the stationary measure, the currents can be explicitly written as functions of the densities

J • = -ρ • (ρ • + β(1 -ρ • -ρ • )), J • = ρ • (ρ • + α(1 -ρ • -ρ • )).
One can consider conserved quantities ρ and v, defined by

ρ v = -α(α+2β) 3 -β(α+2β) 3 α -β ρ • ρ • + (2α+β)(2β+α) 9 β-α 3 (23)
The associated currents are (up to irrelevant additive constants)

J ρ = ρv, J v = ρ + v 2 .
These are the currents of the Leroux system [START_REF] Leroux | Analyse et approximation de problèmes hyperboliques non linéaires[END_REF][START_REF] Serre | Existence globale de solutions faibles sous une hypothèse unilaterale pour un système hyperbolique non linéaire[END_REF] (the particular case α = β = 1 2 is the one considered in [START_REF] Fritz | Derivation of the leroux system as the hydrodynamic limit of a twocomponent lattice gas[END_REF]), which is known to be a Temple class system.

The general case: Riemann variables

For the problem under investigation one could expect that in addition to the generic complexity of the analysis of a coupled system of conservation equations, one has to face the further complication due to the implicit dependence of the currents on the densities, which goes through the auxiliary variables z. Actually, quite unexpectedly, what seems a drawback of the equations turns out to be the main feature which allows to solve them. Indeed the variables z happen to be Riemann variables for our conservation laws, i.e. they diagonalize the system of eqs. [START_REF] Van Beijeren | Exact results for anomalous transport in one-dimensional hamiltonian systems[END_REF] and simplify substantially their analysis. From now on we want to think of both ρ and J as functions of z .

We first notice that, solving eqs. [START_REF] Kipnis | Scaling limits of interacting particle systems[END_REF][START_REF] Derrida | Exact solution of a 1-D asymmetric exclusion model using a matrix formulation[END_REF] for the densities ρ in terms of z and J and replacing them into eqs. [START_REF] Benassi | Hydrodynamical limit for the asymmetric simple exclusion process[END_REF][START_REF] Rezakhanlou | Hydrodynamic limit for attractive particle systems on Z d[END_REF], the currents are the solution of the following linear system of equations

J • z α + J • z α -1 - J • + J • z α -α + 1 = 0 ( 24
)
J • z β + J • z β -1 - J • + J • z β -β + 1 = 0. ( 25
)
Now differentiate the l.h.s. of eq.( 9) with respect to t, differentiate the l.h.s. of eq.( 24) with respect to x and sum the obtained results. Thanks to the conservation laws eqs. [START_REF] Van Beijeren | Exact results for anomalous transport in one-dimensional hamiltonian systems[END_REF] the derivatives ∂ t ρ and ∂ x J cancel and one remains with

∂ t z α + v α (z)∂ x z α = 0, v α (z) = J• z 2 α + J• (zα-1) 2 -J•+J• (zα-α) 2 ρ• z 2 α + ρ• (zα-1) 2 + 1-ρ•-ρ• (zα-α) 2 . ( 26
)
In the same way one obtains the equation for z β

∂ t z β + v β (z)∂ x z β = 0, v β (z) = J• z 2 β + J• (z β -1) 2 -J•+J• (z β -β) 2 ρ• z 2 β + ρ• (z β -1) 2 + 1-ρ•-ρ• (z β -β) 2 . ( 27
)
The speeds v α and v β are the eigenvalues of the linearization matrix ∂ ρ j J i , and on general grounds they can also be written as

v α = ∂ ρ i J i | z β = ∂ zα J i (z) ∂ zα ρ i (z) , v β = ∂ ρ i J i | zα = ∂ z β J i (z) ∂ z β ρ i (z) . (28)
A close inspection of their expression allows to conclude that

v β (z) ≥ v α (z), (29)
with the equality holding for 1-z α -z β = 0, which is a non empty set only for α+β ≥ 1.

We conclude that for α + β < 1, the system in ( 22) is strictly hyperbolic on the whole physical domain D, whereas for for α + β ≥ 1 it is degenerate hyperbolic on the locus 1 -z α -z β = 0, i.e. ρ * = 0, ρ • ≤ α, ρ • ≤ β (green segments in Figs. 2a,2b) and strictly hyperbolic on the rest of the physical domain. Using eqs. [START_REF] Tóth | Perturbation of singular equilibria of hyperbolic two-component systems: a universal hydrodynamic limit[END_REF][START_REF] Leroux | Analyse et approximation de problèmes hyperboliques non linéaires[END_REF] we can easily work out the rarefaction fans. These are continuous solutions of eqs. [START_REF] Van Beijeren | Exact results for anomalous transport in one-dimensional hamiltonian systems[END_REF], which depend only on the self-similarity variable ξ = x/t, z(x, t) = z(ξ = x/t). They are given by the solutions of the equations

(v α -ξ)∂ ξ z α = 0 (v β -ξ)∂ ξ z β = 0. ( 30
)
Locally we have four possibilities.

(i) The trivial solution, namely both z α and z β are constant.

(ii) Both ∂ ξ z α = 0, ∂ ξ z β = 0. In this case we must have v α -ξ = v β -ξ = 0, and in particular v α = v β . As mentioned above, this is possible only if 1 -z α -z β = 0. In this case we get a rarefaction fan of equation

ρ • (ξ) = z α (ξ) = 1 + ξ 2 , ρ • (ξ) = z • (ξ) = 1 -ξ 2 . ( 31
)
Notice that the condition 1 -z α -z β = 0 implies absence of * -particles, and the solution (31) corresponds to the fan solution of the single species TASEP (upon identification of •-particles with empty sites).

(iii) z α constant, ∂ ξ z β = 0. In this case the rarefaction fan (β-fan) is given by

v β (z α , z β (ξ, z α )) = ξ. (32)
This equation, combined with the expression ρ(z α , z β ), allows to write a β-fan in parametrized form (z α is kept constant while z β varies), see Fig. 4b. One can show that (at fixed z α ) v β is a decreasing function of z β , while ρ • is an increasing function of z β . This means that z β (ξ, z α ) and ρ • (ξ, z α ) are decreasing functions of ξ.

(iv) z β constant, ∂ ξ z α = 0. This case is similar to the previous one, but the role of α and β variables is exchanged. So we speak of an α-fan, given by

v α (z α (ξ, z β ), z β ) = ξ. (33)
One can show that (at fixed z β ) v α and ρ • are increasing function of z α . This means that z α (ξ, z β ) and ρ • (ξ, z β ) are increasing functions of ξ.

Projection in the ρ-plane of the three types of fans as well as an example of a β-fan are represented in Fig. 4.

ρ• ρ• 1 1 (a) v - β v + β ρ• ρ * ρ• ρ ξ (b)
Figure 4: (a) Projection in the ρ-plane of the three types of rarefactions fans in the case α, β < 1: in green a TASEP-like rarefaction fan, in blue an α-fan, in red a β-fan. (b) A plot of a β-fan. The local densities correspond to the width of the corresponding colored regions. In particular the line separating the yellow region from the violet region represents the plot of the density ρ • (ξ). The extremes of the fan are given by :

v - β = v β (z α , 1 -z α ), v + β = v β (z α , 0).

Shocks

It is well known that a smooth solution of a nonlinear conservation laws like eqs.( 22) may develop a shock discontinuity in finite time. One needs therefore to admit the notion of weak solution, i.e. solution in the sense of distributions, which need not even be continuous. The discontinuity associated to a shock with trajectory x s (t), has to satisfy the Rankine-Hugoniot jump relations . If we denote by [ρ] and [J] respectively the discontinuities of the densities and of the currents across the shock, i.e.

[ρ] = ρ(x

+ s ) -ρ(x - s ), [J] = J(x + s ) -J(x - s )
then the Rankine-Hugoniot jump relations read

[ρ] -v s [J] = 0, v s = shock's speed. ( 34
)
The Rankine-Hugoniot jump relations, allow to express the speed of the shock in terms of the discontinuity and at the same time put a constraint on the admissible discontinuities, the Hugoniot condition:

det [ρ • ] [J • ] [ρ • ] [J • ] = 0. ( 35
)
In order to analyze the Hugoniot condition, we consider ρ(x - s ) = ρ -as fixed. For strictly hyperbolic systems of N conservation laws, it is known that the set of ρ + = ρ(x + s ) satisfying the Hugoniot condition passes through ρ -and locally decompose around ρ - in N different branches [START_REF] Lefloch | Hyperbolic Systems of Conservation Laws: The theory of classical and nonclassical shock waves[END_REF], each one called a shock curve. In our case we expect two shock curves passing through any point ρ -, which correspond to two different kinds of shocks. Using eq.( 9) and eq.( 24) we see that if

z + α = z - α = z α , then we have [ρ • ] z α + [ρ • ] z α -1 - [ρ • ] + [ρ • ] z α -α = 0 [J • ] z α + [J • ] z α -1 - [J • ] + [J • ] z α -α = 0.
This means that the matrix in the eq.( 35) has the left null vector ( 1 zα -1 zα-α , 1 zα-1 -1 zα-α ) and therefore its determinant vanishes. We conclude that a straight line at z α constant is a shock curves. In the same way we find that the lines at z β constant are also shock curves. In conclusion we have found that we have two kind of admissible shocks

• β-shocks: z α (ρ -) = z α (ρ + ) = z α with speed v s,β (z α ; z - β , z + β ); • α-shocks: z β (ρ -) = z β (ρ + ) = z β , with speed v s,α (z β ; z - α , z + α ). The corresponding shocks speeds v s,β (z α ; z - β , z + β ) and v s,α (z β ; z - α , z + α )
do not have particularly transparent expressions except for some particular cases. In the case α = β = 1, a β-shock is a discontinuity of the density ρ • , with ρ • constant across the discontinuity, while an α-shock is the other way round, i.e. a discontinuity of the density ρ • , with ρ • constant across the discontinuity. Shock curves coincide with rarefaction curves so this is an example of a conservation system of Temple class [START_REF] Temple | Systems of conservation laws with invariant submanifolds[END_REF].

A further analysis of the Hugoniot condition allows to conclude that in the bulk of the physical domain D, the only possible shocks are α and β-shocks. There exists however one more class of shocks when both sides of the discontinuity lie on the boundary line ρ * = 0. These shocks have speed

v s = [J•] [ρ•]
, where the current J • is given by eq.( 16). It is possible to show that at fixed z α the current J • is a concave function of the density ρ • (see fig. 5 ). This implies that, for a fixed value of the densities on one side of the shock, say ρ -, the speed of a β-shock is a decreasing function of ρ +

• . This property can be conveniently reformulated in terms of the z variables. At constant z α , ρ • is an increasing function of z β , hence at fixed z α and z - β , the speed of β-shock v s,β (z α ; z - β , z + β ) is a decreasing function of z + β . In particular it takes its minimum for the largest z β allowed, i.e.

z β,max = 1 -z α v s,β (z α ; z β , z + β ) ≥ v s,β (z α ; z β , 1 -z α ) 0 ≤ z + β ≤ 1 -z α . (36)
The current J • is a convex function of ρ • at fixed z β . A similar reasoning as the one presented above allows to conclude

v s,α (z β ; z α , z + α ) ≤ v s,α (z β ; z α , 1 -z β ) 0 ≤ z + α ≤ 1 -z β . (37)
Now, from an explicit computation, we notice that

v s,α (z β ; z α , 1 -z β ) = v s,β (z α ; z β , 1 -z α ) = z α -z β .
This allows to conclude that

v s,α (z β ; z α , z α ) ≤ v s,β (z α ; z β , z β ). (38)
In words this means that, for a fixed value of the densities on one side of the shock, the speed of any β-shock is larger than the speed of any α-shock.

Since lim zα→zα v s,α (z β ; z α , z α ) = v α (z α , z β ) and lim z β →z β v s,β (z α ; z β , z β ) = v β (z α , z β ), we get also v α (z α , z β ) ≤ v s,β (z α ; z β , z β ) (39) v β (z α , z β ) ≥ v s,α (z β ; z α , z α ). ( 40
)
The Hugoniot condition is not sufficient to select the physical shocks. Indeed eqs. [START_REF] Van Beijeren | Exact results for anomalous transport in one-dimensional hamiltonian systems[END_REF] has to be thought of as the zero viscosity limit of a set of conservation laws which contains a diffusive term, a term which comes from the microscopic corrections to the currents and depends on the derivatives of the densities. Inviscid limits of viscous solutions are typically characterized by entropy conditions, which for shocks take the form of the Liu entropy criterion [START_REF] Lefloch | Hyperbolic Systems of Conservation Laws: The theory of classical and nonclassical shock waves[END_REF][START_REF] Liu | The Riemann problem for general 2× 2 conservation laws[END_REF]. Let the right densities ρ + lie on a Hugoniot curve emanating from ρ -, then for all densities ρ lying on the same Hugoniot curve in between ρ -and ρ + the Liu condition states that

v s (ρ -, ρ + ) ≤ v s (ρ -, ρ). (41)
Physically this condition can be understood as a stability condition: if under a perturbation the shock were to split by inserting an intermediate state ρ, then a violation of condition (41) would imply that the shock between ρ -and ρ would move away from the original shock between ρ -and ρ + . In the case of β-shocks the Liu condition reads

v s,β (z α ; z - β , z + β ) ≤ v s,β (z α ; z - β , z β ), ∀ z β ∈ [min(z - β , z + β ), max(z - β , z + β )] (42)
Since, as mentioned above, at fixed z α the current J β is a concave function of the density

ρ • , ∂ 2 ρ• J β | zα < 0, we conclude that the Liu constraint means ρ - • < ρ + • .
This can be also formulated in terms of the z variables. Indeed, since ∂ z β ρ • | zα > 0, we must have z - β < z + β . A similar analysis can be performed for α-shocks.

Here is a summary of our results. They are schematized in the figure on the right, where we have reported the direction of the possible shock-discontinuities.

• For an α-shocks (blue oriented line), we need ρ -

• > ρ + • or equivalently z - α > z + α . • For a β-shocks (red oriented line), we need ρ - • < ρ + • or equivalently z - β < z + β . ρ• ρ• 1 1

Riemann's problem

With the result for the rarefaction curves and the shock curves at our disposal, it is rather simple to describe the general solution of the Riemann problem, i.e. the solution of eqs. [START_REF] Van Beijeren | Exact results for anomalous transport in one-dimensional hamiltonian systems[END_REF] with domain wall initial conditions

ρ(x, t = 0) = ρ L = (ρ L • , ρ L • ) x < 0 ρ R = (ρ R • , ρ R • ) x > 0. ( 43
)
By uniqueness, the solution of the Riemann problem has to take the form ρ(x, t) = ρ(ξ) with ξ = x/t and is given by a sequence of rarefaction waves and/or shocks. It is best described in terms of the variables

z L = (z L α , z L β ), z R = (z R α , z R β ).
We have four possible situations (these are schematically summarized in figure 6).

• z L α > z R α , z L β < z R β .
The solution is composed of two shocks: an α-shock with z -= (z L α , z L β ) and

z + = (z R α , z L β ) at position ξ α = v s,α (z L β ; z L α , z R α ), followed by a β- shock with z -= (z R α , z L β ) and z + = (z R α , z R β ) at position ξ β = v s,β (z R α ; z L β , z R β )
. This result follows from the inequality (38)

ξ α = v s,α (z L β ; z L α , z R α ) ≤ v s,β (z R α ; z L β , z R β ) = ξ β . • z L α > z R α , z L β > z R β .
The solution is composed of an α-shock and a β-fan. The α-shock has z -= (z L α , z L β ) and z + = (z R α , z L β ) and it is located at position

ξ α = v s,α (z L β ; z L α , z R α ). The β-fan starts with value (z R α , z L β ) at ξ β,1 = v β (z R α , z L β ) and ends with value (z R α , z R β ) at ξ β,2 = v β (z R α , z R β )
. This result follows from the inequality

ξ α = v s,α (z L β ; z L α , z R α ) ≤ v β (z R α , z L β ) = ξ β,1 ≤ v β (z R α , z R β ) = ξ β,2 .
The first inequality is just inequality (40), while the second one follows from the fact that v β is a decreasing function of z β .

• z L α < z R α , z L β < z R β .
The solution is composed of an α-fan and a β-shock. The α-fan starts with value

(z L α , z L β ) at ξ α,1 = v α (z L α , z L β ) and ends with value (z R α , z L β ) at ξ α,2 = v α (z R α , z L β ). The β shock has z -= (z R α , z L β ) and z + = (z R α , z R β ) and is located at position ξ β = v s,β (z R α ; z L β , z R β ). • z L α < z R α , z L β > z R β .
The solution is composed of fans. One has to distinguish two cases depending whether z R α + z L β is larger or smaller than 1. If z R α + z L β < 1, the solution consists of an α-fan starting with value (z L α , z L β ) at ξ α,1 = v α (z L α , z L β ) and ending with value (z R α , z L β ) at ξ α,2 = v α (z R α , z L β ) followed by a β-fan starting at

ξ β,1 = v β (z R α , z L β ) and ending at ξ β,2 = v β (z R α , z R β ). If α + β > 1 then it is possible to have z R α + z L β > 1.
In this case the α-fan cannot reach the value (z R α , z L β ), which lies outside the physical domain. It ends at the value (1 -z L β , z L β ), followed by a degenerate fan till (z R α , 1 -z R α ) and then by a β-fan till (z R α , z R β ).

Monte Carlo simulations

We come back to the original microscopic stochastic model and compare the predictions of the hydrodynamic equations with numerical simulations. We have simulated our model on a finite lattice of integer coordinates, running on the interval [-L, L], with L = 2100. The system is initialized in a random configuration sampled from a product measure of local densities ρ L on sites of coordinate i < 0 and ρ R on sites of coordinates i ≥ 0. At the left and right boundaries the particles are chosen neither to leave nor to enter the system. This means that we expect to find three distinct regions: two kinetic waves coming from the boundaries and a kinetic wave originating from the discontinuity at the origin. Whereas we make no prediction on the boundary waves, we expect that as long as they don't meet the bulk one, they do not influence the latter. Let us introduce the height function, which is defined up to an arbitrary additive constant by

h i n t , t -h i (-1, t) := 1 t -t<k≤n χ i (k, t) -t < n < t, (44)
From our assumption eq.( 21), it follows that at large time and for each sample, the height function should converge to the deterministic shape

h i (ξ) -h i (-1) = ξ -1 ρ i (ξ )dξ , (45)
where ρ i (ξ) is the solution of the Riemann problem found in Section 3.4. In Figs. 7,8 we illustrate some representative results of our simulations. We have run each simulation up to a time t indicated at the top of each plot. At that time we find that the left boundary wave has not yet reached the site i = -t and the right boundary wave has not reached the site i = t, so we can safely use eq.( 44) as a definition of the height function and plot it for -1 ≤ ξ ≤ 1 (for convenience we chose h i (-1) = 1). First we illustrate the simulation of the system for α = 0.3, β = 0.8 (Figs. 7a-7c). In Fig. 7a the left densities are (ρ L • = 0.45, ρ L • = 0.11) and the right densities are (ρ R • = 0.23, ρ R • = 0.67), which in terms of the z variables correspond to a uniform z L α = z R α = 0.15 and z L 1 = 0.1 < z R 1 = 0.6. This choice of boundary densities is expected to lead to the formation of a single β-shock. Indeed the plot of the height functions for the •-particles and •-particles shows two linear regions, corresponding to the two regions of constant densities, separated at the predicted shock speed, situated at ξ = v s,β = 0.217. In Fig. 7b the boundary densities are inverted w.r.t. Fig. 7a. In this case we don't expect any shock to persist, and indeed the result is compatible with a single fan. In Fig. 7c we explore the separation of two shocks. The boundary densities are (ρ L

= 0.3 = 0.8 t =1000 h num h num h h (c) (ρ L • = 0.8, ρ L • = 0.1), (ρ R • = 0.2, ρ R • = 0.7)
= 0.7 = 1.2 t =1100 h num h num h h (d) (ρ L • = 0.2, ρ L • = 0.1), (ρ R • = 0.5, ρ R • = 0.4)
• = 0.8, ρ L • = 0.1), (ρ R • = 0.2, ρ R • = 0.7), which correspond to z L α = 0.265 > z R α = 0.
139 and z L β = 0.097 < z R β = 0.627. In this case the solution of the Riemann problem predicts an α-shock of speed v s,α = 0.008 followed by a β-shock of speed v s,β = 0.186. The β shock is clearly visible at the bend of h • (green line), the α shock is less visible because of the small bends in both height functions. In Fig. 7d, we have α = 0.7, β = 1.2 and we show the result of the simulation where the theoretical analysis predicts an α-fan and a β-shock, the shock being visible at the bend of h • (green line), located at ξ = v s,β = 0.49. Finally in Fig. 8 we explore the formation of a TASEP-like fan corresponding to a region where the density of * -particles vanishes. We chose α = 1.7, β = 2, (ρ L • = 0.1, ρ L • = 0.5, ρ L * = 0.4), (ρ R • = 0.5, ρ R • = 0.05, ρ R * = 0.45) so that initially the * -particles are present everywhere in the system. The plot of h * (purple line) shows a central flat region, which corresponds to a region of vanishing * -particles density.

Conclusion

In this article we have investigated the hydrodynamic behavior of an exactly solvable two species exclusion process, consisting of two kinds of particles moving in opposite directions on a one dimensional lattice and swapping their position when adjacent. By making the assumption of local stationarity and exploiting the knowledge of the particle currents on periodic geometry, we have written down the hydrodynamic conservation laws of the model and investigated their solutions: rarefaction fans, shocks and the solution of the Riemann problem. Then such predictions have been shown to be in agreement with numerical simulations of the microscopic model. The macroscopic conservation laws are shown to belong to a class of conservation laws called Temple systems, which possess coinciding rarefaction and shock curves.

While this property has important implications for the mathematical analysis of the conservation equations, it is not clear to us what is its physical underpinning. We believe it would be interesting to investigate how generic this property is among exclusion processes, and whether it is somehow related to the integrability of the underlying microscopic model. On this line of thoughts, natural candidates to investigate seem to be the multispecies integrable generalizations (with more than 2 species) introduced in [START_REF] Cantini | Inhomogenous Multispecies TASEP on a ring with spectral parameters[END_REF].

Another interesting issue would be to understand the behavior of our model when restricted to a finite lattice in contact with boundary particle reservoirs. On general grounds, novel feautures are expected to emerge like for example the phenomenon of boundary induced phase transitions [START_REF] Krug | Boundary-induced phase transitions in driven diffusive systems[END_REF]. For the case of a single conserved quantity the well-known max-min principle of Krug [START_REF] Krug | Boundary-induced phase transitions in driven diffusive systems[END_REF], later generalized by Popkov and Schütz [START_REF] Popkov | Steady-state selection in driven diffusive systems with open boundaries[END_REF] allows to determine the dependence of the current on the boundary couplings and hence the phase diagram of the model. However, it is not known how to generalize the maxmin principle in the case of more than one conserved species. Preliminary numerical investigations based on the model presented here seems to indicate that the Riemann problem may play a relevant role in order to tackle this problem [38].

  ) ‡ Here we are supposing a, b, c, γ > 0.

Figure 1 :

 1 Figure 1: The shaded triangle is the physical region D. Given z α ∈ [0, min(1, α)], the line α (red) intersects the boundary of D at the bottom side and at the diagonal side at the point C α of coordinates ρ • = z α , ρ • = 1 -z α . Similarly given z β ∈ [0, min(1, β)], the line β (blue) intersects the boundary of D at the left side and at the diagonal side at the point C β of coordinates ρ • = 1 -z β , ρ • = z β . For the two lines to cross we need the point C α to lie above the point C β , which is the case iff z α + z β ≤ 1.

Figure 2 :

 2 Figure 2: On the left the physical domain D in the densities plane, on the right the corresponding domain D z in the z variables plane. On both sides we have drawn in red the lines α at constant z α and in blue the lines β at constant z β . The lines α (red) have slope α(1-zα) (α-1)zα and in particular they have positive slope for α > 1 as in figure (a) and negative slope for α < 1 as in figures (b) and (c). Similarly, the lines β (blue) have slope (β-1)z β β(1-z β ) , they have positive slope for β > 1 as in figure (a) and negative slope for β < 1 as in figures (b) and (c).

Figure 3 :

 3 Figure 3: (a) Current J • for β > 1 and ρ • = 0. (b) Current J • for α > 1 and ρ • = 0. (c) Current J • for α + β < 1 and ρ * = 0. The dashed lines correspond to the current for a strict absence of •-particles (a), •-particles (b), * -particles (c).

Figure 5 :

 5 Figure 5: Current J • as function of ρ • at constant z α = .2. (a) Fixed α = 4 and decreasing β (blue line β = 5, green line β = 1, red line β = 0.3). (b) Fixed β = 0.3 and decreasing α (blue line α = 3, green line α = 1, red line α = 0.3).

Figure 6 :

 6 Figure 6: Projection in the z-plane of the different type of solutions of the Riemann problem. The left values of the z variables are represented by the point L. The points R represent the different distinct possibilities for the right values of the z variables. Continuous lines represent fans, while dashed lines represent shocks. In green is indicated the possible TASEP-like fan, in red either an α-shock or an α-fan and in blue either a β-shock or a β-fan.

Figure 7 :

 7 Figure 7: Height function of the •-particles (green lines) and of the •-particles (red lines). Dashed lines represent the numerical values, continuous lines represent the theoretical prediction. At the bottom of each plot we have reported the predicted densities.

1 .Figure 8 :

 18 Figure 8: Situation diplaying an α and a β-fan with a depletion region in the middle, (ρ L • = 0.1, ρ L • = 0.5, ρ L * = 0.4), (ρ R • = 0.5ρ R • = 0.05, ρ R * = 0.45). On the left: height function of •-particles (green lines) and of * -particles. Continuous lines represent the numerical values, dashed lines represent the theoretical prediction. At the bottom we have reported the predicted densities. Notice the flat central region of h * corresponding to a region of vanishing ρ * . On the right: the projection of this configuration in the densities plane: in blue the β fan, in green the TASEP-like fan with vanishing ρ * and in red the α fan.