Solutions of the Aw-Rascle-Zhang system with point constraints

We revisit the entropy formulation and the wave-front tracking construction of physically admissible solutions of the Aw-Rascle and Zhang (ARZ) “second-order” model for vehicular trafﬁc. A Kruzhkov-like family of entropies is introduced to select the admissible shocks. This tool allows to deﬁne rigorously the appropriate notion of admissible weak solution and to approximate the solutions of the ARZ model with point constraint. Stability of solutions w.r.t. strong convergence is justiﬁed. We propose a ﬁnite volumes numerical scheme for the constrained ARZ, and we show that it can correctly locate contact discontinuities and take the constraint into account.


Introduction
Any macroscopic vehicular traffic model expresses the conservation of the number of vehicles along a highway without entrances or exists with the PDE ρ t + (ρ v) x = 0, where ρ is the density and v is the velocity of the vehicles.In addition to this PDE, the Lighthill-Whitham and Richards (LWR) model, [26,28], assumes that v = V (ρ).Recently, a particular attention was paid to incorporating point constraints of the kind (ρ v)| x=x i ≤ q i (t ) into the LWR model.Here, x i are the locations of the road where a traffic light, a sharply localized construction site or a toll gate, etc., are situated; and for each x i , q i (t ) is the maximal value of the flux q .= ρ v allowed through the point x i at time t .In this case, non-classical (w.r.t. the classical Kruzhkov theory [25]) stationary shocks may appear at x = x i , however, the problem is well-posed both in the BV and the L ∞ frameworks.As suggested in [13], the unilateral constraint is taken into account via an additional singular "compensation term" in the Kruzhkov entropy inequalities supported at x = x i .
We refer to [5,13,15,18,29,30] for the introduction, the fundamental analytic results and a simple approximation strategy for the locally constrained LWR model.Natural extensions of this model are proposed in [2,3,11,14,16,17,31] to reproduce the phenomenon called capacity drop at the exit in pedestrian flows.
These works focus on the description of the crowd behaviour when high density situations occur in a narrow corridor.In particular, in [3] non-local constraints of the form F i (t ) = Q i [ρ(t , •)] are considered.Let us stress that the typical paradoxical features related to capacity drop (Faster is Slower and Braess paradox) can be reproduced within this very simple model, see [1].

The objective of the paper
Although LWR can be used for rough modeling of vehicular or pedestrian traffic, it presents several nonrealistic features.Some of these issues are fixed in the "second-order" macroscopic model of Aw-Rascle and Zhang (ARZ) for vehicular traffic, [7,32], which writes as the hyperbolic 2 × 2 system of conservation laws where p : R + → R + is a function whose properties will be specified later on (typically, p(ρ) = ρ γ , γ > 0).The additional unknown w is a phenomenological quantity interpreted as a Lagrangian marker (indeed, ( implies that formally w t + v w x = 0, i.e. w is transported at the velocity v of the vehicles).It is of interest to understand how the point constraints can be incorporated into ARZ.A first step in this direction was made in [22], where the constrained Riemann problem was described.The goal of the present note is to summarize some recent analytical and numerical results obtained by the authors for the locally constrained ARZ subject to point constraint on the flux (ρ v)| x=x 0 ≤ q 0 (t ), where q 0 is a constant or piecewise constant function of time.In addition to the interest of the result in itself, this is a first step towards the analysis and approximation of ARZ-based models with capacity drop.
The paper is organized as follows: In Section 2 we recall that ARZ has a family of entropies that are conserved by the evolution; we make the latter statement precise for all weak solutions of ARZ by applying Panov's theory [27] of renormalization for scalar continuity equations.We also describe how one can select the admissible shocks and forbid the non-admissible ones by using a specific family of ρ-convex, Kruzhkov-like entropies (E k ) k>0 .We recall that even if contact discontinuities adjacent to vacuum states can be mathematically correct solutions to the ARZ system, they must be disregarded as they do not have physical meaning.Of course, entropy criteria cannot rule out these non-admissible discontinuities, as entropies do not dissipate along discontinuities of the linearly degenerate family.In Section 3 we rely on Riemann invariants to describe a BV-stable Riemann solver able to cope with vacuum states, [6,23].Sections 4 and 5 deal with the constrained ARZ.We give a rigorous definition of admissible solution to ARZ with a local constraint in terms of entropy inequalities with the entropies (E k ) k>0 and suitable singular "compensation terms".Using the aforementioned renormalization property, we prove that this definition is stable w.r.t. the L 1 loc convergence of sequences of solutions.Then we describe the associated Riemann solver and its discretized version suitable for the Wave-Front Tracking (WFT) analysis.Section 6 is devoted to numerical approximation of constrained ARZ.We describe a finite volume numerical scheme for admissible solutions of ARZ, incorporating the point constraint by a simple flux-limiting procedure in the spirit of [5].In the last section we compare an exact solution with the simulation obtained by our finite volume scheme from the same initial condition.This example constitutes an important part of the validation of the numerical scheme.

Notation and structure of ARZ
We assume that p : R + → R + is a continuous function, twice continuously differentiable on ]0, +∞[, such that p(0) = 0, For technical reasons we also assume All these assumptions are satisfied, in particular, for p(ρ) .
Throughout the paper, we choose to represent the unknown states in Recall that the characteristic speed corresponding to the waves of the first (resp., of the second) family is Away from the vacuum the first family is genuinely nonlinear, while the second is linearly degenerate; this motivates the distinction in the terminology used for jumps in the solution (jumps of the first family are called shocks, while jumps of the second family are called contact discontinuities).The system is strictly hyperbolic and a Temple system away from the vacuum; at the vacuum, it is natural to extend the flux by zero, so that the two eigenvalues coincide: λ 1 = 0 = λ 2 , yet the matrix of the system still admits a full basis of eigenvectors also in vacuum.

Entropies and renormalization property for ARZ revisited
In this section, we recall what is known on entropies for ARZ.There are no convex entropies for ARZ when vacuum is allowed, see [23].However, for our result it is sufficient to consider entropies that are preserved on smooth solutions of ARZ and dissipated (in the non-strict sense) by the solutions of the Riemann problems as constructed by Aw and Rascle [7].Indeed, on the one hand, the linearity in w of the second equation of (1.1) implies that fully general functions E that depend solely on the variable w are entropies for the system (and these entropies are conserved along trajectories of ARZ).On the other hand, the very convenient in practice "Kruzhkov-like" entropies E k , introduced below, select the right shocks as they are convex in the variable ρ, while convexity in the conservative variables (ρ, y) is not needed.We also recall that viscous regularization (underlying the convexity requirement for the entropies) has no physical meaning in the context of ARZ.Indeed, it is incorrect to apply a viscous regularization to ARZ because vehicular traffic is anisotropic: unlike gas molecules that respond to stimulus from all sides, the drivers react only to the distance to the vehicle ahead and ignore the location of all other vehicles.
The entropy-entropy flux pairs (E , Q) for ARZ are pairs of real-valued functions of the state of the system T , such that the companion conservation laws E ( W ) t + Q( W ) x = 0 hold for all smooth solutions of (1.1).This leads to the compatibility conditions under the integrability condition, see [20,Eq. (7.4.13)], Accurate calculations show that general solutions of the above system are where f and g are arbitrary sufficiently regular functions, b and c are arbitrary constants and The smoothness assumption for f and g is relaxed by straightforward approximation arguments, indeed, it is enough to have E Lipschitz continuous.Actually, the class of f in the above formula can be even more general, due to where D denotes the space of distributions on ]0, +∞[ × R and B denotes the space of all Borel measurable functions.
Proof.The claim of the lemma is a straightforward application of the result of [27].Indeed, given an L ∞ weak solution W to ARZ, one can consider the L ∞ functions A : (t , x) → ρ( W (t , x)) and B : (t , x) → q( W (t , x)) ≡ ρ( W (t , x)) v(t , x).With this notation, the first equation of (1.1) means that the field (A, B ) in the coordinates (t , x) is divergence-free (in the sense of distributions), while the second equation means that w is a solution (in the sense of distributions) of the continuity equation Under these assumptions, Panov [27] proves that for every Borel function f Remark 1. Due to the result of Lemma 2.1, the terms in (2.4) involving the function f (and obviously, those involving the constants b,c) cannot contribute to select admissible weak solutions.Therefore, they can be set to zero each time entropies are used as admissibility criterion.However, entropies corresponding to a special family of Lipschitz continuous functions (f k ) k will play an important role in the definition of solutions to ARZ with point constraint.
It is not feasible to use entropy dissipation criteria to rule out contact discontinuities (in ARZ as it is described in [7], some contact discontinuities adjacent to vacuum states are declared not admissible on a phenomenological basis).However entropy dissipation can be used to enforce the right notion of admissibility of shocks (in ARZ as described in [7], a shock is admissible if and only if v decreases across the shock).
Indeed, if we fix f ≡ 0, b = c = 0, we have Lemma 2.2.The weak solution to system (1.1) taking the form with w − = w + (i.e., the shock propagating at speed σ . dν are the entropy-entropy flux pairs corresponding to (2.4) with zero f, b,c.
Proof.The claim follows by a straightforward explicit calculation, the inequality sign coming from (1.2).
It was observed in [23] that E g , g ≥ 0, are not convex but they are ρ-convex; we do not insist on this feature, focusing on the fact that these entropies are exactly those that select the admissible shocks in ARZ.However, it is not necessary to use the whole family of entropies {E g : g ≥ 0} in Lemma 2.2: one easily sees that it is enough to consider only the entropies E k .
The explicit expressions of these entropies and of the associated entropy fluxes (2.7) The family (E k ) k>0 spans the cone of all entropies E g , g ≥ 0, in the same sense in which Kruzhkov entropies span the cone of convex entropies in the scalar case.As in the case of Burgers equation, the nonadmissibility of shocks with v − < v + can be witnessed via the anti-dissipation of say that E k are Kruzhkov-like entropies for ARZ.We adopt the following definition holds for all φ ∈ D(R + × R), and moreover, for all k > 0 the entropy inequalities

Riemann solver, L ∞ and BV stability for ARZ
The delicate issue with ARZ concerns the vacuum: already at the level of the Riemann problem, if the initial condition contains a vacuum state there may co-exist two weak solutions satisfying the entropy admissibility criterion of Definition 2.3, one of which is non-physical, [7].Therefore, we adopt the viewpoint that ARZ is most appropriately determined by the Riemann solver rather than by the underlying system (1.1).
Contrarily to [7], we define the Riemann solver in terms of (v, w) and not in terms of (ρ, y) (this has already been the setting of [23] and [6]), which implies making particular choices in presence of vacuum states (v, w) T ∈ W 0 .The Riemann solver RS for ARZ in coordinates (v, w) contains the following elementary waves: with speed of propagation σ( = v and w = w r , then we have a contact discontinuity (3.11) for ARZ is defined as follows: Naturally, we also define RS [ W * , W * ] ≡ W * for any W * ∈ W .
Observe that the cases 1.-7. of Definition 3.1 do not cover all pairs ( W , W r ) ∈ W 2 : namely, we eluded the question of how to define RS for ( W , W r ) ∉ G, where Indeed, first, any initial datum given in physical variables (ρ, y) can be represented in variables (v, w) and discretized with approximate data belonging to G. Second, one can check that the domain G is invariant for the WFT algorithm based on the Riemann solver RS .
On the basis of the Riemann solver, sequences of approximate solutions can be constructed by finite volume schemes, see [6], or by the WFT method, see [19,24] and [21,23].Justification of their convergence relies on the following bounds.

Lemma 3.2. For any given v
is an invariant domain for RS .
The following BV stability property of ARZ is essential for the analysis of the convergence, see [6,23].For a function W = (v, w) T on R, define the total variation TV( W ) as the sum of the total variations of its components v and w.
Lemma 3.3.For all Riemann data ( W , W r ) ∈ G, the total variation of the function given in Definition 3.1 is equal to the total variation of the data.
With these basic ingredients at hand, we will now depart from the pure ARZ system.Notice that the corresponding existence and WFT results are a simple particular case (corresponding to

Locally constrained ARZ: definition and stability of solutions
The goal of this section is to define entropy admissible weak solutions to the constrained ARZ and q| x=0 ≤ q 0 , (4.12) where q 0 ∈ L ∞ (R + ).We start with the following lemma, which can be bypassed in the BV framework, cf.[4], where the existence of strong L 1 loc traces is obvious.
holds for all φ ∈ D(R + ×R).Then for any bounded Borel function f : R + → R, the function x → q( W )f(w) (•, x) admits a trace at x = 0 in the sense of the weak-* convergence in L ∞ (R + ).
Proof.The existence of right and left traces of q W f(w) is a consequence of Lemma 2.1 and of the general result in [12] on the existence of weak traces of L ∞ divergence-measure fields (one considers separately the domains {x < 0} and {x > 0}).The coincidence of the one-sided traces is due to the fact that W is a weak solution also in a neighbourhood of {x = 0}: It expresses the Rankine-Hugoniot condition.
In the context of Lemma 4.1, we will denote the coinciding right and left traces at {x = 0} by γ| x=0 ± [q f(w)] (in the sequel, sometimes they will also be abusively denoted by [q f(w)](t,0 ± ) in order to lighten the notation).For q 0 ≥ 0, k > 0 and with ζ defined in (2.5), let us introduce the auxiliary Lipschitz continuous functions Now we are ready to define solutions to (4.12) for the case of constant q 0 .Definition 4.2.Let W 0 ∈ L ∞ (R; W ). We say that W ∈ L ∞ (R + × R; W ) is an entropy weak solution of the constrained ARZ (4.12) with initial datum W 0 if, firstly, it is a weak solution of the pure ARZ (1.1) with initial datum W 0 (see Definition 2.3); secondly, and finally, the following Kruzhkov-like entropy inequalities with compensation terms hold in D (]0, +∞[ × R): (with the convention q( W )/q 0 .= 1 if q 0 = 0), for all k > 0 and E k , Q k given in (2.7).
Observe that all the trace terms in the above definitions are well defined due to Lemma 4.1.Because the properties required in Definitions 2.3 and 4.2 differ only at the interface {x = 0}, and because the notion of weak entropy solution can be localized by restricting the supports of the test functions, solutions of constrained ARZ (4.12) away from {x = 0} have the same properties as those of pure ARZ (1.1).The constraint at {x = 0} is enforced explicitly by (4.13).Moreover, an essential role is played by the second term in Kruzhkov-like entropy inequalities (4.14) that we call compensation term (cf.[11,13] for the introduction of analogous terms for constrained LWR model).Indeed, we have the following observation (stated, for the sake of simplicity, for solutions having traces at {x = 0} in the strong sense): Proposition 1.Let W be a function which is jump discontinuous at {x = 0}, and is constant in the regions {x < 0} and {x > 0}.Then W is a constrained entropy solution of (4.12) in the sense of Definition 4.2 if and only if the Rankine-Hugoniot jump condition at {x = 0} is satisfied, and one of the following cases occurs: (i) either [q( W )](•, 0 ± ) ≤ q 0 and the shock at {x = 0} is classical in the sense that it is admissible for pure ARZ system (1.1) (cf.Lemma 2.2); (ii) or [q( W )](•, 0 ± ) = q 0 and the shock is not admissible for pure ARZ system (1.1).
As for the constrained LWR case (see [5,13]), from item (ii) in the proposition above we see that an entropy solution to the constrained ARZ problem in the sense of Definition 4.2 may admit non-classical shocks located at {x = 0}, however they are conservative and they occur precisely at the level q 0 of the constraint.
We refer to [4] for a proof of Proposition 1.
The fundamental examples of entropy weak solutions to (4.12) can be found in [22], where a conservative constrained Riemann solver was constructed.Indeed, in a vicinity of {x = 0}, these solutions fall into one of the two cases (i), (ii) of Proposition 1; moreover, their definition implies that they are admissible for pure ARZ away from {x = 0}.Due to the self-similar structure of the Riemann solver, it is enough to check (4.14).
We will use the conservative Riemann solver of [22], written in the coordinates W as the RS of the previous section, as the building block of the WFT algorithm exploited in the next section.
Let us underline that there are simpler choices of "compensation term" in the Kruzhkov-like entropy inequalities (4.14).However, the interest of our definition resides in the following stability result, which strongly relies upon Lemma 2.1.
Proposition 2. Assume that (q i 0 ) i ⊂ R + converges to q 0 and that ( W i 0 ) i ⊂ L ∞ converges a.e. to W 0 as i → +∞.Let W i be the entropy weak solution of constrained ARZ (4.12) associated with the initial datum W i 0 and the constraint q i 0 .Assume that the sequence ( W i ) i is uniformly bounded in L ∞ and a.e.convergent to a limit function W . Then W is an entropy weak solution of constrained ARZ (4.12) associated with the initial datum W 0 and the constraint q 0 .
In contrast to the constrained LWR model where the "compensation term" added to Kruzhkov entropy inequalities depends only on the parameter k of the entropy considered, the last term of (4.14) makes appear a trace of some solution-dependent quantities.For this reason, the proof of Proposition 2 (and actually, the choice of the compensation term that makes this proof possible) is non-trivial.
Proof of Proposition 2. The proof consists in passing to the limit, as i → +∞, in the different relations of Definition 4.2.It is trivial to pass to the limit in all volume integrals, in the initial datum term of the weak solution identity and in the entropy inequalities (4.14).One has to care of the last term of (4.14), and check that the constraint (4.13) with W i , q i 0 , passes to the limit.Both delicate points follow from the Green-Gauss theorem, which reduces the weak convergence of interface traces to the convergence of some volume integrals (cf.[3, Lemma 7.1]).At this point, the renormalization property (2.6) of Lemma 2.1 and the fact that q 0 does not depend on t are essential ingredients of the proof.The details are given in [4].
Let us stress that the stability feature of Proposition 2 is a cornerstone of the existence proof we give in the next section.To conclude this section, let us briefly discuss the case of ARZ with variable in time point constraint q 0 (t ).Remark 2. If q 0 (•) is piecewise constant in time, it is enough to use the same definition, "gluing" in the last term of (4.14) the different terms γ| x=0 ± q( W ) q 0 f k (q 0 ; w) on the respective time intervals [t −1 , t [ where q 0 (t ) takes the value q 0 .
If q 0 (•) varies continuously, it becomes delicate to use the same approach, in particular because the renormalization property of Lemma 2.1 does not extend easily to functions of the form f k (q 0 (t ); •).In this case, one can approximate q 0 (•) in the L 1 loc sense by a family of piecewise constant functions (q m 0 ) m .Then one can use these approximations q m 0 (•) of q 0 (•) in the last term of inequalities (4.14) (in which case the traces are well defined), provided one adds to these inequalities an error term ε(m) R + φ(t , 0) dt .This term should vanish as m → +∞ while controlling, at each fixed m, the approximation error which can be estimated by With such definition, stability of entropy weak solutions under the a.e.convergence of solutions and constraints can be established.A simpler definition and stability proof can be given if we restrict our attention to BV loc constraints q 0 (•), indeed, in this case it is possible to approximate q 0 by a monotone sequence of piecewise constant functions q m 0 (•) m and thus avoid the introduction of error terms.
Right: the curve w = v + p(q 0 /v) in the (v, w)-plane.

The constrained Riemann solver and wave-front tracking
Our first objective is to adapt the Riemann solver RS for pure ARZ (1.1) to the constrained ARZ (4.12).At this point additional notation is needed, see Figure 1.
Definition 5.1.The Riemann solver RS q 0 for ( W , W r ) ∈ G is defined as follows: Thus, in the second case RS q 0 [ W , W r ] introduces non-classical shocks which can be characterized as stationary jumps in v across {x = 0} (w remaining constant across the non-classical shocks) joining two states with the flux equal to q 0 : the state Ŵ (q 0 , w ) on the left and the state W (q 0 , w ) on the right from x = 0, being v(q 0 , w ) < v(q 0 , w ) (that's why this discontinuity is not admissible in the sense of Definition 2.3; cf.Lemma 2.2).However, since the flux q at x = 0 of this solution equals to the constraint level q 0 , in accordance with Proposition 1 this non-classical shock is admissible in the sense of Definition 4.2 (while the classical solution RS [ W , W r ] becomes non-admissible with this new definition).More generally, one easily checks that RS q 0 provides entropy weak solutions to constrained ARZ (4.12) in the sense of Definition 4.2.One can also check that G remains an invariant domain in the WFT algorithm based on RS q 0 , so that there is no need to define RS q 0 on W 2 \ G. Now, we investigate the bounds that can be obtained with the constrained Riemann solver RS q 0 of Definition 5.1.Lemma 3.2 can be substituted by Then the domain (v, w) T ∈ W : v ≥ v * and w * ≤ w ≤ w * is invariant for RS q 0 .
At the same time, it is obvious from Definition 5.1 that Lemma 3.3 fails if RS is replaced with RS q 0 .E.g., any constant datum W c = (v c , w c ) T ∈ W such that q( W c ) > q 0 gives rise to a non-constant solution As a preliminary observation, note that the amplitude of the non-classical shock created at x = 0 in this case is given by J (q 0 , w c ) .= v(q 0 , w c ) − v(q 0 , w c ).In fact, the map J (q 0 , •) is a basic ingredient of the "interaction potential" needed to control the increase of the BV norm of the WFT approximations (such increase may occur at times where a non-classical shock appears or disappears at x = 0); we defer to [4].
The WFT algorithm of [19], further developed in [9,24], has become the classical tool for construction of solutions to hyperbolic systems of conservation laws, particularly in the situations where the solutions are defined by a (non-classical) Riemann solver.We pursue this strategy for approximation of the Cauchy problem for (4.12).
Let us briefly describe the preparatory work needed to implement the WFT for (4.12).One picks a sufficiently large h ∈ N and sets ε h .= 2 −h V 0 , where V 0 is an a priori L ∞ bound for approximate solutions with given initial datum W 0 (cf.Lemma 5.2).Then W is discretized by . The initial datum will be discretized with values in W h ; the Riemann solver RS (to be used for interactions at x = 0) and the constrained Riemann solver RS q 0 (to be used at x = 0) for ARZ will be substituted with the approximate Riemann solvers RS h , RS h q 0 such that for all ( Figure 2: Left: W h (dots), C h (highlighted) in the (ρ, q)-plane.
Right: sets W h and C h pictured in the (v, w)-plane.
To this end, the rarefactions (3.10) are discretized in a classical way, carefully choosing the jump velocities; in particular, the resulting fans of non-admissible (in the sense of Lemma 2.2 and Definition 2.3) shocks provides a weak solution to ARZ.Next, the classical shocks (3.9) and the contact discontinuities (3.11) in RS , RS q 0 are not modified.The last modification that only concerns RS q 0 : one has to approximate the non-classical shocks so that their end states Ŵ , Wr belong to W h .This forces us to relax the equality q( Ŵ ) = q 0 = q( Wr ) which holds true for non-classical shocks in RS h q 0 : actually we may even loose the desirable conservation property at x = 0, whenever a non-classical shock is present.Here, let us give the details.The region of states W ∈ W that do not satisfy the constraint q( W ) We approximate it with the set (pictured in Figure 2) Let wh (q 0 ) be the minimal value of w ∈ ε h N for which the set C h (q 0 ) has two elements with the second coordinate equal to w.We introduce the approximations in C h (q 0 ) of v and v given in (5.15), (5.16): vh (q 0 , w) .= max ε h N ∩ 0, v q 0 , w , vh (q 0 , w) .= max ε h N ∩ 0, v q 0 , w .Now, by following [24] it is easy to introduce the WFT algorithm to construct approximate solutions of the Cauchy problem for (4.12).The painstaking analysis developed in [4] aims at showing that the total number of wave-fronts appearing in the construction remains finite; moreover, the total variation of W h (t , •) remains bounded uniformly in time, in spite of the fact that the constrained Riemann solver RS h q 0 is not variation-diminishing.With BV bounds in hand, one derives compactness of solutions; then, convergence (up to a subsequence) of the algorithm can be justified combining classical ideas of analysis of WFT schemes and the idea of Proposition 2. This material can be found in [4].

Numerical approximation of constrained ARZ model
In this section we describe a finite volumes numerical scheme for the constrained ARZ model in the conservative form.The same scheme is to be employed in our ongoing research project to investigate the ARZ system subject to more general point constraints, possibly non-local.Our scheme is obtained combining the techniques presented in [22] to treat the constrained Riemann problem for ARZ, together with the Transport-Equilibrium scheme introduced in [10] to deal with the contact discontinuities.The main differences between our scheme and the aforementioned ones is that we do not use an exact Riemann solver in the implementation of the scheme as we consider general initial conditions and we wish to keep our scheme as fast as possible.
In order to fix the notation we recall that the conservative variables for ARZ are Y = (ρ, y) T , and that the Cauchy problem for the system writes where F ( Y ) .= (ρ v, y v) T , v .= y ρ − p(ρ), q = ρ v as above, and Y .= ( ρ, ȳ) is the initial condition.
Let us denote by [x , x r ]×[0, T ] the computational domain and let ∆x and ∆t be the constant space and time increments.We introduce the points x j +1/2 .= x + j ∆x for 0 ≤ j ≤ M , where M is the natural number such that x M +1/2 = x r .Then for any 1 ≤ j ≤ M , we define the cells K j .= [x j −1/2 , x j +1/2 [ with cell centres We introduce the index j c such that x j c +1/2 is the location of the constraint (a tollgate or a traffic light).We define the time discretization t n = n ∆t for 0 ≤ n ≤ N , where N is the natural number such that t N = T .For 1 ≤ j ≤ M and 0 ≤ n ≤ N , we denote by Y n j = (ρ n j , y n j ) T the approximation of the average of Y (t n , • ) on the cell K j , namely Observe that Y n j is a vector: we omit the upper arrow in order to keep the notation as light as possible.To each state Y n j we can associate a unique pair of "Riemann invariants" (v n j , w n j ) as follows: • If ρ n j = 0, then w n j = y n j /ρ n j and v n j = w n j − p(ρ n j ).
• If ρ n j = 0 and j = 1, then w n j = v n j = w n j −1 .
Note that the definition above is not motivated by any physical insight, but just designed to optimize computations, as we already pointed out in Section 1.2.
Assume for the moment that we do not enforce the constraint to our problem.Even in this simpler case we cannot apply directly a classical conservative scheme because it may generate important non-physical oscillations near contact discontinuities.Note that contact discontinuities are unavoidable in solutions of the ARZ model with a generic initial condition.In [10] the authors propose an efficient numerical strategy, inspired by the random sampling technique, to address this issue.The main idea is to separate the treatment of the contact discontinuity from the computation of other waves in the solution.
Assume (Y n j ) 1≤ j ≤M are given: Instead of computing immediately (Y n+1 j ) 1≤ j ≤M Chalons and Goatin [10] introduce an intermediate step n +1/2 so that from n to n +1/2 they only update the position of the possibly present contact discontinuities, while from n + 1/2 to n they use a (essentially) standard scheme to update the values in all the cells.
In practice to go from n to n +1/2 one generates a random number a n+1 and, for each value of j , defines Here Y * (Y Then, to go from n + 1/2 to n + 1 one uses the standard scheme and In this paper, the numerical flux F .= ( 1 , 2 ) T is the HLL flux, see [8], defined as follows Now we are ready to consider the Cauchy problem (6.17) subject to a constraint located at x j c +1/2 .We proceed exactly as in [22] (see also [5,11] for the LWR case) and define the new numerical flux F n j +1/2 .= (ˆ1,ˆ2) T as follows: • • If j = j c , then the components of the vector F n j c +1/2 are ˆn 1, j c +1/2 .= min n 1, j c +1/2 , q n and ˆn 2, , where q n 0 is an approximation of q 0 (t n ).

Elements of validation of the scheme
In this section we compare an exact solution of a constrained Cauchy problem for the ARZ model obtained by WFT method, with the numerical simulations obtained from the same initial condition by the numerical scheme described in Section 6.
This first comparison allows us to check that our scheme correctly locates contact discontinuities and enforces the constraint.This example also shows that the total variation of the solution can increase abruptly due to the non-classical Riemann solver at the constraint.A complete validation of the scheme together with its convergence and more detailed examples will be presented in a forthcoming paper.
Consider the domain of computation [−30, 30] × [0, 3], take the anticipation function p(ρ) = ρ 3 , and assume that the initial condition takes the form and is such that V the common value of v A and v B .The constraint is located at x = 0 and its constant value, q 0 , satisfies The solution of this problem consists of a contact discontinuity originated at (t , x) = (0, −10), which reaches the constraint location at time T 1 .Then, as the state on the left to the contact discontinuity does not satisfy the constraint, the solution changes its profile and new waves are created, according to the nonclassical Riemann solver RS q 0 : A shock wave of negative speed, a stationary non-classical shock and a shock wave of positive speed after which we finally recover the contact discontinuity.
The exact solution, corresponding to the values Figure 3: The (exact) solution described in Section 7 represented in the (t , x)-plane (left) and in the (ρ, q)plane (right).
is represented in the (t , x)-plane and in the (ρ, q)-plane in Figure 3.The time at which the contact discontinuity reaches the position of the constraint is T 1 = 5/3 and the solution is computed up to the final time In Figure 4, we compare the exact solution with the numerically computed solution at time T = 3.The parameters of discretization for the numerically computed solution are ∆x = 5×10 −3 and ∆t = 10 −4 .We observe a good agreement between the two solutions, in particular the contact discontinuity is well captured.
In Figure 5, we show the comparison between the two solutions when focusing on the shock discontinuity at (3, −19.5869).The approximate solution is computed for different values of the space step while the time step is kept to ∆t = 10 −4 .We observe that the accuracy of the approximate solution increases as the space step decreases.This numerical convergence is also observed in Table 1, where the relative L 1 errors between the two solutions are computed for different space step at time T = 3 and for the fixed ∆t = 10 −4 .Using the logarithm scale we can deduce from Table 1 the orders of convergence which are approximately 0.906 and 0.869 for ρ and y respectively.
Finally, as due to the random sampling the scheme presented in Section 6 is non-conservative, we check in        We can easily see from Table 2 that these conservative errors are small and converge to zero with the mesh size.

( 1 . 1 )
in variables W . = (v, w) T , which are the Riemann invariants for ARZ.Then the conservative variables for (1.1) are (ρ, y) T , with ρ .= p −1 (w − v) and y .= w p −1 (w − v), defined in the domain W . = {(v, w) T ∈ R 2 + : v ≤ w}.We denote the density flux ρ v by q.Since ρ = 0 if and only if v = w, in the (v, w) coordinates we can distinguish different vacuum states which turns out to be very convenient from the mathematical viewpoint, even if it has no physical meaning.Hence, we denote by W 0 .= {(v, w) T ∈ R 2 + : v = w} the set of all vacuum states and by W c 0 .= {(v, w) T ∈ R 2 + : v < w} the set of non-vacuum states.

Figure 4 :
Figure 4: The components of the exact solution and the numerically computed solution at time T = 3 corresponding to the values (7.19) for the initial datum (7.18).

Figure 5 :
Figure 5: The shock discontinuity at (t , x) = (3, −19.5869) and its numerical approximations for different values of ∆x corresponding to the values (7.19) for the initial datum (7.18).

Table 1 :
Relative L 1 errors at time T = 3.

Table 2 :
Measurement of the relative total mass errors on both ρ and y as defined by(7.20)