Applied Mathematics and Computation 266 (2015) 390–403
Contents lists available at ScienceDirect
Applied Mathematics and Computation journal homepage: www.elsevier.com/locate/amc
A cutting plane projection method for bi-level area traffic control optimization with uncertain travel demand Suh-Wen Chiou∗ Department of Information Management, National Dong Hwa University, Da Hsueh Rd., Shou-Feng, Hualien, 97401, Taiwan
a r t i c l e
i n f o
Keywords: Congestion pricing Bi-level min–max model Cutting plane method Nash–Stackelberg game Equilibrium flow
a b s t r a c t A robust congestion pricing (ROPRICE) of signal-controlled road network is considered for equilibrium network flow with uncertain travel demand. A min–max bi-level program is presented to mitigate vulnerability of area traffic control road network against growing travel demand and reduce traffic congestion. A cutting plane projection approach (CPP) is proposed to effectively solve the ROPRICE problem with global convergence. Numerical computations are performed using various test road networks and comparisons are made with recently proposed heuristics. Computational results indicate that the proposed solution scheme can substantially achieve greater system performance as compared to other alternatives. © 2015 Elsevier Inc. All rights reserved.
1. Introduction For most urban road networks, severe traffic congestion would be incurred at signal-controlled junctions by road travelers as a result of inappropriate design of signal settings when uncertain travel demand is prevailing. Congestion pricing for system optimum flow has long been regarded as the most efficient instrument to mitigate traffic congestion and achieve efficient use of road network [1–7]. Among them, the marginal social cost pricing (MSCP) recognized as the first-best congestion pricing has been theoretically considered a most economically efficient tool to deal with congestion and make best usage of road network [8–14]. However, for some road networks, a second-best congestion pricing is pursued instead when only a subset of links can be tolled due to public acceptance and social and political restrictions. Ways of using mathematical programming to solve a constrained optimization problem of a second-best congestion pricing in a subset of road network links with certain travel demand have received ample attention [15–20]. Yang and Lam [15] proposed a sensitivity-analysis based algorithm (SAB) for simple congestion pricing problem in a general road network. Unfortunately, the relation between equilibrium flow and congestion pricing variables may not always be differentiable, which has been theoretically investigated in [21] and [22] and numerically illustrated in [23] and [24]. The combined problem of congestion pricing and signal settings for equilibrium flow can be formulated as a case of mathematical program with equilibrium constraints (MPEC) or MPCC (mathematical program with complementarity constraints) [25–27]. Lawphongpanich and Hearn [18] proposed a cutting constraint algorithm (CCA) approach to solving a second-best toll pricing problem. Using the concept of MPEC, constrained non-linear programming problems for second-best tolling optimization with fixed and variable travel demand were addressed. As it has been recognized widely from literature [28–30], the MPEC problem is generally a non-convex problem due to the implicit form in constraints. Chiou [12] introduced a projected subgradient approach to effectively solve a first-best congestion pricing problem for a general road network. More recently, Chiou [13,14] explored several heuristics respectively for variable travel demand in a general road network and for ∗
Tel.: +886 3 8633108; fax: +886 3 8633100. E-mail address:
[email protected]
http://dx.doi.org/10.1016/j.amc.2015.05.009 0096-3003/© 2015 Elsevier Inc. All rights reserved.
S.-W. Chiou / Applied Mathematics and Computation 266 (2015) 390–403
391
signal-controlled road networks with regard to certain travel demand. On the other hand, Sumalee [20] presented a GA-based heuristic to optimize a charging cordon in a general traffic network and gave promising results on finding optimal charging cordon design from empirical studies. The effectiveness of proposed GA-based heuristics has been demonstrated numerically at example test networks with considerable travel cost savings. In order to mitigate vulnerability of system against uncertainty, there are basically two popular approaches to optimization under uncertainty. One is the approach of stochastic programming [31–33] modeling uncertain parameters as random variables with a priori probability distribution and optimizing the expected value of the objective function. For example, Gardner et al. [5] presented a stochastic mathematical program with equilibrium constraints for robust congestion pricing of transportation networks under uncertain demand and gave solution methods for robust prices which are resilient to demand variations. Li et al. [34] proposed a dynamic pricing for perishable products involving hybrid uncertainty in demand. Ban et al. [35] proposed a risk-neutral second-best toll pricing model to account for the possible non-uniqueness of user equilibrium solutions. Instead of knowing the probability distribution of uncertain parameters, the other approach, robust optimization, where one optimizes a worst possible case of a considered problem, has long been considered an effective tool to hedge against uncertainty [36–38]. In this regard, Lou et al. [39] presented a robust congestion pricing for a traffic road network under boundedly rational user equilibrium (BRUE). Chung et al. [2] also presented a dynamic congestion pricing with uncertain travel demand. Because of the non-linearity of the equilibrium constraint with respect to congestion pricing variables, most solution heuristics mentioned above can simply solve a non-convex and non-linear inequality constrained problem with congestion pricing only locally. In this paper, a robust congestion pricing (ROPRICE) of a signal-controlled road network is presented following earlier work in Chiou [12]. In order to hedge against uncertain travel demand, a min–max model is proposed. A set of robust signals for a ROPRICE problem in the presence of uncertain travel demand can be determined against a worst-case that might be taken by the opponent in anticipation of road users’ responding strategy via route choice. In this regard, the worst-case performance measure serves as an upper bound estimate. A cutting plane projection (CPP) is presented and numerical computations are performed using road networks of realistic size. The contributions made from this paper are summarized as follows. First, a ROPRICE problem is performed optimally to determine robust signal settings and link congestion pricing in the presence of uncertain travel demand. A Nash–Stackelberg equilibrium is established. The performance measure (PM), maximized with respect to a travel demand growth factor, is minimized with respect to signal-setting and link congestion pricing. A cutting plane projection is presented for a ROPRICE problem with uncertain demand. Computational results obtained indicate that proposed solution scheme can solve the ROPRICE problem successfully. As compared to deterministic solution for nominal travel demands using the MSCP, the proposed scheme achieved greater improvement while incurring a relatively slighter loss of optimality. The rest of the paper is organized as follows. Section 2 introduces a bi-level min–max model for a ROPRICE problem with uncertain travel demand. A Nash–Stackelberg solution can be effectively characterized by a tractable computation scheme proposed in Section 3. Numerical computations of proposed scheme are performed in Section 4. Conclusions for this paper and extensions of the proposed approach to topics of interest are briefly summarized in Section 5. 2. A ROPRICE problem with uncertain travel demand In the presence of uncertain travel demand, a bi-level min–max model is introduced. Notation used for a ROPRICE problem with link congestion pricing is introduced first. 2.1. Notation G(N, L ) W Rw = (ζ , θ , φ )
β λa λmin τ jlm m ( j, l ) ρa
sa q
μ
f h
δ c ( , β , f ) π
a road network with node set N and link set L. a set of origin–destination (OD) pairs. a set of routes between OD pair w, ∀w ∈ W . set of signal setting variables. a set of link tolls. duration of effective green for link a, ∀a ∈ L. minimum green. clearance time between the end of green for group j and the start of green for incompatible group l at junction m. collection of numbers 0 and 1 for each pair of incompatible signal groups at junction m; where m ( j, l ) = 0 if the start of green for signal group j proceeds that of l and m ( j, l ) = 1, otherwise. maximum degree of saturation for link a, ∀a ∈ L. saturation flow on link a, ∀a ∈ L. a matrix of travel demand for OD pairs. OD demand growth factor. a vector of link flow, i.e. f = [ fa ], ∀a ∈ L. a vector of route flow between OD trips, i.e. h = [h p ],∀ p ∈ Rw ,∀w ∈ W . a link-route incidence matrix. a OD-route incidence matrix. a vector of link flow travel cost, i.e. c ( , β , μ ) = [ca ( , β , f )], ∀a ∈ L. a vector of minimum travel cost between OD, i.e. π = [πw ], ∀w ∈ W .
392
S.-W. Chiou / Applied Mathematics and Computation 266 (2015) 390–403
a vector of route flow travel cost, i.e. C = [C p ], ∀ p ∈ Rw , ∀w ∈ W . rate of delay on link a, ∀a ∈ L. number of stops per unit time on link a, ∀a ∈ L. link-specific weighting factor for rate of delay on link a, ∀a ∈ L. link-specific weighting factor for number of stops on link a, ∀a ∈ L. monetary factor associated with Da , ∀a ∈ L. monetary factor associated with Sa , ∀a ∈ L.
C Da Sa WaD WaS MD MS
2.2. A user equilibrium traffic assignment A user equilibrium flow following Wardrop’s principle [40] can be identified by a variational inequality. Let denote a feasible set for network flow, i.e.
= { f : f = δ h, h = q, h ≥ 0}
(1)
A user equilibrium flow can be characterized if and only if for every f¯ ∈ there exists a f ∈ such that
c ( f )( f¯ − f ) ≥ 0
(2)
2.3. A parametric user equilibrium flow with uncertain travel demand For a ROPRICE problem with a set of link congestion pricing β and signal settings in the presence of a realization taken by unknown demand growth factor μ, a responding user equilibrium flow f ( , β , μ ) can be characterized in this section. Let (μ ) denote a feasible set for a parametric user equilibrium flow with respect to some realization taken by unknown demand growth factor μ, we have
(μ ) = { f : f = δ h, h = μq, h ≥ 0}
(3)
Therefore a parametric user equilibrium flow with OD demand growth factor μ can be characterized if and only if for every f¯ ∈ (μ ) there exists a f ( , β , μ ) ∈ (μ ) such that
c ( , β , f )( f¯ − f ) ≥ 0
(4)
∗ ( , β , μ )
Let denote a solution set determined by (4) consisting of responding flow f ( , β , μ ) to a set of link congestion pricing β and signal settings in the presence of some realization taken by unknown demand growth factor μ, i.e.
f ( , β , μ ) ∈ ∗ ( , β , μ )
(5)
2.4. A bi-level min–max model A ROPRICE problem in the presence of uncertain travel demand can be expressed as a bi-level min–max model. In the bilevel min–max model, a link congestion pricing β ∗ and signal settings ∗ in the presence of uncertain travel demand can be determined with a worst-case realization of a travel demand growth μ∗ and equilibrium flow f ∗ . In this regard, the objective function evaluated at the upper level constrained with a responding flow f ∗ , maximized with respect to a realization taken by unknown demand growth factor μ∗ , is minimized with respect to link congestion pricing β ∗ and signal settings ∗ .
Min Max
P ( , β , f ( , β , μ ))
(6)
Sub ject to
ζ ≤ζ ≤ζ
(7)
,β , f
μ
λmin ζ ≤ φ jm ≤ 1, fa ( , β , μ ) ≤ ρa sa λa ,
∀ j, m
(8)
∀a ∈ L
θ jm + φ jm + τ jlm ζ ≤ θlm + m ( j, l ),
(9)
j = l, ∀ j, l, m
(10)
According to [51], the performance measure (PM) of a ROPRICE problem can be effectively evaluated in terms of traffic delays and queues, i.e.
P ( , β , f ) =
a∈L
Da ( , f ( , β , μ ))WaD MD + Sa ( , f ( , β , μ ))WaS MS +
βa fa ( , β , μ )
(11)
a
In constraints (7)–(10): constraint (7) is a bound for common cycle time, and the constraints (8)–(10) are respectively for green phase, link capacity and clearance time for signal setting variables. The performance measure (PM), P ( , β , f ( , β , μ )) in (11), is defined as a monetary sum of a weighted linear combination of the rate of delay and number of stops per unit time for all traffic
S.-W. Chiou / Applied Mathematics and Computation 266 (2015) 390–403
393
streams and selected link tolls. A responding flow f ∗ in a bi-level min–max problem (6)–(10) is solved by a parametric variational inequality (4). Both of link congestion pricing β ∗ and signal settings ∗ together with a travel demand growth factor μ∗ can be optimally determined by solving a bi-level min–max model. A triple of saddle points ( ∗ , β ∗ , μ∗ ) exists with equilibrium network flow f ∗ if the following condition holds:
P ( ∗ , β ∗ , f ( ∗ , β ∗ , μ )) ≤ P ( ∗ , β ∗ , f ( ∗ , β ∗ , μ∗ )) ≤ P ( , β , f ( , β , μ∗ )) According to (5), a responding user equilibrium flow variational inequality (4). Therefore, we have
f∗
to a triple of saddle points
(12)
( ∗ , β ∗ , μ∗ )
is a solution of a parametric
P ( ∗ , β ∗ , f ( ∗ , β ∗ , μ )) ≤ P ( ∗ , β ∗ , f ∗ ) ≤ P ( , β , f ( , β , μ∗ ))
(13)
μ∗
is maximized under link In (13), at the first inequality a worst-case realization taken by unknown demand growth factor congestion pricing β ∗ and signal settings ∗ . Analogously, at the second inequality of (13) a link congestion pricing β ∗ and signal settings ∗ are minimized under a worst-case scenario of uncertain demand with a growth factor μ∗ . 2.5. Directional derivatives for responding flow Regarding a link congestion pricing β ∗ and signal settings ∗ for a worst-case scenario of uncertain demand with growth factor μ∗ , the direction of change that occurs in responding flow f ∗ can be evaluated by the first-order directional derivatives. According to [21–24], a responding flow f ∗ is not always differentiable. Therefore it would be preventive from direct use of sensitivity analysis results in [41] for responding flow. However, according to Rademacher’s theorem [42], the first-order directional derivatives exist almost everywhere for general equilibrium network flows. According to results given in [21,22] and [43], the perturbation ∇ f of responding flow f ∗ with respect to change in ,β and μ can be characterized by solving a linearized variational inequality in the following way. Introduce
(μ ) = {∇ f : ∇ f = δ (h ), (h ) = μ, ∃h ∈ K0 } and
h∗p > 0, (i )h p f ree, if K0 = h : (ii )h p ≥ 0, i f h∗p = 0, Cp = π w , (iii )h p = 0, i f h∗p = 0, Cp > π w
∀ p ∈ Rw , ∀w ∈ W
(14)
For every fˆ ∈ (μ ), a directional derivative ∇ f ∈ (μ ) along a direction β and can be determined in the following way.
(∇β c( ∗ , β ∗ , f ∗ )β + ∇ c( ∗ , β ∗ , f ∗ ) + ∇ f c( ∗ , β ∗ , f ∗ )∇ f )( fˆ − ∇ f ) ≥ 0
(15)
The gradients ∇β c ( ∗ , β ∗ , f ∗ ),∇ c ( ∗ , β ∗ , f ∗ )and ∇ f c ( ∗ , β ∗ , f ∗ ) in (15) are evaluated at ( ∗ , β ∗ , f ∗ ) when perturbations in ,β and μ are specified. The directional derivatives ∇ f ,∇β f and ∇μ f are piecewise linear and differentiable almost everywhere according to Rademacher’s theorem in [42]. According to (5), the generalized gradients for responding flow f ∗ can be characterized in the following way. Let co denote a convex hull, it implies
∂ ∗ ( ∗ , β ∗ , μ∗ ) = co n→∞ lim ∇ f ( (n) , β (n) , μ(n) ) : ( (n) , β (n) , μ(n) ) → ( ∗ , β ∗ , μ∗ ),∇ f ( (n) , β (n) , μ(n) )exist
(16)
∂β ∗ ( ∗ , β ∗ , μ∗ ) = co n→∞ lim ∇β f ( (n) , β (n) , μ(n) ) : ( (n) , β (n) , μ(n) ) → ( ∗ , β ∗ , μ∗ ),∇β f ( (n) , β (n) , μ(n) )exist (17)
∂μ ∗ ( ∗ , β ∗ , μ∗ ) = co n→∞ lim ∇μ f ( (n) , β (n) , μ(n) ) : ( (n) , β (n) , μ(n) ) → ( ∗ , β ∗ , μ∗ ),∇μ f ( (n) , β (n) , μ(n) )exist (18) 2.6. A single-level problem According to results given in (16) and (17), the generalized gradients for performance measure in (11) can be calculated as follows.
∇ P = (∇ D + ∇ f D( ∗ , f ∗ )∇ f )WD MD + (∇ S + ∇ f S( ∗ , f ∗ )∇ f )WS MS + β∇ f ( , β , μ )
(19)
∇β P = (∇ f D( ∗ , f ∗ )∇β f )WD MD + (∇ f S( ∗ , f ∗ )∇β f )WS MS + f ( , β , μ ) + β∇β f ( , β , μ )
(20)
394
S.-W. Chiou / Applied Mathematics and Computation 266 (2015) 390–403
The performance measure in (11) can be re-expressed as
P˜1 ( , β ) =
Da ( )WaD MD + Sa ( )WaS MS
βa fa
(21)
a
a∈L
And the generalized gradients can be expressed as
∇ P˜1 = ∇ P
(22)
∇β P˜1 = ∇β P
(23)
Now, a bi-level min–max model for (6)–(10) can be reduced into a single-level optimization problem according to (18) and (21).
P1 ( , β , μ ) =
Min Max μ
,β
Da ( )WaD MD + Sa ( )WaS MS +
βa fa
a
a∈L
ζ ≤ζ ≤ζ ∀ j, m λmin ζ ≤ φ jm ≤ 1, fa ( , β , μ ) + ∇ fa + β∇β fa + μ∇μ fa ≤ ρa sa λa , θ jm + φ jm + τ jlm ζ ≤ θlm + m ( j, l ), j = l, ∀ j, l, m
Sub ject to
Analogous to (13), a triple of saddle points
( ∗ , β ∗ , μ∗ )
∀a ∈ L (24)
exists for a min–max problem (24) if the following condition holds:
P1 ( , β , μ ) ≤ P1 ( , β , μ ) ≤ P1 ( , β , μ ) ∗
∗
∗
∗
∗
∗
(25)
μ∗
At the first inequality of (25) a worst-case realization taken by unknown demand growth factor is determined such that the performance measure in (25) is maximized. Similarly, at the second inequality of (25) the link congestion pricing β ∗ and signal settings ∗ are determined such that performance measure in (25) is minimized. Let
P1 (μ ) = P1 ( ∗ , β ∗ , μ ) = Min P1 ( , β , μ )
(26)
,β
where ( ∗ , β ∗ ) = Arg min P1 ( , β , μ ), and ,β
P1 ( , β ) = P1 ( , β , μ∗ ) = Max P1 ( , β , μ )
(27)
μ
where μ∗ = Arg max P1 ( , β , μ ). The inequalities in (25) can be re-expressed. μ
P1 (μ ) ≤ P1 ( ∗ , β ∗ , μ∗ ) ≤ P1 ( , β )
(28)
( ∗ , β ∗ , μ∗ )
for a min–max problem (24) can be optimally determined by alternately A triple of Nash–Stackelberg solutions solving two following sub-problems. First, under initial demand with a factor μ, according to the generalized gradient ∇μ f for responding flow determined in (15) and characterized in (18), it is to find a μ at ( ∗ , β ∗ ) such that μ∗ = μ + μ; i.e.
Max μ
P1 (μ ) =
Sub ject to
a∈L fa∗ +
Da ( ∗ )WaD MD + Sa ( ∗ )WaS MS +
βa fa∗
a
μ∇μ fa ≤ ρa sa λ∗a ,
∀a ∈ L
(29)
Next, under link capacity β and signal-setting , according to the generalized gradients ∇ f and ∇β f for responding flow determined in (15) and characterized in (16) and (17), it is to determine a pair of ( , β ) at a worst-case scenario of uncertain demand with a growth factor μ∗ such that ∗ = + and β ∗ = β + β ; i.e.
Min
( ,β )
P¯1 ( , β ) =
Da ( )WaD MD + Sa ( )WaS MS +
a∈L
ζ ≤ζ ≤ζ ∀ j, m λmin ζ ≤ φ jm ≤ 1, ∀a ∈ L fa∗ + ∇ fa + β∇β fa ≤ ρa sa λa , θ jm + φ jm + τ jlm ζ ≤ θlm + m ( j, l ), j = l, ∀ j, l, m
βa fa
a
Sub ject to
(30)
According to a maximization problem (29) with respect to μ and a minimization problem (30) with respect to ( , β ), a Nash–Stackelberg solution for (24) can be defined as follows. Definition 1. (A Nash–Stackelberg solution) We say that ( ∗ , β ∗ , μ∗ ) is a triple of Nash–Stackelberg solutions for (24) if and only if there exists a link congestion pricing β ∗ and signal settings ∗ for a minimization problem (30) and a worst-case scenario
S.-W. Chiou / Applied Mathematics and Computation 266 (2015) 390–403
395
of demand with a growth factor μ∗ for a maximization problem (29). In other words, there exists a link congestion pricing β ∗ and signal settings ∗ under a worst-case realization of a demand growth factor μ∗ such that the equalities in (28) hold, i.e.
P1 (μ∗ ) = P1 ( ∗ , β ∗ , μ∗ ) = P1 ( ∗ , β ∗ )
(31)
In (31) we have
μ∗ = Arg max P1 (μ )
(32)
( ∗ , β ∗ ) = Arg min P1 ( , β )
(33)
According to (26) and (27), it concludes at a Nash–Stackelberg solution we have a link congestion pricing β ∗
∗ in the presence of a worst-case realization of a demand growth factor μ∗ such that
and signal settings
P1 ( ∗ , β ∗ , μ∗ ) = Max Min P1 ( , β , μ ) = Min Max P1 ( , β , μ ) μ
,β
,β
(34)
μ
Furthermore, by means of (13), we have
P ( ∗ , β ∗ , f ∗ ) = Max Min P ( , β , f ( , β , μ )) = Min Max P ( , β , f ( , β , μ )) μ
,β
,β
μ
(35)
3. A tractably computational solution scheme In this section, a tractable computation solution scheme is proposed to solve two sub-problems (29) and (30). In the presence of a realization taken by unknown demand growth factor that is most unfavorable, perturbations of parametric flow with respect to link congestion pricing and signal settings can be evaluated using piecewise linear directional derivatives given in (16)–(18). For a maximization problem (29), a demand growth factor μ∗ can be determined by solving a linear programming problem. For a minimization problem (30), a cutting plane projection (CPP) approach is introduced. Regarding the CPP approach, two relevant properties for a minimization problem (30) are introduced first. Definition 2.
We say an objective function P1 ( , β ) in (30) is semi-smooth on the domain set if P1 ( , β ) is locally Lipschitz. There exists a vector (v1 , v2 ) ∈ ∂ P¯1 ( + l , β + l β ), → and β → β such that the limit
(v1 , v2 β )
lim
(v1 ,v2 )∈∂ P1 ( +l ,β +β ), → ,β →β ,l→0
(36)
exists for every direction and β . Lemma 1. Suppose P1 (y, ) in (30) is a locally Lipschitzian function and the directional derivative ∇ P1 ( ∗ , β ∗ ; , β ) exists for every direction and β in the change of ∗ and β ∗ . Then (1). ∇ P1 (·, ·; , β ) is Lipschitzian. (2). For any direction and β , there exists a subgradient (v1 , v2 ) ∈ ∂ P1 ( ∗ , β ∗ ) such that
∇ P1 ( ∗ , β ∗ ; , β ) = (v1 , v2 β )
(37)
3.1. A cutting plane model for a minimization problem In order to effectively solve a minimization problem (30), in this section, a cutting plane method is employed. The cutting plane method has been widely used in the field of non-smooth optimization [46–48]. For the objective function P1 (·, · ) in (30), a linear approximation can be constructed using a bundle of past computed subgradients {P1 ( (i ) , β (i ) ), ∇ P1 ( (i ) , β (i ) ), 1 ≤ i ≤ k}. Let co denote a convex hull, thus a convex hull of all points of the form lim ∇ P1 ( (k ) , β (k ) ) can be constructed below for the generalized gradients of (30).
k→∞
∂ P1 ( ∗ , β ∗ ) = co lim ∇ P1 ( (k) , β (k) ) : ( (k) , β (k) ) → ( ∗ , β ∗ ), ∇ P1 ( (k) , β (k) )exists k→∞
Let P1
P1
(k )
(k )
(38)
denote a linear approximation of P1 ( , β ) close to ( (k ) , β (k ) ) for (30) at iteration i, ∀1 ≤ i ≤ k, we have
≈ Max
1≤i≤k
Supposing that P1 (k )
∇ P1 ( (i) , β (i) ) (k )
εi,k = P1 ( , β
− (i ) β − β (i )
+ P1 ( (i) , β (i) )
is convex, let
(k )
(i )
(i )
(i )
(i )
) − P1 ( , β ) + ∇ P1 ( , β )
(k ) − (i ) β (k ) − β (i )
(39)
(40)
396
S.-W. Chiou / Applied Mathematics and Computation 266 (2015) 390–403
denote an error bound for a linear approximation of P1 ( (k ) , β (k ) ) and εi,k ≥ 0 due to local convexity of P1 (β ). Let ( (k) , β (k) ) be a search direction from current iterate ( (k) , β (k) ) to the next one and denoted by
( (k) , β (k) ) = ( − (k) , β − β (k) ) A linear approximation of P1
(k )
(41)
in (39) can be expressed in terms of bundle gradients in the following way:
(k ) (k) (i ) (i ) P1 ≈ Max ∇ P1 ( , β ) − εi,k + P1 ( (k) , β (k) ) β ( k ) 1≤i≤k
(42)
Since P1 ( (k ) , β (k ) ) is constant in (42), for convenience let us skip this, a cutting plane model for a linear approximation of P1 ( (k ) , β (k ) ) in terms of bundle gradients, along a perturbed direction ( (k ) , β (k ) ) can be constructed.
(k )
P1 cp ( (k) , β (k) ; (k) , β (k) ) = Max
1≤i≤k
∇ P1 ( (i) , β (i) )
(k) β ( k )
− εi,k
(43)
The general concept of cutting plane model for (30) is based on a piecewise linear lower approximation of P1 ( (k ) , β (k ) ) in the neighborhood of iterate ( (k ) , β (k ) ). Such approximations are built by employing the supporting hyper-planes to the (k )
corresponds to one step of graph of minimization problem, received at the iterates {( (i ) , β (i ) )}i≤k . The minimization of P1 cutting plane method which may give us very slow convergence. In order to improve slow convergence in (43) for a minimization problem (30), a projection approach is employed.
3.2. A cutting plane projection (CPP) approach In this section, we propose a cutting plane projection (CPP) solution scheme to effectively solve a minimization problem (30) with global convergence. The bundle subgradients of the objective function in (30) are projected onto null space of active constraints in search for local optimum. By consecutive implementing CPP, successive projections of the bundle of past computed subgradients of the objective function in (30) will help us to dilate the direction provided by the negativity of the subgradients, which may effectively improve the local solutions obtained for (30). Theorem 1. A cutting plane projection (CPP) approach: Let d (k ) = ( (k ) , β (k ) ) denote a solution for problem (30) which can be determined in the following way:
Min d (k )
2 1 v + d ( k )
sub ject to
2
∇ P1 ( (i) , β (i) )d (k) − εi,k ≤ ν,
1≤i≤k
(44)
Let x(k ) = ( (k ) , β (k ) ) and Pr (x ) denote the projection of x on a domain set for a minimization problem (30) such that
x − Pr (x ) = inf x − z
(45)
z∈
A sequence of iterates {x(k ) } for a minimization problem (30) can be determined in accordance with
x(k+1) = Pr (x(k) + ld (k) ),
k = 1, 2, . . .
(46)
where l ∈ (0, 2 ) is the step length which minimizes P¯1(k ) in (30) and d (k ) solves problem (44). Let x∗ be a minimum point, then x(k ) − x∗ is monotonically decreasing and x(k ) − x∗ ≤ x(1 ) − x∗ , i.e. the sequence of points {x(k ) } generated by CPP is bounded. Proof. For any x1 and x2 in the domain set for (30), by definition of the projection, we have
Pr (x1 ) − Pr (x2 ) ≤ x1 − x2
(47)
thus for x(k+1 ) we have
(k+1) ∗ 2 2 x − x = Pr (x(k) + ld (k) ) − x∗ 2 ≤ x(k) + ld (k) − x∗ 2 2 = x(k) − x∗ + l 2 d (k) + 2l (x(k) − x∗ )d (k)
(48)
S.-W. Chiou / Applied Mathematics and Computation 266 (2015) 390–403
397
let d (k ) = x∗ − x(k ) , thus in (48)
(k+1) ∗ 2 (k) ∗ 2 2 (k) 2 x − x ≤ x − x + l d − 2l (x(k) − x∗ )(x(k) − x∗ ) 2 2 = x(k) − x∗ + l 2 (x∗ − x(k) ) − 2l (x(k) − x∗ )(x(k) − x∗ ) 2 2 = x(k) − x∗ + l (l − 2 )(x∗ − x(k) ) 2 2 = x(k) − x∗ + l (l − 2 )d (k)
(49)
since 0 < l < 2, we have
(k+1) ∗ 2 (k) ∗ 2 x − x < x − x
(50)
for k = 1, 2, 3... It implies x(k ) − x∗ is monotonically decreasing and x(k ) − x∗ ≤ x(1 ) − x∗ . Theorem 2. (Convergence of CPP) If the solution set ∗ for a minimization problem (30) is nonempty then the proposed CPP is globally convergent, i.e.
2
lim x(k) − x∗ = 0
∀x∗ ∈ ∗
k→∞
(51)
Proof. Since x∗ ∈ ∗ , following the results in Theorem 1, it is easy to check that every sequence {x(k ) } generated by the proposed CPP is bounded. We prove this theorem by contradiction. Suppose
2
lim x(k) − x∗ = η > 0
(52)
k→∞
then
x(k) ⊂ Y = x ∈ : η ≤ ||x − x∗ ||2 , ||x − x∗ ||2 ≤ ||x(1) − x∗ ||2 , ∀x∗ ∈ ∗
(53)
and Y is a closed bounded set. From (49), we have
(k+1) ∗ 2 (k) ∗ 2 2 x − x ≤ x − x − l (2 − l )d (k)
where 0 < l < 2. Because
Y ∩ ∗
2 X ( x ) = l ( 2 − l ) d ( k )
(54)
= φ , then on Y we let
(k )
(55)
Since X (x(k ) ) is Lipschitz continuous on Y, there exists a ε ≥ 0 such that
inf X (x(k) ) = ε
(56)
x(k) ∈Y
From (52), there is a k0 > 0 such that for all k > k0 ,
( k ) ∗ 2 x − x < η + ε
(57)
2
On the other hand, from (54) and (56),
(k+1) ∗ 2 (k) ∗ 2 ε x − x ≤ x − x < η −
2
(58)
This contradicts (52) and completes this proof. 3.3. A tractable computation scheme According to Definition 1, a Nash–Stackelberg solution for a ROPRICE problem can be characterized. The corresponding link congestion pricing and signal settings are simultaneously determined by solving two sub-problems (29)–(30) alternately. The corresponding performance measure in (30) can be continuously improved from iteration to iteration when link congestion pricing and signal settings are updated in the following way until a search direction (k ) and β (k ) defined in (41) vanishes.
( (k+1) , β (k+1) ) = ( (k) , β (k) ) + ( (k) , β (k) )
(59)
In this section, a tractable computation scheme solving for two sub-problems (29)–(30) is proposed. Step 1. Start with initial link congestion pricing β (k ) together with signal settings (k ) and demand growth factor μ(k ) and set index k = 0. Set stopping threshold ε . Step 2. Solve a maximization problem (29) to obtain a demand growth factor μ(k+1 ) in the presence of link congestion pricing β (k) and signal settings (k) . Update an upper bound estimate P¯1 = P1 ( (k) , β (k) , μ(k+1) ) of performance measure in (29) with a worst-case realization taken by a demand growth factor μ(k+1 ) . Details about setting an upper bound estimate about P¯1 with respect to μ are described in the following sub-steps.
398
S.-W. Chiou / Applied Mathematics and Computation 266 (2015) 390–403 Table 1 Initial data for Sioux Falls city network. 1st set data
2nd set data
φ1,10 /ζ φ2,10 /ζ φ1,15 /ζ φ2,15 /ζ φ1,16 /ζ φ2,16 /ζ φ1,17 /ζ φ2,17 /ζ φ1,19 /ζ φ2,19 /ζ β(1,2) β(2,1) β(12,13) β(13,12) β(18,20) β(20,18) 1/ζ μ
25 25 25 25 25 25 25 25 25 25 0.0 0.0 0.0 0.0 0.0 0.0 60 0.91
55 55 55 55 55 55 55 55 55 55 5.0 5.0 5.0 5.0 5.0 5.0 120 0.92
PM (in $)
1126
1633
where 1/ζ and φ jm /ζ respectively denote the common cycle time and green duration for junction m with signal group j measured in s.
Step 2-1. Characterize a responding flow f ( (k ) , β (k ) , μ(k ) ) by solving an inequality (4) subject to link congestion pricing β (k ) and signal settings (k ) under a worst-case realization taken by a demand growth factor μ(k ) . Step 2-2. Perform a first-order sensitivity analysis for flow f ( (k ) , β (k ) , μ(k ) ). Characterize directional derivative ∇μ f ( (k) , β (k) , μ(k) ) along direction of change μ(k) by solving an inequality (15). Step 2-3. Determine an optimal direction of change μ(k ) by solving a maximization problem (29) such that a new demand growth factor μ(k+1 ) along a determined direction μ(k ) can be determined in the following way: μ(k+1 ) = μ(k ) + μ(k ) . Step. 2-4. Update an upper bound estimate P¯1 = P1 ( (k ) , β (k ) , μ(k+1 ) ) for a new demand growth factor μ(k+1 ) . Step 3. Solve a minimization problem (30) to obtain a link congestion pricing β (k+1 ) and signal settings (k+1 ) under a worstcase of realization of a demand growth factor μ(k+1 ) . Update a lower bound estimate P1 = P1 ( (k+1 ) , β (k+1 ) , μ(k+1 ) ) of performance measure in (30). Details about setting a lower bound estimate about P1 with respect to link congestion pricing β and signal settings are described in the following sub-steps. Step 3-1. Determine a responding flow f ( (k ) , β (k ) , μ(k+1 ) ) by solving an inequality (4) with link congestion pricing β (k ) and signal settings (k ) in the presence of a worst-case of a realization of a demand growth factor μ(k+1 ) . Step 3-2. Perform a first-order sensitivity analysis and characterize directional derivatives ∇β f ( (k ) , β (k ) , μ(k+1 ) ) and
∇ f ( (k) , β (k) , μ(k+1) )along direction of change β (k) and (k) by solving an inequality (15). Step 3-3. Perform a first-order sensitivity analysis for an upper bound estimate P1 ( , β ) in (30) with respect to perturbation (k) and β (k) . Characterize generalized gradient ∇ P˜1 and ∇β P˜1 via Eqs. (22) and (23). Step 3-4. Determine optimal directions of change and β by constructing a cutting plane model for a linear approximation of lower bound estimate P1 ( (k ) , β (k ) ) via (43). Step 3-5. Conduct a CPP approach given in Theorem 1 and determine signal settings (k+1 ) and link congestion pricing β (k+1 ) along a determined direction β and such that ( (k+1 ) , β (k+1 ) ) = ( (k ) , β (k ) ) + ( (k ) , β (k ) ). Step. 3-6. Update a lower bound estimate P1 = P1 ( (k+1 ) , β (k+1 ) , μ(k+1 ) ) for new link congestion pricing β (k+1 ) and signal settings (k+1 ) in the presence of a worst-case of a realization of demand growth factor μ(k+1 ) . Step 4. If the distance of performance measure between upper bound estimate P1 and lower bound estimate P1 is within a threshold ε , i.e. |P1 − P1 | ≤ ε then stop. According to (31), we have a Nash–Stackelberg equilibrium solution for a ROPRICE
problem in (24), i.e. ( ∗ , β ∗ ) = ( (k+1 ) , β (k+1 ) ). Otherwise, increase index k by 1 and go to Step 2. 4. Numerical computations and comparisons
In this section, numerical computations were performed at a benchmark problem of aggregated Sioux Falls city network [49] and [50] as shown in Figs. 1 and 2. The performance measure (PM) used in problem (6)–(12) for a ROPRICE problem with signal-controlled links was calculated on the basis of results given in [51]. The travel time and link investment cost functions for a ROPRICE problem with link congestion pricing in problem (6)–(12) were adopted from [49] and [50]. This numerical test includes five signal-controlled junctions and six candidate links for congestion pricing. In order to illustrate effectiveness of the proposed CPP for second-best congestion pricing, two sets of initial data with travel demand growth factors were given in Table 1. The PM values calculated for problem (6)–(12) were also given. Three solution alternatives for ROPRICE with second-best congestion
S.-W. Chiou / Applied Mathematics and Computation 266 (2015) 390–403
Fig. 1. Sioux Falls real-data test network [49].
Fig. 2. 25-node grid network.
399
400
S.-W. Chiou / Applied Mathematics and Computation 266 (2015) 390–403 Table 2 Computational results for Sioux Falls city network. 1st set data
2nd set data
SAB
CCA
GA
CPP
SAB
CCA
GA
CPP
PM (in $) CPU time (in s)
40 40 42 38 37 43 41 39 35 45 5.6 10.2 9.8 8.7 6.4 5.5 90 1.02 982 107
55 55 54 56 49 61 60 50 52 58 6.2 11.5 12.5 10.5 9.2 8.5 120 1.03 890 55
44 46 45 45 50 40 42 48 51 39 4.5 10.5 12.3 11.6 9.1 7.5 100 1.05 887 923
62 60 59 63 58 64 65 57 61 61 5.2 9.9 11.5 7.4 11.9 8.1 132 1.09 875 19
39 51 51 39 45 45 45 45 44 46 9.6 9.2 11.2 12.7 10.2 9.7 100 1.01 1017 112
50 60 60 50 52 58 62 48 49 61 7.2 10.5 8.5 8.5 7.8 9.3 120 1.02 912 76
38 42 40 40 35 45 45 35 45 35 10.1 10.5 8.6 11.2 6.8 9.6 90 1.05 905 987
56 56 57 55 52 60 61 51 59 53 11.5 9.8 9.2 7.5 10.1 5.9 122 1.10 879 17
Improvement rate (in %)
12.79
20.96
21.23
22.29
37.72
44.15
44.58
46.17
φ1,10 /ζ φ2,10 /ζ φ1,15 /ζ φ2,15 /ζ φ1,16 /ζ φ2,16 /ζ φ1,17 /ζ φ2,17 /ζ φ1,19 /ζ φ2,19 /ζ β(1,2) β(2,1) β(12,13) β(13,12) β(18,20) β(20,18) 1/ζ μ
pricing: the sensitivity-based method (SAB) [15], the GA approach employed in [20] and recently employed CCA in [52] for continuous network design problem and in [18] are chosen for comparison with CPP for ROPRICE problem. Implementations for carrying out the following computations were made on SUN SPARC SUNW, 900 MHZ processor with 4 GB RAM under Unix SunOS 5.8 using C++ compiler. The stopping criterion is set when the relative difference in the PM value is less than 0.15%. Computational results for CPP and three other alternatives were summarized in Table 2. As it shown in Table 2, the improvement rates in the value of PM of heuristics ranged from 12% to 45% over initial data set. The proposed CPP achieved the greatest improvement rate in the value of PM of all heuristics at two sets of initial data. Also the proposed CPP took the least CPU time to solve a ROPRICE problem when compared to other solution alternatives. As for other heuristics, the GA approach took the greatest CPU time to solve a ROPRICE problem. The SAB proposed in [15] achieved the least improvement rate of all heuristics for a ROPRICE problem. As shown in Table 2, the recently proposed CCA achieved a comparable improvement rate in the value of PM for a ROPRICE problem as that did the CPP. However, as it observed in Table 2, the CCA took more CPU time to find a Nash–Stackelberg for a min–max problem (24) than that did the proposed CPP. A robust solution is often regarded as a solution for a worst-case scenario where the feasibility of solutions is pursued at the expense of solution optimality against parameter ambiguity. A link congestion pricing and signal settings for a ROPRICE problem can be optimally determined under a worst-case scenario of uncertain travel demand. For a ROPRICE problem with uncertain travel demand, it becomes a trade-off between system performance and mitigation against uncertain travel demand. As it noted in literature [8–11], the MSCP congestion pricing achieves the most efficient usage of road networked system based on the design of system optimal flow. In order to understand the effectiveness and robustness of the proposed heuristic over MSCP (the firstbest congestion pricing), two computational indices for the PM value of a ROPRICE problem are introduced: the infeasibility gain and optimality loss. First, the infeasibility gain of a proposed solution for a ROPRICE problem can be defined as a percent gain for the difference between the PM value of the upper bound estimate and that of a proposed solution in the following manner. Let P˜ denote an upper bound estimate for a MSCP congestion pricing under a worst-case travel demand. Let P ∗ represent a Nash–Stackelberg solution for a min–max problem (24). We define the infeasibility gain of proposed link congestion pricing and signal-setting in the following way:
P˜ − P ∗ P∗
× 100%
(60)
Second, let Pˆ denote a sub-optimal value for a proposed heuristic at nominal condition and let Pl denote a lower bound estimate of MSCP at a nominal condition for a ROPRICE problem. The optimality loss of a proposed solution for a ROPRICE problem can be defined as a percent loss for the difference between the PM value of a lower bound estimate and that of a proposed solution at a nominal condition in the following manner:
Pˆ − Pl Pl
× 100%
(61)
S.-W. Chiou / Applied Mathematics and Computation 266 (2015) 390–403
401
Fig. 3. Comparison for infeasibility gain to MSCP between CPP and CCA.
Fig. 4. Comparison for optimality loss to MSCP between CPP and CCA.
Numerical comparisons for CPP and CCA alternative in terms of the infeasibility gain and optimality loss defined in (60) and (61) with respect to the performance of MSCP were also performed using 10 sets of initial data and computational results were numerically plotted in Figs. 3 and 4. As seen from Figs. 3 and 4, the infeasibility gains of proposed CPP over MSCP under the worst-case scenario of uncertain travel demand were over 30% on average while those of CCA alternative were around 10%. The optimality loss of CPP at nominal condition with regard to MSCP was within 5% while most of CCA alternative was more than 10%. Furthermore, in order to investigate the effectiveness and robustness of the proposed CPP in a general road network, extensive numerical calculations were performed at a grid-size road network as shown in Fig. 2. Numerical comparisons were also made with CCA. In this numerical example, 18 sets of initial data were taken into account. According to (60) and (61), the infeasibility gain and sub-optimality loss of CPP and CCA were calculated and the results were numerically plotted in Figs. 5 and 6 respectively. As expected, the infeasibility gains of CPP in Fig. 5 were between 30% and 50% while those of CCA were no more than 10%. Again, the sub-optimality loss of CPP in Fig. 6 was below 4% while most of CCA was not less than 9%.
402
S.-W. Chiou / Applied Mathematics and Computation 266 (2015) 390–403
Fig. 5. Comparison for infeasibility gain between CPP and CCA at grid-size network.
Fig. 6. Comparison for optimality loss between CPP and CCA at grid-size network.
5. Conclusions and further issues In this paper, we presented a bi-level min–max model for a ROPRICE problem with uncertain travel demand. A proposed cutting plane projection heuristic was presented to effectively reduce overall travel delay to all road users. Using the first-order directional derivatives, a bi-level min–max problem can be reduced to a single-level optimization problem and solved by a computationally tractable solution scheme. A Nash–Stackelberg solution can be effectively characterized. Numerical computations were performed using a realistic road network and good results were obtained. In order to understand the effectiveness and robustness of the proposed heuristic, further numerical tests were also performed and comparisons were made with the first-best congestion pricing. Computational results obtained indicate that the proposed CPP has considerably reduced the worst-case PM value of a ROPRICE problem while incurring a relatively slight loss of optimality at nominal travel demand when comparing those with MSCP. Moreover, in comparison with a recently proposed heuristic CCA, particularly, our computation results showed that the proposed CPP enjoyed greater advantage in gain of infeasibility against the MSCP at worst-case scenario of uncertain travel demand while incurring very modest loss of sub-optimality. Considering a general road network, the growing complexity of signal-controlled road junctions has increased a formidable risk of sudden loss of road network capacity. The vulnerability of road networks as a result of capacity reduction at links downstream is of growing concern. The method investigated in this paper can be applied to a general urban road network design problem with signal-controlled junctions in the presence of uncertain capacity. The present work on robust signal-setting with equilibrium flow for a ROPRICE problem indicated that a bi-level min–max model is attractive to deal with instances of network design under uncertainty. We will discuss these issues of interest in subsequent papers.
S.-W. Chiou / Applied Mathematics and Computation 266 (2015) 390–403
403
Acknowledgments The author was very grateful to anonymous reviewers for their valuable comments on earlier version of this paper. The author is also thankful for helpful arrangement given by Editor-in-Chief: Prof. T.E. Simos. The work reported in this paper has been supported by grants NSC 98-2410-H-259-009-MY3 and NSC-101-2628-H-259-001-MY2 from Taiwan National Science Council. References [1] A. de Palma, R. Lindsey, Traffic congestion pricing methodologies and technologies, Transport. Res. Part C 19 (2011) 1377–1399. [2] B.D. Chung, T. Yao, T.L. Friesz, H. Liu, Dynamic congestion pricing with demand uncertainty: a robust optimization approach, Transport. Res. Part B 46 (2012) 1504–1518. [3] J. Ekstrom, A. Sumalee, H.K. Lo, Optimizing toll locations and levels using a mixed integer linear approximation approach, Transport. Res. Part B 46 (2012) 834–854. [4] T. Tsekeris, S. Voβ , Design and evaluation of road pricing: state-of-the-art and methodological advances, Netnomics 10 (2009) 5–52. [5] L.M. Gardner, A. Unnikrishnan, S. Waller, Solution methods for robust pricing of transportation networks under uncertain demand, Transport. Res. Part C 18 (2010) 656–667. [6] E. Verhoef, A. Koh, S. Shepherd, Pricing, capacity and long-run cost functions for first-best and second-best network problems, Transport. Res. Part B 44 (2010) 870–885. [7] S. Lawphongpanich, D.W. Hearn, M. Smith, Mathematical and Computational Models for Congestion Charging, Springer, New York, 2006. [8] P. Bergendorff, D.W. Hearn, M.V. Ramana, Congestion toll pricing of traffic networks, in: P.M. Pardalos, D.W. Hearn, W.W. Hager (Eds.), Network Optimization, Springer, Berlin, 1997, pp. 51–71. [9] D.W. Hearn, M.V. Ramana, Solving congestion toll pricing models, in: P. Marcotte, S. Nguyen (Eds.) Equilibrium and Advanced Transportation Modeling, Kluwer Academic Publishers, Boston, pp.109–124. [10] M. Yildirim, D.W. Hearn, A first best toll pricing framework for variable demand traffic assignment problems, Transport. Res. Part B 39 (2005) 659–678. [11] H. Yang, H-J. Huang, Mathematical and Economic Theory of Road Pricing, Elsevier, Amsterdam, 2005. [12] S-W. Chiou, An optimal scheme for toll pricing problem, Appl. Math. Comput. 182 (2006) 1127–1136. [13] S-W. Chiou, An optimization model for area traffic control with link tolls, Appl. Math. Comput. 192 (2007) 520–532. [14] S-W. Chiou, Optimization of congestion pricing road network with variable demands, Appl. Math. Comput. 195 (2008) 382–391. [15] H. Yang, W.H.K. Lam, Optimal road tolls under conditions of queuing and congestion, Transport. Res. Part A 30 (5) (1996) 319–332. [16] M. Labbe, P. Marcotte, G. Savard, A bilevel model of taxation and its application to optimal highway pricing, Manage. Sci. 44 (12) (1998) 1608–1622. [17] E. Verhoef, Second-best congestion pricing in general static transportation networks with elastic demands, Reg. Sci. Urban Econ. 32 (3) (2002) 281–310. [18] S. Lawphongpanich, D.W. Hearn, An MPEC approach to second-best toll pricing, Math. Program. 101 (1) (2004) 33–55. [19] X. Zhang, H. Yang, The optimal cordon-based network congestion pricing problem, Transport. Res. Part B 38 (2004) 517–537. [20] A. Sumalee, Optimal road user charging cordon design: a heuristic optimization approach, Computer-Aided Civil Infrastruct. Eng. 19 (5) (2004) 377–392. [21] M. Josefsson, M. Patriksson, Sensitivity analysis of separable traffic equilibrium equilibria with application to bilevel optimization in network design, Transport. Res. Part B 41 (2007) 4–31. [22] S. Lu, Sensitivity of static traffic user equilibria with perturbations in arc cost function and travel demand, Transport. Sci. 42 (1) (2008) 105–123. [23] M.G.H. Bell, Y. Iida, Transportation Network Analysis, Wiley, Chichester, 1997. [24] H. Yang, M.G.H. Bell, Sensitivity analysis of network traffic equilibria revisited: the corrected approach, in: B. Heydecker (Ed.), Mathematics in Transport, Selected Proceedings of the Fourth IMA International Conference on Mathematics in Transport in Honor of Richard Allsop, Elsevier, Oxford, 2007, pp. 373– 395. [25] Y. Li, T. Yan, X. Li, A log-exponential smoothing method for mathematical programs with complementarity constraints, Appl. Math. Comput. 218 (2012) 5900–5909. [26] T. Yan, A class of smoothing methods for mathematical programs with complementarity constraints, Appl. Math. Comput. 186 (2007) 1–9. [27] J.X. Ban, H.X. Liu, M.C. Ferris, B. Ran, A general MPCC model and its solution algorithm for continuous network design problem, Math. Comput. Model. 43 (2006) 493–505. [28] Z.Q. Luo, J.S. Pang, D. Ralph, Mathematical Programs with Equilibrium Constraints, Cambridge University Press, Cambridge, UK, 1996. [29] K. Shimizu, Y. Ishizuka, J.F. Bard, Nondifferentiable and Two-Level Mathematical Programming, Kluwer Academic Press, Boston, 1997. [30] J. Outrata, M. Kocvara, J. Zowe, Nonsmooth Approach to Optimization Problems with Equilibrium Constraints: Theory Applications and Numerical Results, Kluwer Academic Publishers, Boston, MA, 1998. [31] J.R. Birge, F. Louveaux, Introduction to Stochastic Programming, Springer, New York, 1997. [32] A. Ruszczynski, A. Shapiro, Stochastic Programming, Handbooks in Operations Research and Management Science, Elsevier, 2003. [33] A. Shapiro, D. Dentcheva, A. Ruszczynski, Lectures on Stochastic Programming Modeling and Theory, SIAM, Philadelphia, 2009. [34] G. Li, Z. Xiong, Y. Zhou, Y. Xiong, Dynamic pricing for perishable products with hybrid uncertainty in demand, Appl. Math. Comput. 219 (2013) 10366–10377. [35] X. Ban, M. Ferris, L. Tang, S. Lu, Risk-neutral second best toll pricing, Transport. Res. Part B 48 (2013) 67–87. [36] J.M. Mulvey, R.J. Vanderbei, S.A. Zenios, Robust optimization of large-scale systems, Oper. Res. 43 (2) (1995) 264–281. [37] A. Ben-Tal, L. El Ghaoui, A. Nemirovski, Robustness Optimization, Princeton University Press, Princeton New Jersey, 2009. [38] D. Bertsimas, D.B. Brown, C. Caramanis, Theory and applications of robust optimization, SIAM Rev. 53 (3) (2011) 464–501. [39] Y. Lou, Y. Yin, S. Lawphongpanich, Robust congestion pricing under boundedly rational user equilibrium, Transport. Res. Part B 44 (1) (2010) 15–28. [40] J. Wardrop, Some theoretical aspects of road traffic approach, in: Proc. Inst. Civil Eng. II, 1952, pp. 325–378. [41] R.L. Tobin, T.L. Friesz, Sensitivity analysis for equilibrium network flow, Transport. Sci. 22 (1988) 242–250. [42] F.F. Clarke, Optimization and Nonsmooth Analysis, John Wiley & Sons, New York, 1983. [43] Y. Qiu, T.L. Magnanti, Sensitivity analysis for variational inequalities defined on polyhedral sets, Math. Oper. Res. 14 (3) (1989) 410–432. [44] R. Mifflin, Semismooth and semiconvex functions in constrained optimisation, SIAM Control Optim. 15 (1977) 959–972. [45] L. Qi, J. Sun, A nonsmooth version of Newton’s method, Math. Program. 58 (1993) 353–368. [46] C. Lemarechal, Bundle methods in nonsmooth optimization, in: C. Lemarechal, R. Mifflin (Eds.), Nonsmooth Optimization, Pergamon Press, Oxford, 1978, pp. 79–102. [47] K.C. Kiwiel, Methods of descent for nondifferentiable optimization, Lecture Notes in Mathematics, Springer-Verlag, Berlin, 1985, p. 1133. (362). [48] N. Shor, Nondifferentiable Optimization and Polynomial Problems, Kluwer Academic Publishers, Boston, 1998. [49] M. Abdulaal, L.J. LeBlanc, Continuous equilibrium network design models, Transportation Res. Part B 13 (1) (1979) 19–32. [50] C. Suwansirikul, T.L. Friesz, R.L. Tobin, Equilibrium decomposed optimization: a heuristic for continuous equilibrium network design problem, Transport. Sci. 21 (1987) 254–263. [51] S-W. Chiou, TRANSYT derivatives for area traffic control optimisation with network equilibrium flows, Transport. Res. Part B 37 (3) (2003) 263–290. [52] G. Wang, Z. Gao, M. Xu, An MPEC formulation and its cutting constraint algorithm for continuous network design with multi-user classes, Appl. Math. Model. 38 (2014) 1846–1858.