Liudmila Tumash

Carlos Canudas De Wit

Maria Laura

Delle Monache

M L Delle Monache

Boundary and In-Domain Control for Large-Scale Urban Traffic Networks

come

Boundary and In-Domain Control for Large-Scale Urban Traffic Networks

Liudmila Tumash, Carlos Canudas-de-Wit, Fellow, IEEE and Maria Laura Delle Monache, Member, IEEE Abstract-This paper presents a novel approach to design control for traffic on large-scale urban networks by analyzing the structure of one single partial differential equation (PDE). In particular, we elaborate a method that represents a curvilinear coordinate transformation translating a 2D conservation law into a parametrized set of equations each having a structure of inhomogeneous 1D LWR equation. This resulting system can be explicitly analyzed for traffic evolving on urban networks of arbitrary size. As an application example, we demonstrate how the 2D model in curvilinear coordinates can be used to design two different controllers for urban traffic. First, a boundary controller is designed to track any desired space-and timedependent vehicle density profile with the help of Hamilton-Jacobi formalism. Second, we design an in-domain variable speed limit (VSL) controller that steers traffic flow such that any space-varying equilibrium can be achieved. We validate the control results numerically using the structure of Grenoble downtown. Index Terms-boundary control, Hamilton-Jacobi, conservation law, urban traffic control, variable speed limit.

I. INTRODUCTION

C ONTINUING urbanization caused by ever-growing pop- ulation of the planet implies a growing demand for transportation. This entails formation of severe traffic congestions that cost people hundreds of hours per year and that also have a significant negative impact on the environment, see Urban Mobility Report for 2019 in USA [START_REF] Schrank | Urban Mobility Report 2019[END_REF]. This requires development of scalable models able to predict congestion formations and control techniques to resolve them.

The most common and simple model to describe traffic behaviour is the LWR model presented by Lighthill, Whitham [START_REF] Lighthill | On kinematic waves, II: A theory of traffic flow on long crowded roads[END_REF] and Richards [START_REF] Richards | Shock waves on the highway[END_REF] in the fifties. This macroscopic model describes temporal evolution of aggregated quantities (traffic density and kinematic wave speed) as fluids. In particular, the LWR model is based on the conservation principle, where the conserved quantity is the density of vehicles. Mathematically, this model is a first-order hyperbolic PDE with a concave flux function that represents an empirical relation between flow and density, see [START_REF] Li | A generic characterization of equilibrium speed-flow curves[END_REF] for a review on flow-density functions.

However, the LWR model was originally designed to describe traffic flow on a single infinite road. Thus, additional conditions and constraints had to be imposed in order to model traffic flow on networks that consist of links (roads) and nodes (junctions). A methodology for intersection modelling within LWR framework was proposed, e.g., in [START_REF] Lebacque | First-order macroscopic traffic flow models: Intersection modeling, network modeling[END_REF]. For a general theory of traffic flow on networks see [START_REF] Garavello | Traffic flow on networks[END_REF]. The Cauchy problem for complex networks (with more than two incoming and outgoing roads at junctions) was considered in [START_REF] Garavello | Conservation laws on complex networks[END_REF].

The main challenge in this link-level representation of traffic networks is a large computational time that significantly exaggerates optimization of large networks consisting of thousands of roads [START_REF] Ziliaskopoulos | Large-scale dynamic traffic assignment: implementation issues and computational analysis[END_REF]. Moreover, these models are often not scalable and mathematically intractable. Alternatively, for modelling of urban traffic one can use two-dimensional continuous models. These models describe urban traffic as if it were a fluid evolving on a continuum plane in R 2 , and the urban network (collection of roads and junctions) is embedded into this plane. The advantages of using such continuous models is that they do not require a high computational effort and require much less data to tune than road-by-road models.

The first works describing transportation networks in terms of aggregated variables appeared several decades ago, see [START_REF] Smeed | The road capacity of city centers[END_REF], [START_REF] Herman | A two-fluid approach to town traffic[END_REF]. These early models, however, failed in capturing traffic flow dynamics during a rush hour due to the lack of any knowledge on a flow-density function. The existence of a Macroscopic Fundamental Diagram (MFD) in congested urban regions has been observed empirically [START_REF] Geroliminis | Existence of urban-scale macroscopic fundamental diagrams: Some experimental findings[END_REF], and was generalized in [START_REF] Daganzo | An analytical approximation for the macroscopic fundamental diagram of urban traffic[END_REF]. This discovery led to appearance of reservoir models, which describe traffic as evolution of the total number of cars in some zone (reservoir) depending on its inflows and outflows. For a detailed review on several MFD-based models see [START_REF] Aghamohammadi | Dynamic traffic assignment using the macroscopic fundamental diagram: A Review of vehicular and pedestrian flow models[END_REF]. Their main drawback is that in case of heterogeneously congested links, MFD might not be well-defined, see [START_REF] Geroliminis | Properties of a well-defined macroscopic fundamental diagram for urban traffic[END_REF] for properties of well-defined MFDs. In general, there must exist only one flow value for a given number of vehicles. This feature is preserved only in regions that consist of links characterized by similar congestion levels. To overcome this limitation, [START_REF] Leclercq | Macroscopic traffic dynamics with heterogeneous route patterns[END_REF] presented partitioning algorithms to divide a urban area into multiple zones each having its own well-defined MFD. However, these reservoir models are non-adaptive to changing traffic conditions, as partitioning depends on the current level of congestion in each zone, see also [START_REF] Tumash | Equilibrium Manifolds in 2D Fluid Traffic Models[END_REF] for more details.

An alternative way to describe the propagation of urban traffic is to use models based on 2D conservation laws that describe traffic state by its density, see [START_REF] Aghamohammadi | Dynamic traffic assignment using the macroscopic fundamental diagram: A Review of vehicular and pedestrian flow models[END_REF] for a review. Recently, [START_REF] Mollier | Twodimensional macroscopic model for large scale traffic networks[END_REF] proposed an extension of the original LWR model to two dimensions assuming that the direction of traffic flow is determined by geometry and infrastructure of the underlying urban network. This is captured by introducing space-dependency into the flow-density function.

Transportation networks are controlled in order to improve their overall efficiency, e.g., decrease total travel time. There exists a vast choice on literature devoted to practical network control techniques, such as routing of traffic [START_REF] Gugat | Optimal control for traffic flow networks[END_REF], traffic light control [START_REF] Goettlich | Traffic light control: a case study[END_REF], ramp metering [START_REF] Papageorgiou | ALINEA: a local feedback control law for on-ramp metering[END_REF], variable speed limits [START_REF] Delle Monache | Traffic Regulation via Controlled Speed Limit[END_REF], see also [START_REF] Papageorgiou | Review of Road Traffic Control Strategies[END_REF] for a review. However, most of them are based on the first-discretize-then-optimize procedure or require using much of traffic data. On the contrary, being able to derive a controller for a PDE-like model without discretizing it, allows to obtain analytical results in a closed form.

In this paper we propose a new technique that brings a 2D conservation law model into a form that can be analyzed in a similar way as the well-known LWR model. Information on boundary flow data and network (its infrastructure and geometry) is incorporated as model parameters. This modified model represents a single equation that can be used for a pure analytic control design. Hence, two different control designs are performed: boundary control for tracking any desired space-and time-varying profile, and in-domain VSL control to steer traffic to any desired space-varying equilibrium.

Our main contributions are the following: 1) We present a method of transforming the 2D LWR model [START_REF] Mollier | Twodimensional macroscopic model for large scale traffic networks[END_REF] into a parametrized set of 1D systems, which enables an explicit elaboration of control strategies for various tasks to solve on arbitrarily large urban networks. The only limitation thereby is that this 2D conservation law is applicable only to urban traffic with a preferred direction of motion. 2) We explicitly derive a boundary controller for this 2D LWR equation to track a space-and time-dependent trajectory that admits discontinuities. This is the first time such a general control result is obtained in the context of large-scale traffic. For this we use the Hamilton-Jacobi framework, as it was done in [START_REF] Tumash | Boundary Control Design for Traffic with Nonlinear Dynamics[END_REF], but extending it to 2D and handling explicit space-dependency of the fundamental diagram by applying the viability theory for the case of space-dependent Hamiltonians explained in [START_REF] Désilles | Optimal dynamic route guidance: A model predictive approach using the macroscopic fundamental diagram[END_REF], [START_REF] Aubin | Viability Theory: New Directions[END_REF]. There is also another boundary controller for 2D LWR [START_REF] Tumash | Topologybased control design for congested areas in urban networks[END_REF] for a special case of congested traffic regime with a stationary desired state. 3) We explicitly derive a VSL controller for a 2D conservation law that drives traffic to any space-dependent desired equilibrium: we extend the work of [START_REF] Karafyllis | Feedback control of scalar conservation laws with application to density control in freeways by means of variable speed limits[END_REF], who treated 1D LWR, to two dimensions that implies spacevarying equilibria. These are steady states that are non trivial due to space-dependency of the flux function in 2D. Compared to [START_REF] Karafyllis | Feedback control of scalar conservation laws with application to density control in freeways by means of variable speed limits[END_REF], we also consider more general fundamental diagrams: space-dependent and nonmonotonically depending on speed limits, which are more realistic (see [START_REF] Carlson | Optimal Motorway Traffic Flow Control Involving Variable Speed Limits and Ramp Metering[END_REF]). We also design a desired equilibrium, at which urban networks are used at maximal capacity by the maximal possible number of drivers.

II. 2D LWR MODEL

We use a macroscopic model described in [START_REF] Mollier | Twodimensional macroscopic model for large scale traffic networks[END_REF] that predicts traffic evolution in a urban network represented by a 2D continuum plane (x, y) ∈ Ω ⊂ R 2 that is a bounded rectangular domain, i.e., Ω = [x min , x max ] × [y min , y max ]. The size of this domain is determined by the size of urban network, i.e., x min is related to an intersection with the minimal x space coordinate among all intersections (the rest is defined similarly). This 2D model can be seen as an extension of the classical LWR model to two dimensions, and it describes traffic over a continuum plane ∀(x, y, t) ∈ Ω × R + as:

             ∂ρ(x, y, t) ∂t + ∇ • Φ(x, y, ρ(x, y, t)) = 0, Φ(x, y, t) = φ in (x, y, t) d θ (x, y), ∀(x, y) ∈ Γ in φ out (x, y, t) d θ (x, y), ∀(x, y) ∈ Γ out ρ(x, y, 0) = ρ 0 (x, y), (1)
where ρ(x, y, t) : Ω × R + → R + is a 2D density that aggregates the number of vehicles per square meter, ρ 0 (x, y) is its value at initial time (a function of bounded variation). The flux function Φ(x, y, ρ) in ( 1) is a space-dependent vector function with magnitude Φ(x, y, ρ) : E → R + , and its set of departure is

E = {(x, y, ρ) : (x, y) ∈ Ω, ρ ∈ [0, ρ max (x, y)]}.
The flux magnitude Φ(x, y, ρ) is a concave Lipschitz continuous function with the maximum φ max (x, y) ∀(x, y) ∈ Ω (capacity) achieved at the critical density ρ c (x, y), and the minimum is achieved at Φ(x, y, 0) = Φ(x, y, ρ max (x, y)) = 0 with ρ max (x, y) being the space-dependent traffic jam density. Function Φ(x, y, ρ) is known as a fundamental diagram (FD) that reflects an empirically established law relating the average flow φ with the average density ρ in a 2D plane, i.e., Φ (x, y, ρ) = φ(x, y, t). There exist a lot of shapes for FDs, see [START_REF] Li | A generic characterization of equilibrium speed-flow curves[END_REF] for a review. The most simple FD has a triangular shape suggested in [START_REF] Daganzo | The cell transmission model: A dynamic representation of highway traffic consistent with the hydrodynamic theory[END_REF], and it is defined as:

Φ(x, y, ρ) = vρ ρ ∈ [0, ρ c ], -ω(ρ -ρ max ) ρ ∈ (ρ c , ρ max ], (2)
where ρ c , ρ max , v and ω are space-dependent functions. For triangular FD, the critical density and the maximal flow are given ∀(x, y) ∈ Ω as

ρ c = ω v + ω ρ max , φ max = vω v + ω ρ max . ( 3
)
Triangular FD is characterized by having only two slopes: Φ (ρ) = v > 0 for traffic in the free-flow regime achieved for ρ(x, y) ∈ [0, ρ c (x, y)], and Φ (ρ) = -ω < 0 for the congested traffic regime when ρ(x, y) ∈ (ρ c (x, y), ρ max (x, y)]. These slopes v and -ω are kinematic wave speeds.

The flux vector function Φ in (1) is then defined as a product of its magnitude Φ (2) and direction vector d θ (unit vector):

Φ(x, y, ρ) = Φ(x, y, ρ) d θ (x, y), (4)
where

d θ = cos(θ(x, y)) sin(θ(x, y)) (5)
is a vector that depends on network geometry given by angle θ(x, y) : Ω → [0, 2π) that is a differentiable function, i.e., θ ∈ C 1 . Angle θ(x, y) is related to the orientation of roads in a network, thus, it determines the direction of traffic flow. We call d θ the direction field to stress its physical meaning.

Let us now explain the boundary conditions in [START_REF] Lighthill | On kinematic waves, II: A theory of traffic flow on long crowded roads[END_REF]. These are defined on a set Γ ⊂ Ω that is the boundary of domain Ω. It consists of two subsets Γ = Γ in ∪ Γ out . Thereby, Γ in is a set of boundary points (x, y) for which n(x, y) • d θ (x, y) > 0, where n(x, y) is a unit normal vector to the boundary oriented inside the domain. In a similar way, we also define Γ out such that ∀(x, y) ∈ Γ out : n(x, y) • d θ (x, y) < 0. Further, inflows φ in (x, y, t) and outflows φ out (x, y, t) are respectively defined ∀(x, y) ∈ Γ in and ∀(x, y) ∈ Γ out as

φ in (x, y, t) = min {D(ρ in (x, y, t)), S (ρ(x, y, t))} , φ out (x, y, t) = min {D (ρ(x, y, t)) , S (ρ out (x, y, t))} , (6)
where D(ρ) and S(ρ) are demand and supply functions:

D(ρ) = Φ(x, y, ρ), if 0 ≤ ρ ≤ ρ c (x, y), φ max (x, y), if ρ c (x, y) < ρ ≤ ρ max (x, y), (7)
S(ρ) = φ max (x, y), if 0 ≤ ρ ≤ ρ c (x, y), Φ(x, y, ρ), if ρ c (x, y) < ρ ≤ ρ max (x, y). (8)
Note that boundary control for system (1) is usually implemented by actuating D(ρ in (x, y, t)) or S(ρ out (x, y, t)). We will do this later, for now these functions should be seen as some given exterior signals.

The existence and uniqueness of solutions for multi-variable conservation laws like (1) were shown in [START_REF] Kruzhkov | First order quasilinear equations in several independent variables[END_REF] (see p.223 for the conditions of uniqueness, and existence is discussed on p.230). The boundary conditions [START_REF] Herman | A two-fluid approach to town traffic[END_REF] are given in a weak form [START_REF] Strub | Weak formulation of boundary conditions for scalar conservation laws: An application to highway traffic modelling[END_REF] that is needed to guarantee the uniqueness of (weak) solution, as will be shown later.

A continuous model like (1) requires applying some interpolation techniques to define its parameters everywhere from the corresponding values at real roads. Thus, we define ∀(x, y) ∈ Ω the direction field d θ (x, y) and the kinematic wave speed v(x, y) by Inverse Distance Weighting, see [START_REF] Tumash | Equilibrium Manifolds in 2D Fluid Traffic Models[END_REF], [START_REF] Mollier | Twodimensional macroscopic model for large scale traffic networks[END_REF]. Thereby, we assume that d θ (x, y) depends on v(x, y) and on a weighted distance to physical roads, and the weights can be tuned depending on the desired sensitivity of flux to follow precisely the location of roads. We also define the maximal density ρ max (x, y) ∀(x, y) ∈ Ω by placing vehicles on every road of a network as densily as possible (here we assume that the headway distance is 6 m). Then, density is reconstructed from vehicle positions, while treating every vehicle as a Gaussian kernel with standard deviation of 50 m centred at its position (see [START_REF] Mollier | Twodimensional macroscopic model for large scale traffic networks[END_REF]).

The main limitation of model (1): it is applicable only to networks without loops. In particular, all roads need to be uni-directional, and there exists a preferred direction of motion. Otherwise there would be no possibility to define a differentiable vector field d θ on the whole domain Ω. Indeed, if there would be a loop, then there would be a point inside of every loop where θ is undefined, since integral lines cannot cross each other (these should be seen as unique solutions to a differential equation governed by d θ ). Thus, we require that any integral line of d θ begins and ends at the boundary of the domain. Moreover, any loop would have no boundary, thus cars following this path would never be created nor destroyed.

This model can indeed predict real traffic evolution in some particular scenarios, which was validated in [START_REF] Mollier | Twodimensional macroscopic model for large scale traffic networks[END_REF], [START_REF] Tumash | Equilibrium Manifolds in 2D Fluid Traffic Models[END_REF]. It has a structure of a conservation law with state being the vehicle density. The maximum density is constructed to approximate the maximal possible number of cars in every area of the corresponding urban network, and velocity vectors are reconstructed from real roads. In the free-flow regime vehicles move freely at speeds approximating speed limits of the roads, while in the congested regime a traffic jam starts to propagate according to the approximate maximal car density. The model restriction to networks with no loops limits its applicability for general traffic although being still useful in several situations, e.g., imagine a rush hour, when many people are driving to a business district located in some particular point.

III. COORDINATE TRANSFORMATION

A. General idea

The structure of 2D LWR model [START_REF] Lighthill | On kinematic waves, II: A theory of traffic flow on long crowded roads[END_REF] implies that the direction field d θ (5) depends only on network geometry and not on state. This enables us to describe traffic flow trajectories that do not change with time. These trajectories are integral lines obtained by building tangents to d θ . In the following, we will perform a curvilinear coordinate transformation that translates these integral curves into a set of straight parallel lines as illustrated in Fig. 1. Afterwards, state evolution along a straight line can be treated as a 1D system, which is significantly easier to analyze than the 2D system (1). This idea is related to rectification of vector fields from [START_REF] Arnold | Ordinary Differential Equations[END_REF] with a difference that here we want to define a global transformation valid for traffic in the whole domain.

We introduce new coordinates (ξ, η) in a differential form:

dξ dη = C(x, y)R θ (x, y) dx dy , (9)
where R θ (x, y) is a rotation matrix given by

R θ (x, y) = cos (θ(x, y)) sin (θ(x, y)) -sin (θ(x, y)) cos (θ(x, y)) , (10)
and C(x, y) is a diagonal scaling matrix given by

C(x, y) = α(x, y) 0 0 β(x, y) , (11)
where α(x, y) and β(x, y) are positive and bounded scaling parameters needed for the existence of the coordinate transformation (see Lemma 1). In Fig. 1a) we have used the geometry of Grenoble downtown (grey arrows). The direction at each road is set such that loops are impossible. On a global scale, motion is oriented towards North-East of the city. Matrix R θ provides rotation of the integral lines (in green) in (x, y)-plane, and scaling matrix C acts such that these lines have the same metric.

B. Intuition: straight lines

In case of straight lines depicted in Fig. 1b) we would have θ = 0 ∀(x, y) ∈ Ω, which implies that rotation [START_REF] Arfken | Mathematical Methods for Physicists[END_REF] and scaling matrices [START_REF] Papageorgiou | Review of Road Traffic Control Strategies[END_REF] become identity matrices, i.e., C = R θ = I. By (9) the new coordinates (ξ, η) would completely coincide with (x, y) up to a constant shift. In this case, the direction field in (5) becomes d θ (ξ, η) = (1, 0), and by ( 4) we obtain:

Φ = Φ(ξ, η, ρ) 1 0 , (12)
which can be inserted into the divergence term in (1):

∂ ∂ξ , ∂ ∂η 1 0 Φ(ξ, η, ρ) = ∂Φ(ξ, η, ρ) ∂ξ . (13)
Thus, in case of straight lines the divergence (13) contains only one term instead of two as it was in [START_REF] Lighthill | On kinematic waves, II: A theory of traffic flow on long crowded roads[END_REF]. Flow evolves only along ξ coordinates, which are tangent to the flow motion, and there is no motion in the orthogonal direction of η, which can be treated as a parameter (a label that numbers flow path). Afterwards, we can treat each such line of constant η as a 1D equation, for which we will solve different control tasks.

C. Coordinate transformation

Let us provide the necessary and sufficient condition on the existence of this curvilinear coordinate transformation.

Lemma 1. Assume θ ∈ C 1 (Ω) and α, β ∈ C 1 (Ω). Then there exists a bijective transformation (ξ, η) in C 2 (Ω) satisfying (9)
if and only if the following PDEs hold ∀(x, y) ∈ Ω:

-sin θ ∂ (ln α) ∂x + cos θ ∂ (ln α) ∂y = cos θ ∂θ ∂x + sin θ ∂θ ∂y , (14)
cos θ ∂ (ln β) ∂x + sin θ ∂ (ln β) ∂y = sin θ ∂θ ∂x -cos θ ∂θ ∂y . (15)
Proof. For any function in C 2 mixed partial derivatives must be equal by the Schwarz theorem. In our case this is equivalent to the invariance in the order of taking partial derivatives of ξ and η w.r.t. x and y, i.e.,

∂ ∂y ∂ξ(x, y) ∂x = ∂ ∂x ∂ξ(x, y) ∂y , (16)
and

∂ ∂y ∂η(x, y) ∂x = ∂ ∂x ∂η(x, y) ∂y . ( 17
)
By applying ( 16) and ( 17) to ( 9) using the definitions ( 10) and ( 11) we obtain ( 14)- [START_REF] Gugat | Optimal control for traffic flow networks[END_REF]. Finally, ξ and η can be obtained by integrating [START_REF] Daganzo | The cell transmission model: A dynamic representation of highway traffic consistent with the hydrodynamic theory[END_REF]. Bijectivity follows since the determinant of the Jacobian ( 9) is given by α(x, y)β(x, y), and by ( 14)-( 15) both α(x, y) and β(x, y) are strictly positive.

Note that α(x, y) and β(x, y) being functions of d θ (x, y) only can be computed from the network geometry.

D. Model in (ξ, η)-space

Let us now perform this coordinate transformation to the original system (1) to turn it into a parametrized set of 1D LWR equations. According to Chapter 2 of [START_REF] Arfken | Mathematical Methods for Physicists[END_REF], we can apply the divergence formula to calculate ∇ • Φ in (ξ, η)-space:

∇ • Φ(ξ, η, ρ) = 1 h ξ h η   ∂ Φ ξ h η ∂ξ + ∂ Φ η h ξ ∂η   , (18)
where h ξ and h η are known as Lamé coefficients, which correspond to the lengths of basis vectors in (ξ, η)-space:

h ξ = ∂x ∂ξ , ∂y ∂ξ T and h η = ∂x ∂η , ∂y ∂η T . (19)
A detailed calculation of ∇ • Φ(ξ, η, ρ) is given in [START_REF] Tumash | Equilibrium Manifolds in 2D Fluid Traffic Models[END_REF], thus, we will directly state the result that reads:

∇ • Φ(ξ, η, ρ) = α(ξ, η)β(ξ, η) ∂(Φ(ξ, η, ρ)/β) ∂ξ . (20)
Thus, we have shown that the curvilinear coordinate transformation ( 9) makes the divergence operator uni-dimensional. This means that the temporal change of density in a 2D plane is caused only by the change of flow along η-lines, as we were showing by [START_REF] Daganzo | A variational formulation of kinematic waves: basic theory and complex boundary conditions[END_REF] for the case of straight lines. We rescale all density-, flow-and velocity-related functions:

ρ = ρ αβ , Φ = Φ β , v = αv. (21)
Finally, we define a spatial domain, on which the system in new coordinates will evolve as:

Ω = {(ξ, η) : ∃ (x, y) ∈ Ω, ξ = ξ(x, y), η = η(x, y)} .
Then, the domain boundary in (x, y)-space can be uniquely projected into the boundary in (ξ, η)-space, i.e., Γ Ω → ΓΩ. In particular, ΓΩ consists of points (ξ min (η), ξ max (η)) such that

ξ min (η) = min (x,y)∈Ω, η(x,y)=η ξ(x, y), ξ max (η) = max (x,y)∈Ω, η(x,y)=η ξ(x, y).
and we can also define maximal and minimal values of η as

η min = min{η | ∃ξ : (ξ, η) ∈ Ω}, η max = max{η | ∃ξ : (ξ, η) ∈ Ω}.
Using the divergence term in (ξ, η)-space (20), we can now rewrite the 2D model (1) that reads ∀(ξ, η, t) ∈ Ω × R + :

             ∂ ρ(ξ, η, t) ∂t + ∂ Φ(ξ, η, ρ) ∂ξ = 0, φin (η, t) = min D (ρ in (η, t)) , S (ρ (ξ min (η), η, t)) , φout (η, t) = min D (ρ (ξ max (η), η, t)) , S (ρ out (η, t)) , ρ(ξ, η, 0) = ρ0 (ξ, η), (22)
where Φ(ξ, η, ρ) is now a scalar function that preserves all the FD properties such as being Lipschitz continuous and concave. Traffic flow evolves only along lines of constant η in (ξ, η)space. Thus, ( 22) is a continuous set of inhomogeneous 1D LWR equations each following a path parametrized by η. This means that one can also analyze its solution in the same way as in case of 1D LWR. Notice that boundary conditions φin (η, t) and φout (η, t) are formulated using the demandsupply concept (weak formulation). Thus, the initial boundary value problem ( 22) is well-posed, see more details in [START_REF] Bürger | Difference schemes, entropy solutions, and speedup impulse for an inhomogeneous kinematic traffic flow model[END_REF], [START_REF] Wu | Steadystate traffic flow on a ring road with up-and down-slopes[END_REF] for entropy conditions for inhomogeneous LWR model. In the following, we will always refer to the system in (ξ, η)space [START_REF] Geroliminis | Existence of urban-scale macroscopic fundamental diagrams: Some experimental findings[END_REF]. For simplicity of notations, bars are omitted for the rest of the paper, it is left only for domain Ω.

IV. BOUNDARY CONTROL

The system in new coordinates [START_REF] Geroliminis | Existence of urban-scale macroscopic fundamental diagrams: Some experimental findings[END_REF] offers a variety of possibilities for explicit control design for urban networks by analyzing only its structure. Let us now formulate the following boundary control problem as an application example.

A. Problem Statement

Boundary control for the system in (ξ, η)-space [START_REF] Geroliminis | Existence of urban-scale macroscopic fundamental diagrams: Some experimental findings[END_REF] implies setting actuators at entry and exit of η-lines, along which the solution propagates. Thus, we control traffic by modifying demand and supply functions, i.e., ∀(η, t)

∈ [η min , η max ]×R + u in (η, t) = D (ρ in (η, t)) , u out (η, t) = S (ρ out (η, t)) . (23) Problem 1. Design ∀(η, t) ∈ [η min , η max ] × R + boundary controllers u in (η, t
) and u out (η, t) such that density ρ(ξ, η, t) governed by [START_REF] Geroliminis | Existence of urban-scale macroscopic fundamental diagrams: Some experimental findings[END_REF] tracks a desired trajectory as t → ∞.

In [START_REF] Tumash | Boundary Control Design for Traffic with Nonlinear Dynamics[END_REF] a similar problem was posed for one road with a homogeneous FD. Here we extend this result by solving the tracking problem for a large urban area whose infrastructure is captured by an explicit space-dependency of FD.

B. Assumptions

Let us introduce capacity of the strongest bottleneck along the η-line ∀η ∈ [η min , η max ] as

φ min max (η) = min ξ∈[ξmin(η),ξmax(η)] φ max (ξ, η). (24)
To solve Problem 1, we need to assume the following:

Assumption 1. Inflows φ in (η, t) and outflows φ out (η, t) from (22) must satisfy ∀(η, t) ∈ [η min , η max ] × R + φ in (η, t) ≤ φ min max (η), φ out (η, t) ≤ φ min max (η), (25)
where φ min max (η) is defined in [START_REF] Garavello | Conservation laws on complex networks[END_REF]. Moreover, there exists ε > 0 such that φ in (η, t) and φ out (η, t) additionally satisfy

t+tc(η) t φ in (η, τ )dτ ≤ t c (η)φ min max (η) -ε and t+tc(η) t φout (η, τ )dτ ≤ t c (η)φ min max (η) -ε, (26)
where t c (η) is time needed for a solution evolving from one end of η-line to reach the opposite end:

t c (η) = min      ξmax(η) ξmin(η) 1 v( ξ) d ξ, ξmax(η) ξmin(η) 1 ω( ξ) d ξ     . ( 27
)
It means that inflows and outflows at each η-line are not allowed to pass the capacity of the strongest bottleneck of the corresponding line instantly (see [START_REF] Carlson | Optimal Motorway Traffic Flow Control Involving Variable Speed Limits and Ramp Metering[END_REF]). Moreover, they must be strictly lower during the time interval given by t c (η) (see [START_REF] Claudel | Lax-Hopf Based Incorporation of Internal Boundary Conditions Into Hamilton-Jacobi Equation. Part I: Theory[END_REF]), i.e., a road is never filled with the maximal number of vehicles, which gives more possibilities for control.

Assumption 2. Solution of ( 22) is determined only by boundary data, i.e., the influence of initial conditions left the system.

Influence of initial condition is meant in sense of episolution evolving from initial data. For example, if a road gets blocked at a boundary, then cars from initial condition propagate until the boundary and become a part of boundary data.

Remark 1.

If Assumption 1 is satisfied, then by taking t ≥ t min , where t min is derived in (84), Assumption 2 holds trivially, as it is shown in Appendix I-E.

C. Boundary Control Design

We are now going to solve Problem 1, where the goal is formulated in terms of tracking a desired trajectory ρ d . This is a time-and space-dependent function of bounded variation ρ d (ξ, η, t) governed by the same system ( 22) ∀(ξ, η, t) ∈ Ω × R + , and its maximal value is also bounded by ρ max that is given by network geometry. Boundary flows φ in d (η, t) and φ out d (η, t) correspond to weak boundary conditions in [START_REF] Geroliminis | Existence of urban-scale macroscopic fundamental diagrams: Some experimental findings[END_REF].

Theorem 1. Consider density ρ(ξ, η, t) governed by system (22) ∀(ξ, η, t) ∈ Ω×R + , for which Assumptions 1 and 2 hold. The desired density ρ d (ξ, η, t) is also governed by [START_REF] Geroliminis | Existence of urban-scale macroscopic fundamental diagrams: Some experimental findings[END_REF], for which Assumption 2 holds. Then if

∀(η, t) ∈ [η min , η max ] × R + the controls (23) are set to (1) u in (η, t) = φ in d (η, t) -ke(η, t), (2) u out (η, t) = φ out d (η, t) + ke(η, t), with k > 0 and e(η, t) = ξmax(η) ξmin(η) ρ( ξ, η, t) -ρ d ( ξ, η, t) d ξ, ( 28
) then ∀a, b: ξ min (η) ≤ a < b ≤ ξ max (η) we obtain ∀η lim t→∞ b a ρ( ξ, η, t) -ρ d ( ξ, η, t) d ξ = 0.
Remark 2. Here we consider vehicle density governed by system [START_REF] Geroliminis | Existence of urban-scale macroscopic fundamental diagrams: Some experimental findings[END_REF], which is analyzed in Hamilton-Jacobi formulation (H-J) given in Appendix I. Thereby, traffic is described in terms of cumulative vehicle number rather than flow and density. The solution of [START_REF] Geroliminis | Existence of urban-scale macroscopic fundamental diagrams: Some experimental findings[END_REF] with triangular FD (2) is obtained explicitly in H-J formulation, and it is given by (85). Its step-by-step derivation being quite technical is presented in Appendix I. It is then used to analyze properties of system [START_REF] Geroliminis | Existence of urban-scale macroscopic fundamental diagrams: Some experimental findings[END_REF] in a cumulative (integral) form. In particular, it enables treatment of weak boundary conditions that are given by the minimum between demand and supply in [START_REF] Geroliminis | Existence of urban-scale macroscopic fundamental diagrams: Some experimental findings[END_REF]. These imply that control [START_REF] Geroliminis | Properties of a well-defined macroscopic fundamental diagram for urban traffic[END_REF] may not be accepted at boundaries at any time. Thus, (85) is used to analyze the time during which controls are not accepted by the system in terms of control restriction functions, similarly as it was done in [START_REF] Tumash | Boundary Control Design for Traffic with Nonlinear Dynamics[END_REF]. Remark 3. Control gain k determines the convergence rate towards the desired trajectory, i.e., the larger it gets the faster is convergence. However, due to weak boundary conditions, for high k values the convergence rate remains almost constant.

Remark 4. Note that the integral convergence of densities stated in Theorem 1 implies that the state approximates the desired trajectory as time goes to infinity, since a and b can be arbitrarily close in space, i.e., ρ ≈ ρ d as t → ∞.

Proof of Theorem 1. The proof is quite technical and would take several pages. It follows the same major steps as in 1D case in [START_REF] Tumash | Boundary Control Design for Traffic with Nonlinear Dynamics[END_REF], and it shows that the convergence is exponential for controller [START_REF] Geroliminis | Properties of a well-defined macroscopic fundamental diagram for urban traffic[END_REF]. However, since we deal with a multidimensional inhomogeneous system, we list the differences that need to be taken into account.

1) Length L of a 1D road varies as a function of line number

η, i.e., [0, L] → [ξ min (η), ξ max (η)]. This implies that L v → T v (ξ max (η), η) and L ω → T ω (ξ min (η), η) where T v (ξ max (η), η) should be taken from (70) for ξ = ξ max (η) and T ω (ξ min (η), η) from (74) for ξ = ξ min (η).
2) Every occurrence of Lρ max should be substituted by

ξmax(η) ξmin(η) ρ max ( ξ, η)d ξ.
3) Equation ( 28) of [START_REF] Tumash | Boundary Control Design for Traffic with Nonlinear Dynamics[END_REF] should be rewritten as:

g in (η, t) = 0 ⇒ ∀t ∈ [t -Tω(ξ min (η), η), t] : R(η, t ) ≥ ξmax(η) ξ min (η) ρc( ξ, η)d ξ, g out (η, t) = 0 ⇒ ∀t ∈ [t -Tv(ξmax(η), η), t] : R(η, t ) ≤ ξmax(η) ξ min (η) ρc( ξ, η)d ξ. (29)
We obtain (29) by using the following upper bound:

t t-Tω (ξ min (η),η) φ out (η, τ )dτ + t t φ in (η, τ )dτ ≤ Tω(ξ min (η), η)φ min max (η) ≤ ξmax(η) ξ min (η) φmax( ξ, η) ω( ξ, η) d ξ = ξmax(η) ξ min (η) ρmax( ξ, η) -ρc( ξ, η) d ξ.

D. Numerical Example

The efficiency of boundary controller ( 28) is demonstrated on a numerical example, where the controller is used to track a desired density profile that is space-dependent and periodic in time. The controller is set to boundaries of a urban area that has the same structure as an area of Grenoble downtown of size 1.4 × 1 km 2 (the selected area is shown in grey in Fig. 3). Road directions were (numerically) modified such that there is a preferred direction of motion in this area (here towards North-East). Speed limits are set as in the selected area: some roads are 30 km/h and other roads are 50 km/h.

We define a numerical grid in Ω × R + and deploy the Godunov scheme in 2D for [START_REF] Geroliminis | Existence of urban-scale macroscopic fundamental diagrams: Some experimental findings[END_REF]. First, discretize the η dimension into m = 180 cells. Then, we use the 2D Godunov scheme for every j ∈ {1, . . . , m} with a discretization step ∆ξ = 5 m (space cell size in ξ dimension). We also set the time cell size ∆t = 0.1 s, which provides that the CFL condition is satisfied. In order to compute the integral in [START_REF] Geroliminis | Properties of a well-defined macroscopic fundamental diagram for urban traffic[END_REF] we perform the Riemann summation for every j ∈ {1, . . . , m} over all ξ cells, i.e., i ∈ {1, . . . , n j }, where n j is the number of ξ cells contained in each cell j.

For triangular FD (2) we use ρ c = ρ max /3. The system to be controlled is initially given as a traffic jam, see Fig. 2a):

ρ 0 (ξ, η) = ρ max (ξ, η), ∀(ξ, η) ∈ Ω.
We set the desired inflow demand D(ρ in d )(η, t) and outflow supply S(ρ out d )(η, t) to be time-periodic functions:

D(ρ in d )(η, t) =φ min max (η) 0.6+ 0.4 sin 2π t 1200 + 2 η -η min η max -η min , S(ρ out d )(η, t) =φ min max (η) 0.6+ 0.4 sin 2π t 2400 + 2 η -η min η max -η min .
Hence, these boundary flow functions will never exceed the minimal capacity on the corresponding η-line. They were chosen such to generate a desired trajectory ρ d with a period of τ = 2400 seconds, as drawn on the right column in Fig. 3. We demonstrate here, how the boundary control law from Theorem 1 enhances traffic state if there is a feedback, i.e., k > 0 in [START_REF] Geroliminis | Properties of a well-defined macroscopic fundamental diagram for urban traffic[END_REF]. Control is applied at domain boundaries, and it physically corresponds to demand at the entry and supply of the exit, as illustrated in Fig. 2a). Fig. 3 shows evolution of traffic density within the time interval of 2τ = 4800 sec. The top row is related to desired trajectory with the boundary data described above. The middle and bottom row correspond to density evolution for k = 5 • 10 -5 and k = 0, respectively. We observe convergence to the desired profile for the case with feedback (middle row) that becomes visible already at t = 2τ , while this does not happen for the case without feedback (bottom row).

In Fig. 2b) the L 1 error norm is depicted as a function of time for different control gains. It can be computed as follows: We see that a higher control gain k = 10 -3 provides a higher convergence speed in comparison to k = 5 • 10 -5 . Notice that we do not obtain higher convergence speed than for k = 10 -3 due to weak boundary conditions (control is sometimes not imposed). In general, we achieve exponential convergence with control [START_REF] Li | A generic characterization of equilibrium speed-flow curves[END_REF]. On the contrary, k = 0 does not work even if we start from an empty network unless the initial data is the same as in the desired profile (hardly ever possible).

ρ -ρ d 1 = ηmax ηmin ξmax(η) ξmin(η) |ρ(ξ, η, t) -ρ d (ξ, η, t)|dξdη. (30)

V. VARIABLE SPEED LIMIT CONTROL

Let us now demonstrate how to solve control tasks using variable speed limit in a 2D-plane by stating a new problem in (ξ, η)-space. We consider the following initial-boundary value problem given ∀(ξ, η, t)

∈ Ω × R + as ∂ρ(ξ, η, t) ∂t + ∂Φ(ξ, η, ρ(ξ, η, t), u) ∂ξ = 0, (31)
where the initial and boundary data are the same as in [START_REF] Geroliminis | Existence of urban-scale macroscopic fundamental diagrams: Some experimental findings[END_REF]. The only difference is that the flux function Φ now depends also on a control parameter u ∈ [0, 1] that represents a variable speed limit ratio. Applying VSL should be understood as a flexible restriction on speed at which vehicles can drive on a given road stretch. Its value varies according to the current environmental and road conditions and is displayed on electronic traffic signs. Setting u = 1 implies that vehicles can drive at speeds bounded by the legal maximum (e.g., 130 km/h on French highways). Moreover, if u = 0, then no movement is allowed and Φ(ξ, η, ρ, 0) = 0. One should see u as an indomain controller that affects the flux. It is applied in the whole domain including its boundaries. Note that Φ is still a concave function wrt ρ, and Φ is continuous in u.

A. Contribution

The material of this chapter was inspired by a previous work [START_REF] Karafyllis | Feedback control of scalar conservation laws with application to density control in freeways by means of variable speed limits[END_REF]. However, there are some major points that were not considered in [START_REF] Karafyllis | Feedback control of scalar conservation laws with application to density control in freeways by means of variable speed limits[END_REF], and thus will be addressed here:

1) 2D systems: this is the first time that VSL control is applied on a large-scale network directly using the intrinsic properties of the model only (no discretization). 2) Space-dependent FD: we extend the result of [START_REF] Karafyllis | Feedback control of scalar conservation laws with application to density control in freeways by means of variable speed limits[END_REF] by considering space-dependent FD that imply obtaining non trivial desired equilibria (with space-dependency). This is an essential extension, since space-dependent equilibria necessarily arise in urban networks that in general do not have a homogeneous structure. 3) Realistic FDs: in [START_REF] Karafyllis | Feedback control of scalar conservation laws with application to density control in freeways by means of variable speed limits[END_REF] ∂Φ(ξ, η, ρ, u)/∂u > 0 holds, see Fig. 4a). This assumption was made for simplicity to avoid multi-valued functions (always only one value of u for each flow φ). Here we omit this condition by allowing more general forms of FD. In general, applying speed limits (u < 1) can shift the critical density towards larger values in realistic FDs. This is schematically depicted in Fig. 4b), see red FD for u = 0.7 and green FD for u = 0.5 and compare ρ c3 and ρ c2 with ρ c1 achieved with u = 1. This means that applying speed limits can increase the range of vehicle density, for which the free-flow regime is preserved. There it is also shown how VSL can enhance flow for some densities in the congested regime, e.g., compare flows φ 2 with φ 1 that can be achieved with different speed limits for the same density ρ c2 . These VSL effects on FD were validated by real data collected from European VSLequipped motorway, see [START_REF] Carlson | Optimal Motorway Traffic Flow Control Involving Variable Speed Limits and Ramp Metering[END_REF]. 4) Study the controller smoothness: considering a general class of fundamental diagrams may lead to irregular control policies. We investigate whether any conditions must be imposed on the functional dependence of FD on VSL in order to provide smoothness.

B. Problem Statement

Let us also introduce the following notations: ,

and now we can formulate the stabilization problem as follows.

Problem 2. Given ∀(ξ, η) ∈ Ω the fundamental diagram Φ(ξ, η, ρ, u) and initial density ρ 0 (ξ, η) with dynamics governed by [START_REF] Wu | Steadystate traffic flow on a ring road with up-and down-slopes[END_REF], find a VSL controller u such that

lim t→∞ ρ(ξ, η, t) = 0, ∀(ξ, η) ∈ Ω, (32)
where ρ = ρ -ρ d is the deviation from a desired equilibrium

ρ d (ξ, η) ∈ (0, ρ max (ξ, η)). ρ Φ(ρ) a) 0 ρ c ρ max 0 ρ Φ(ρ) b) ρ c1 ρ c2 ρ c3 ρ max f (ρ) φ 1 φ 2
Fig. 4. FDs and their dependence on VSL: a) monotonic dependence ∂Φ(ξ, η, ρ, u)/∂u > 0 used in [START_REF] Karafyllis | Feedback control of scalar conservation laws with application to density control in freeways by means of variable speed limits[END_REF]; b) dependence we assume here, i.e., possible increase of ρc for stronger speed limits (from real data, see [START_REF] Carlson | Optimal Motorway Traffic Flow Control Involving Variable Speed Limits and Ramp Metering[END_REF]). Blue line: u = 1. Red line: u = 0.7. Green line: u = 0.5. Bold dashed line: maximal flow function defined in [START_REF] Delle Monache | Traffic Regulation via Controlled Speed Limit[END_REF].

C. VSL Control Design

Let us define ∀(ξ, η) ∈ Ω and ∀ρ ∈ [0, ρ max (ξ, η)] a maximal flow function f (ξ, η, ρ), which is the maximum possible flow that can be achieved at a given point for a given density over all the VSL values (see the thick dashed line in Fig. 4):

f (ξ, η, ρ) = max u∈[0,1] Φ(ξ, η, ρ, u). (33)
We also introduce a multi-valued function G(ξ, η, ρ, φ), which is the inverse image of FD with respect to speed limit:

G(ξ, η, ρ, φ) = {u ∈ (0, 1] : Φ(ξ, η, ρ, u) = φ}. (34)
In general, it is possible that several values of speed limits u provide the same flow value, see the black dot in Fig. 4. Therefore, G(ξ, η, ρ, φ) for a fixed set of parameters represents a set, not a single value.

Theorem 2. Assume ρ d ∈ C 1 ( Ω). Let the controller u(ξ, η, ρ) be given ∀(ξ, η) ∈ Ω and for ρ = ρ(ξ, η, t) by the following inclusion u(ξ, η, ρ) ∈ G(ξ, η, ρ, φ d (ξ, η, ρ)), with

φ d (ξ, η, ρ) = B(ξ, η, ρ) min ξ f (ξ , η, ρ(ξ , η, t)) B(ξ , η, ρ) and B(ξ, η, ρ) = 1 + γ ξ ξmin(η) ρ( ξ, η, t)d ξ, (35)
where control gain γ is a positive constant defined as

0 < γ < min η    ξmax(η) ξmin(η) ρ max ( ξ, η)d ξ   -1
.

Then there exists c = c(γ, ρ 0 ) > 0 such that ∀ρ 0 ∈ C 1 ( Ω) the system (31) with initial condition ρ(ξ, η, 0) = ρ 0 (ξ, η) has a unique solution ρ ∈ C 1 ( Ω × R + ) which satisfies

max (ξ,η)∈ Ω |ρ(ξ, η, t)| ≤ e -ct max (ξ,η)∈ Ω |ρ(ξ, η, 0)| ∀t ∈ R + , (36)
and, moreover,

∀(ξ, η) ∈ Ω lim t→∞ Φ(ξ, η, ρ(ξ, η, t), u(ξ, η, ρ)) = min ξ f (ξ , η, ρ d (ξ , η)). ( 37
)
Remark 5. Note that u(ξ, η, ρ) depends on state (feedback).

Let us give several comments on its structure:

1) Such a controller choice [START_REF] Mollier | Twodimensional macroscopic model for large scale traffic networks[END_REF] assures that the system flow is immediately set to desired flow φ d , i.e.

Φ (ξ, η, ρ, u) = φ d (ξ, η, ρ), ∀(ξ, η, t) ∈ Ω × R + .
2) The desired flow φ d (ξ, η, ρ) is designed such that it does not exceed the maximal flow function f (ξ, η, ρ) at any point of space-time. The space-dependency of the desired flow is incorporated into function B(ξ, η, ρ). 3) Function B(ξ, η, ρ) acts as a feedback linearization for system [START_REF] Wu | Steadystate traffic flow on a ring road with up-and down-slopes[END_REF], which loses the conservation law structure, and we do not have to handle discontinuities in the solution. This will be shown in the proof of Theorem 2. 4) The lower and upper bound on control gain γ are set such to guarantee that B(ξ, η, ρ) is positive, i.e., B : Ω × R + → R + . The upper bound on γ is required for situations when the density error ρ is negative, which can appear since we design a general controller that drives any state to any desired equilibrium.

Proof of Theorem 2. First of all, we need to prove that controller ( 35) is well-defined. Namely, we will show that the set G(ξ, η, ρ, φ d (ξ, η, ρ)) is not empty, i.e., the desired flow takes values in a bounded range that can be achieved by the VSL control. Indeed, for all (ξ, η) ∈ Ω we obtain from (35) that

φ d (ξ, η, ρ) B(ξ, η, ρ) = min ξ f (ξ , η, ρ) B(ξ , η, ρ) ≤ f (ξ, η, ρ) B(ξ, η, ρ) , (38)
and by positivity of B(ξ, η, ρ) we get φ d (ξ, η, ρ) ∈ [0, f (ξ, η, ρ)] ∀(ξ, η, t) ∈ Ω × R + . This interval exactly corresponds to the range of function Φ(ξ, η, ρ, u(ξ, η, ρ)) w.r.t. u, therefore the set function G(ξ, η, ρ, φ d (ξ, η, ρ)) is not empty.

Then, we substitute the constructed flux function

Φ(ξ, η, ρ, u(ξ, η, ρ)) = B(ξ, η, ρ) min ξ f (ξ , η, ρ) B(ξ, η, ρ) (39)
into system (31) and obtain

∂ ρ(ξ, η, t) ∂t + min ξ f (ξ , η, ρ) B(ξ , η, ρ) ∂B(ξ, η, ρ) ∂ξ = 0.
By inserting the definition of B(ξ, η, ρ) from ( 35), this equation can be further simplified as

∂ ρ(ξ, η, t) ∂t = -γ ρ(ξ, η, t) min ξ f (ξ , η, ρ) B(ξ , η, ρ) . (40)
Equation [START_REF] Tumash | Boundary Control Design for Traffic with Nonlinear Dynamics[END_REF] does not contain any partial space derivatives, and thus the conservation law structure is lost in the closedloop system. This dynamic equation has a stable equilibrium at zero. By [START_REF] Karafyllis | Feedback control of scalar conservation laws with application to density control in freeways by means of variable speed limits[END_REF] we obtain smoothness and an exponential convergence to a desired equilibrium with rate c(γ, ρ 0 ). Finally, we see that density convergence ∀(ξ, η) ∈ Ω ρ(ξ, η, t) → ρ d (ξ, η) as t → +∞ implies that B(ξ, η, ρ) → 1, and thus (39) results into

Φ(ξ, η, ρ(ξ, η, t), u(ξ, η, ρ)) → min ξ f (ξ , η, ρ d (ξ , η)),
which coincides with [START_REF] Aghamohammadi | Dynamic traffic assignment using the macroscopic fundamental diagram: A Review of vehicular and pedestrian flow models[END_REF], and thus concludes the proof. Remark 6. Property [START_REF] Aghamohammadi | Dynamic traffic assignment using the macroscopic fundamental diagram: A Review of vehicular and pedestrian flow models[END_REF] means that the highest possible equilibrium flow is achieved for a given ρ d (ξ, η). With definition [START_REF] Delle Monache | Traffic Regulation via Controlled Speed Limit[END_REF], the following double inequality holds ∀η ∈ [η min , η max ]

min ξ Φ(ξ, η, ρ d , 1) ≤ min ξ f (ξ, η, ρ d ) ≤ φ min max (η), (41)
where the left inequality implies that the same or higher traffic flow can be achieved with lower speed limits than for u = 1.

D. Smoothness of VSL Controller

VSL controller ( 35) is defined via inclusion, and in general it can result in a discontinuous function in space. For example, imagine there are two different speed limits able to provide the desired flow. In this case, one should fear that speed limits jump from one value to another along a road infinitely many times. However, if we assume additional properties on how the flux function should depend on speed limit, we will obtain that u(ξ, η, ρ) is differentiable almost everywhere. Theorem 3. Assume that ∀(ξ, η) ∈ Ω, ∀ρ ∈ [0, ρmax (ξ, η)] and ∀u ∈ [0, 1] flux Φ(ξ, η, ρ, u) is differentiable. Moreover, assume that it is either strictly concave in u (congested regime) or monotonic in u and reaches its maximum at u = 1 (freeflow regime). Then using controller provided in Theorem 2 and assuming ρ ∈ C 1 ( Ω) ∀t ∈ R + , we can choose u(ξ, η, ρ) such that it is differentiable almost everywhere w.r.t. ξ.

Remark 7. This additional assumption on functional dependence of Φ(ξ, η, ρ, u) on u can be interpreted as follows. In the congested regime when speed limit decreases, the flow can first increase for a fixed density as in Fig. 5b), and then it drops to zero as the speed limit approaches zero. In the free-flow regime, flow is maximal for u = 1 and decreases monotonically as u decreases as in Fig. 5a).

Remark 8. Notice that by Theorem 2, density is a differentiable function ρ ∈ C 1 ( Ω) ∀t ∈ R + if the initial condition function of system (31) is differentiable, i.e., ρ 0 ∈ C 1 ( Ω).

Proof of Theorem 3. Let us fix time t and line η, and consider an interval of all possible ξ values and split it in two subsets H 1 and H 2 as [ξ min (η), ξ max (η)] = H 1 ∪ H 2 , where

H 1 = ξ ∈ [ξ min (η), ξ max (η)] ∂Φ(ξ, ρ(ξ), u(ξ, ρ)) ∂u = 0 , H 2 = ξ ∈ [ξ min (η), ξ max (η)] ∂Φ(ξ, ρ(ξ), u(ξ, ρ)) ∂u = 0 .
We introduce also interiors of H 1 and H 2 as:

E 1 = int(H 1 ), E 2 = int(H 2 ).
Moreover, we introduce a complementary subset E 0 as

E 0 = (H 1 \ E 1 ) ∪ (H 2 \ E 2 ), such that E 0 ∪ E 1 ∪ E 2 = H 1 ∪ H 2 = [ξ min (η), ξ max (η)].
It is clear that sets E 1 and E 2 have the same Lebesgue measure as sets H 1 and H 2 , respectively. This implies that set E 0 is of measure zero. Thus, showing that controller u(ξ) = u(ξ, ρ(ξ)) is differentiable on E 1 and E 2 would imply that it is differentiable almost everywhere. Let us first consider set E 1 with the following function defined from (39):

F 1 (ξ, u) = Φ(ξ, ρ(ξ), u) -min ξ f (ξ , ρ(ξ )) B(ξ , ρ(ξ )) B(ξ, ρ(ξ)).
This function is differentiable by the assumptions made in Theorem 3 and is equal to zero by [START_REF] Tumash | Topologybased control design for congested areas in urban networks[END_REF]. Moreover, the derivative of Φ(ξ, ρ(ξ), u) w.r.t. u is non-zero on set E 1 by its definition. This immediately implies that the derivative of F 1 (ξ, u) with respect to u is also non-zero. Therefore, we can use the Implicit Function Theorem, which assures that there exists a differentiable function u(ξ) on this set satisfying [START_REF] Tumash | Topologybased control design for congested areas in urban networks[END_REF]. On the second set E 2 we define another function as

F 2 (ξ, u) = ∂Φ(ξ, ρ(ξ), u) ∂u .
Notice that F 2 (ξ, u) is zero by the definition of set E 2 , and it has a negative derivative w.r.t. u, since we assumed concavity of the flux function for the congested traffic regime (in a pure free-flow regime set E 2 would be empty). This means that we can use the Implicit Function Theorem again, thus a differentiable function u(ξ) exists on set E 2 as well.

Combining these results, we obtain that function u(ξ) is differentiable on E 1 ∪ E 2 , i.e., almost everywhere. Proposition 1. In case of concave dependence of FD on speed limits, u(ξ, η, ρ) can sometimes be chosen from two values G(ξ, η, ρ, φ d ) for ρ being in congested regime, see Fig. 5 b). Then the most appropriate choice is the minimal value, since it provides the free-flow traffic regime:

u(ξ, η, ρ) := min{G(ξ, η, ρ, φ d )}.
As an example, consider the intersection point (black dot in Fig. 4b)) corresponding to a flow-density pair that can be achieved using either u = 1 or u = 0.7. In this case we choose u = 0.7, since this provides the free-flow regime and, thus, a more smooth traffic motion.

E. Parametrization of Fundamental Diagram

Before applying the designed VSL controller [START_REF] Mollier | Twodimensional macroscopic model for large scale traffic networks[END_REF] in practice (or in our case it will be a numerical example), we should first discuss flux functions depending on u by suggesting an explicit relation satisfying assumptions made in Theorem 3.

Let us consider triangular FD as in [START_REF] Richards | Shock waves on the highway[END_REF], which should be modified due to the dependence on speed limits. We assume a linear dependence of kinematic wave speeds on speed limits:

v(ξ, η, u) = u v 1 (ξ, η), ω(ξ, η, u) = ω 1 (ξ, η) + (1 -u)ω add (ξ, η), ( 42
)
where v 1 (ξ, η) and ω 1 (ξ, η) are kinematic wave speeds for u = 1, and ω add (ξ, η) is an additional value expressing the effect of speed limit on the kinematic wave speed in the congested regime. Thus, if speed limits are high (u 1), drivers are moving slowly, and therefore start braking late (larger safety distance for lower speeds). Let us estimate a range of reasonable values for ω add (ξ, η) such that ∀(ξ, η) ∈ Ω

∂φ max (ξ, η, u) ∂u ≥ 0. ( 43
)
Condition ( 43) means that it is not possible to enhance the transportation capacity by applying speed limits, see (41).

Transportation capacity is a property of network geometry, i.e., φ max is determined by the number of lanes and free-flow kinematic wave speed, and thus it should not be changed with a VSL. In the following we skip the dependence on (ξ, η) in notations. We insert ω(u) and v(u) from ( 42) into (3) and get

φ max (u) = v 1 ρ max u (ω 1 + (1 -u)ω add ) ω 1 + v 1 u + (1 -u)ω add . ( 44
)
We take the partial derivative of (44) w.r.t. u and obtain

∂φmax(u) ∂u = v 1 ρmax (ω 1 + (1 -u) ω add ) 2 -u 2 v 1 ω add (ω 1 + v 1 u + (1 -u)ω add ) 2 . ( 45
)
In accordance with (43), we need to find such range of ω add that (45) is positive. We distinguish two different cases for which nominator of (45) takes non-negative values ∀u ∈ [0, 1]: 1) ω add ≤ 0: then ∂φ max (u)/∂u > 0 holds always.

2) ω add > 0: then we must provide that

ω 1 + (1 -u) ω add ≥ u √ ω add v 1 ⇒ ω 1 + ω add ≥ u (ω add + √ ω add v 1 ) .
In the worst case this inequality must be satisfied for u = 1, which results into ω add ≤ ω 2 1 /v 1 (upper bound for ω add ). By definition (42) and the fact that ω(u) should be non-negative, the lower bound is -ω 1 . Thus, a reasonable range reads

ω add ∈ -ω 1 , ω 2 1 v 1 .
For a numerical example, we will pick the largest value ω add = ω 2 1 /v 1 , since by ( 45) it provides ∂φ max (u)/∂u = 0 at u = 1. This choice implies the biggest influence of VSL on FD in congested regime (the largest surface enclosed by the blue line in congested regime and the thick dashed line in Fig. 4).

F. Optimal Equilibrium

The controller given by ( 35) can be applied to achieve any type of desired equilibrium ρ d (ξ, η) ∈ (0, ρ max (ξ, η)) ∀(ξ, η) ∈ Ω. However, for the upcoming numerical example, we seek to design an optimal equilibrium ρ opt d that corresponds to throughput maximization and, at the same time, to density maximization, i.e., the highest possible number of cars should be able to pass a network at maximal flow. Thereby, the number of cars in a urban area is directly related to its density in this area that can be increased due to a change in the shape of FD caused by u(ξ, η, ρ), as it is shown in Fig. 6.

A method to find equilibrium profiles providing the maximal flow in the system was presented in [START_REF] Tumash | Topologybased control design for congested areas in urban networks[END_REF]. However, there it was done for u = 1 (no speed limits). With the help of speed limits, we are now able to extend the result of [START_REF] Tumash | Topologybased control design for congested areas in urban networks[END_REF] by maximizing also the number of vehicles. In particular, we seek to find ∀(ξ, η) ∈ Ω speed limits u opt (ξ, η) such that φ max (ξ, η, u opt ) = φ min max (η, 1). The VSL controller must provide that this steady state flow is achieved, and at the same time

ρ opt d (ξ, η) = ρ c (ξ, η, u opt
). Thus, desired density corresponds to the critical density achieved for u opt . In terms of Fig. 6, if φ min max (u = 1) = φ max (u opt ) for some (ξ, η) ∈ Ω, then u opt is such that ρ opt d = ρ c (u opt ). In terms of Theorem 2, the desired flow φ d = φ max (u opt ). Hence, the controller should provide the same maximal possible flow, while the density is increased, since ρ opt d > ρ 1 . Notice that due to the change of FD shape, at the desired equilibrium traffic operates only at critical density, i.e., there are no congestions in the whole area. Let us again skip (ξ, η) in the notations for simplicity. In order to find u opt ∀(ξ, η) ∈ Ω, we use (44) and (3), and obtain

φmax(u opt ) = v 1 v 1 + ω 1 ω 1 ρc u opt ω 1 + (1 -u opt )ω add ω 1 + v 1 u opt + (1 -u opt )ω add , (46)
where ρ c corresponds to the critical density as in (3) for v = v 1 and ω = ω 1 . Further, we use ρ c v 1 = φ max1 with φ max1 being the highest possible flow for some (ξ, η) ∈ Ω reached with u = 1, and

ω add = ω 2 1 /v 1 to rewrite (46) as φmax(u opt ) = φmax 1 u opt v 1 + (1 -u opt )ω 1 ω 1 + (v 1 -ω 1 ) u opt . (47)
Let us now introduce a coefficient κ ∈ (0, 1] to denote the ratio of the flow at the strongest bottleneck along η-line to the maximal possible flow at space point (ξ, η) for u = 1:

κ(ξ, η) = φ min max (η, 1) φ max (ξ, η, 1)
.

From (47) we get the following equation ∀(ξ, η) ∈ Ω to be solved for u opt :

κ = u opt (v 1 + (1 -u opt )ω 1 ) ω 1 + (v 1 -ω 1 ) u opt ,
which can be further expanded as

u opt 2 + u opt κ v 1 ω 1 -1 - v 1 ω 1 -1 + κ = 0.
This is a quadratic equation with respect to u opt , which yields two solutions. We pick the one with the minus sign, since this guarantees that u opt remains below 1:

u opt = µ + 1 -κ(ν -1) -(ν + 1 -κ(ν -1)) 2 -4κ 2 , (48)
with ν = v 1 /ω 1 . Finally, the optimal equilibrium is the critical density defined in (3) obtained for u opt from (48):

ρ opt d = ω(u opt ) v(u opt ) + ω(u opt ) ρ max , (49)
where v(u opt ) and ω(u opt ) can be taken from (42) for u = u opt and ω add = ω 2 1 /v 1 .

G. Numerical Setup

As a network we again take the structure of Grenoble downtown with its real infrastructure parameters. The critical density in triangular FD is again ρ c = ρ max /3. The initial datum is given ∀(ξ, η) ∈ Ω by ρ 0 (ξ, η) = 3ρ max (ξ, η)/4, thus, it is in the congested traffic regime. Inflow demand and outflow supply are set to the maximal possible steady state flows for u = 1, that is

D (ξ min , η, ρ in (η), u) = S (ξmax, η, ρ out (η), u) = φ min max (η, 1),
which should be chosen to maximize system's throughput. The desired optimal steady state (49) is constructed following the steps described above, and it is depicted in Fig. 7b). This state is characterized by the maximal possible flow through the system achieved for the maximal possible number of vehicles. The numerical scheme needed to discretize the PDE system is again the Godunov scheme in 2D as in the previous numerical example. The only difference is that the for every grid point in space-time ∀(i, j, k) ∈ {1, . . . , m} × {1, . . . , n j } × Z + , the flux function must include dependence on VSL controller as in (42) for u = u opt from (48).

In [START_REF] Mollier | Twodimensional macroscopic model for large scale traffic networks[END_REF] there exists an upper bound for the controller gain γ that guarantees that B(ξ, η, ρ) > 0 ∀(ξ, η, t) ∈ Ω × R + . However, one can accelerate the convergence rate by choosing the maximal possible γ(η, t) for each line of constant η and for each time. Thus, we will compare control results obtained with two different control gains:

1) A constant control gain γ = 0.14 that is the largest possible value for a given urban network (Grenoble) that matches the bounds stated in Theorem 2. 2) A time-and space-varying control gain γ(η, t):

γ(η, t) = 1 - max -min η ξ ξmin(η) ρ( ξ, η, t)d ξ, δ , (50)
where δ > 0 is chosen to get γ > 0 even if the minimum is positive (and in this case an arbitrarily large γ can be used), and > 0 provides the lower bound for B(ξ, η, ρ). Notice that Theorem 2 was proved for the case of constant γ (as in item 1). However, convergence can be accelerated also with γ that depends on η and t as in (50). The only issue is that function B must be always positive, and also that γ can not depend on dimension ξ, since in this case the feedback linearization would not work such that the dynamic equation turns into [START_REF] Tumash | Boundary Control Design for Traffic with Nonlinear Dynamics[END_REF] due to an additional derivative term w.r.t. ξ. Fig. 7c) -f) illustrates the temporal evolution of traffic density under the VSL control (48) with a time-varying gain given by (50) with = 0.01 and δ = 0.1. Thereby, at every time step, demand and supply functions at domain boundaries are set to the desired (maximal possible) flow of every η-line. State converges to the desired equilibrium, which becomes visible already after t = 2 hours of simulation time, (compare plots b) and f)).

Remark 9. At the desired equilibrium the critical density at each point of space will be higher than at initial time, since the VSL control changes the FD shape and affects the desired density as in (49). Therefore, results presented in Fig. 7 may look like driving traffic towards more congested regime, although it is still in free-flow (recall that at the desired equilibrium traffic operates at critical density that becomes higher under VSL control). Traffic flow corresponds to the maximal possible steady state flow determined by the network geometry (capacities at strongest bottlenecks).

Further, we compute L 1 error norm as in [START_REF] Goettlich | Traffic light control: a case study[END_REF] with ρ opt d (ξ, η) from (49) being the desired state. Its evolution within under control within 3 hours is shown in Fig. 7a) for two different control gains. As in the previous example, we again observe that a larger control gain (50) provides a higher convergence speed in comparison to constant γ = 0.14. Recall that as soon as we start applying control, the traffic system is completely set to the free-flow regime, since we always choose the minimal VSL value (see Proposition 1).

VI. CONCLUSIONS

We have elaborated control of traffic on urban networks of any size. We found an approach to analyze the 2D conservation law model such that one gets information about vehicle trajectories in urban area. Such analysis became possible since by its structure this model is applicable only for networks that consist of uni-directional roads and have no loops. The direction field depends only on network geometry and not on state. Thus, we were able to define a curvilinear coordinate transformation that transforms the 2D traffic system into a parametrized set of 1D systems with space-dependent FD, which is similar to inhomogeneous 1D LWR. Although this coordinate transformation could be defined due to specific model restrictions, it can still be used to predict traffic evolution in several frequently occurring situations, e.g., when during a morning rush hour all vehicles stream to the business district.

Further, we have elaborated two control strategies to demonstrate the benefits of rewriting the system in curvilinear coordinates. First, we have posed a boundary control design problem for this system to approximate a desired vehicle trajectory in asymptotic time. The problem was solved using the Hamilton-Jacobi formalism that enabled us to handle discontinuities that always occur in mixed-regime traffic. An additional difficulty was caused by the explicit space-dependency in the FD such that we had to apply the viability theory for space-dependent Hamiltonians. For a numerical example, we took the structure of Grenoble downtown as a urban network. The simulation results confirmed the effectiveness of the boundary controller in achieving convergence towards a desired density profile.

Second, we analyzed the system in curvilinear coordinates to design an in-domain controller using VSL that affects traffic flow by imposing temporary restrictions on allowed speed. It acts as a feedback linearization such that the state equation loses its conservation law structure. This VSL controller is able to stabilize the 2D system to any desired space-varying equilibrium. We have shown that the controller is differentiable almost everywhere in space if FD depends on VSL in a special (physically reasonable) way. We have also designed an optimal steady state that corresponds to throughput maximization achieved for the maximal possible number of cars. In a numerical example, it was shown how vehicle density convergences to the desired equilibrium under the VSL controller.

An appealing direction for future studies might be to extend this analysis for truly multi-directional traffic that includes flow crossings and loops.

APPENDIX I HAMILTON-JACOBI PDE FOR SPACE-DEPENDENT FD

Here we consider the system in (ξ, η)-space [START_REF] Geroliminis | Existence of urban-scale macroscopic fundamental diagrams: Some experimental findings[END_REF] in Hamilton-Jacobi formulation (H-J). It has been known [START_REF] Daganzo | A variational formulation of kinematic waves: basic theory and complex boundary conditions[END_REF], [START_REF] Daganzo | A variational formulation of kinematic waves: solution methods[END_REF], [START_REF] Daganzo | On the variational theory of traffic flow: Well-posedness, duality and applications[END_REF] that an LWR-type equation can be solved exactly in H-J formalism, which is an integral formulation of LWR. Its solution does not contain shocks, which facilitates analysis of its properties (see [START_REF] Tumash | Boundary Control Design for Traffic with Nonlinear Dynamics[END_REF]). H-J formalism for [START_REF] Geroliminis | Existence of urban-scale macroscopic fundamental diagrams: Some experimental findings[END_REF] is not trivial due to explicit space-dependency of FD and additional space parameter η ∈ [η min , η max ] used as a flow path label.

In H-J formulation, traffic is described in terms of cumulative number of vehicles M (ξ, η, t) counted at a given position after a given time. This function is called the Moskowitz function named after an engineer who first used it to investigate traffic, although it was first mentioned only some decades later in [START_REF] Newell | A simplified theory of kinematic waves in highway traffic, part i,ii and iii[END_REF]. In the following, we skip writing η in the arguments to make the notations less heavy. Since all functions in this section depend on η in the same way, all the following steps hold for each line of constant η. Derivatives of M (ξ, η, t) w.r.t. time and space correspond to flow φ (value of flux function Φ(ρ) for a fixed ρ) and density ρ, respectively:

ρ(ξ, t) = - ∂M (ξ, t) ∂ξ , φ(ξ, t) = ∂M (ξ, t) ∂t . (51)
Let us now express M (ξ, t) through the boundary values (inflows φ in (t), outflows φ out (t) and initial condition ρ 0 (ξ)), as well as through the current state ρ(ξ, t). This can be simply done by using the definitions from (51). Namely, we can define a conservative field (-ρ(ξ, t), φ(ξ, t)), which is a gradient of M (ξ, t) ∀(ξ, t) ∈ [ξ min , ξ max ]×R + . By the gradient theorem, it follows that the value of line integral of this field does not depend on a particular chosen path, and equals only to the difference in the values of the Moskowitz function between ending and starting points of the path in space-time. Since M (ξ, t) is an integral function that is defined up to a constant, we are free to assign a reference value to this function at some particular point in space-time. Let us choose a starting point (ξ max , 0) corresponding to the end of traffic line at inital time. Then, we also set M (ξ max , 0) = 0, since this is a decreasing function of position and increasing function of time, see Chapter 14 of [START_REF] Aubin | Viability Theory: New Directions[END_REF]. Thus, taking the ending point of the path as (ξ, t), one possible integration path is

(ξ max , 0) → (ξ max , t) → (ξ, t), which yields ∀(ξ, t) ∈ [ξ min , ξ max ] × R + M (ξ, t) = t 0 φ out (τ )dτ + ξmax ξ ρ( ξ, t)d ξ. (52)
or taking another integration path (ξ max , 0) → (ξ min , 0) → (ξ min , t) → (ξ, t) results into

M (ξ, t) = ξmax ξmin ρ 0 ( ξ)d ξ + t 0 φ in (τ )dτ - ξ ξmin ρ( ξ, t)d ξ. (53)
Let us use (51) to rewrite the fundamental flow-density relation Φ (ξ, ρ) = φ(ξ, t) as

∂M (ξ, t) ∂t -Φ ξ, - ∂M (ξ, t) ∂ξ = 0. (54)
This is a Hamilton-Jacobi PDE (or Moskowitz PDE in the context of traffic). In terms of viability theory, its state M (ξ, t) can also be called the congestion function (see [START_REF] Aubin | Viability Theory: New Directions[END_REF]), since (54) can be viewed as an optimal control problem minimizing a congestion functional M (ξ, t), i.e., vehicles tend to minimize congestion by adapting their individual (microscopic) velocities to the kinematic wave velocity (a macroscopic quantity). In (54) Φ plays the role of a Hamiltonian.

A. General Solution of H-J

The Moskowitz PDE (54) can be solved analytically in accordance with the variational principle using the boundary data. This requires to specify the initial M Ini (ξ) and the boundary condition functions M Up (t), M Down (t).

For convenience of notation, let us introduce the value condition function c(ξ, t) : Dom(c) → R + , where Dom(c) = ({ξ min , ξ max } × R + ) ∪ ([ξ min , ξ max ] × {0}), which aggregates the initial and boundary conditions of (54):

c(ξ, t) =      M Ini (ξ, ) t = 0, M Up (t) ξ = ξ min , M Down (t) ξ = ξ max . (55)
Note that due to (52) and (53), function c is well-defined and continuous on its domain. The IBVP (54)-(55) represents an integral form of [START_REF] Geroliminis | Existence of urban-scale macroscopic fundamental diagrams: Some experimental findings[END_REF]. The boundary conditions are meant to be given in a weak sense (demand-supply problem), and H-J IBVP is well-posed. Let us solve it explicitly.

First, we specify the value condition function (55) by calculating M Up (t), M Down (t) and M Ini (ξ). Thus, the upstream boundary condition M Up (t) can be obtained by considering (53) for ξ = ξ min , which results ∀t ∈ R + into

M Up (t) = c(ξ min , t) = t 0 φ in (τ )dτ + ξmax ξmin ρ 0 ( ξ)d ξ. (56)
Then, M Down (t) can be expressed from (52) for ξ = ξ max :

M Down (t) = c(ξ max , t) = t 0 φ out (τ )dτ, ∀t ∈ ×R + . (57)
Finally, M Ini (ξ) can be expressed from (53) or (52) for t = 0:

M Ini (ξ) = c(ξ, 0) = ξmax ξ ρ 0 ( ξ)d ξ. (58)
Second, we introduce a Legendre-Fenchel transform of the flux function Φ(ξ, ρ):

∀v ∈ [-ω(ξ), v(ξ)] : L(ξ, v ) = sup ρ∈[0,ρmax(ξ)] (Φ(ξ, ρ) -v ρ), (59)
where v and -ω are kinematic wave speeds (slopes of FD) in case of ρ = 0 and ρ = ρ max , respectively. From Chapter 14.3 of [START_REF] Aubin | Viability Theory: New Directions[END_REF] we get the Legendre transform (59) of a triangular FD:

L(ξ, v ) = φ max (ξ) -ρ c (ξ)v , ∀v ∈ [-ω(ξ), v(ξ)]. (60)
The closed-form unique solution of (54)-(55) in Barron-Jensen/Frankowska sense (see Theorem 13.10.3 of [START_REF] Aubin | Viability Theory: New Directions[END_REF]) corresponds to the infimum among all viable evolutions starting from a boundary at t -T and arriving at ξ at terminal time t:

M (ξ, t) = inf (T,v )∈S   c ξ(0), t -T + T 0 L ξ(τ ), v (τ ) dτ   , (61)
where the infimum is taken over domain S defined as:

S = (T, v ) T ∈ R + , v (•) ∈ L 1 (0, T ), ξ(τ ) = v (τ ), ξ(T ) = ξ, v (τ ) ∈ -ω ξ(τ ) , v ξ(τ ) , ξ(0), t -T ∈ Dom(c) . (62)
Here ξ(τ ) denotes trajectory of an observer moving along a traffic stream at speed v (τ ) that might be non-constant due to inhomogeneity of network infrastructure. Trajectory ξ(τ ) originates at τ = 0 on a boundary of domain of c and arrives at point ξ at terminal time τ = T .

We introduce two-argument functions M Up (ξ, t), M Down (ξ, t) and M Ini (ξ, t) as viability episolutions [START_REF] Claudel | Lax-Hopf Based Incorporation of Internal Boundary Conditions Into Hamilton-Jacobi Equation. Part I: Theory[END_REF] to (61) for corresponding domains of function c, which are M Up (t) (56), M Down (t) (57) and M Ini (ξ) (58). For example, M Up (ξ, t) is a "solution candidate" that arrives at space-point (ξ, t) from the upstream boundary with a given "initial cost" M Up (t), while M Down (ξ, t) and M Ini (ξ, t) are defined similarly. This enables us to restate the unique solution of (54)-(55) as a minimum of three functions

∀(ξ, t) ∈ [ξ min , ξ max ] × R + M (ξ, t) = min{M Up (ξ, t), M Down (ξ, t), M Ini (ξ, t)}. (63)
Notice that in the following we will consider only solutions for large enough time

t ≥ max      ξmax ξmin 1 v( ξ)d ξ , ξmax ξmin 1 ω( ξ)d ξ      . (64)
In Sections I-B-I-D we calculate episolutions M Up (ξ, t), M Down (ξ, t) and M Ini (ξ, t), and then the unique solution is the minimum of them (63).

B. Upstream Boundary Condition

M Up (ξ, t) is the minimal cumulative vehicle number that originates from the upstream boundary ξ min at initial time.

By definition of the value condition function (55) we get c( ξ(0), t -T ) = M Up (t -T ) in (61). The upstream boundary condition is assigned to ξ min , which implies the following start and end points of observer trajectory that starts traveling with non-constant speed v (τ ) ∈ [-ω (ξ(τ )) , v (ξ(τ ))]:

ξ(0) = ξ min , ξ(t) = ξ min + t 0 v (τ )dτ. (65)
Using ( 61) with (60) for c ξ(0), t -T = M Up (t -T ) from (56), we get the following infimum problem:

M Up (ξ,t) = inf (T,v )∈S Up t-T 0 φ in (τ )dτ + ξmax ξ min ρ 0 ( ξ)d ξ + T 0 φmax ξ(τ ) dτ - T 0 ρc ξ(τ ) v (τ )dτ . ( 66
)
where infimum is taken over domain S Up defined as in (62) but with ξ(0), t -T ∈ Dom(c Up ), where c Up = M Up (t) as in (55). Let us consider the last term in (66). By definition d ξ = v (τ )dτ , and we perform a change of variables:

T 0 ρ c ξ(τ ) v (τ )dτ = ξ ξmin ρ c ( ξ)d ξ =: R c (ξ), (67)
where R c (ξ) is a new variable that denotes cumulative critical density. Further, we can decompose the integrals in (66) as and recall that v (τ ) is related to ξ(τ ) by (65). By Assumption 1 we have φ in (t) ≤ φ max (ξ) ∀(ξ, t) ∈ [ξ min , ξ max ] × R + . Hence, the infimum is achieved when the traveling time T is minimized, i.e., the solution is assigned to a traveler that moves with the maximal speed v everywhere, thus (65) becomes

ξ(t) = ξ min + t 0 v( ξ(τ ))dτ. (69)
Thus, in the infimum, T is the solution to (69) for t = T :

We rewrite the first term on the right-hand side of (71) as Using Assumption 1, we obtain that the infimum is achieved for the minimal traveling time T , which corresponds to:

T ω (ξ) =

D. Initial Condition

Finally, we calculate M Ini (ξ, t) related to vehicle with a known label at initial time (58) that follows the path of viable evolution. We can establish that T = t, since the viability evolution starts its path at initial time. Thus, using (61) with initial condition (58), we can state the infimum problem as

where φ min max is as in [START_REF] Garavello | Conservation laws on complex networks[END_REF]. In the following, we will show that starting from t min the effect of initial conditions will have left the system and thus can be excluded from the minimum.

E. Time When Initial Conditions Leave The System

We seek to estimate t min such that ∀(ξ, t) ∈ [ξ min , ξ max ] × [t min , +∞): M Ini (ξ, t) ≥ M Up (ξ, t) or M Ini (ξ, t) ≥ M Down (ξ, t). First, we estimate the time after which M Ini (ξ, t) ≥ M Up (ξ, t), then we do the same for M Ini (ξ, t) ≥ M Down (ξ, t), and t min is the smallest of two values. We combine the result for M Up (ξ, t) (72) with the lower bound for M Ini (ξ, t) (79), and write Finally, t min is the minimum between (82) and (83):

t min = t c   1 +      1 ξmax ξmin ρ max ( ξ) + ρ c ( ξ) d ξ        . ( 84
)
Notice that t min in (84) is different for every η-line. The common t min is then just determined as the maximal value.

F. Unique Solution

The final solution M (ξ, t) of ( 54)-( 55) can be obtained as a minimum of (72) and (75) ∀t ∈ [t min , +∞), thus, the effect of initial conditions is negligible:

M (ξ, t) = min

Fig. 1 .

 1 Fig. 1. Coordinate transformation mapping: a) curved trajectories in Grenoble downtown in (x, y)-plane into b) straight lines in (ξ, η)-plane.

Fig. 2 .Fig. 3 .

 23 Fig. 2. a) Initial state to be controlled from the boundaries (black arrows). b) L 1 error as a function of time for different control gains.

1 Fig. 5 .

 15 Fig. 5. FD as a function of u: a) monotonic dependence for a fixed ρ in free-flow; b) concave dependence for a fixed ρ in congested regime.

Fig. 6 .

 6 Fig. 6. Blue line: FD for u = 1. Red line: FD for u = u opt .

a) L 1

 1 Fig. 7. a) L 1 norm of density error as a function of time for different control gains, b) the desired optimal equilibrium as in (49). Traffic flow control by VSL in Grenoble downtown. Density ρ(x, y, t) at: c) t = 0; d) t = 10 min; e) t = 30 min; f) t = 2 hours.

φ

  max ( ξ(τ )) -φ in (t -T + τ ) dτ.(68) Thus, using (67) and (68) we can rewrite (66) asM Up (ξ, t) = inf (T,v )∈S Up   T 0 φmax( ξ(τ )) -φ in (t -T + τ ) dτ ξ)d ξ -Rc(ξ),

  the viability episolution yields from (66):M Up (ξ, t) = Tv(ξ) 0 φ max ( ξ(τ ))dτ + t-Tv(ξ) 0 φ in (τ )dτ + ξmax ξmin ρ 0 ( ξ)d ξ -R c (ξ).

  ξ(τ ))d ξ = R c (ξ).where (69) and (67) were used. With this result, two R c (ξ) terms with opposite signs in (71) cancel each other, and the upstream boundary solution reads:M Up (ξ, t) = t-Tv(ξ) 0 φ in (τ )dτ + ξmax ξmin ρ 0 ( ξ)d ξ.(72)C. Downstream Boundary ConditionNow we need to obtain M Down (ξ, t) related to the downstream boundary ξ max . Here, viable evolutions are characterized by the following start and end points of traveling:ξ(0) = ξ max , ξ(t) = ξ max + t 0 v (τ )dτ,where v (τ ) ∈ [-ω (ξ(τ )) , v (ξ(τ ))].

0 φ

 0 previous case, we use M Down (t) from (57) and the result from (68) to write the following infimum problemM Down (ξ, t) = inf (T,v )∈S Down -T 0 ρc( ξ(τ ))v (τ )dτ + T 0 φmax( ξ(τ )) -φ out (t -T + τ ) dτ + t out (τ )dτ,where domain S Down is defined as in (62) with ξ(0), t -T ∈ Dom(c Down ), where c Down = M Down (t) as in (55).

ξmax ξ 1 ω

 1 ( ξ) d ξ and v = -ω.(74)From the definition of triangular FD (3) we getρ c = ρ max ω v + ω ⇒ ρ max ω = ρ c (v + ω),which is then used together with (74), φ max = ρ c v and the change of variables d ξ = ω(τ )dτ to solve the infimum: M Down (ξ, t) = t-Tω(ξ) 0 φ out (τ )dτ + ξmax ξ ρ max ( ξ)d ξ. (75)

M 0 v

 0 Ini (ξ, t) the domain S Ini is defined as in (62) for T = t:S Ini = v v (•) ∈ L 1 (0, t), ξ(τ ) = v (τ ), ξ(t) = ξ, v (τ ) ∈ -ω ξ(τ ) , v ξ(τ ) , ξ(0) ∈ [ξ min , ξ max ] .(77)Notice that now the viable evolution starts its path from ξ0 at initial time with v(τ ) ∈ [-ω (ξ(τ )) , v (ξ(τ ))] ξ(0) = ξ0 , ξ(t) = ξ0 + t (τ )dτ.(78)Using the change of variables v (τ )dτ = d ξ we rewrite (76) and estimate the lower bound term by term:M Ini (ξ, t) ≥ 0 -

MM

  Ini (ξ, t) -M Up (ξ, t) ≥ -Ini (ξ, t) -M Up (ξ, t) ≥ -ξmax ξ min ρmax( ξ) + ρc( ξ) d ξ + t-Tv (ξ) 0 φ min max -φ in (τ ) dτ.(80)Using Assumption 1, we are able to estimate the following lower bound for the second term on the right-hand side (80):t-Tv(ξ) 0 φ min max -φ in (τ ) dτ ≥ t -T v (ξ) t c ,with > 0 and t c is from[START_REF] Aubin | Viability Theory: New Directions[END_REF], which lets us rewrite (80) asM Ini (ξ, t)-M Up (ξ, t) ≥ t -T v (ξ max ) t c ξmax ξmin ρ max ( ξ) + ρ c ( ξ) d ξ.(81)The right-hand side of (81) becomes always non-negative after same steps are performed to obtain the minimal time, for which M Ini (ξ, t) -M Down (ξ, t) ≥ 0 holds:

  where T v (ξ) and T ω (ξ) are given by (70) and (74).

03167733v2

This work is a part of Scale-FreeBack project that received funding from the European Research Council (ERC) under the European Unions Horizon 2020 research and innovation programme (grant agreement N 694209). L. Tumash and C. Canudas-de-Wit are with Univ.