A Godunov type scheme for a class of LWR traffic flow models with non-local flux

We present a Godunov type numerical scheme for a class of scalar conservation laws with non-local flux arising for example in traffic flow models. The proposed scheme delivers more accurate solutions than the widely used Lax-Friedrichs type scheme. In contrast to other approaches, we consider a non-local mean velocity instead of a mean density and provide $L^\infty$ and bounded variation estimates for the sequence of approximate solutions. Together with a discrete entropy inequality, we also show the well-posedness of the considered class of scalar conservation laws. The better accuracy of the Godunov type scheme in comparison to Lax-Friedrichs is proved by a variety of numerical examples.


Introduction
Over recent years, non-local conservation laws gained growing interest for a wide field of applications such as supply chains [3], sedimentation [4], conveyor belts [12], crowd motion [8] or traffic flow [5,10]. In the latter case, the well-known Lighthill-Whitham-Richards (LWR) model [15,16] has been extended by considering non-local velocity terms depending on the downstream traffic so that drivers adapt their velocity to the mean traffic in front, see [5,10].
The well-posedness of special non-local flux problems has been investigated in for example [1,2,7,9]. However, only a few numerical schemes have been applied so far to solve these type of equations. The most common approach are first order Lax-Friedrichs (LxF) type schemes [2,4,5,7,10], while recently second-and higher-order schemes have been introduced [6,11]. We remark that also these higher-order methods rely on LxF type numerical flux functions, which imply the same drawbacks known from local conservation laws. Certainly, the LxF type scheme offers a powerful tool to numerically analyze non-local flux problems but typically leads to approximate solutions with strong diffusive behavior. As we are interested in a more accurate approach, we present a Godunov type scheme for a class of scalar conservation laws with non-local flux. In addition, by deriving several properties of the scheme, we prove the well-posedness for these special non-local conservation laws, which in contrast to other models [5,7,10] focus on a non-local mean velocity of the downstream traffic. Furthermore, the Godunov type scheme approach allows for physically reasonable solutions meaning that a maximum principle is satisfied and negative velocities as well as negative fluxes are avoided. This work is organized as follows: In Section 2, we present the considered class of non-local conservation laws for traffic flow. Afterwards, we introduce the Godunov type scheme and derive important properties of the scheme such as L ∞ and bounded variation (BV) estimates in Section 3. Those are also used to show the well-posedness of the proposed traffic model. In Section 4, we present numerical examples, which demonstrate the better accuracy of the Godunov type scheme in comparison to the widely used LxF type scheme and also provide a comparison to the model in [7].

Modeling
We briefly present an already existing traffic flow model with non-local flux originally introduced in [5,7,10]. Based on the modeling ideas therein, we propose an adapted model and show its well-posedness.
The key difference appears in the flux function, where instead of a mean downstream density a mean downstream velocity is considered.

An existing model with mean downstream density
The model considered in [7] is given by a scalar conservation law of the form where For initial conditions the existence and uniqueness of weak entropy solutions is stated in [7,Theorem 1] if the following hypotheses are satisfied: Note that a whole family of kernel functions w η is considered to also analyze the limit behaviour of the model (η → 0 and η → ∞).

Model considering a mean downstream velocity
The non-local model (1) to (3) can be applied in the context of traffic flow and chooses the velocity based on a mean downstream traffic density. In contrast to this approach, it is also reasonable to assume that drivers adapt their speed based on a mean downstream velocity, anticipating the future space in front of them, see Figure 1. Therefore, we consider a slightly different model compared to (1) to (3), namely where and we have given initial conditions For simplicity, let us also define In (4) and (5), we assume the same hypotheses as for (1) and (2) and one additional restriction: As we will see in Section 3, the reformulation of the original model (1) to (3) keeps the main properties and allows for a straightforward application of a Godunov type scheme.
The weak entropy solutions of problem (4) to (6) are intended in the following sense: Our main result concerning the new model is given by the following theorem, which states the well-posedness of problem (4) to (6). Theorem 2.3. Let ρ 0 ∈ BV(R; [0, ρ max ]) and hypotheses (H2) hold. Then, the Cauchy prob- admits a unique weak entropy solution in the sense of Definition 2.2 and The proof consists of two parts: existence and uniqueness of entropy solutions. While the uniqueness proof follows from the Lipschitz continuous dependence of weak entropy solutions on the initial data, the existence proof is based on a construction of a converging sequence of approximate solutions defined by a Godunov type scheme.

Uniqueness of entropy solutions
One part of the proof to Theorem 2.3 is to show uniqueness of entropy solutions for the model (4) to (6). Therefore, we prove the Lipschitz continuous dependence of weak entropy solutions with respect to the initial data. Here, we follow [5,7,10] and use Kruzkov's doubling of variables technique [13]. Note that in the following · denotes · L ∞ . Theorem 2.4. Under hypotheses (H2), let ρ and σ be two entropy solutions of (4) to (6) with initial data ρ 0 and σ 0 , respectively. Then, for any T > 0, there holds with K given by (12).
Proof. The functions ρ and σ are weak entropy solutions of respectively, and V , U are bounded measurable functions and Lipschitz continuous w.r.t. x since ρ, σ ∈ (L 1 ∩ L ∞ ∩ BV )(R + × R; R). Using the classical doubling of variables technique, we get the following inequality: where ρ x must be understood in the sense of measures. Applying the mean value theorem and using the properties of the kernel function, we deduce Using the Leibniz integral rule and again the mean value theorem, we can also obtain for a.e.
The existence of weak entropy solutions is now proved in Section 3. We therefore introduce a Godunov type scheme used to construct approximate solutions.

A Godunov type scheme
The main new contribution of this work is to develop a suitable Godunov type numerical scheme for the non-local model (4) to (6). We derive L ∞ and BV bounds for the approximate solutions and further provide, due to a discrete entropy inequality, all ingredients to finally prove the well-posedness result Theorem 2.3.

Numerical scheme
We take a space step h such that for the size of the (downstream) kernel η = N h with N ∈ N holds. The time step is denoted by τ and we set λ = τ /h. We define the space grid by x j+ 1 2 = (j + 1 2 )h being the cell interfaces and x j = jh the cell centers for j ∈ Z, see Figure 2. The finite volume approximate ]. Within the proposed scheme, we intend to mimic the numerical flux function of the Godunov scheme for local conservation laws, i.e., minimizing or maximizing the flux within [ρ n j , ρ n j+1 ] or [ρ n j+1 , ρ n j ] depending on whether ρ n j ≤ ρ n j+1 or ρ n j+1 ≤ ρ n j , respectively. As the convolution term V (t, x) is a non-local velocity function (see (7)), there is no straightforward way of adapting this result to the non-local case. Therefore, we first examine the flux at the interface x j+ 1 2 , where the convolution term for t ∈ [t n , t n+1 [ is then approximated by V n with which is motivated by Remark 3.1. For (14) we follow [4, Equation (3. 2)] to satisfy if W 0 = 1 and γ k is computed exactly or by an appropriate quadrature formula (such that all γ k are non-negative and N −1 k=0 γ k = W 0 = 1 holds). An example for the computation of V n Summarizing, the entire Godunov type scheme is initialized by the initial data ρ 0 as and can be computed using (13) and (16) by the finite volume scheme Remark 3.2. Analogously to the scheme (17) and (18) with (13), a Godunov type scheme for the model (1) to (3) considered in [7] can be derived. For this, similar properties can be shown analogously to the following sections. One major advantage of our Godunov type scheme unlike the LxF type scheme used in [5,7,10] is that the numerical fluxes F n

Maximum principle
The approximate solutions constructed by the Godunov type scheme (18) satisfy a strict maximum principle: Let hypotheses (H2) hold. For a given initial datum ρ 0 j , j ∈ Z with ρ M = sup j∈Z ρ 0 j and ρ m = inf j∈Z ρ 0 j , the approximate solutions constructed by the scheme (18) satisfy the bounds if the following Courant-Friedrichs-Levy (CFL) condition holds: Proof. We prove the claim per induction. For n = 0 the claim is obvious, so we suppose ρ m ≤ ρ n j ≤ ρ M , ∀j ∈ Z holds for a fixed n ∈ N. Before considering ρ n+1 j we show some general inequalities: where we used the monotonicity of w η and v, and the mean value theorem. Analogously, we obtain from (20) By multiplying inequality (21) by g(ρ M ), subtracting V n j+ 1 2 g(ρ n j ) and applying the mean value theorem, we obtain Therefore, under the CFL condition (19), we have Analogously, we obtain which gives us the claim.
The maximum principle ensures that the numerical solution to (4) to (6) is bounded from above by ρ M = sup j∈Z ρ 0 j ∈ [0, ρ max ] and hence does not exceed the maximal density ρ max . In addition, the scheme is positivity preserving as the solution stays non-negative, since ρ m = inf j∈Z ρ 0 j ≥ 0.

BV estimates
Next, we derive a BV estimate for the approximate solutions constructed by the Godunov type scheme (18). Similar to the LxF type scheme analyzed in [5,7,10], BV estimates cannot be derived using the standard general approaches. In particular, the Godunov type scheme also does not fit into the classical assumptions of total variation diminishing (TVD) schemes, as the total variation may slightly increase (as it is the same for the analytical solution). Nevertheless, the numerical scheme has a bounded total variation. Further, to finally also prove the existence of solutions to the model (4) to (6), we also need to provide a bound on the (discrete) variation in space and time. We begin with the BV estimate in space: Theorem 3.4. Let hypotheses (H2) hold, ρ 0 ∈ BV (R; [0, ρ max ]) and letρ be given by (18). If the CFL condition (19) holds, then for every T > 0 the following discrete space BV estimate is satisfied: Proof. Let us define ∆ n j+k− 1 2 := ρ n j+k − ρ n j+k−1 .
Then, we obtain ∆ n+1 where ξ j+ 1 2 is between ρ n j and ρ n j+1 . With (20) we derive where ζ j+ 1 2 is again between ρ n j and ρ n j+1 . Thus, we have ∆ n+1 Due to the CFL condition (19) and hypotheses (H2), all terms before the differences are positive and we get j |∆ n+1 Rearranging the indices we obtain j |∆ n+1 Using inequality (21), for which (21) ≤ γ 0 v ρ max holds, and with γ 0 ≤ hw η (0), we obtain j |∆ n+1 Therefore, we recover the following estimate for the total variation We are now able to also provide an estimate for the discrete total variation in space and time: Theorem 3.5. Let hypotheses (H2) hold, ρ 0 ∈ BV (R; [0, ρ max ]) and letρ be given by (18). If the CFL condition (19) holds, then for every T > 0 the following discrete space and time total variation estimate is satisfied: If we consider the scheme (18), we obtain Taking absolute values yields Summing over j and rearranging the indices gives us Therefore, we recover as desired.

Discrete entropy inequality
As another desirable property and final ingredient regarding the proof of Theorem 2.3, we next show that the approximate solutions obtained by the Godunov type scheme (18) fulfill a discrete entropy inequality. Therefore, we follow [2,5,7,10] and define max(a, b) and a ∨ b = min(a, b).
Proof. The proof closely follows [2,5,7]. We set which is a monotone non-decreasing function with respect to each variable under the CFL condition (19) since we have Moreover, we have the identitỹ Then, by monotonicity, the definition of the scheme (18) and by using (for the last inequality) the non-negativity of (a, b) → (sgn(a + b) − sgn(a))(a + b), we get (24): ).

Proof of Theorem 2.3
Since we have already shown uniqueness of weak entropy solutions to the model (4) to (6), it remains to finalize the existence proof. Similar to [5, Section 4] and [7, Proof of Theorem 1], the convergence of the approximate solutions constructed by the Godunov type scheme (18) towards the unique weak entropy solution can be proven by applying Helly's theorem. The latter can be applied due to Theorems 3.3 and 3.5 and states that there exists a sub-sequence of the constructedρ that converges to some ρ ∈ (L 1 ∩ L ∞ ∩ BV )(R + × R; I). Following a Lax-Wendroff type argument, see [14,Theorem 12.1], one can show that the limit function ρ is a weak entropy solution of (4) to (6)

Numerical examples
In this section, we present some numerical examples demonstrating the advantages of the Godunov type scheme in comparison to the widely used LxF type scheme. The latter will be briefly introduced in the following section. In addition, we also comment on the differences between the model considered in this work (4) to (6) and the earlier one (1) to (3).

A LxF type scheme
The common scheme used so far for the problem (1) to (3) is a LxF type scheme, where the downstream velocity of the convolution term is computed by with w k η = w η (kh). This discretization of w η implies (for W 0 = 1) Thus, the approximation (25) of the convolution term slightly overestimates the traffic density and therewith underestimates the velocity. Further, unphysical densities beyond ρ max or negative velocities are possible. This can be both avoided by discretizing the kernel function as proposed in (14). In the following subsections we will use this discretization to avoid that the accuracy studies are biased by different quadrature rule errors. The numerical flux function of the LxF scheme is given by with α ≥ 0 being the viscosity coefficient. This leads to the scheme For the corresponding CFL condition and restrictions on α, we refer to [7, Proposition 2].
Remark 4.1. Note that the LxF type scheme can be adapted to the model (4) and vice versa the Godunov type scheme to model (1), where in comparison the LxF type scheme adds more diffusion to the numerical solution (see also Section 4.2).

Accuracy of the Godunov type scheme
In the following sections, we consider the solution of model (4) to (6) and model (1) to (3) in the case of traffic flow. So we have g(ρ) = ρ, a velocity function v(ρ) specified below, and a road of length L = 1. We are interested in the solution at a certain final time T for the initial conditions For simplicity, we use periodic boundary conditions in all our examples. To compute the L 1 error of an approximate solution ρ h with step size h compared to a reference solution ρh with step sizeh at time T , we apply Now, in order to compare the accuracy of the Godunov scheme and the LxF scheme, we first use the linear velocity function v(ρ) = 1 − ρ, as for this choice both models coincide (see Remark 2.1) and the different discretization schemes both have been well analyzed. We consider the final time T = 0.1 and the kernel function w η (x) = 3(η 2 − x 2 )/(2η 3 ) with η = 0.1. We compare the L 1 distances of the LxF and the Godunov type scheme to a reference solution computed with the LxF type scheme andh = 0.02 · 2 −9 . The spatial step size is given by h = 0.02 · 2 −n with n ∈ {0, . . . , 6}. The time step parameter τ is given by the minimum of the CFL condition (19) and the CFL condition of the LxF type scheme.
The results for this first test case are given in Table 1. Obviously, the L 1 errors of the Godunov type scheme are significantly smaller than the ones of the LxF scheme. Table 1: n Godunov LxF 0 9.38e-03 1.99e-02 1 6.97e-03 1.30e-02 2 4.29e-03 9.31e-03 3 3.00e-03 6.41e-03 4 1.96e-03 4.27e-03 5 1.33e-03 2.71e-03 6 9.05e-04 1.64e-03 The better accuracy of the Godunov type scheme can also be seen directly in Figure 3. We notice that in the presence of the two discontinuities, the Godunov type scheme in particular shows a better resolution of the solution structure than the LxF scheme at the left-hand side, while the resolution for the rarefaction wave close to the jump on the right-hand side is quite similar for the short time period T = 0.1.
If we consider a longer time period, i.e. T = 1, the good performance of the Godunov type scheme can be observed in Figure 4. Here, the reference solution is computed with the LxF type scheme with a spatial step size ofh = 0.02 · 2 −7 .
In addition, we consider the accuracy for a non-linear velocity function, v(ρ) = 1 − ρ 5 .
We choose the constant kernel w η (x) = 1 η with η = 0.1 and consider again the solution of the initial conditions (27) but at time T = 0.05 on a road of length L = 1. The spatial step size is given as above by h = 0.02 · 2 −n with n ∈ {0, . . . , 6}. The reference solution is computed by the LxF scheme adapted to the model (4) and withh = 0.02 · 2 −9 . The results for the non-linear test case can be seen in Table 2. Similar to the linear test case, the L 1 errors of the Godunov type scheme are significantly smaller than the ones of the LxF scheme. In addition, the better resolution of the Godunov type scheme can be seen in Figure 5. Again, similar to the linear case, the Godunov type scheme in particular shows a better resolution of the solution structure than the LxF scheme at the left-hand side, while the resolution for the rarefaction wave close to the jump on the right-hand side is quite similar.

Comparison of the models
Next, we aim to discuss the differences between the models (4) to (6) and (1) to (3). To see the different dynamics within the two models, we have to choose a non-linear velocity function and we choose the same non-linear velocity v(ρ) = 1 − ρ 5 as before with the same parameter η = 0.1 for the constant kernel function and final time T = 0.05. For a fair comparison of the evolution of densities, we apply the Godunov type scheme to both models. Figure 6 shows the results for a spatial step size h = 0.01 and a time step size τ determined by the CFL condition. It can be observed that the approximate solution of the earlier model (1) to (3) contains rather large oscillations while resolving the traffic jam. In contrast to that, the model (4) to (6) resolves the traffic jam in a more monotone way.

Limit η → 0
Finally, we take a look at the limit case η → 0 and investigate numerically whether the approximate solutions constructed by the proposed Godunov type scheme converge towards the solution of the (local) LWR traffic model. Note that this property is far from obvious since the constants in Theorems 2.4 and 3.5 blow up (see also [7]). For a numerical investigation, we consider the same (non-linear) scenario as above with a fixed space step size h = 0.5 · 10 −4 and vary η ∈ {10 −1 , 10 −2 , 10 −3 , 10 −4 }. As final time we take T = 0.05.
To evaluate the convergence, we compute the L 1 distances between the approximate solutions obtained for the proposed Godunov type scheme applied to (4) to (6) and the result of a classical Godunov scheme for the corresponding local LWR problem. Obviously, the corresponding L 1 distances shown in Table 3 demonstrate the convergence towards the solution of the local traffic model. The results are further illustrated in Figure 7. Table 3: L 1 distances between the approximate solutions to the local LWR model and the non-local model for different η at T = 0.05. η 10 −1 10 −2 10 −3 10 −4 L 1 distance 4.46e-02 6.85e-03 9.90e-04 1.60e-04

Conclusion
In this work, we have presented a Godunov type scheme for a class of non-local conservation laws. For this novel scheme we provide L ∞ and BV bounds as well as a discrete entropy inequality. Based on these results, we also proved the well-posedness, i.e., existence and uniqueness of weak entropy solutions. The proposed Godunov type scheme can be adapted to other classes of non-local conservation laws and is very promising as it adds less numerical diffusion to the solution compared to the LxF type scheme. Certainly, this advantage can also be exploited by using the underlying numerical flux function within higher-order methods. In future work we aim at constructing several higher order methods based on the presented Godunov type scheme. In addition, the considered model with mean downstream velocity may be advantageous in the context of networks. Here we aim at investigating appropriate coupling conditions and suitable discretization schemes.