Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations

ARTICLE IN PRESS

JID: TRB

[m3Gsc;June 16, 2017;12:13]

Transportation Research Part B 0 0 0 (2017) 1–24

Contents lists available at ScienceDirect

Transportation Research Part B journal homepage: www.elsevier.com/locate/trb

Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations Edward S. Canepa a,∗, Christian G. Claudel b a b

Electrical Engineering Department, King Abdullah University of Science and Technology, Thuwal, Saudi Arabia Department of Civil, Architectural and Environmental Engineering, University of Texas Austin, USA

a r t i c l e

i n f o

Article history: Received 23 May 2016 Revised 30 May 2017 Accepted 31 May 2017 Available online xxx Keywords: Traffic estimation Mixed integer linear programming Optimization

a b s t r a c t Nowadays, traffic management has become a challenge for urban areas, which are covering larger geographic spaces and facing the generation of different kinds of traffic data. This article presents a robust traffic estimation framework for highways modeled by a system of Lighthill Whitham Richards equations that is able to assimilate different sensor data available. We first present an equivalent formulation of the problem using a Hamilton– Jacobi equation. Then, using a semi-analytic formula, we show that the model constraints resulting from the Hamilton–Jacobi equation are linear ones. We then pose the problem of estimating the traffic density given incomplete and inaccurate traffic data as a Mixed Integer Program. We then extend the density estimation framework to highway networks with any available data constraint and modeling junctions. Finally, we present a travel estimation application for a small network using real traffic measurements obtained obtained during Mobile Century traffic experiment, and comparing the results with ground truth data. © 2017 Elsevier Ltd. All rights reserved.

1. Introduction Transportation research is currently at a tipping point; the emergence of new transformative technologies and systems, such as vehicle connectivity, automation, shared-mobility, and advanced sensing is rapidly changing the individual mobility and accessibility. This will fundamentally transform how transportation planning and operations should be conducted to enable smart and connected communities. The transport systems can be highly beneficiated and become safer, more efficient and reliable. Nowadays, dynamic routing and traffic-dependent navigation services are available for users. Such applications need to estimate the present traffic situation and that of the near future at a forecasting horizon based on data that are available in real-time. Traffic state estimation for a road network refers to estimate all the traffic variables (e.g. cars density, speed) of the network at an instant of time based of traffic measurements. This is, for a limited amount of traffic data the estimator obtains a complete view of the traffic scenario. This estimation requires the fusion or traffic data and traffic models, the latter are typically formulated as partial differential equations (PDEs). For this framework, we will use the Lighthill–Whitham–Richards (LWR) partial differential equation (Lighthill and Whitham, 1956; Richards, 1956) which is commonly used to model highway traffic; derivating the model constraints is a complex problem. Other estimation tech-

∗

Corresponding author. E-mail addresses: [email protected] (E.S. Canepa), [email protected] (C.G. Claudel).

http://dx.doi.org/10.1016/j.trb.2017.05.016 0191-2615/© 2017 Elsevier Ltd. All rights reserved.

Please cite this article as: E.S. Canepa, C.G. Claudel, Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations, Transportation Research Part B (2017), http://dx.doi.org/10.1016/j.trb.2017.05.016

JID: TRB 2

ARTICLE IN PRESS

[m3Gsc;June 16, 2017;12:13]

E.S. Canepa, C.G. Claudel / Transportation Research Part B 000 (2017) 1–24

niques such as Extended Kalman Filtering (Alvarez-Icaza et al., 2004) (EKF), Ensemble Kalman Filtering (Work et al., 2010) or Particle Filtering (PF) rely on approximations to determine the model constraints, either through linearization or sampling. Recent estimation traffic network techniques involve the Link Transmission Model (LTM) (Yperman et al., 2005; Jin, 2015) in which the Lax-Hopf formula is used to compute the demand and supply of the links and invariant junction models are used to calculate the boundary flows. Traffic estimation on a freeway segment using heterogeneous data set has also been addressed on Deng et al. (2013), is based on the establishment of a cumulative vehicle count-based state estimation models for using AVI and GPS data. Along this line, in Nantes et al. (2016) the fusion of three heterogeneous data sources on a single EKF-based estimator was proposed,they provided a comprehensive analysis of the estimation accuracy using data from loop detectors, GPS and Bluetooth scanners. No approximation of the model is required by the framework presented on this article. An example of the usage of this framework is determining the ranges input flows (or any convex function of the boundary data) compatible with the traffic model and measurement data. The exact estimation technique presented in this paper is based on the Moskowitz function (Moskowitz, 1965; Newell, 1993); it is used here as an intermediate computational abstraction. The Moskowitz function can be understood as the integral form of the density function, and solves an Hamilton–Jacobi (HJ) PDE, whereas the density function itself solves the LWR PDE. An advantage of using the HJ PDE is that its solutions can be expressed semi analytically (Claudel and Bayen, 2010), which enables the derivation of the model constraints explicitly. We will now summarize the key differences between our estimation approach and the works presented earlier on this section: •

•

•

•

•

Kalman filtering considers Gaussian model noise, which is not considered in this approach. It is also considering a Gaussian error model (probabilistic), whereas the present algorithm considers a deterministic error model (for example L1, L2 or L infinity error under some threshold). Note that the present algorithm can also adjust the confidence we put on the data, which can be a term in the objective function In the context of the triangular diagram, the Extended Kalman filter requires mode identification (as in Munoz et al., 2003). Other approaches can be used, such as the EnKF, which do not require mode identification (Piccoli and Bayen, 2010). Virtually all approaches are based on the CTM (the discretized LWR model), which introduces discretization errors, unlike the approach used in the present article Incorporating travel time constraints requires the integration of velocities over multiple time steps, which can slow down Kalman Filter approaches significantly. In addition, the integration errors cause an additional amount of uncontrolled errors when incorporating the travel time constraints, unlike the proposed approach. The present approach is not a minimum mean square error estimator: the objective of the optimization problem can be chosen freely, allowing different types of problem to be solved, for instance to solve L1 norm minimization problems, to find traffic state estimates that are as sparse as possible.

1.1. Contributions of the article The present article builds on Canepa and Claudel (2012), Canepa et al. (2013), Li et al. (2014a), Li et al. (2014b) and Anderson et al. (2013) which introduced a Mixed Integer Linear Programming framework for solving data assimilation and data reconciliation problems, for specific objective functions. In the present article, the framework initially described in Canepa and Claudel (2012) is extended to network traffic density estimation. The present article has the following contributions over the previous work presented in Canepa and Claudel (2012): •

•

•

The integration of internal traffic density data, or arbitrary travel time data (not necessarily defined as the travel time required to cross the entire physical domain), which was not considered in earlier articles (Canepa and Claudel, 2012; Canepa et al., 2013). The extension of the traffic state estimation framework defined in Canepa and Claudel (2012) and Canepa et al. (2013) to transportation networks, which require the proper modeling of junctions, and the integration of the entropy condition to junction flows. The formulation of estimation problems that do not involve a minimum variance estimation, unlike classical estimation schemes derived from the Kalman Filter. Examples of non minimum variance estimation include compressed sensing (L1 norm minimization), shown in Section 6.

The outline of this article is the following. In Section 2 we define the solution to the LWR PDE and its equivalent formulation as a HJ PDE. In Section 3, we recall the analytical expressions of the solutions to HJ PDEs for the triangular flux functions investigated in this article, and show that the LWR PDE constraints correspond to convex constraints in the unknown initial, boundary and internal condition parameters. A first estimation example is shown in Section 4, using boundary and internal conditions from measurement data the unknown initial conditions are estimated. The framework is extended to Highway Networks in Section 6, where we also validated it using experimental traffic flow data (e.g. density, point velocity and travel time) collected during the Mobile Century traffic experiment. Please cite this article as: E.S. Canepa, C.G. Claudel, Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations, Transportation Research Part B (2017), http://dx.doi.org/10.1016/j.trb.2017.05.016

ARTICLE IN PRESS

JID: TRB

E.S. Canepa, C.G. Claudel / Transportation Research Part B 000 (2017) 1–24

[m3Gsc;June 16, 2017;12:13] 3

2. Background 2.1. The Lighthill–Whitham–Richards traffic flow model For the remainder of the article, we will assume that the spatial domain representing the highway section is [ξ , χ ], where ξ and χ respectively represent the upstream and downstream boundaries of the domain. Traffic flow on this section can be described by the density function, denoted as ρ ( ·, ·). The density function represents an aggregated number of vehicles per space unit, and can is modeled by the Lighthill–Whitham–Richards (LWR) PDE:

∂ρ (t, x ) ∂ψ (ρ (t, x )) + = 0 ∂t ∂x

(1)

The function ψ ( · ) is named flux function. It depends on several empirical parameters (e.g. number of lanes, the drivers habits). Different models have been proposed for ψ , in particular the triangular model defined below, this model is widely used in the literature (Daganzo, 1994; 2005; 2006).

ψ (ρ ) =

 vρ w ( ρ − ρm )

if ρ ≤ ρc otherwise

(2)

In the remainder of this article, we assume that the flux function is triangular and given by (2). While other concave flux functions could be used and would also yield convex constraints, the instantiation of model constraints as linear inequalities requires piecewise linear flux functions, such as (2). 2.2. Hamilton–Jacobi equation Equivalently, the state of traffic can be described by a scalar function M( ·, ·) of both time and space, known as Moskowitz function (Moskowitz, 1965; Newell, 1993). The Moskowitz function is a macroscopic description of traffic flow, which appears naturally in the context of traffic. It can be interpreted as follows: let consecutive integer labels be assigned to vehicles entering the highway at location x = ξ . The Moskowitz function M( ·, ·) satisfies M(t, x ) = n where n is the label of the vehicle located in x at time t (Daganzo, 2005; 2006), and is assumed to be continuous. The density function ρ ( ·, ·) is related (Newell, 1993) to the spatial derivative of the Moskowitz function M( ·, ·) as follows:

ρ (t, x ) = −

∂ M(t, x ) ∂x

(3)

If the density function is to be modeled by the LWR PDE, the Moskowitz function satisfies an Hamilton–Jacobi (HJ) PDE obtained (Aubin et al., 2008; C.G.Claudel and Bayen, 2010) by integration of the LWR PDE:

  ∂ M(t, x ) ∂ M(t, x ) −ψ − = 0 ∂t ∂x

(4)

Several classes of weak solutions to Eq. (4) exist, such as viscosity solutions (Crandall and Lions, 1983; Bardi and CapuzzoDolcetta, 1997) or Barron–Jensen/Frankowska (B-J/F) solutions (Barron and Jensen, 1990; Frankowska, 1993). For the problem investigated in this article, these solutions are equivalent, and can be computed implicitly using a Lax-Hopf formula. 2.3. Barron–Jensen/Frankowska solutions to Hamilton Jacobi equations In order to characterize the B-J/F solutions, we first need to define the Legendre–Fenchel transform of the Hamiltonian

ψ ( · ) as follows.

Definition 2.1 (Legendre–Fenchel transform). For an upper semicontinuous Hamiltonian ψ ( · ), the Legendre–Fenchel transform ϕ ∗ ( · ) is given by:

ϕ ∗ (u ) := sup p∈Dom(ψ ) [ p · u + ψ ( p)]

(5)

Solving the HJ PDE (4) requires the definition of value conditions, which encode the traditional concepts of initial, boundary and internal conditions. Definition 2.2 (Value condition). A value condition c( ·, ·) is a lower semicontinuous function defined on a subset of [0, tmax ] × [ξ , χ ]. In the remainder of the article, a value condition can encode an initial condition, an upstream boundary condition or a downstream boundary condition. Each of these functions is defined on a subset of R+ × [ξ , χ ]. For each value condition c( ·, ·), we define the partial solution (Mazare et al., 2011) to the HJ PDE (4) using the Lax-Hopf formula (Aubin et al., 2008; C.G.Claudel and Bayen, 2010). Please cite this article as: E.S. Canepa, C.G. Claudel, Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations, Transportation Research Part B (2017), http://dx.doi.org/10.1016/j.trb.2017.05.016

ARTICLE IN PRESS

JID: TRB 4

[m3Gsc;June 16, 2017;12:13]

E.S. Canepa, C.G. Claudel / Transportation Research Part B 000 (2017) 1–24

Proposition 2.3 (Lax-Hopf formula). Let ψ ( · ) be a concave Hamiltonian, and let ϕ ∗ ( · ) be its Legendre-Fenchel transform (5). Let c( ·, ·) be a lower semicontinuous value condition, as in Definition 2.2. The B-J/F solution Mc ( ·, ·) to (4) associated with c(·, ·) can be algebraically represented (Aubin et al., 2008; C.G.Claudel and Bayen, 2010) by:

Mc (t, x ) =

inf

(u,T )∈Dom(ϕ ∗ )×R+

(c(t − T , x + T u ) + T ϕ ∗ (u ) )

(6)

Eq. (6) implies the existence of a B-J/F solution Mc ( ·, ·) for any value condition function c( ·, ·). However, the solution itself may be incompatible with the value condition that we imposed on it, i.e. we do not necessarily have ∀(t, x ) ∈ Dom(c ), Mc (t, x ) = c(t, x ). 2.4. Properties of the solutions to scalar Hamilton–Jacobi equations The structure of the Lax-Hopf formula (6), implies the following important property, known as inf-morphism property. The inf-morphism property can be formally derived through capture basins, such as in Aubin et al. (2008). Proposition 2.4 (Inf-morphism property). Let the value condition c( ·, ·) be minimum of a finite number of lower semicontinuous functions:

∀(t, x ) ∈ [0, tmax ] × [ξ , χ ], c(t, x ) := min c j (t, x ) j∈J

(7)

The solution Mc ( ·, ·) associated with the above value condition can be decomposed (Aubin et al., 2008; C.G.Claudel and Bayen, 2010; Claudel and Bayen, 2010) as:

∀(t, x ) ∈ [0, tmax ] × [ξ , χ ], Mc (t, x ) = min Mc j (t, x ) j∈J

(8)

In the present work, the value conditions cj ( ·, ·) are not known exactly, either because of measurement uncertainty (case of the upstream and downstream boundary condition) or because of the lack of measurements (case of the initial condition). However, even if the real values of cj ( ·, ·) are not known exactly, they cannot be arbitrary as they have to apply in the strong sense (see Strub and Bayen, 2006 for a mathematical definition) to be compatible with the LWR model. In the remainder of this article, we define the model constraints as the set of constraints that applies on the value conditions cj ( ·, ·) to ensure that all value conditions apply in the strong sense. The inf-morphism property is critical for the derivation of the LWR PDE model constraints, allowing us to instantiate these constraints as inequalities. Proposition 2.5 (Model compatibility constraints for block value conditions). Let c(·, · ) = min j∈J c j (·, · ) be given, and let Mc ( ·, ·) be defined as in (6). The value condition c( ·, ·) satisfies ∀(t, x ) ∈ Dom(c ), Mc (t, x ) = c(t, x ) if and only if the following inequality constraints are satisfied:

Mc j (t, x ) ≥ ci (t , x )

∀(t , x ) ∈ Dom(ci ), ∀(i, j ) ∈ J2

(9)

The proof of this proposition is available in Claudel and Bayen (2011). Note that equation (9) represents an important improvement, as the model constraints are now semi-explicit. In order to solve the problem completely, we still need to evaluate the functions Mc j (·, · ) explicitly. These explicit solutions were derived in Claudel and Bayen (2010) for affine initial, boundary and internal conditions blocks which are presented in the next section. 3. Explicit solutions to piecewise affine initial, internal and boundary conditions As mentioned before, different types of value conditions can be incorporated into the estimation problem. In the present article we will handle initial, internal, upstream and downstream conditions. These value conditions are typically measured (with some error) using fixed sensors, such as inductive loop detectors, magnetometers or GPS devices. 3.1. Definition of affine initial, upstream/downstream boundary and internal density conditions The formal definition of initial, upstream/downstream boundary and internal conditions associated with the HJ PDE (4) is the subject of the following definition. Definition 3.1 (Affine initial, upstream/downstream boundary and internal conditions). Let us define K = {0, . . . , kmax }, N = {0, . . . , nmax }, M = {0, . . . , mmax } and U = {0, . . . , umax }. For all k ∈ K, n ∈ N, m ∈ M and u ∈ U, we define the following functions, respectively called initial, upstream, downstream internal flow and internal density conditions:

⎧ k−1 ⎪ ⎨− i=0 ρini (i )X −ρini (k )(x − kX ) Mk (t, x ) = ⎪ ⎩ +∞

if t = 0 and x ∈ [kX, (k + 1 )X ] otherwise

(10)

Please cite this article as: E.S. Canepa, C.G. Claudel, Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations, Transportation Research Part B (2017), http://dx.doi.org/10.1016/j.trb.2017.05.016

ARTICLE IN PRESS

JID: TRB

E.S. Canepa, C.G. Claudel / Transportation Research Part B 000 (2017) 1–24

[m3Gsc;June 16, 2017;12:13] 5

Fig. 1. Domains of the initial, upstream/downstream boundary and internal conditions. The block upstream and downstream boundary conditions respectively denoted by γ n ( ·, ·) and β n ( ·, ·) are defined on line segments corresponding to the upstream and downstream boundaries of the physical domain. In contrast, the block initial conditions Mk ( ·, ·) are defined on line segments corresponding to the initial time. Note that the actual problem involves block initial conditions covering the entire physical domain [ξ , χ ], and block boundary conditions covering the temporal domain [0, tmax ]. The block of internal flow and density conditions (μm ( ·, ·) and ϒu ( ·, ·) respectively) are defined and any position inside the computational domain. All these functions are unknown (or only partially known).

⎧ n−1 ⎪ ⎨ i=0 qin (i )T γn (t, x ) = +qin (n )(t − nT ) ⎪ ⎩ +∞

⎧ n−1 ⎪ i=0 qout (i )T ⎪ ⎪ ⎨+qout (n )(t − nT ) max βn (t, x ) = − kk=0 ρ ( k )X ⎪ ⎪ ⎪ ⎩ +∞

μm (t, x ) =

if x = ξ and t ∈ [nT , (n + 1 )T ] otherwise

if x = χ and t ∈ [nT , (n + 1 )T ] otherwise

⎧ ⎪ ⎨L(m ) + r (m )(t − tmin (m )) ⎪ ⎩

+∞

ϒu (t, x ) = where vmeas (m ) =

L(u ) − ρ (u )(x − xminρ (u )) +∞

(11)

(12)

if x = xmin (m ) +vmeas (m )(t − tmin (m )) and t ∈ [tmin (m ), tmax (m )] otherwise

(13)

if x ∈ [xminρ (u ), xmaxρ (u )] and t = tρ (u ) otherwise

(14)

xmax (m )−xmin (m ) tmax (m )−tmin (m )

In the above definition, internal density conditions (14) are specific to model density sensors that are inside our computational domain. Regarding flow sensors also located inside the computational domain, they can be thought as an internal flow condition (13) associated with a zero velocity (vmeas (m ) = 0). Note that the affine initial, upstream/downstream boundary and internal conditions defined above for the HJ PDE (4) are equivalent to constant initial, upstream/downstream boundary and internal conditions for the LWR PDE (1). The domains of definitions of these functions are illustrated in Fig. 1.

3.2. Analytical solutions to affine initial, upstream/downstream boundary and internal conditions Given the affine initial, upstream/downstream boundary and internal conditions defined above, the corresponding solutions MMk (·, · ), Mγn (·, · ), Mβn (·, · ) and and Mϒu (·, · ) are given (Claudel and Bayen, 2011; Mazare et al., 2011) by the following Please cite this article as: E.S. Canepa, C.G. Claudel, Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations, Transportation Research Part B (2017), http://dx.doi.org/10.1016/j.trb.2017.05.016

ARTICLE IN PRESS

JID: TRB 6

[m3Gsc;June 16, 2017;12:13]

E.S. Canepa, C.G. Claudel / Transportation Research Part B 000 (2017) 1–24

formulas:

MMk (t, x ) =

⎧ +∞ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪− k−1 ρ (i )X ⎪ ⎪ i=0 ini ⎪ ⎪ +ρini (k )(t v + kX − x ) ⎪ ⎪ ⎪ ⎪ and ⎪ ⎪ ⎪ and ⎪ ⎪ k−1 ⎪ ⎪ − ρ ( i ) X ⎪ ini i =0 ⎪ ⎪ ⎪ +ρc (t v + kX − x ) ⎪ ⎪ ⎪ and ⎪ ⎨

if x ≤ kX + wt or x ≥ (k + 1 )X + vt if kX + t v ≤ x ( k + 1 )X + t v ≥ x ρini (k ) ≤ ρc

and

−1 − ki=0 ρini (i )X ⎪ ⎪ ⎪ ⎪ +ρini (k )(tw + kX − x ) ⎪ ⎪ ⎪ −ρm tw ⎪ ⎪ ⎪ ⎪ and ⎪ ⎪ ⎪ and ⎪ ⎪ ⎪ ⎪ − ki=0 ρini (i )X ⎪ ⎪ ⎪ ρc (tw + (k + 1 )X − x ) ⎪ ⎪ ⎪ ⎪ −ρm tw ⎪ ⎪ ⎪ and ⎪ ⎩

and

Mγn (t, x ) =

⎧ +∞ ⎪ ⎪ n−1 ⎪ ⎪ qin (i )T ⎪ ⎪ ⎪+qi=0(n )(t − x−ξ − nT ) ⎨ in v ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ni=0 qin (i )T ⎪ ⎩ +ρc v(t − (n + 1 )T −

if kX + t v ≥ x kX + tw ≤ x ρini (k ) ≤ ρc

(15)

if kX + tw ≤ x (k + 1 )X + tw ≥ x ρini (k ) ≥ ρc if (k + 1 )X + t v ≥ x (k + 1 )X + tw ≤ x ρini (k ) ≥ ρc if t ≤ nT + x−v ξ if nT + x−v ξ ≤ t and t ≤ (n + 1 )T + x−v ξ

x −ξ

v

(16)

) otherwise

⎧ ⎪ ⎪+∞ ⎪ ⎪ ⎪ ⎪ ⎪− kmax ρini (k )X + n−1 qout (i )T ⎪ i=0 ⎪ k=0 ⎪ ⎪+qout (n )(t − x−wχ − nT ) ⎪ ⎨ −ρm (x − χ ) Mβn (t, x ) = ⎪ ⎪ ⎪ and ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ kmax ⎪ ⎪ − k=0 ρ (k )X + ni=0 qout (i )T ⎪ ⎩ +ρc v(t − (n + 1 )T − x−vχ )

if t ≤ nT + x−wχ

if nT + x−wχ ≤ t t ≤ ( n + 1 )T + x−wχ

(17)

otherwise

 ⎧ vmeas m )(t−tmin (m )) Lm + rm t − x−xmin (m )− − tmin (m ) meas m ) ⎪ v − v ⎪ ⎪ if x ≥ xmin (m ) + vmeas m )(t − tmin (m )) ⎪ ⎪ ⎪ ⎪and x ≥ xmax (m ) + v(t − tmax (m )) ⎪ ⎪ ⎪ ⎪and x ≤xmin (m ) + v(t − tmin (m )) ⎪ ⎪ −vmeas m )(t−tmin (m )) ⎪ Lm + rm t − x−xmin (m )w − tmin (m ) ⎪ −vmeas m ) ⎪ ⎪ x−xmin (m )−vmeas m )(t−tmin (m )) ⎪ ⎪ w−vmeas m ) ⎨+kc (v − w ) meas m )(t − tmin (m )) Mμm (t, x ) = if x ≤ xmin (m ) + v ⎪ ⎪and x ≤ xmax (m ) + w(t − tmax (m )) ⎪ ⎪ and x ≥ xmin (m ) + w(t − tmin (m )) ⎪ ⎪ ⎪Lm + rm (tmax (m ) − tmin (m ) )+ ⎪ ⎪ ⎪ ⎪(t − tmax (m ) )kc v − x−xmax (m) ⎪ t−tmax (m ) ⎪ ⎪ ⎪ if x ≤ xmax (m ) + v(t − tmax (m )) ⎪ ⎪ ⎪ ⎪ ⎩and x ≥ xmax (m ) + w(t − tmax (m ))

(18)

+∞ otherwise

Please cite this article as: E.S. Canepa, C.G. Claudel, Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations, Transportation Research Part B (2017), http://dx.doi.org/10.1016/j.trb.2017.05.016

ARTICLE IN PRESS

JID: TRB

[m3Gsc;June 16, 2017;12:13]

E.S. Canepa, C.G. Claudel / Transportation Research Part B 000 (2017) 1–24

Mϒu (t, x ) =

⎧ +∞ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ L (u ) ⎪ ⎪ ⎪+ρ (u )(v(t − tρ (u )) + xminρ (u ) − x ) ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ and ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ and ⎪ ⎪ ⎪ L ( u ) ⎪ ⎪ ⎪ ⎪+ρc (v(t − tρ (u )) + xminρ (u ) − x ) ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎨ and ⎪ ⎪ and ⎪ ⎪ ⎪ L ( u ) ⎪ ⎪ ⎪ ⎪ +ρ (u )(w(t − tρ (u )) + xminρ (u ) − x ) ⎪ ⎪ ⎪ −ρm w(t − tρ (u )) ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ and ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ and ⎪ ⎪ ⎪ ⎪ L (u ) ⎪ ⎪ ⎪+ρc (w(t − tρ (u )) + xmaxρ − x ) ⎪ ⎪ ⎪ −ρm w(t − tρ (u )) ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ ⎪ and ⎪ ⎪ ⎪ ⎪ ⎩ and

7

if x ≤ xminρ (u ) +w(t − tρ (u )) or x ≥ xmaxρ (u ) +v(t − tρ (u )) or t ≤ tρ (u ) if xminρ (u ) +v(t − tρ (u )) ≤ x xmaxρ (u ) +v(t − tρ (u )) ≥ x ρ ( u ) ≤ ρc if xminρ (u ) +w(t − tρ (u )) ≤ x xminρ (u ) +v(t − tρ (u )) ≥ x ρ ( u ) ≤ ρc

(19)

if xminρ (u ) +w(t − tρ (u )) ≤ x xmaxρ (u ) +w(t − tρ (u )) ≥ x ρ ( u ) ≥ ρc if xmaxρ (u ) +w(t − tρ (u )) ≤ x xmaxρ (u ) +v(t − tρ (u )) ≥ x ρ ( u ) ≥ ρc

3.3. Properties of the affine initial, upstream/downstream boundary and internal conditions In

this

section,

we

show

that

the

LWR

model

constraints

(ρ (1 ), ρ (2 ), . . . , ρ (kmax , qin (1 ), . . . , qin (nmax ), qout (1 ), . . . , qout (nmax )).

(9)

are

convex

in

the

variable

Proposition 3.2 (Linearity property of the initial, upstream/downstream boundary and internal conditions). Let us fix (t, x ) ∈ R+ × [ξ , χ ]. The initial, upstream and downstream boundary condition functions Mk (t, x), γ n (t, x) and β n (t, x) are linear functions of the coefficients (ρ (1 ), . . . , ρ (kmax ), qin (1 ), . . . , qin (nmax ), qout (1 ), . . . , qout (nmax )). The proof of this proposition is straightforward and follows directly from Eqs. (10), (14) and (12). Proposition 3.3 (Concavity property of the solution associated with the initial condition). Let us fix (t, x ) ∈ R+ ∈ [ξ , χ ]. The solution MMk (t, x ) associated with the initial condition (10) is a concave function of the coefficients ρ ( · ). Proof. The Lax-Hopf formula (6) associated with the solution MMk (·, · ) can be written (C.G.Claudel and Bayen, 2010; Claudel and Bayen, 2010) as:

MMk (t, x ) =

inf

u∈Dom(ϕ ∗ ) s. t. (x+tu )∈[kX,(k+1 )X]

−

k−1 

 ρ (i )X − ρ (k )(x − kX ) + t ϕ ∗ (u )

(20)

i=0

−1 Let us fix u × Dom(ϕ ∗ ). The function f defined as f (ρ (1 ), . . . , ρ (kmax )) = − ki=0 ρ (i )X − ρ (k )(x − kX ) + t ϕ ∗ (u ) is concave (indeed, affine). Hence, the function MMk (t, x ) is a concave function of (ρ (1 ), . . . , ρ (kmax )), since it is the infimum of concave functions (Boyd and Vandenberghe, 2004; Rockafellar, 1970).  Proposition 3.4 (Concavity property of the solutions associated with upstream and downstream boundary conditions). Let us fix (t, x ) ∈ R+ × [ξ , χ ]. The solutions Mγn (t, x ) and Mβn (t, x ) respectively associated with the upstream and downstream boundary conditions (14) and (12) are concave functions of the coefficients ρ ( · ), qin ( · ) and qout ( · ). Please cite this article as: E.S. Canepa, C.G. Claudel, Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations, Transportation Research Part B (2017), http://dx.doi.org/10.1016/j.trb.2017.05.016

ARTICLE IN PRESS

JID: TRB 8

[m3Gsc;June 16, 2017;12:13]

E.S. Canepa, C.G. Claudel / Transportation Research Part B 000 (2017) 1–24

Proof. The Lax-Hopf formula (6) associated with the solution Mγn (·, · ) can be written (C.G.Claudel and Bayen, 2010; Claudel and Bayen, 2010) as:

Mγn (t, x ) =

inf

s∈R+ ∩[t −(n+1 )T,t −nT ]

n−1 

 qin (i )T + qin (n )(t − s ) + sϕ ∗

i=0

ξ −x

 (21)

s

Let us fix s ∈ R+ ∩ [t − (n + 1 )T , t − nT ]. The function d defined as d (qin (1 ), . . . , qin (nmax )) = sϕ ∗ ( ξ −x s )

n−1 i=0

qin (i )T + qin (n )(t − s ) +

is concave (indeed, affine). Hence, the solution Mγn (t, x ) is a concave function of (qin (1 ), . . . , qin (nmax )), since it is the infimum of concave functions (Boyd and Vandenberghe, 2004; Rockafellar, 1970). The same property applies for Mβn (t, x ), which is a concave function of (ρ (1 ), . . . , ρ (kmax ), qout (1 ), . . . , qout (nmax ))  Propositions 3.2, 3.3 and 3.4 thus imply the following convexity property: Proposition 3.5 (Convexity property of model constraints). The model constraints (9) are convex functions of (ρ (1 ), ρ (2 ), . . . , ρ (kmax ), qin (1 ), . . . , qin (nmax ), qout (1 ), . . . , qout (nmax )). Proof. The set of inequality constraints (9) can be written as:

Mci (t, x ) ≥ cj (t, x ), ∀(t, x ) ∈ Dom(cj ) ∀ j ∈ I such that(t, x ) ∈ Dom(c j ), ∀i ∈ I

(22)

Note that Proposition 3.2 implies that the term cj (t, x) in (22) is a linear function (labeled lj, t, x ( · )) of (ρ (1 ), ρ (2 ), . . . , ρ (kmax ), qin (1 ), . . . , qin (nmax ), qout (1 ), . . . , qout (nmax )). In addition, by Propositions 3.3 and 3.4, the term Mci (t, x ) is a concave function (labeled ci, t, x ( · )) of (ρ (1 ), ρ (2 ), . . . , ρ (kmax ), qin (1 ), . . . , qin (nmax ), qout (1 ), . . . , qout (nmax )). Hence, the equality (22) can be written as:

−ci,t,x (ρ (1 ), ρ (2 ), . . . , ρ (kmax ), qin (1 ), . . . , qin (nmax ), qout (1 ), . . . , qout (nmax ) ) +l j,t,x (ρ (1 ), ρ (2 ), . . . , ρ (kmax ), qin (1 ), . . . , qin (nmax ), qout (1 ), . . . , qout (nmax ) ) ≤ 0,

∀ j ∈ I, ∀(t, x ) ∈ Dom(c j ), ∀i ∈ I

(23)

This last inequality is a convex inequality (Boyd and Vandenberghe, 2004) in (ρ (1 ), ρ (2 ), . . . , ρ (kmax ), qin (1 ), . . . , qin (nmax ), qout (1 ), . . . , qout (nmax )), that is, an inequality of the form f( · ) ≤ 0 where f( · ) is a convex function.  The above property is very important, and can be thought of as follows. Consider the vector space V of all parameters of the initial, upstream and downstream boundary conditions. Each point of this vector space corresponds to a known value condition (encompassing initial, upstream and downstream boundary conditions). However, the solution to the LWR PDE (1) associated with this arbitrary value condition will satisfy the value condition itself on its boundaries if and only if the model constraints (9) are satisfied. Proposition 3.5 essentially states that the set of value conditions compatible with the LWR PDE model is convex.

4. Formulation of the density estimation problem as a mixed integer linear program 4.0.1. Decision variable As outlined in Proposition 3.5, the variable (ρ (1 ), ρ (2 ), . . . , ρ (kmax ), qin (1 ), . . . , qin (nmax ), qout (1 ), . . . , qout (nmax )) plays an important role in our estimation problem, and will be defined as the decision variable of our optimization framework. Definition 4.1 (Decision variable). Let us consider a finite set of, initial, upstream and downstream boundary conditions to be defined as in (10) (14) and (12). The decision variable v associated with this finite set of value conditions is defined by:

v := (ρ (1 ), ρ (2 ), . . . , ρ (kmax ), qin (1 ), . . . , qin (nmax ), qout (1 ), . . . , qout (nmax ) )

(24)

We denote by V the vector space of the decision variables v defined by Eq. (24).

4.0.2. Model and data constraints Let v denote the value of the decision variable associated with the true state of the system (which is not known in practice, and can only be estimated). Because of model and data constraints, v must satisfy the set of constraints outlined in Propositions 4.2 and 4.4 below. Proposition 4.2 (Model constraints). The model constraints (9) can be expressed as the following finite set of convex inequality constraints: Please cite this article as: E.S. Canepa, C.G. Claudel, Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations, Transportation Research Part B (2017), http://dx.doi.org/10.1016/j.trb.2017.05.016

ARTICLE IN PRESS

JID: TRB

E.S. Canepa, C.G. Claudel / Transportation Research Part B 000 (2017) 1–24

[m3Gsc;June 16, 2017;12:13] 9

⎧ MMk (t0 , x p ) ≥ M p (t0 , x p ) ⎪ ⎪ ⎪ ⎪ ∀(k, p) ∈ K2 (i ) ⎪ ⎪ ⎪ M ( pT , χ ) ≥ β ( pT , χ ) p Mk ⎪ ⎪ ⎪ ∀ ( k, p ) ∈ K × N ( ii )( a ) ⎪ ⎪ ⎨M (t + χ −xk+1 , χ ) ≥ β (t + χ −xk+1 , χ ) p 0 Mk 0 v v ∀(k, p) ∈ K × N s. t. t0 + χ −xv k+1 ∈ [ pT , ( p + 1 )T ] (ii )(b) ⎪ ⎪ ⎪ MMk ( pT , ξ ) ≥ γ p ( pT , ξ ) ⎪ ⎪ ⎪ ⎪∀(k, p) ∈ K × N (iii )(a ) ⎪ ⎪ ⎪ k k M (t + ξ −x , ξ ) ≥ γ p (t0 + ξ −x ,ξ) ⎪ w w ⎪ ⎩ Mk 0 k ∀(k, p) ∈ K × N s. t. t0 + ξ −x ∈ [ pT , ( p + 1 )T ] (iii )(b) w

(25)

⎧ MMk (tmin (m ), xmin (m )) ≥ μm (tmin (m ), xmin (m )) ⎪ ⎪ ⎪∀k ∈ K, ∀m ∈ M (iv )(a ) ⎪ ⎪ ⎪ ⎪MMk (tmax (m ), xmax (m )) ≥ μm (tmax (m ), xmax (m )) ⎪ ⎪ ⎪ ⎪∀k ∈ K, ∀m ∈ M (iv )(b) ⎪ ⎪ ⎪ ⎨MMk (t1 (m, k ), x1 (m, k )) ≥ μm (t1 (m, k ), x1 (m, k )) ∀k ∈ K, ∀m ∈ M s. t. ∀t1 (m, k ) ∈ [tmin (m ), tmax (m )] MMk (t2 (m, k ), x2 (m, k )) ≥ μm (t2 (m, k ), x2 (m, k )) ⎪ ⎪ ⎪ ∀k ∈ K, ∀m ∈ M s. t. ∀t2 (m, k ) ∈ [tmin (m ), tmax (m )] ⎪ ⎪ ⎪ ⎪ M ( t3 (m, k ), x3 (m, k )) ≥ μm (t3 (m, k ), x3 (m, k )) M ⎪ k ⎪ ⎪ ∀ k ∈ K, ∀m ∈ M s. t. ∀t3 (m, k ) ∈ [tmin (m ), tmax (m )] ⎪ ⎪ ⎪ ⎪ ⎩MMk (t4 (m, k ), x4 (m, k )) ≥ μm (t4 (m, k ), x4 (m, k )) ∀k ∈ K, ∀m ∈ M s. t. ∀t4 (m, k ) ∈ [tmin (m ), tmax (m )]

(26)

⎧ MMk (tρ (u ), xminρ (u )) ≥ ϒu (tρ (u ), xminρ (u )) ⎪ ⎪ ⎪ ⎪ ∀k ∈ K, ∀u ∈ U (v )(a ) ⎪ ⎪ ⎪ MMk (tρ (u ), xmaxρ (u )) ≥ ϒu (tρ (u ), xmaxρ (u )) ⎪ ⎪ ⎪ ∀k ∈ K, ∀u ∈ U (v )(b) ⎪ ⎪ ⎪ ⎪ ⎪ ⎨MMk (tρ (u ), x5 (u, k )) ≥ ϒu (tρ (u ), x5 (u, k )) ∀k ∈ K, ∀u ∈ U s. t. x5 (u, k ) ∈ [xminρ (u ), xmaxρ (u )] MMk (tρ (u ), x6 (u, k )) ≥ ϒu (tρ (u ), x6 (u, k )) ⎪ ⎪ ⎪ ∀k ∈ K, ∀u ∈ U s. t. x6 (u, k ) ∈ [xminρ (u ), xmaxρ (u )] ⎪ ⎪ ⎪ ⎪ M ( Mk tρ (u ), x7 (u, k )) ≥ ϒu (tρ (u ), x7 (u, k )) ⎪ ⎪ ⎪ ∀k ∈ K, ∀u ∈ U s. t. x7 (u, k ) ∈ [xminρ (u ), xmaxρ (u )] ⎪ ⎪ ⎪ ⎪ MMk (tρ (u ), x8 (u, k )) ≥ ϒu (tρ (u ), x8 (u, k )) ⎪ ⎩ ∀k ∈ K, ∀u ∈ U s. t. x8 (u, k ) ∈ [xminρ (u ), xmaxρ (u )] ⎧ Mγn ( pT , ξ ) ≥ γ p ( pT , ξ ) ⎪ ⎪ ⎪ ⎨Mγn ( pT , χ ) ≥ β p ( pT , χ ) ξ Mγn (nT + χ − v , χ ) ≥ β p (nT + ⎪ ⎪ ⎪ ⎩

(iv )(c ) (iv )(d ) (iv )(e ) (iv )( f )

(v )(c )

(27)

(v )(d ) (v )(e ) (v )( f )

∀(n, p) ∈ N2 (vi ) ∀(n, p) ∈ N2 (vii )(a ) χ −ξ 2 v , χ ) ∀ (n, p) ∈ N s. t. nT + χ −ξ v ∈ [ pT , ( p + 1 )T ] (vii )(b)

(28)

⎧ Mγn (tmin (m ), xmin (m )) ≥ μm (tmin (m ), xmin (m )) ⎪ ⎪ ⎪∀n ∈ N, ∀m ∈ M (viii )(a ) ⎪ ⎨ Mγn (tmax (m ), xmax (m )) ≥ μm (tmax (m ), xmax (m )) ∀n ∈ N, ∀m ∈ M (viii )(b) ⎪ ⎪ ⎪ ⎪ ⎩Mγn (t9 (m, n ), x9 (m, n )) ≥ μm (t9 (m, n ), x9 (m, n )) ∀n ∈ N, ∀m ∈ M s. t. t9 (m, n ) ∈ [tmin (m ); tmax (m )] (viii )(c )

(29)

⎧ Mγn (tρ (u ), xminρ (u )) ≥ ϒu (tρ (u ), xminρ (u )) ⎪ ⎪ ⎪ ⎪ ⎨∀n ∈ N, ∀u ∈ U (ix )(a ) Mγn (tρ (u ), xmaxρ (u )) ≥ ϒu (tρ (u ), xmaxρ (u )) ∀n ∈ N, ∀u ∈ U (ix )(b) ⎪ ⎪ ⎪ ⎪ ⎩Mγn (tρ (u ), x10 (u, n )) ≥ ϒu (tρ (u ), x10 (u, n )) ∀n ∈ N, ∀u ∈ U s. t. x10 (u, n ) ∈ [xminρ (u ); xmaxρ (u )] (ix )(c )

(30)

Please cite this article as: E.S. Canepa, C.G. Claudel, Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations, Transportation Research Part B (2017), http://dx.doi.org/10.1016/j.trb.2017.05.016

JID: TRB 10

ARTICLE IN PRESS

[m3Gsc;June 16, 2017;12:13]

E.S. Canepa, C.G. Claudel / Transportation Research Part B 000 (2017) 1–24

⎧ Mβn ( pT , ξ ) ≥ γ p ( pT , ξ ) ⎪ ⎪ ⎪ ∀(n, p) ∈ N2 (x )(a ) ⎪ ⎪ ⎨ Mβn (nT + ξ −wχ , ξ ) ≥ γ p (nT + ξ −wχ , ξ ) ∀(n, p) ∈ N2 s. t. nT + ξ −wχ ∈ [ pT , ( p + 1 )T ] (x )(b) ⎪ ⎪ ⎪ ⎪ ⎪ ⎩Mβn ( pT , χ )2≥ β p ( pT , χ ) ∀(n, p) ∈ N (xi ) ⎧ Mβn (tmin (m ), xmin (m )) ≥ μm (tmin (m ), xmin (m )) ⎪ ⎪ ⎪∀n ∈ N, ∀m ∈ M (xii )(a ) ⎪ ⎪ ⎨ Mβn (tmax (m ), xmax (m )) ≥ μm (tmax (m ), xmax (m )) ∀n ∈ N, ∀m ∈ M (xii )(b) ⎪ ⎪ ⎪ ⎪ ⎪ ⎩Mβn (t11 (m, n ), x11 (m, n )) ≥ μm (t11 (m, n ), x11 (m, n )) ∀n ∈ N, ∀m ∈ M s. t. t11 (m, n ) ∈ [tmin (m ); tmax (m )] (xii )(c ) ⎧ Mβn (tρ (u ), xminρ (u )) ≥ ϒu (tρ (u ), xminρ (u )) ⎪ ⎪ ⎪ ⎪ ∀n ∈ N, ∀u ∈ U (xiii )(a ) ⎪ ⎨ Mβn (tρ (u ), xmaxρ (u )) ≥ ϒu (tρ (u ), xmaxρ (u )) ∀n ∈ N, ∀u ∈ U (xiii )(b) ⎪ ⎪ ⎪ ⎪ ⎪ ⎩Mβn (tρ (u ), x12 (u, n )) ≥ ϒu (tρ (u ), x12 (u, n )) ∀n ∈ N, ∀m ∈ M s. t. x12 (m, n ) ∈ [xminρ (u ); xmaxρ (u )] (xiii )(c ) ⎧ Mμm ( pT , ξ ) ≥ γ p ( pT , ξ ) ∀(m, p) ∈ M × N(xiv )(a ) ⎪ ⎪ ⎪Mμ (t13 (m ), ξ ) ≥ γ p (t13 (m ), ξ ) ∀(m, p) ∈ M × N s. t. ⎨ m t13 (m ) ∈ [ pT , ( p + 1 )T ] (xiv )(b) ⎪ ⎪ Mμm (t14 (m ), ξ ) ≥ γ p (t14 (m ), ξ ) ∀(m, p) ∈ M × N s. t. ⎪ ⎩ t14 (m ) ∈ [ pT , ( p + 1 )T ] (xiv )(c ) ⎧ Mμm (pT, χ ) ≥ βp (pT, χ ) ∀(m, p) ∈ M × N (xv )(a ) ⎪ ⎪ ⎪ ⎨Mμm (t15 (m ), χ ) ≥ β p (t15 (m ), χ ) ∀(m, p) ∈ M × N s. t. t9 (m ) ∈ [ pT , ( p + 1 )T ] (xv )(b) ⎪ ⎪ ⎪Mμm (t16 (m ), χ ) ≥ β p (t16 (m ), χ ) ∀(m, p) ∈ M × N s. t. ⎩ t10 (m ) ∈ [ pT , ( p + 1 )T ] (xv )(c ) ⎧ Mμm (tmin (p ), xmin (p )) ≥ μp (tmin (p ), xmin (p )) ⎪ ⎪ ⎪ ⎪ ∀(m, p) ∈ M2 (xvi )(a ) ⎪ ⎪ ⎪ ⎪Mμm (tmax ( p), xmax ( p)) ≥ μ p (tmax ( p), xmax ( p)) ⎪ ⎪ ⎪ ∀(m, p) ∈ M2 (xvi )(b) ⎪ ⎪ ⎪ ⎪ Mμm (t17 (m, p), x17 (m, p)) ≥ μ p (t17 (m, p), x17 (m, p)) ⎪ ⎪ ⎪ ⎪ ∀(m, p) ∈ M2 s. t. t17 (m, p) ∈ [tmin ( p), tmax ( p)] (xvi )(c ) ⎪ ⎨ Mμm (t18 (m, p), x18 (m, p)) ≥ μ p (t18 (m, p), x18 (m, p)) ∀(m, p) ∈ M2 s. t. t18 (m, p) ∈ [tmin ( p), tmax ( p)] (xvi )(d ) ⎪ ⎪ ⎪ ⎪ ⎪Mμm (t19 (m, p), x19 (m, p)) ≥ μ p (t19 (m, p), x19 (m, p)) ⎪ ⎪ 2 ⎪ ⎪∀(m, p) ∈ M s. t. t19 (m, p) ∈ [tmin ( p), tmax ( p)] (xvi )(e ) ⎪ ⎪ M ( t ( m, p ), x20 (m, p)) ≥ μ p (t20 (m, p), x20 (m, p)) ⎪ μm 20 ⎪ ⎪ ⎪∀(m, p) ∈ M2 s. t. t20 (m, p) ∈ [tmin ( p), tmax ( p)] (xvi )( f ) ⎪ ⎪ ⎪ ⎪Mμm (t21 (m, p), x21 (m, p)) ≥ μ p (t21 (m, p), x21 (m, p)) ⎪ ⎩ ∀(m, p) ∈ M2 s. t. t21 (m, p) ∈ [tmin ( p), tmax ( p)] (xvi )(g) ⎧ Mμm (tρ (u ), xminρ (u )) ≥ ϒu (tρ (u ), xminρ (u )) ⎪ ⎪ ⎪∀m ∈ M, ∀u ∈ U (xvii )(a ) ⎪ ⎪ ⎪ ⎪ Mμm (tρ (u ), xmaxρ (u )) ≥ ϒu (tρ (u ), xmaxρ (u )) ⎪ ⎪ ⎪ ⎪ ∀m ∈ M, ∀u ∈ U (xvii )(b) ⎪ ⎪ ⎪ ⎪ ⎪Mμm (tρ (u ), x22 (m, u )) ≥ ϒu (tρ (u ), x22 (m, u )) ⎪ ⎪ ⎪ ∀m ∈ M, ∀u ∈ U s. t. x22 (m, u ) ∈ [xminρ (u ), xmaxρ (u )] (xvii )(c ) ⎪ ⎨ Mμm (tρ (u ), x23 (m, u )) ≥ ϒu (tρ (u ), x23 (m, u )) ∀m ∈ M, ∀u ∈ U s. t. x23 (m, u ) ∈ [xminρ (u ), xmaxρ (u )] (xvii )(d ) ⎪ ⎪ ⎪ ⎪ Mμm (tρ (u ), x24 (m, u )) ≥ ϒu (tρ (u ), x24 (m, u )) ⎪ ⎪ ⎪ ⎪ ∀m ∈ M, ∀u ∈ U s. t. x24 (m, u ) ∈ [xminρ (u ), xmaxρ (u )] (xvii )(e ) ⎪ ⎪ ⎪ Mμm (tρ (u ), x25 (m, u )) ≥ ϒu (tρ (u ), x25 (m, u )) ⎪ ⎪ ⎪ ⎪ ∀m ∈ M, ∀u ∈ U s. t. x25 (m, u ) ∈ [xminρ (u ), xmaxρ (u )] (xvii )( f ) ⎪ ⎪ ⎪ ⎪ ⎪Mμm (tρ (u ), x26 (m, u )) ≥ ϒu (tρ (u ), x26 (m, u )) ⎩ ∀m ∈ M, ∀u ∈ U s. t. x26 (m, u ) ∈ [xminρ (u ), xmaxρ (u )] (xvii )(g)

(31)

(32)

(33)

(34)

(35)

(36)

(37)

Please cite this article as: E.S. Canepa, C.G. Claudel, Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations, Transportation Research Part B (2017), http://dx.doi.org/10.1016/j.trb.2017.05.016

JID: TRB

ARTICLE IN PRESS E.S. Canepa, C.G. Claudel / Transportation Research Part B 000 (2017) 1–24

⎧ Mϒu ( pT , χ ) ≥ β p ( pT , χ ) ⎪ ⎪ ⎪ ∀(u, p) ∈ U × N (xviii )(a ) ⎪ ⎪ ⎪ χ −xmaxρ (u ) ⎪ ρ (u ) ⎪Mϒu (tρ (u ) + χ −xmax , χ ) ≥ β p (tρ (u ) + ,χ) ⎪ v v ⎪ ⎨ χ −xmaxρ (u ) ∀(u, p) ∈ U × N s. t. tρ (u ) + ∈ [ pT , ( p + 1 )T ] (xviii )(b) v Mϒu ( pT , ξ ) ≥ γ p ( pT , ξ ) ⎪ ⎪ ⎪ ⎪ ⎪∀(u, p) ∈ U × N (xviii )(c ) ⎪ ξ −xminρ (u ) ξ −xminρ (u ) ⎪ ⎪ Mϒ (tρ (u ) + , ξ ) ≥ γ p (tρ (u ) + ,ξ) ⎪ w w ⎪ ⎩ u ξ −xminρ (u ) ∀(u, p) ∈ K × N s. t. tρ (u ) + ∈ [ pT , ( p + 1 )T ] (xviii )(d ) w ⎧ Mϒu (tmin ( p), xmin ( p)) ≥ μ p (tmin ( p), xmin ( p)) ⎪ ⎪ ⎪ ⎪ ∀(u, p) ∈ U × M (xix )(a ) ⎪ ⎪ ⎪ Mϒu (tmax ( p), xmax ( p)) ≥ μ p (tmax ( p), xmax ( p)) ⎪ ⎪ ⎪ ⎪ ⎪∀(m, p) ∈ U × M (xix )(b) ⎪ ⎪ ⎪ Mϒu (t27 (u, p), x27 (u, p)) ≥ μ p (t27 (u, p), x27 (u, p)) ⎪ ⎪ ⎪ ⎪ ∀(u, p) ∈ U × M s. t. t27 (u, p) ∈ [tmin ( p), tmax ( p)] (xix )(c ) ⎪ ⎨ Mϒu (t28 (u, p), x28 (u, p)) ≥ μ p (t28 (u, p), x28 (u, p)) ∀(u, p) ∈ U × M s. t. t28 (u, p) ∈ [tmin ( p), tmax ( p)] (xix )(d ) ⎪ ⎪ ⎪ ⎪ Mϒu (t29 (u, p), x29 (u, p)) ≥ μ p (t29 (u, p), x29 (u, p)) ⎪ ⎪ ⎪ ⎪ ∀(u, p) ∈ U × M s. t. t29 (u, p) ∈ [tmin ( p), tmax ( p)] (xix )(e ) ⎪ ⎪ ⎪ Mϒu (t30 (u, p), x30 (u, p)) ≥ μ p (t30 (u, p), x30 (u, p)) ⎪ ⎪ ⎪ ⎪ ⎪∀(u, p) ∈ U × M s. t. t30 (u, p) ∈ [tmin ( p), tmax ( p)] (xix )( f ) ⎪ ⎪ ⎪Mϒu (tρ (u ), x26 ( p, u )) ≥ μ p (tρ (u ), x26 ( p, u )) ⎪ ⎩ ∀(u, p) ∈ U × M s. t. x26 ( p, u ) ∈ [xminρ (u ), tmaxρ (u )] (xix )(g) ⎧ Mϒu (tρ ( p), xminρ ( p)) ≥ ϒ p (tρ ( p), xminρ ( p)) ⎪ ⎪ ⎪∀(u, p) ∈ U2 (xx )(a ) ⎪ ⎪ ⎪ ⎪ Mϒu (tρ ( p), xmaxρ ( p)) ≥ ϒ p (tρ ( p), xmaxρ ( p)) ⎪ ⎪ ⎪ ⎪ ∀(u, p) ∈ U2 (xx )(b) ⎪ ⎪ ⎪ ⎪ ⎪Mϒu (tρ ( p), x31 (u, p)) ≥ μ p (tρ ( p), x31 (u, p)) ⎨ ∀(u, p) ∈ U2 s. t. x31 (u, p) ∈ [xminρ ( p), tmaxρ ( p)] (xx )(c ) ⎪Mϒu (tρ ( p), x32 (u, p)) ≥ μ p (tρ ( p), x32 (u, p)) ⎪ ⎪ ⎪ ∀(u, p) ∈ U2 s. t. x32 (u, p) ∈ [xminρ ( p), tmaxρ ( p)] (xx )(d ) ⎪ ⎪ ⎪ ⎪ Mϒu (tρ ( p), x33 (u, p)) ≥ μ p (tρ ( p), x33 (u, p)) ⎪ ⎪ ⎪∀(u, p) ∈ U2 s. t. x33 (u, p) ∈ [xmin ( p), tmaxρ ( p)] (xx )(e ) ⎪ ρ ⎪ ⎪ ⎪ ⎪ ⎩Mϒu (tρ ( p),2x34 (u, p)) ≥ μ p (tρ ( p), x34 (u, p)) ∀(u, p) ∈ U s. t. x34 (u, p) ∈ [xminρ ( p), tmaxρ ( p)] (xx )( f )

[m3Gsc;June 16, 2017;12:13] 11

(38)

(39)

(40)

where the coefficients t1 (m, k), x1 (m, k), t2 (m, k), x2 (m, k), t3 (m, k), x3 (m, k), t4 (m, k), x4 (m, k), x5 (u, k), x6 (u, k), x7 (u, k), x8 (u, k), t9 (m, n), x9 (m, n), x10 (u, n), t11 (m, n), x11 (m, n), x12 (u, n), t13 (m), t14 (m), t15 (m), t16 (m), t17 (m, p), x17 (m, p), t18 (m, p), x18 (m, p), t19 (m, p), x19 (m, p), t20 (m, p), x20 (m, p), t21 (m, p), x21 (m, p), x22 (m, u), x23 (m, u), x24 (m, u), x25 (m, u), x26 (m, u), t27 (u, p), x27 (u, p), t28 (u, p), x28 (u, p), t29 (u, p), x29 (u, p), t30 (u, p), x30 (u, p), x31 (u, p), x32 (u, p), x33 (u, p) and x34 (u, p) are given by Eqs. (41)–(44) below:

⎧ x (m )−xk+1 −vmeas (m )tmin (m ) ⎪ t1 (m, k ) = min ⎪ v− vmeas (m ) ⎪ ⎪ xmin (m )−xk+1 −vmeas (m )tmin (m ) ⎪ meas ⎪ x ( m, k ) = v ( m ) 1 ⎪ v−vmeas (m ) ⎪ ⎪ ⎪ −t ( m ) + x ( m ) ⎪ ) min min ⎪ meas ⎪ (m )tmin (m ) k −v ⎪ t2 (m, k ) = xmin (m )−x ⎪ w−vmeas (m ) ⎪ ⎪ xmin (m )−xk −vmeas (m )tmin (m ) ⎪ meas ⎪ x2 (m, k ) = v (m ) ⎪ w−vmeas (m ) ⎪ ⎪ −t ( m ) + x ⎪ ) min min (m ) ⎪ ⎪ vmeas (m )tmin (m ) ⎪ ⎨t3 (m, k ) = xmin (m)−xv−k −vmeas (m ) k −vmeas (m )tmin (m ) x3 (m, k ) = vmeas (m ) xmin (m )−x v−vmeas (m ) ⎪ ⎪ −tmin (m ) ) + xmin (m ) ⎪ ⎪ ⎪ vmeas (m )tmin (m ) ⎪t4 (m, k ) = xmin (m)−xk+1 −meas ⎪ ⎪ w−v (m ) ⎪ ⎪x (m, k ) = vmeas (m ) xmin (m)−xk −vmeas (m)tmin (m) ⎪ 4 ⎪ w−vmeas (m ) ⎪ ⎪ −tmin (m ) ) + xmin (m ) ⎪ ⎪ ⎪ ⎪x5 (u, k ) = xk+1 + v(tρ (u ) − t0 ) ⎪ ⎪ ⎪ ⎪ ⎪x6 (u, k ) = xk+1 + w(tρ (u ) − t0 ) ⎪ ⎪ ⎪ ⎩x7 (u, k ) = xk + v(tρ (u ) − t0 ) x8 (u, k ) = xk + w(tρ (u ) − t0 )

(41)

Please cite this article as: E.S. Canepa, C.G. Claudel, Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations, Transportation Research Part B (2017), http://dx.doi.org/10.1016/j.trb.2017.05.016

JID: TRB 12

ARTICLE IN PRESS E.S. Canepa, C.G. Claudel / Transportation Research Part B 000 (2017) 1–24

⎧ nT v−vmeas (m )tmin (m )+xmin (m )−ξ ⎪ v−vmeas (m ) ⎪t9 (m, n ) = ⎪ ⎪ x9 (m, n ) = xmin (m )+ ⎪ ⎪ meas )tmin (m )+xmin (m )−ξ ⎪ ⎪ vmeas (m ) nT v−v (m − tmin (m ) meas (m ) ⎪ v − v ⎪ ⎪ ⎪ x10 (u, n ) = ξ + v(tρ (u ) − nT ) ⎪ ⎪ meas (m )+xmin (m )−χ ⎪ ⎪ t (m, n ) = nT w−v (wm−)tvmin meas (m ) ⎪ ⎨ 11 x11 (m, n )= xmin (m )+ meas (m )+xmin (m )−χ vmeas (m ) nT w−v (wm−)tvmin − tmin (m ) ⎪ meas (m ) ⎪ ⎪ ⎪ x12 (u, n ) = χ + w(tρ (u ) − nT ) ⎪ ⎪ ⎪ ⎪ t13 (m ) = ξ −xmin (m )w+wtmin (m ) ⎪ ⎪ ⎪ ⎪ t14 (m ) = ξ −xmax (m )w+wtmax (m ) ⎪ ⎪ ⎪ χ −xmin (m )+vtmin (m ) ⎪ ⎪ ⎩t15 (m ) = χ −xmax (mv)+vtmax (m) t16 (m ) = v and

[m3Gsc;June 16, 2017;12:13]

⎧ vmeas ( p)tmin ( p)−vmeas (m )tmin (m ) t17 (m, p) = xmin (m )−xmin ( p)v+meas ⎪ ( p)−vmeas (m ) ⎪ ⎪ ⎪ x17 (m, p) = xmin ( p) + vmeas ( p)(−tmin ( p)+ ⎪ ⎪ xmin (m )−xmin ( p)+vmeas ( p)tmin ( p)−vmeas (m )tmin (m ) ⎪ ⎪ vmeas ( p)−vmeas (m ) ⎪ ⎪ xmax (m )−xmin ( p)+vmeas ( p)tmin ( p)−vtmax (m ) ⎪ ⎪ t ( m, p ) = 18 ⎪ vmeas ( p)−v ⎪ ⎪ meas ⎪ x ( m, p ) = x ( p ) + v ( p)(−t min ( p)+ 18 min ⎪ ⎪ xmax (m )−xmin ( p)+vmeas ( p)tmin ( p)−vtmax (m ) ⎪ ⎪ meas v ( p ) −v ⎪ ⎪ x (m )−xmin ( p)+vmeas ( p)tmin ( p)−vtmin (m ) ⎪ ⎪ ⎪t19 (m, p) = min vmeas ( p)−v ⎪ ⎪ meas ⎪ x ( m, p ) = x ( p ) + v 19 min ⎪ x (m)−x ( p)+vmeas ( p)t ( p)−vt( p()m(−t min ( p)+ ⎪ ) ⎪ min min min min ⎪ meas v ( p ) −v ⎨ vmeas ( p)tmin ( p)−vtmax (m ) t20 (m, p) = xmax (m )−xmin ( p)v+meas ( p)−w ⎪ meas ⎪ x20 (m, p) = xmin ( p) + v ( p)(−t ( p )+ ⎪ ⎪ xmax (m)−xmin ( p)+vmeas ( p)tmin ( p)−vtmax (m) min ⎪ ⎪ meas v ( p)−w ⎪ ⎪ vmeas ( p)tmin ( p)−vtmin (m ) ⎪ ⎪ t21 (m, p) = xmin (m )−xmin ( p)v+meas ⎪ ( p)−w ⎪ ⎪ meas ⎪ x21 (m, p) = xmin ( p) + v ( p)(−t ( p)+ ⎪ ⎪ xmin (m)−xmin ( p)+vmeas ( p)tmin ( p)−vtmin (m) min ⎪ ⎪ ⎪ vmeas ( p)−w ⎪ ⎪ x22 (m, u ) = xmin (m ) + w(tρ (u ) − tmin (m )) ⎪ ⎪ ⎪ ⎪x23 (m, u ) = xmax (m ) + w(tρ (u ) − tmax (m )) ⎪ ⎪ ⎪ ⎪x24 (m, u ) = xmin (m ) + v(tρ (u ) − tmin (m )) ⎪ ⎪ ⎪ ⎩x25 (m, u ) = xmax (m ) + v(tρ (u ) − tmax (m )) x26 (m, u ) = xmin (m ) + vmeas (m )(tρ (u ) − tmin (m )) ⎧ x ( p)−xmaxρ (u )−vmeas ( p)tmin ( p)+vtρ (u ) t27 (u, p) = min ⎪ v−vmeas ( p) ⎪ ⎪ ⎪ x27 (u, p) = xmin ( p) + vmeas ( p)(−t ⎪  min ( p)+ ⎪ ⎪ ⎪ xmin ( p)−xmaxρ (u)−vmeas ( p)tmin ( p)+vtρ (u) ⎪ ⎪ v−vmeas ( p) ⎪ ⎪ ⎪ xmin ( p)−xminρ (u )−vmeas ( p)tmin ( p)+wtρ (u ) ⎪ ⎪ t28 (u, p) = ⎪ w−vmeas ( p) ⎪ meas ⎪ x ( u, p ) = x ( p ) + v ( p)(−t ( p)+ ⎪ 28 min ⎪ ⎪ xmin ( p)−xminρ (u)−vmeas ( p)tmin ( p)+wtρ (u)  min ⎪ ⎪ ⎪ w−vmeas ( p) ⎪ ⎪ ⎪ xmin ( p)−xminρ (u )−vmeas ( p)tmin ( p)+vtρ (u ) ⎪ ⎪t29 (u, p) = ⎨ v−vmeas ( p) x29 (u, p) = xmin ( p) + vmeas ( p)(−t  min ( p)+ ⎪ xmin ( p)−xminρ (u )−vmeas ( p)tmin ( p)+vtρ (u ) ⎪ ⎪ ⎪ v−vmeas ( p) ⎪ ⎪ xmin ( p)−xmaxρ (u )−vmeas ( p)tmin ( p)+wtρ (u ) ⎪ ⎪ t30 (u, p) = ⎪ w−vmeas ( p) ⎪ ⎪ ⎪ x30 (u, p) = xmin ( p) + vmeas ( p)(−t ⎪  min ( p)+ ⎪ ⎪ xmin ( p)−xmaxρ (u )−vmeas ( p)tmin ( p)+wtρ (u ) ⎪ ⎪ ⎪ w−vmeas ( p) ⎪ ⎪ ⎪ x ( u, p ) = x ⎪ 31 minρ (u ) + v (tρ ( p) − tρ (u )) ⎪ ⎪ ⎪ x32 (u, p) = xmaxρ (u ) + v(tρ ( p) − tρ (u )) ⎪ ⎪ ⎪ ⎪ ⎩x33 (u, p) = xminρ (u ) + w(tρ ( p) − tρ (u )) x34 (u, p) = xmaxρ (u ) + w(tρ ( p) − tρ (u ))

(42)

(43)

(44)

Proof. Note that ∀(k, n ) ∈ [0, kmax ] × [0, nmax ], Dom(Mk ) ∩ Dom(Mγn ) = ∅ and that ∀(k, n ) ∈ [0, kmax ] × [0, nmax ], Dom(Mk ) ∩ Dom(Mβn ) = ∅. Thus, the set of inequality constraints (9) can be written in the case of initial, upstream, downstream and Please cite this article as: E.S. Canepa, C.G. Claudel, Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations, Transportation Research Part B (2017), http://dx.doi.org/10.1016/j.trb.2017.05.016

JID: TRB

ARTICLE IN PRESS E.S. Canepa, C.G. Claudel / Transportation Research Part B 000 (2017) 1–24

[m3Gsc;June 16, 2017;12:13] 13

internal conditions as:

⎧ MMk (t, ξ ) ≥ γp (t, ξ ) ⎪ ⎪ ⎪ MMk (t, χ ) ≥ β p (t, χ ) ⎪ ⎪ ⎪ MMk (t, x ) ≥ μm (t, x ) ⎪ ⎪ ⎪ ⎪ MMk (t, x ) ≥ ϒu (t, x ) ⎪ ⎪ ⎪ Mγn (t, ξ ) ≥ γ p (t, ξ ) ⎪ ⎪ ⎪ ⎪ Mγn (t, χ ) ≥ β p (t, χ ) ⎪ ⎪ ⎪ Mγn (t, x ) ≥ μm (t, x ) ⎪ ⎪ ⎪ ⎪ Mγn (t, x ) ≥ ϒu (t, x ) ⎪ ⎪ ⎪ ⎪ ⎨Mβn (t, ξ ) ≥ γ p (t, ξ ) Mβn (t, ξ ) ≥ β p (t, χ ) Mβn (t, x ) ≥ μm (t, x ) ⎪ ⎪ ⎪ Mβn (t, x ) ≥ ϒu (t, x ) ⎪ ⎪ ⎪ ⎪ Mμk (t, ξ ) ≥ γ p (t, ξ ) ⎪ ⎪ ⎪ Mμk (t, χ ) ≥ β p (t, χ ) ⎪ ⎪ ⎪ ⎪ Mμk (t, x ) ≥ μm (t, x ) ⎪ ⎪ ⎪ Mμk (t, x ) ≥ ϒu (t, x ) ⎪ ⎪ ⎪ ⎪ Mϒk (t, ξ ) ≥ γ p (t, ξ ) ⎪ ⎪ ⎪ Mϒk (t, χ ) ≥ β p (t, χ ) ⎪ ⎪ ⎪ ⎪ ⎩Mϒk (t, x ) ≥ μm (t, x ) Mϒk (t, x ) ≥ ϒu (t, x )

∀t ∈ [ pT , ( p + 1 )T ], ∀k ∈ K, ∀ p ∈ N ∀t ∈ [ pT , ( p + 1 )T ], ∀k ∈ K, ∀ p ∈ N ∀(t, x ) ∈ Dom(μm ), ∀(k, m ) ∈ K × M ∀(t, x ) ∈ Dom(ϒu ), ∀(k, u ) ∈ K × U ∀t ∈ [ pT , ( p + 1 )T ], ∀(n, p) ∈ N2 ∀t ∈ [ pT , ( p + 1 )T ], ∀(n, p) ∈ N2 ∀(t, x ) ∈ Dom(μm ), ∀(n, m ) ∈ N × M ∀(t, x ) ∈ Dom(ϒu ), ∀(k, u ) ∈ K × U ∀t ∈ [ pT , ( p + 1 )T ], ∀(n, p) ∈ N2 ∀t ∈ [ pT , ( p + 1 )T ], ∀(n, p) ∈ N2 ∀(t, x ) ∈ Dom(μm ), ∀(n, m ) ∈ N × M ∀(t, x ) ∈ Dom(ϒu ), ∀(n, u ) ∈ N × U ∀t ∈ [ pT , ( p + 1 )T ], ∀k ∈ M, ∀ p ∈ N ∀t ∈ [ pT , ( p + 1 )T ], ∀k ∈ M, ∀ p ∈ N ∀(t, x ) ∈ Dom(μm ), ∀(k, m ) ∈ M2 ∀(t, x ) ∈ Dom(ϒu ), ∀k ∈ M, ∀u ∈ U ∀t ∈ [ pT , ( p + 1 )T ], ∀k ∈ U, ∀ p ∈ N ∀t ∈ [ pT , ( p + 1 )T ], ∀k ∈ U, ∀ p ∈ N ∀(t, x ) ∈ Dom(μm ), ∀k ∈ U, ∀m ∈ M ∀(t, x ) ∈ Dom(ϒu ), ∀(k, u ) ∈ U2

(45)

The conditions (45) all involve checking that a function of (t, x) is greater than another function of (t, x) on a line segment of R+ × [ξ , χ ]. Yet, because of the affine structure of the initial and boundary conditions (10), (14) and (12) as well as the piecewise affine structure of their solutions (15), (16) and (17), the inequalities of the form ∀(t, x ) ∈ Dom(ci ), Mc j (t, x ) ≥ ci (t, x ) are equivalent to a finite number of inequalities of the form ∀ p ∈ {0, . . . , pmax }, Mc j (t p , x p ) ≥ ci (t p , x p ). This arises from the fact that a piecewise affine function is positive on all points of a segment if and only if it is positive on each extremity of the segment, and on the finite number of points of the segment on which the function is not differentiable. In the present case, this property implies the equivalence of (45) and of (25)–(40). The equality constraints (47) arise from continuity conditions of the Moskowitz function (Newell, 1993).  Proposition 4.3 (Continuity constraints). Let a set of initial, boundary and internal conditions be defined as in (10), and let the corresponding partial solutions be defined as MMk (·, · ), Mγn (·, · ), Mβn (·, · ) and Mμm (·, · ). Let us also assume that the model constraints (9) are satisfied. Let Mp ( ·, ·) be defined as M p (·, · ) = mink,n,u,m|m = p (MMk (·, · ), Mγn (·, · ), Mβn (·, · ), Mϒu (·, · ), Mμm (·, · )) and Mo ( ·, ·) be defined as Mo (·, · ) = mink,n,m,u|u =o (MMk (·, · ), Mγn (·, · ), Mβn (·, · ), Mμm (·, · )), Mϒu (·, · ) . The solution M( ·, ·) to the HJ PDE (7) defined by M(·, · ) = mink,n,m,u (MMk (·, · ), Mγn (·, · ), Mβn (·, · ), Mμm (·, · ), Mϒu (·, · )) is continuous if and only if the following conditions are satisfied:

∀ p ∈ M, M p (tmin ( p), xmin ( p)) = μ p (tmin ( p), xmin ( p)) ⎧  μm (tmin (m ), xmin (m )) = min MMk (tmin (m ), xmin (m )), ⎪ ⎪ ⎪ ⎪Mγn (tmin (m ), xmin (m )), Mβn (tmin (m ), xmin (m )), ⎪ ⎪ ⎪ ⎪ Mϒu (tmin (m ), xmin (m )), Mμ p (tmin (m ), xmin (m )) ∀k ∈ K, ⎪ ⎨ ∀n ∈ N, ∀u ∈ U, ∀(m, p) ∈ M2  ⎪μm (tmax (m ), xmax (m )) = min MMk (tmax (m ), xmax (m )), ⎪ ⎪ ⎪ ⎪Mγn (tmax (m ), xmax (m )), Mβn (tmax (m ), xmax (m )), ⎪ ⎪ ⎪ M (t (m ), xmax (m )), Mμ p (tmax (m ), xmax (m )) ∀k ∈ K, ⎪ ⎩ ϒu max ∀n ∈ N, ∀u ∈ U, ∀(m, p) ∈ M2 ∀ p ∈ M, Mo (tρ (o), xminρ (o)) = ϒo (tρ (o), xminρ (o)) ⎧  ϒu (tρ (u ), xminρ (u )) = min MMk (tρ (u ), xminρ (u )), ⎪ ⎪ ⎪ ⎪ Mγn (tρ (u ), xminρ (u )), Mβn (tρ (u ), xminρ (u )), ⎪ ⎪ ⎪ M (t (u ), x (u )), M (t (u ), x (u )) ∀k ∈ K, ⎪ ϒo ρ minρ minρ ⎪ ⎨ μp ρ ∀n ∈ N, ∀ p ∈ M, ∀(o, u ) ∈ U2 ϒu (tρ (u ), xmaxρ (u )) = min MMk (tρ (u ), xmaxρ (u )), ⎪ ⎪ ⎪ ⎪Mγn (tρ (u ), xmaxρ (u )), Mβn (tρ (u ), xmaxρ (u )), ⎪ ⎪ ⎪ ⎪ Mμ p (tρ (u ), xmaxρ (u )), Mϒo (tρ (u ), xmaxρ (u )) ∀k ∈ K, ⎪ ⎩ ∀n ∈ N, ∀ p ∈ M, ∀(o, u ) ∈ U2

(46)

(47)

(48)

(49)

Please cite this article as: E.S. Canepa, C.G. Claudel, Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations, Transportation Research Part B (2017), http://dx.doi.org/10.1016/j.trb.2017.05.016

ARTICLE IN PRESS

JID: TRB 14

[m3Gsc;June 16, 2017;12:13]

E.S. Canepa, C.G. Claudel / Transportation Research Part B 000 (2017) 1–24

Furthermore, the equality constraints (46) and (48) can be written as a set of mixed integer linear inequalities involving the continuous variables ρ (1 ), ρ (2 ), . . . , ρ (kmax ), qin (1 ), . . . , qin (nmax ), qout (1 ), . . . , qout (nmax ), L1 , . . . , Lmmax and r1 , . . . , rmmax , as well as auxiliary integer variables. The proof of (46) is straightforward, and follows directly (Claudel and Bayen, 2011) from the piecewise affine structure of the partial solutions MMk (·, · ), Mγn (·, · ), Mβn (·, · ) and Mμm (·, · ). The fact that (46) can be written as a set of mixed integer linear inequalities is more involved. It can be shown that since M p (·, · ) = mink,n,m|m = p (MMk (·, · ), Mγn (·, · ), Mβn (·, · ), Mμm (·, · )), equation (46) can be written as a set of inequalities involving the continuous variables ρ (1 ), ρ (2 ), . . . , ρ (kmax ), qin (1 ), . . . , qin (nmax ), qout (1 ), . . . , qout (nmax ), L1 , . . . , Lmmax and r1 , . . . , rmmax , as well as boolean variables. An example of such derivation is shown in Canepa and Claudel (2012) for the case in which mmax = 1. These inequalities can be further rewritten as mixed integer linear inequalities using the piecewise affine dependency of the partial solutions with respect to the variables ρ (1 ), ρ (2 ), . . . , ρ (kmax ), qin (1 ), . . . , qin (nmax ), qout (1 ), . . . , qout (nmax ), L1 , . . . , Lmmax and r1 , . . . , rmmax . Proposition 4.4 (Data constraints). In the remainder of our article, we assume that the data constraints are convex in the decision variable v. Different choices of error models yield convex data constraints, such as the two examples outlined below. Example of convex data constraints (1) — Consider a sensor measuring the boundary flows (qin (0 ), . . . qin (nmax )) with 5% relative uncertainty, a loop detector measuring the initial density ρ (3) with 10% absolute uncertainty, and no downstream sensor. In this situation, the constraints are convex inequalities (indeed, linear inequalities) in the decision variable:



0.95qmeasured (n ) ≤ qin (n ) ≤ 1.05qmeasured (n ) ∀n ∈ [0, nmax ] in in measured ρ (3 ) − 0.1ρm ≤ ρ (3 ) ≤ ρ (3 )measured + 0.1ρm

(50)

Example of convex data constraints (2) — Consider two identical sensors measuring the boundary flows

(qin (0 ), . . . qin (nmax )), and (qout (1 ), . . . qout (nmax )) which are characterized by a RMS relative error of 3%. In this situation, the constraints are convex inequalities (quadratic convex inequalities) in the decision variable:



nmax n=0 nmax n=0

 2 max  measured 2 qin (n ) − qmeasured ≤ 0.03 nn=0 qin in  nmax  measured 2 measured 2 qout (n ) − qout ≤ 0.03 n=0 qout

(51)

In this situation the estimation problem becomes a Quadratic Program. Another option to represent the data constraints are the so-called chance constraints which could be used to enforce constraints with the knowledge of practical error distribution.However, we think that the computational time may be severely affected by using chance constraints, and thus we kept deterministic constraints. Proposition 4.5 (Travel time constraints). Travel time data can be used in the estimation process, in order to properly define this information we define the travel time as:

Ttime = t ftravel − t0travel

(52)

The travel time constraints are specified as the following equality:



Mγn (t0travel , ξ ) = Mβ p (t ftravel , χ )

∀(n, p) ∈ N2

(53)

5. Traffic estimation on a single highway link We now present an implementation of the estimation framework presented earlier on an experimental dataset. The dataset includes fixed sensor data (obtained from inductive loop detectors in the present case) travel time and mobile sensor data. 5.1. Experimental setup In the following sections, the effectiveness of the method is illustrated on different traffic flow estimation problems which are formulated as LPs and MIPs, using the Mobile Century (http://traffic.berkeley.edu/ ), (Work et al., 2010) dataset. The Mobile Century field experiment demonstrated the use of Nokia N-95 cellphones as mobile traffic sensors in 2008, and was a joint UC-Berkeley/Nokia project. For the numerical applications, a spatial domain of 1.2 km is considered, located between the PeMS (http://pems.eecs. berkeley.edu ) VDS stations 400,536 and 401,529 on the Highway I - 880 N around Hayward, California. The data used in this implementation was generated on February 8th, 2008, between times 18: 30 and 18: 55. In our scenario, we only consider inflow and outflow data qmeasured (· ) and qmeasured (· ) generated by the above PeMS stations, i.e. we do not assume out in to know any density data. Of course the framework presented in this article allows incorporation of density data, see for instance Example 1 in the previous section. The layout of the spatial domain is illustrated in Fig. 2. Please cite this article as: E.S. Canepa, C.G. Claudel, Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations, Transportation Research Part B (2017), http://dx.doi.org/10.1016/j.trb.2017.05.016

ARTICLE IN PRESS

JID: TRB

E.S. Canepa, C.G. Claudel / Transportation Research Part B 000 (2017) 1–24

[m3Gsc;June 16, 2017;12:13] 15

Fig. 2. Spatial domain considered for the numerical implementation. The upstream and downstream PeMS stations are delimiting a 1.2 km spatial domain, outlined by a solid line. The direction of traffic flow is represented by an arrow.

5.2. Initial density estimation on systems modeled by the Lighthill–Whitham–Richards PDE In this first scenario, our objective is to find the minimal or maximal values of a function of the decision variables, assuming that boundary flow data is available from the PeMS sensors. For this specific application the objective function max is chosen as the total number of vehicles at initial time, defined by ki=0 ρ (i ), though any convex piecewise affine function of the decision variable would be acceptable. The constraints are linear inequalities in (24), and comprise both model (25), (28) and (30) as well as data constraints. For this specific application, the data constraints are (1 − e )qmeasured (n ) ≤ in/out

(n ) ∀n ∈ [0, nmax ], where e = 0.01 = 1% is chosen the worst-case relative error of the sensors. qin/out (n ) ≤ (1 + e )qmeasured in/out The maximal densities are solutions to the following linear program: kmax

ρ (i )

Minimize

(respectively Maximize )

such that

⎧ (25 ) ⎪ ⎪ ⎪ (28 ) ⎪ ⎪ ⎪ ⎨ (30 ) (1 − e )qmeasured (n ) ≤ qin (n ) ∀n ∈ [0, nmax ] in ⎪ ⎪ q ( n ) ≤ ( 1 + e ) qmeasured (n ) ∀n ∈ [0, nmax ] in ⎪ in ⎪ measured ⎪ ( 1 − e ) q ( n ) ≤ q (n ) ∀n ∈ [0, nmax ] ⎪ out out ⎩ qout (n ) ≤ (1 + e )qmeasured (n ) ∀n ∈ [0, nmax ] out

i=0

(54)

For this implementation, we chose 20 pieces of upstream and downstream flow data corresponding to 10 min of data. The parameters of the triangular flux function are set with standard values: v = 65mph, w = −10mph, ρc = 30veh/(l ane.mil e ). The optimal solutions to (54) are illustrated in Fig. 3.

6. Extension to highway networks 6.1. Junction model In this section, we generalize the above framework to traffic state estimation on road networks. For this, we first need to derive the boundary flows occurring at the junctions. This process is done through a junction model, which we now outline. Proposition 6.1 (Conservation of vehicles). We consider a general junction (without storage capacity) integrating on-ramps, off-ramps and incoming as well as outgoing links, and illustrated in Fig. 4. The conservation of vehicles at the junction imposes the following equality constraint.



qouti +



qramp = in k



qin j +



qramp outn

(55)

Please cite this article as: E.S. Canepa, C.G. Claudel, Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations, Transportation Research Part B (2017), http://dx.doi.org/10.1016/j.trb.2017.05.016

ARTICLE IN PRESS

JID: TRB 16

[m3Gsc;June 16, 2017;12:13]

E.S. Canepa, C.G. Claudel / Transportation Research Part B 000 (2017) 1–24

Fig. 3. Single road traffic state estimation. In this problem, we want to reconstruct the state of traffic using three different types of data: GPS data generated by a probe vehicle, a travel time measurement generated by a vehicle entering and exiting the physical domain, and a radar generating a density measurement. The present objective was to maximize the initial average density (worst-case average density), and the problem involves 49 variables and 929 constraints.

Fig. 4. Junction convention and flow conservation.

qouti are the flows leaving their respective road (to cross the junction) and qin j represent the flows entering ramp the roads (after crossing the junction). qin are the external flows entering the network through the intersection, conversely, k ramp qoutn are the flows exiting the network. The quantities



Note that the above equality constraint is insufficient in practice to derive the actual incoming and outgoing flows. Following Garavello and Piccoli (2005) and Piccoli and Bayen (2010), we assume that the flow at the junctions is the solution to a LP. This LP maximizes the total inflow through the junction under the constraints of an allocation matrix (which splits the flow according to the driver preferences), and physical constraints of demand and supply occurring at the boundaries of the incoming and outgoing links respectively. The equation of conservation of flows (55) is linear in the decision variable. We assume that the volume of traffic en ramp tering the junction (Ui s ) through qouti and qin is distributed among the exit options (O j s ) according to an allocation k

parameter α Outs, Ins ≥ 0. Since the junction has no storage capacity the sum of all the allocation parameters from a fixed incoming option among all the output options must be equal to one:



αO j ,Ui = 1

(56)

j

The relation of the incoming-outgoing flows in terms of the allocation parameters is encoded by the allocation matrix:







Oj



=A

j





Ui

(57)

i

The well-posedness constraints imposed by the LWR model on the junctions can be written as

⎧ q i ≤ di ⎪ ⎨qout ramp ≤d ink

≤ sj ⎪ j ⎩qin ramp

u,k

(58)

qoutn ≤ sd,n

where di correspond to the demand of link i, du, k corresponds to the demand of the on-ramp k, sj correspond to the supply of link j and sd, n correspond to the supply of off-ramp n. Since the demand and supplies of incoming and outgoing links is Please cite this article as: E.S. Canepa, C.G. Claudel, Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations, Transportation Research Part B (2017), http://dx.doi.org/10.1016/j.trb.2017.05.016

ARTICLE IN PRESS

JID: TRB

E.S. Canepa, C.G. Claudel / Transportation Research Part B 000 (2017) 1–24

[m3Gsc;June 16, 2017;12:13] 17

a mixed integer linear function of the decision variable (Han et al., 2016; 2014) (while the flows are linear in the decision variable), these constraints are mixed integer linear. Finally, the maximization of the sum of the flows through the junction imposes that at least one of the inequalities in (58) becomes an equality, adding additional integer variables to the problem. The network topology used is based on the Supply and Demand approach (as proposed by Daganzo, 1995 and Costeseque and Lebacque, 2014) in which based on the traffic conditions of the links involved in the junction, the maximum allowable flow is chosen. Hence, by combining all of the above constraints the junction constraints can be written as mixed integer linear inequalities in the decision variable, while the model constraints (on each links) are also mixed integer linear. If the objective to be minimized (as part of the estimation) is a linear function of the decision variable (which is typically the case), then the problem of estimating the state of traffic on a general highway network remains a MILP. Hence, the estimation of the state of traffic on a general network can be posed as a mixed integer linear program, in which the integer variables are a consequence of the junction models and of internal data (either internal boundary conditions or internal density conditions). 6.2. Merge junction model Junction models involving more outgoing links than incoming links do not yield a unique solution using allocation matrices only (Daganzo, 1995; Han et al., 2014). However, incorporating priority constraints in the merge model can be done in a mixed integer linear framework. For example, in a merge model proposed by Lebacque (1996), the problem states:

q = qin1 + qin2 qin1 = min(D1 , S1 ) sub ject to qin2 = min(D2 , S2 ) qin1 , qin2 ≥ 0

(59)

Where Di is the Demand of each link and Si = αi S includes a non-negative distribution factor α i for the Supply that satisfies α1 + α2 = 1. The supply of the downstream cell is first distributed to the two upstream segments with two fractions α 1 and α 2 , the flows q1 and q2 reaches their individual maximums since their respective Demands are computed by the model. The maximization of the total flow through the junction with fixed distribution factors imposes that:

q = min{D1 + D2 , α1 S + α2 S} which can be rewritten as a Mixed Integer Linear Programming (MILP) problem with an arbitrary objective function and appending the following constraints:

⎧ ⎪ ⎨q ≥= D1 + D2 q ≥= α1 S + α2 S ⎪q ≤ D1 + D2 + b1C ⎩ q ≤ α1 S + α2 S + (1 − b1 )C

(60)

Where bi is binary variable and C is a constant with a high value; with this framework is guaranteed that only one of the options in (60) is selected for the total flow in the merge. This junction model is based on Lebacque (1996) with a MILP approach. 6.3. Implementation We illustrate the above results by implementing a MILP to solve two travel time estimation problems (with distinct network structures) involving different segments of the I-880 highway located in the San Francisco bay area. The first problem will involve three roads, two junctions with a ramp on and off respectively as can be seen in Fig. 5. For this problem, we consider data generated by the PeMS (http://pems.eecs.berkeley.edu ) VDS stations 400,674 and 400,640 for the upstream and downstream position respectively (Fig. 5). The stations are located on the Highway I - 880 N around Hayward, California, the space domain of the entire network is 6 miles. We assumed that we do not have data in the junctions except for the allocation matrix which is predefined by the user, we are modeling the junctions with the approach described earlier. The traffic flow allocation matrix for both junctions in this example is the following:



qin1 rampout1





=

0.9 0.1



1 qout1 0 rampin1

(61)

An example of the travel time estimation of the structure mentioned above is show in Fig. 6 below, the result obtained by the estimation toolbox is compared with the ground truth, that is, GPS data obtained during the Mobile Century experiment. This particular example has 447 variables and 3127 constraints and took 4.57 s to solve on an iMac with a processor Intel Core i5-2400 @2.5 GHz. A series of estimations with different traffic conditions were analyzed, comparing the estimated travel time with the ground truth obtained from GPS data. After 30 estimations the RMS error obtained is 11 s. Please cite this article as: E.S. Canepa, C.G. Claudel, Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations, Transportation Research Part B (2017), http://dx.doi.org/10.1016/j.trb.2017.05.016

ARTICLE IN PRESS

JID: TRB 18

[m3Gsc;June 16, 2017;12:13]

E.S. Canepa, C.G. Claudel / Transportation Research Part B 000 (2017) 1–24

Fig. 5. Traffic Network Structure 1. In this problem we want to estimate the travel time of the network from the upstream to downstream position; the network consists of three roads connected through two junctions with a ramp on and off respectively, the latter is used to consider the flow leaving and entering the higway on the intersections.

The second structure is a merge and will involve three roads and one junction with a ramp on and off as can be seen in Fig. 7. For this problem, we consider data generated by the PeMS (http://pems.eecs.berkeley.edu ) VDS stations 400,490 and 400,685 for the upstream and downstream position respectively (Fig. 5). The stations are located on the Highway I 880 N around Hayward, California, the space domain of the entire network is 10.1 miles. Also, travel time data obtained on the Mobile Century experiment was added as data constraints to the problem. The road merging to the I-880 is Decoto Road, from which there is no data available; a random generator of reasonable flows in the upstream position of Decoto Rd. was included. It is important to mention that for a merge structure, some priority parameters must be defined in order to have a well-posed problem (Daganzo, 1995), these parameters differ from the allocation matrix and are included in the estimation framework. We assumed that we don’t have data in the junction except for the allocation matrix which is again predefined by the user. The traffic flow allocation matrix in this problem is the following:



qin1 rampout1





0.9 = 0.1

1 0

q 1 out1 q 0 out2 rampin1

(62)

Two different objective functions were selected to estimate the travel time. The first objective function is the minimiza max tion of the initial densities: ki=0 ρ (i ), with this objective objective the most optimistic initial number or vehicles in the highway network is obtained; an example of this objective function is shown in Fig. 8. The second objective function is the minimization of the L1 norm among all the variables of the decision vector (24); this is done to obtain a more uniform density map (Fig. 9). The results obtained by the estimation toolbox are compared with the ground truth, that is, travel time data validated during the Mobile Century experiment. The example in Fig. 8 has 546 variables and 4336 constraints and took 5.12 s to solve on an iMac with a processor Intel Core i5-2400 @2.5GHz. The MATLAB toolbox can be downloaded from https://www.dropbox.com/sh/852eu9xisrehtcv/AAA61cb2-vmKnFNS8UNTlXEea?dl=0. The travel time estimation was evaluated using different traffic conditions, a comparison of the results obtained by the estimator with both objective functions can be observed in Fig. 10. From the figure we can conclude that with the L1 norm as an objective function the travel time is always overestimated, this is due to estimated initial densities which (in contrast to the previous objective function) are not minimized, allowing more vehicles in the road. The L1 norm optimization can be understood as a worst case scenario estimation which is an important for the highway management and users. 6.4. Model parameter evaluation The inclusion of heterogeneous data during estimation is highly important. In order to evaluate the robustness of the framework, the following data were added extensively for the estimation: • • •

GPS data from probe vehicles Travel time data obtained from cameras Boundary flow data

An example of the data added to the network is shown in figure (11). Considerable amount of data was added until the problem became unfeasible; once unfeasibility was reached a relaxation of the continuity constraints (47) was introduced. The optimization problem is defined now as the minimization of the relaxation term with following modification on the Please cite this article as: E.S. Canepa, C.G. Claudel, Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations, Transportation Research Part B (2017), http://dx.doi.org/10.1016/j.trb.2017.05.016

ARTICLE IN PRESS

JID: TRB

E.S. Canepa, C.G. Claudel / Transportation Research Part B 000 (2017) 1–24

[m3Gsc;June 16, 2017;12:13] 19

Fig. 6. Traffic Time Estimation Example Structure 1. In all subfigures, we compute the density map for which the initial number of vehicles is the lowest, given the boundary data or the junction model, the sample car path is highlighted in red. For this example the estimated travel time is 443 s and the ground truth is 445 s. The objective function is the minimization of the initial densities. Top: Density map of the last highway segment, this segment has data in the downstream end. Center: Density map of the middle highway segment, this segment has a junction in both ends. Bottom: Density map of the initial segment of the highway, this segment has data on the upstream end. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

continuity constraints:

minimize such

that



||

∀ p ∈ M, M p (tmin ( p), xmin ( p)) = μ p (tmin ( p), xmin ( p)) ±

(63)

One of the disadvantages of the model is the uncertainty of the parameters; in order to select values that better describe the traffic (according to both the model and data constraints) the above optimization problem was evaluated several times with modifications on two parameters of the model: the free flow (vf ) and congested (w) velocity. As can be seen in figure Please cite this article as: E.S. Canepa, C.G. Claudel, Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations, Transportation Research Part B (2017), http://dx.doi.org/10.1016/j.trb.2017.05.016

JID: TRB 20

ARTICLE IN PRESS

[m3Gsc;June 16, 2017;12:13]

E.S. Canepa, C.G. Claudel / Transportation Research Part B 000 (2017) 1–24

Fig. 7. Traffic Network Merge Structure. We want to estimate the travel time of the network from the upstream to downstream position; the network consists of three roads connected through one junction with a ramp on and off, the latter is used to consider the flow leaving and entering the highway on the intersection.

Fig. 8. Traffic Time Estimation Example 1 Merge Structure. In these subfigures we compute the density map for which the initial number of vehicles is the lowest, given the boundary data or the junction model, the sample car path is highlighted in red. For this example the estimated travel time is 767 s and the ground truth is 769 s. Top: Density map of the last highway segment, this segment has data in the downstream end and a travel time constraint. Bottom: Density map of the initial segment of the highway, this segment has data on the upstream end and a travel time constraint. (For interpretation of the references to color in this figure legend, the reader is referred to the web version of this article.)

Please cite this article as: E.S. Canepa, C.G. Claudel, Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations, Transportation Research Part B (2017), http://dx.doi.org/10.1016/j.trb.2017.05.016

JID: TRB

ARTICLE IN PRESS E.S. Canepa, C.G. Claudel / Transportation Research Part B 000 (2017) 1–24

[m3Gsc;June 16, 2017;12:13] 21

Fig. 9. Traffic Time Estimation Example 2 Merge Structure. For these subfigures we compute the density map for which the L1 norm of the decision variables is the minimum, creating a more uniform density map. Top: Density map of the last highway segment, this segment has data in the downstream end and a travel time constraint. Bottom: Density map of the initial segment of the highway, this segment has data on the upstream end and a travel time constraint.

Fig. 10. Travel times Comparison.The ground truth are the travel times experienced on February 8th, 2008 during the Mobile Century experiment between Mowry and Winton Avenue (10.1 miles) compared with the results obtained by the travel time estimation using the merge structure and two different objective functions.

Please cite this article as: E.S. Canepa, C.G. Claudel, Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations, Transportation Research Part B (2017), http://dx.doi.org/10.1016/j.trb.2017.05.016

JID: TRB 22

ARTICLE IN PRESS

[m3Gsc;June 16, 2017;12:13]

E.S. Canepa, C.G. Claudel / Transportation Research Part B 000 (2017) 1–24

Fig. 11. Hetereogenous Data. Internal and travel time conditions were added to the estimation problem. This is an example of one link of the network.

Fig. 12. Results of the relaxation term minimization.

(12) low values in both parameters give better results (small value of Delta), with this procedure the parameters can be tune for better results in the long term. 7. Limitations 7.1. Model parameters It is important to point out that the calibration of model parameters does not necessarily lead to a strong convergence of the framework since multiple model parameters can correspond to the same (or nearly the same) value of the objective. In the long term, this issue can be addressed by computing the intersection of most likely model parameter domains associated with each dataset. For example, if n datasets are available (corresponding for example to different time instants), the set Sn defined by that Sn = {(v, w, kc )|O(v, w, kc ) < minν,ω,κ O(v, w, κ ) + } (for some configurable  ) would correspond to the domain of acceptable model parameters. The intersection ∪n∈N Sn would correspond to a more robust convergence, if this intersection is non empty (which can require to choose a higher value of  ). 7.2. Scalability The proposed algorithm is best suited to large networks containing low numbers of nodes. Indeed, since it requires Boolean variables only for internal conditions or junctions, the proposed optimization formulation involves less Boolean variables than discretization-based methods such as the Cell Transmission Model. However, being non-recursive, the comPlease cite this article as: E.S. Canepa, C.G. Claudel, Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations, Transportation Research Part B (2017), http://dx.doi.org/10.1016/j.trb.2017.05.016

JID: TRB

ARTICLE IN PRESS E.S. Canepa, C.G. Claudel / Transportation Research Part B 000 (2017) 1–24

[m3Gsc;June 16, 2017;12:13] 23

putational time required by this algorithm is a function of the chosen time horizon, and iterative Kalman Filter or Particle Filter based methods can be faster if the time horizon corresponding to the dataset is large. 7.3. Junction model The proposed framework assumes that the allocation matrices and demand/supply functions at the boundaries of the network are known. Jointly estimating the state and these parameters would not result in a standard form optimization program, since the constraints and objective are not piecewise linear functions of the state and parameters. Some authors have investigated the joint state/parameter estimation problem, for example in Wang et al. (2009); 2016). 8. Conclusion This article illustrates the travel time estimation application of a new Mixed Integer Programming estimation framework for highway traffic state estimation, in which the state of the system is modeled by the Lighthill Whitham Richards PDE. Using a Lax-Hopf formula, we show that the constraints arising from the model, as well as the measurement data result in linear inequality constraints for a specific decision variable, and that the problem of estimating a linear function of this decision variable is a Mixed Integer Linear Program. We also show that the method can be generalized to highway networks, at the expense of increasing the decision variable size, hence affecting the computational time. A numerical implementation of the estimation on an experimental dataset containing fixed sensor data as well as probe data is performed, and illustrates the ability of the method to quickly and efficiently compute all traffic scenarios compatible both with both the LWR model and given traffic states. This framework has the advantage of being exact and efficient for small-scale networks and gives the ability for the user to select any objective function, to explore the possible state estimates associated with a given dataset. Future work will involve investigating a two phase flow model that has a realistic upper bound on acceleration, to model the state of traffic in urban environments. Preliminary analysis shows that this type of model could be integrated within a similar optimization framework. References Alvarez-Icaza, L., Munoz, L., Sun, X., Horowitz, R., 2004. Adaptive observer for traffic density estimation. Boston, MA, pp. 2705–2710. Anderson, L., Canepa, E., Horowitz, R., Claudel, C., Bayen, A., 2013. Optimization-based queue estimation on an arterial traffic link with measurement uncertainties. Transportation Research Board 93rd Annual Meeting. Paper 14-4570. Aubin, J.-P., Bayen, A.M., Saint-Pierre, P., 2008. Dirichlet problems for some Hamilton-Jacobi equations with inequality constraints. SIAM Journal on Control and Optimization 47 (5), 2348–2380. Bardi, M., Capuzzo-Dolcetta, I., 1997. Optimal Control and Viscosity Solutions of Hamilton-Jacobi-Bellman equations. Birkhäuser, Boston, MA. Barron, E.N., Jensen, R., 1990. Semicontinuous viscosity solutions for Hamilton-Jacobi equations with convex Hamiltonians. Commun. Partial Differ. Equ. 15, 1713–1742. Boyd, S., Vandenberghe, L., 2004. Convex Optimization. Cambridge University Press, Cambridge, U.K.. Canepa, E.S., Bayen, A.M., Claudel, C.G., 2013. Spoofing cyber attack detection in probe-based traffic monitoring systems using mixed integer linear programming. Networks and Heterogeneous Media 8 (3), 783–802. doi:10.3934/nhm.2013.8.783. Canepa, E.S., Claudel, C.G., 2012. Exact solutions to traffic density estimation problems involving the Lighthill-Whitman-Richards traffic flow model using Mixed Integer Linear Programing. In: Proceedings of the 15th International IEEE Conference on Intelligent Transportation Systems. Anchorage, AK. Claudel, C.G., Bayen, A., 2010. Lax-Hopf based incorporation of internal boundary conditions into Hamilton-Jacobi equation. Part I: theory. IEEE Trans. Automat. Control 55 (5), 1142–1157. doi:10.1109/TAC.2010.2041976. Claudel, C., Bayen, A., 2010. Lax-Hopf based incorporation of internal boundary conditions into Hamilton-Jacobi equation. Part II: Computational methods. IEEE Trans. Automat. Control 55 (5), 1158–1174. Doi:10.1109/TAC.2010.2045439 Claudel, C.G., Bayen, A.M., 2011. Convex formulations of data assimilation problems for a class of Hamilton-Jacobi equations. SIAM Journal on Control and Optimization 49, 383–402. Costeseque, G., Lebacque, J.-P., 2014. Discussion about traffic junction modelling: conservation laws vs hamilton-jacobi equations. Discrete Continuous Dynamical Systems - Series S 7 (3), 411–433. doi:10.3934/dcdss.2014.7.411. Crandall, M.G., Lions, P.-L., 1983. Viscosity solutions of Hamilton-Jacobi equations. Trans. Am. Math. Soc. 277 (1), 1–42. Daganzo, C., 1994. The cell transmission model: a dynamic representation of highway traffic consistent with the hydrodynamic theory. Transportation Research 28B (4), 269–287. Daganzo, C., 1995. The cell transmission model, part II: network traffic. Transportation Research 29B (2), 79–93. Daganzo, C., 2006. On the variational theory of traffic flow: well-posedness, duality and applications. Networks and Heterogeneous media 1, 601–619. Daganzo, C.F., 2005. A variational formulation of kinematic waves: basic theory and complex boundary conditions. Transporation Research B 39B (2), 187–196. Deng, W., Lei, H., Zhou, X., 2013. Traffic state estimation and uncertainty quantification based on heterogeneous data sources: a three detector approach. Transportation Research Part B 57, 132–157. http://dx.doi.org/10.1016/j.trb.2013.08.015. Frankowska, H., 1993. Lower semicontinuous solutions of Hamilton-Jacobi-Bellman equations. SIAM Journal of Control and Optimization 31 (1), 257–272. Garavello, G.M.Coclite.M., Piccoli, B., 2005. Traffic flow on a road network. SIAM J. Math. Anal. 36 (6), 1862–1886. Han, K., Gayah, V.V., Piccoli, B., Friesz, T.L., Yao, T., 2014. On the continuum approximation of the on-and-off signal control on dynamic traffic networks. Transportation Research Part B 61, 73–97. http://dx.doi.org/10.1016/j.trb.2014.01.001. Han, K., Liu, H., Gayah, V.V., Friesz, T.L., Yao, T., 2016. A robust optimization approach for dynamic traffic signal control with emission considerations. Transportation Research Part C 70, 3–26. http://dx.doi.org/10.1016/j.trc.2015.04.001. Jin, W.-L., 2015. Continuous formulations and analytical properties of the link transmission model. Transportation Research Part B. 74, 88–103. http://dx. doi.org/10.1016/j.trb.2014.12.006. Lebacque, J.-P., 1996. The godunov scheme and what it means for first order traffic flow models. Internaional symposium on transportation and traffic theory 13th. Lyon, France. Li, Y., Canepa, E., Claudel, C., 2014. Efficient robust control of first order scalar conservation laws using semi-analytical solutions. Discrete and Continuous Dynamical Systems - Series S 7 (3), 525–542. doi:10.3934/dcdss.2014.7.525.

Please cite this article as: E.S. Canepa, C.G. Claudel, Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations, Transportation Research Part B (2017), http://dx.doi.org/10.1016/j.trb.2017.05.016

JID: TRB 24

ARTICLE IN PRESS

[m3Gsc;June 16, 2017;12:13]

E.S. Canepa, C.G. Claudel / Transportation Research Part B 000 (2017) 1–24

Li, Y., Canepa, E., Claudel, C., 2014. Optimal control of scalar conservation laws using linear/quadratic programming: application to transportation networks. IEEE Trans. Control Netw. Syst. 1 (1), 28–39. doi:10.1109/TCNS.2014.2304152. Lighthill, M.J., Whitham, G.B., 1956. On kinematic waves. II. A theory of traffic flow on long crowded roads. Proceedings of the Royal Society of London 229 (1178), 317–345. Mazare, P.E., Dehwah, A., Claudel, C.G., Bayen, A.M., 2011. Analytical and grid-free solutions to the lighthill-whitham-richards traffic flow model. Transportation Research Part B 45 (10), 1727–1748. Moskowitz, K., 1965. Discussion of ‘freeway level of service as influenced by volume and capacity characteristics’ by D.R. Drew and C. J. Keese. Highway Research Record 99, 43–44. Munoz, L., Sun, X., Horowitz, R., Alvarez, L., 2003. Traffic density estimation with the cell transmission model. In: American Control Conference, 2003. Proceedings of the 2003, 5, pp. 3750–3755 vol.5. doi:10.1109/ACC.2003.1240418. Nantes, A., Ngoduy, D., Bhaskar, A., Miska, M., Chung, E., 2016. Real-time traffic state estimation in urban corridors from heterogeneous data. Transportation Research Part C 66, 99–118. http://dx.doi.org/10.1016/j.trc.2015.07.005. Advanced Network Traffic Management: From dynamic state estimation to traffic control Newell, G.F., 1993. A simplified theory of kinematic waves in highway traffic, Part (I), (II) and (III). Transporation Research B 27B (4), 281–313. Piccoli, D., Work, S., Blandin, O., Tossavainen, B., Bayen, A., 2010. A traffic model for velocity data assimilation. Applied Research Mathematics eXpress (ARMX) (1) 1–35. Richards, P.I., 1956. Shock waves on the highway. Oper. Res. 4 (1), 42–51. Rockafellar, R., 1970. Convex Analysis. Princeton University Press, Princeton, NJ. Strub, I.S., Bayen, A.M., 2006. Weak formulation of boundary conditions for scalar conservation laws. International Journal of Robust and Nonlinear Control 16 (16), 733–748. Wang, R., Work, D.B., Sowers, R., 2016. Multiple model particle filter for traffic estimation and incident detection. IEEE Trans. Intell. Transp. Syst. 17 (12), 3461–3470. Wang, Y., Papageorgiou, M., Messmer, A., Coppola, P., Tzimitsi, A., Nuzzolo, A., 2009. An adaptive freeway traffic state estimator. Automatica 45 (1), 10–24. Work, D., Blandin, S., Tossavainen, O., Piccoli, B., Bayen, A., 2010. A distributed highway velocity model for traffic state reconstruction. Applied Research Mathematics eXpress (ARMX) 1, 1–35. Yperman, I., Logghe, S., Immers, B., 2005. The link transmission model: An efficient implementation of the kinematic wave theory in traffic networks. In: In Advanced OR and AI Methods in Transportation, Proceedings of.

Please cite this article as: E.S. Canepa, C.G. Claudel, Networked traffic state estimation involving mixed fixed-mobile sensor data using Hamilton-Jacobi equations, Transportation Research Part B (2017), http://dx.doi.org/10.1016/j.trb.2017.05.016