Transportation Research Part B 35 (2001) 83±105
www.elsevier.com/locate/trb
An equivalent continuously dierentiable model and a locally convergent algorithm for the continuous network design problem Q. Meng a, H. Yang a,*, M.G.H. Bell b a
Department of Civil Engineering, The Hong Kong University of Science & Technology, Clear Water Bay, Kowloon, Hong Kong, People's Republic of China b Department of Civil Engineering, The University of Newcastle upon Tyne, NE1 7RU, UK
Abstract The continuous network design problem (CNDP) is characterized by a bilevel programming model and recognized to be one of the most dicult and challenging problems in transportation. The main diculty stems from the fact that the bilevel formulation for the CNDP is nonconvex and nondierentiable, and indeed only some heuristic methods have been so far proposed. In this paper, the bilevel programming model for CNDPs is transferred into a single level optimization problem by virtue of a marginal function tool. By exploring the inherent nature of the CNDP, the marginal function for the lower-level user equilibrium problem is proved to be continuously dierentiable and its functional value and derivative in link capacity enhancement can be obtained eciently by implementing a user equilibrium assignment subroutine. Thus a continuously dierentiable but still nonconvex optimization formulation of the CNDP is created and a locally convergent augmented Lagrangian method is applied to solve this equivalent problem. The descent direction in each step of the inner loop of the solution method can be found by doing an all or nothing assignment. These favorable characteristics indicate the potential of the algorithm to solve large CNDPs. Numerical examples are presented to compare the proposed method with some existing algorithms. Ó 2000 Elsevier Science Ltd. All rights reserved. Keywords: Network design; Bilevel programming; Optimization; Network equilibrium
*
Corresponding author. Tel.: +852-2358-7178; fax: +852-2358-1534. E-mail address:
[email protected] (H. Yang).
0191-2615/00/$ - see front matter Ó 2000 Elsevier Science Ltd. All rights reserved. PII: S 0 1 9 1 - 2 6 1 5 ( 0 0 ) 0 0 0 1 6 - 3
84
Q. Meng et al. / Transportation Research Part B 35 (2001) 83±105
1. Introduction To accommodate the growing level of trac demand in transportation networks, one way people usually adopt is to expand the capacities of the existing congested links or build new links. In this case, how to select the location of these new links and how much additional capacity is to be added to each of these existing links is an interesting problem, trying to minimize the total system costs under limited expenditure, while accounting for the route choice behavior of network users. This problem is referred to as the network design problem (NDP) and can be roughly classi®ed into three categories: the discrete network design problem (DNDP) that deals with the selection of the optimal locations (expressed by 0±1 integer decision variables) of new links to be added; the continuous network design problem (CNDP) that determines the optimal capacity enhancement (expressed by continuous decision variable) for a subset of existing links; and the mixed network design problem (MNDP) that combines both CNDP and DNDP in a network. The objective for the NDP is toward achieving a system optimum by choosing the optimal decision variables in terms of link capacity enhancement or location for the new link. The system optimum normally represents minimization of the total travel time in the network. Recently maximization of network reserve capacity was also introduced for NDP (Wong and Yang, 1997; Yang and Bell, 1998a, b). The decision variables by the road planner aect the route choice behavior of network users, which is normally described by a network user equilibrium model. Two types of user equilibrium model are often adopted: deterministic user equilibrium (DUE) and stochastic user equilibrium (SUE). For either DUE or SUE, the equilibrium link ¯ow can be expressed by the solution of an optimization problem. Mathematically, the bilevel programming is a good technique to describe this hierarchical property of the NDP with an equilibrium constraint. The upper-level problem is to minimize the total system cost and the lower-level problem is to characterize either the DUE or SUE trac ¯ow pattern. LeBlanc (1975) made a pioneer study of the DNDP. Poorzahedy and Turnquist (1982) presented a typical approximate algorithm for solving the integer programming model of the DNDP. Various variations of the DNDP were also formulated and solved by Boyce and Janson (1980) and Chen and Alfa (1991). Recently, Yang and Bell (1998a) provided a comprehensive review of the models and algorithms for the NDP, in which MNDP was also mentioned. In addition, Yang and Meng (2000) adopted bilevel programming methods for modeling the contemporary Build± Operate±Transfer (BOT) highway design problem that can be considered one of the MNDPs as well. Up to date, studies have been overwhelmingly focused on the CNDP and substantial achievements in algorithmic development have been made. Abdulaal and LeBlanc (1979) formulated the CNDP under DUE as a bilevel programming model and the Hook±Jeeves heuristic algorithm was also introduced. Generally, bilevel programming problem is dicult to solve, designing ecient algorithms for CNDP is recognized to be one of the most challenging problems in transportation. One class of the existing methods tries to derive a set of equivalent dierentiable equations for the DUE assignment problem. The CNDP is then reformulated as a constrained dierentiable optimization problem, which can be solved by the existing convergent methods. Tan et al. (1979) expressed the DUE problem by a set of nonlinear and noncovex, but dierentiable constraints in terms of path ¯ow variables. Friesz (1981) extended this result to the multiclass DUE problem. As
Q. Meng et al. / Transportation Research Part B 35 (2001) 83±105
85
an application, Friesz et al. (1993) used a simulated annealing approach to solve the multiobjective equilibrium network design problem as a single level minimization problem. Since the number of paths in the networks of a realistic size is huge, this approach can only be suitable for small, hypothetical networks. In view of the fact that the DUE problem can be described by a variational inequality (Dafermos, 1980; Smith, 1979), Marcotte (1983) transferred the CNDP into a single level equivalent dierentiable optimization problem. The required constraints involve all the extreme points of the closed convex polyhedron for the feasible acyclic multicommodity ¯ow patterns. It is generally dicult to identify all the extreme points for a polyhedron, and in particular, for a moderately large network problem the constraint set might become huge and intractable. In fact, this technique, like the approach used by Tan et al. (1979), belongs to the type of enumeration methods, it is thus dicult to implement the existing algorithms developed based on the dierentiable optimization formulation of the CNDP. Fisk (1984a, b) developed an alternative single level optimization model for the optimal signal control problem. A max-min technique is used to de®ne a gap function in terms of variational inequality formulation of the DUE problem. The gap function is, however, generally not dierentiable. In view of the solution diculty, various heuristic algorithms for the CNDP are developed to try to produce acceptable solutions for large problems, without necessarily guaranteeing optimality. The relation between DUE link ¯ow and link capacity enhancement is nonlinear and implicit. Under certain strict conditions, the derivative of DUE link ¯ow with respect to link capacity enhancement can be obtained by using the sensitivity analysis method (Tobin and Friesz, 1988; Cho, 1988; Yang, 1995; 1997). Various sensitivity analysis-based heuristic algorithms are proposed for the CNDP and relevant problems using the derivative information (Friesz et al., 1990; Yang and Yagar, 1994; Yang et al., 1994). Unfortunately, the ¯ow-capacity enhancement relation may not always be dierentiable. Another class of heuristic algorithms is the so-called iterative optimization-assignment (IOA) algorithms that iteratively solve the upper and lower level optimization problems of the CNDP. Marcotte (1986), Friesz and Harker (1985) and Marcotte and Marquis (1992) conducted detailed performance analyses of this type of algorithm. Furthermore, Suwansirikul et al. (1987) developed an alternative heuristic method called the equilibrium decomposed optimization (EDO) algorithm by approximating the derivative of the objective function in the upper level problem. This approximation requires that the approximated derivatives should have the same sign as the original true derivatives, which is dicult to verify, in particular, for realistically large network problems. With the increasing availability of powerful computers, computationally demanding, but globally optimal methods such as Monte Carlo optimization method or simulated annealing algorithms are rapidly becoming practical. The simulated annealing algorithm is motivated by an analogy to the statistical mechanics of the annealing of solids and has been applied by Friesz et al. (1992) for the CNDP. The algorithm is potentially globally convergent at the expense of extensive numerical computation and search. Instead of assuming the DUE ¯ow pattern, an alternative would be to describe the CNDP under the route choice behavior of network users satisfying SUE principle. Under the Logit-based SUE problem, Davis (1994) provided a group of continuously dierentiable constraints to represent the SUE condition. In this case, the link ¯ow has an analytical expression associated with the capacity enhancement variables and its derivative can also be obtained by implementing stochastic network loading procedure. Thus locally convergent nonlinear programming
86
Q. Meng et al. / Transportation Research Part B 35 (2001) 83±105
algorithms such as GRG2 (Murtagh and Saunders, 1982) and sequential quadratic programming (SQP) are applied and compared. Recently, Marcotte and Zhu (1996) and Luo et al. (1996) have obtained some interesting results in optimality conditions and provided some algorithms for a class of general bilevel programming problem in which the lower level problem is described by variational inequalities. Certain exact penalty functions and the corresponding algorithms are investigated and established using the theory of exact penalization for mathematical programs with subanalytic constraints under certain regularity conditions. In this approach, the key part is to replace the lower level problem with its equivalent KKT conditions. For the DUE problem with either the minimization model or variational inequality formulation, the KKT condition must include all of the path ¯ows. This means that the approach suers the same limitation of using path ¯ows in the formulation discussed previously for the large-scale problems. Therefore, how to design an ecient convergent algorithm for large-scale CNDPs remains a meaningful research topic. This study is devoted to the development of ecient algorithms for the CNDP. We will reformulate the CNDP under the DUE constraint into an equivalent single level continuously dierentiable problem by virtue of a marginal function tool. The DUE conditions of the lower level problem are represented by a single constraint in terms of the marginal function. Chen (1992) and Shimizu et al. (1997) have used marginal function to investigate the optimality property and algorithms for the general bilevel programming problems, where the constraints of lower-level problem include both the upper- and lower-level decision variables. Basically, the marginal function in this case is directionally dierentiable under certain conditions (DemÕynov and Malzememov, 1974, 1990), but generally nondierentiable (Gauvin and Dubeau, 1982). Due to the nondierentiable property of the marginal function, the relevant results including the models and solution methods are derived based on the analysis of nonsmooth problem, one of the hardest problems to solve in optimization. Nevertheless, the nonlinear bilevel CNDP with the DUE constraint entails a special structure. Namely, the link capacity enhancements do not alter the link-path ¯ow incidence relationship, or the upper-level decision variables do not enter the lower-level constraints and thus do not aect the feasible set of the lower-level problem. This is actually a favorable property, from which we are motivated to develop an ecient algorithm for the CNDP. With this property we prove that the marginal function for the CNDP is a continuously dierentiable function, thereby leading to a single level continuously dierentiable but still a nonconvex formulation of the CNDP. Furthermore, we show that an analytical expression of the gradient of the marginal function can be derived, and readily calculated together with its functional value by implementing a DUE assignment only. Consequently, these characteristics lead us to the application of the existing ecient algorithms for conventional nonlinear dierentiable optimization problems. A convergent augmented Lagrangian method is thus employed to ®nd an exact local solution for the CNDP. The subproblem in the algorithm is solved by the convex combination method, in which the descent direction can be obtained by performing an all or nothing trac assignment. This paper is organized as follows. In the next section, an equivalent single-level model of the CNDP and some propositions for the marginal function are presented. Section 3 designs an augmented Lagrangian algorithm and examines how to ®nd the descent direction by using the decomposition technique for the inner loop of the algorithm. In Section 4, the performances of our algorithm are compared with some existing methods with two test examples. Conclusions are
Q. Meng et al. / Transportation Research Part B 35 (2001) 83±105
87
summarized in Section 5. Finally, the mathematical proof of the propositions on the dierentiability of the marginal function for the CNDP is provided in Appendix A.
2. An equivalent model 2.1. Notations and postulates In this section, we present some notations, de®nitions and postulates used in this paper. These postulates are reasonable and are the same as the key assumptions made in the study by Suwansirikul et al. (1987). Notations A W Rw frw va v ya y qw ta
va ; ya ga
ya dwar h
the set of links in the network the set of origin±destination (OD) pairs the set of routes between OD pair w 2 W the path ¯ow on route r 2 Rw the link ¯ow on link a 2 A the vector of link ¯ow, v
. . . ; va ; . . .T the capacity enhancement of link a 2 A the vector of link capacity enhancement, y
. . . ; ya ; . . .T the ®xed OD demand between OD pair w 2 W the travel cost on link a 2 A described as a function of link ¯ow va and capacity enhancement ya the cost function of capacity enhancement of link a 2 A 1 if route r between OD pair w uses link a, and 0 otherwise the relative weight of total capacity enhancement cost and total travel cost in the system design objective function, or the dual variable for the budget constraint P g
ya 6 B, where B is a predetermined budget for network capacity a a2A enhancement
Postulates (1) The link travel cost function ta
va ; ya ; a 2 A is strictly increasing and continuously dierentiable with respect to link ¯ow va ; a 2 A for any ®xed capacity enhancement ya ; a 2 A. (2) The link travel cost function ta
va ; ya ; a 2 A and ota
va ; ya =oya ; a 2 A both are continuous functions with respect to
va ; ya . (3) The capacity enhancement cost function ga
ya ; a 2 A is continuously dierentiable with respect to ya . 2.2. The CNDP under DUE and the relevant marginal function Mathematically, the CNDP under DUE can be expressed by the following bilevel programming model:
88
Q. Meng et al. / Transportation Research Part B 35 (2001) 83±105
min SC
y y
X
X ta
va ; ya va h ga
ya ;
a2A
1
a2A
subject to la 6 ya 6 ua ;
a 2 A;
2
where fva ; a 2 Ag is the solution of the following convex optimization problem: X Z va ta
x; ya dx; min v
0
a2A
subject to X frw qw ;
w2W;
3
4
r2Rw
va
XX w2W r2Rw
frw dwar ;
frw P 0; r 2 Rw ;
a 2 A;
5
w2 W;
6
where la and ua are the lower and upper bounds for capacity enhancement of link a 2 A, respectively. Let X1 denote the set of the feasible link ¯ows for the above lower level user-equilibrium problem. That is, ( ) X X frw dwar ; frw satisfies
4 and
6 : X1 v
va ; a 2 A va w2W r2R w
The marginal function is de®ned as follows: X Z va ta
x; ya dx: u
y min v2X1
Let F
v; y
0
a2A
XZ a2A
va 0
7
ta
x; ya dx and X2 f y j la 6 ya 6 ua g:
Based on Postulate (1), the DUE problem is a strictly convex optimization problem and thus there is a unique solution of the user equilibrium link ¯ow. This means that Assumption 3 in Appendix A holds. According to the above postulates, the other three conditions in Appendix A are also satis®ed. Thus from Theorems 2 and 3 shown in Appendix A, we can directly obtain the following important properties of this marginal function. Property 1. Let
v
y argmin v2X1
XZ a2A
0
va
ta
x; ya dx
be a user equilibrium link ¯ow, then v
y is a continuous function with respect to y.
Q. Meng et al. / Transportation Research Part B 35 (2001) 83±105
It is clear that the value of this marginal function is X Z va
y u
y ta
x; ya dx: a2A
0
89
8
Although the equilibrium link ¯ow va
y; a 2 A is a function with respect to the capacity enhancement vector y, this marginal function is dierentiable and its partial derivative with respect to individual link capacity enhancement yb ; b 2 A is determined by holding va
y ®xed at the optimal value va
y; a 2 A (refer to Theorem 3 in Appendix A). Thus, we have the following useful proposition. Property 2. u
y is a dierentiable function and its gradient is T ou
y ru
y . . . ; ;... ; oya ou
y oya
Z
va
y 0
ota
x; ya dx; a 2 A: oya
9
10
For any ®xed link capacity enhancement y, the corresponding user equilibrium link ¯ow, v
y can be obtained easily by implementing an ecient DUE trac assignment procedure (e.g., Patriksson, 1994). Then the value and the derivative of the marginal function can be obtained straightforwardly through Eqs. (8)±(10). In addition, from Theorem 1 in the appendix if the travel cost function ta
va ; ya ; a 2 A is convex in ya then the marginal function u
y is also a convex function. A simple example: Here we use a simple example to demonstrate the aforementioned two properties of the marginal function. The example network, depicted in Fig. 1, consists of one OD pair, two nodes, two links and the OD demand q 100. It is straightforward to ®nd the user equilibrium link ¯ows associated with the capacity enhancement y as:
Fig. 1. A simple example.
90
Q. Meng et al. / Transportation Research Part B 35 (2001) 83±105
v1
y v2
y
100 ÿ y; 0; y; 100;
0 < y 6 100; y > 100;
0 < y 6 100; y > 100:
Clearly, v1
y and v2
y are two continuous functions in the interval
0; 1. By the de®nition of marginal function (7), analytical expression of the marginal function for this example is 200 ÿ
y=2; 0 < y 6 100; u
y 2 100
100 =2y; y > 100: It is easy to verify that u
y is dierentiable within any y 2
0; 1 and its derivative directly is given as ÿ1=2; 0 < y 6 100; 0 u
y 2 2 ÿ
100 =2y ; y > 100: On the other hand, Z v
y 1 dt1
x; y dx 0; dy 0 Z
v2
y 0
dt2
x; y
v
y2 dx ÿ 2 2 dy 2y
ÿ
1=2; ÿ
1002 =2y 2 ;
0 < y 6 100; y > 100:
This con®rms that Z v
y Z v
y 1 2 dt1
x; y dt1
x; y 0 u
y dx dx: dy dy 0 0 Namely, Eq. (10) holds. Note that although the link cost function t1
v1 ; y in this example is only a constant instead of a strictly increasing function, Properties (1) and (2) of the marginal function still hold. Actually these two properties are conducted in the appendix under the requirement of unique solution for the lower level minimization problem. In other words, the requirements of link travel cost function in Postulate 1 can be relaxed for a speci®c problem like the above simple example. 2.3. An equivalent continuously dierentiable CNDP formulation For simplicity of notation, we de®ne a function termed as a gap function below: X Z va ta
x; ya dx ÿ u
y: h
v; y a2A
0
11
Obviously, this function is nonnegative, continuously dierentiable for any feasible link ¯ow v and capacity enhancement y. Based on the uniqueness property of link ¯ow in the DUE trac assignment, the equation h
v; y 0 is equivalent to these link ¯ows satisfying the DUE conditions for any given capacity enhancement variable y.
Q. Meng et al. / Transportation Research Part B 35 (2001) 83±105
91
Therefore, an equivalent continuously dierentiable CNDP formulation can be stated below: X X min SC
v; y ta
va ; ya va h ga
ya ;
12 v;y
a2A
a2A
subject to la 6 ya 6 ua ; a 2 A; X frw qw ; w 2 W ;
13
14
r2Rw
va
XX w2W r2Rw
frw dwar ;
a 2 A;
h
v; y 0; frw P 0; r 2 Rw ;
15
16
w2W:
17
In this model, the DUE condition is expressed by the set of Eqs. (14)±(17). Indeed, the complicated constraint in this formulation is Eq. (16), in which the explicit formula of the gap function can not be obtained. Fortunately, the value of the gap function and the value of its gradient can both be simply calculated by implementing a DUE trac assignment. These two values are key parameters in designing an ecient solution algorithm here. It should be pointed out that the gap function h
v; y is generally nonconvex. In other words, the above formulation is still a nonconvex minimization program. 3. A solution algorithm For a continuously dierentiable minimization problem, there exist many methods to ®nd a local minimum. Among these algorithms, the augmented Lagrangian algorithm is one of the ecient and locally convergent methods for the optimization problem with nonlinear constraints (Bertsekas, 1982). Here, we employ this type of algorithm to design a solution method for the aforementioned CNDP. One of the advantages of the method is that the nonlinear, implicit constraint (16) can be incorporated into the objective function as a penalty term. The partially augmented Lagrangian function for problem (12)±(17) is de®ned as L
v; y; l; q SC
v; y lh
v; y 12q
h
v; y2 :
18
3.1. The scheme of augmented Lagrangian algorithm Here, we only describe the framework of the augmented Lagrangian algorithm. A detailed description can be found in Bertsekas (1982) and Bazaraa et al. (1993).
0 Step 0. Initialization: Select a feasible link ¯ow v
0 a ; a 2 A, capacity enhancement ya ; a 2 A
0 and a positive constant b. Choose a multiplier l P 0 and a penalty parameter q
0 > 0. Let n : 0.
92
Q. Meng et al. / Transportation Research Part B 35 (2001) 83±105
Step 1. Inner loop: Solve the following subproblem min L
v; y; l; q;
19
v;y
n1 subject to Constraints and y
n1 denote the solution and ÿ
n1 (13)±(15) and (17). Let v
n1 by executing a DUE trac assignment. ;y then the value h v ÿ can be obtained Step 2. Verify stopping criterion: If h v
n1 ; y
n1 6 e, then terminate. Otherwise, go to Step 3. ÿ ÿ Step 3. Update the Lagrangian multiplier: If h v
n1 ; y
n1 6 0:25h v
n ; y
n , then set ÿ l
n1 l
n q
n h v
n1 ; y
n1 ;
20
let n : n 1 and go to Step 1. Otherwise, go to Step 4. Step 4. Update penalty parameter: Set q
n1 b q
n , l
n1 l
n , n : n 1, and go to Step 1. In this algorithm, various penalty parameter-updating scenarios can be adopted for the problems considered here as long as the following condition is satis®ed: 1 X
q
n 1:
21
n0
3.2. A convex combination method for subproblem (19) In the aforementioned augmented Lagrangian method, the key computational issue is how to solve subproblem (19) in an ecient manner. In this section, we use the convex combination method (Frank±Wolfe algorithm) to ®nd a local solution for this subproblem. Since the constraint set contains only simple bounds of link capacity enhancements and ¯ow conservation equations, the linear programming (LP) problem for ®nding the descent direction in each iteration can be decomposed into two simple LP subproblems. For any given feasible solution denoted by
v
k ; y
k , the multiplier l, the penalty parameter q, and the gradient at this point of the Lagrangian function (18) are ÿ ÿ oL
k ota
k ; a 2 A;
22 va 1 l qh v
k ; y
k ta v
k a ; ya
k ova
v;y
v
k ;y
k ova
va ;ya
v
k a ;ya oL dga
k ota h va
k oya
v;y
v
k ;y
k dya ya ya
k oya
va ;ya
v
k a ;ya ! ÿ
k
k oh ; a 2 A;
23 l qh v ; y oya
v;y
v
k ;y
k where the partial derivative value of the gap function is Z v
k a oh oya
v;y
v
k ;y
k 0
! ota
x; ya ou
k dx ÿ oy
k : oya a yy ya ya
24
Q. Meng et al. / Transportation Research Part B 35 (2001) 83±105
93
De®ne the two partial derivative values (22) and (23) as the generalized link cost and link capacity enhancement denoted by ta
k and ga
k , respectively. Namely, oL
k ta ; a 2 A;
25 ov
k
k a
v;y
v
ga
k
;y
oL ; oya
v;y
v
k ;y
k
a 2 A:
26
In the event, solving the following LP problem can generate the descent direction at the current solution: X X ta
k va va ; ya
27 ga
k ya ; min Zk
va ; ya
a2A
a2A
subject to la 6 ya 6 ua ; a 2 A; X frw qw ; w 2 W ;
28
29
r2Rw
va
XX w2W r2Rw
frw P 0;
frw dwar ;
a 2 A;
r 2 Rw ; w 2 W :
30
31
By inspection of its structure, the LP model can be decomposed into two independent LP subproblems by portioning the coecients in the objective function and the constraints below. LP1: X ta
k va ; min
32 va
a2A
subject to (29)±(31). LP2: X ga
k ya ; min
33
subject to (28). In terms of path ¯ow variables, the model (LP1) can be rewritten as XX cw
k frw ; min r w
34
ya
fr
a2A
w2W r2Rw
; r 2 Rw is the path cost associated with the generalized link subject to (29) and (31), where cw
k r cost between OD pair w, that is, X ta
k dwar ; r 2 Rw ; w 2 W :
35 cw
k r a2A
94
Q. Meng et al. / Transportation Research Part B 35 (2001) 83±105
Therefore, the optimal solution of model (LP1) can be simply obtained by the following allor-nothing assignment procedure: (a) For each OD pair w 2 W ; ®nd the shortest path based on the generalized link cost de®ned in (25) and assign the OD demand qw to this path. Let rw denote this path. (b) From these path ¯ows, the optimal solution called auxiliary link ¯ow is X w qw darw ; a 2 A:
36 v
k a w2W
ya
k ;
Let form of ya
k
a 2 A denote the optimal solution of the model (LP2), which can be expressed in the (
la ; ua ;
if ga
k > 0; otherwise;
a 2 A;
37
The descent direction for subproblem (19) can thus be de®ned by using the optimal solutions of LP1 and LP2, respectively. In summary, we state the solution method for the subproblem (19):
0 Step 0. Initialization: Select a feasible link ¯ow v
0 a ; a 2 A, capacity enhancement ya ; a 2 A and a small positive constant e. Let l and q denote the known multiplier and penalty parameter. Let k 0. Step 1. DUE assignment: For given ya
k ; a 2 A, use the DUE assignment procedure to obtain the user equilibrium link ¯ows denoted by v
k a ; a 2 A. Step 2. Calculation of the marginal function: Use Eqs. (8) and (10) to compute the value and the derivative of the marginal function. Step 3. Calculation of the generalized cost: Use Eqs. (22) and (23) to compute ta
k and ga
k ; a 2 A. Step 4. Calculations of auxiliary variables: Perform all-or-nothing assignment based on the generalized travel cost ta
k ; a 2 A to obtain auxiliary link ¯ow vka ; a 2 A. Use Eq. (37) to obtain the optimal solution of LP2 ya
k ; a 2 A. Step 5. Stopping criterion of inner loop: De®ne the descent direction by
k da
k v
k a ÿ va ;
a 2 A;
38
39 da
k ya
k ÿ ya
k ; a 2 A: P ÿ
k
k If a2A ta da ga
k da
k 6 e, then stop. Otherwise, go to Step 6. Step 6. Line search: According to the descent direction given by (38) and (39), obtain an optimal step size denoted by ak by using a line search procedure. Step 7. Updating of link flow and capacity enhancement: Set
k v
k1 v
k a a ak da ;
a 2 A;
40
ya
k1 ya
k ak da
k ;
a 2 A:
41
Let k : k 1, go to Step 1. In this method, we can adopt several line search methods such as exact line search and ArmijoÕs rule. The convergence of this method by using these rules has already been investigated by
Q. Meng et al. / Transportation Research Part B 35 (2001) 83±105
95
Bazaraa et al. (1993). However, the procedure of an exact line search needs a number of DUE trac assignments to ®nd the step size and therefore is computationally intensive for large-scale problems. In actual applications, a predetermined step size ak satisfying the following conditions can be adopted: X X a2k < 1:
42 ak 1 and This type of step size has been used by Powell and She (1982) in solving the stochastic userequilibrium assignment problem. As the subproblem (19) is also a nonconvex optimization problem, the performance of above algorithm heavily depends on the step size determined in Step 7. For a dierent predetermined step size framework, the method may generate dierent solutions and the speed of the method may be dierent. Friesz et al. (1990) provided and investigated some very useful strategies for choosing step size for CNDP. This may guide us to ®nd a relatively better procedure to determine the step size for a speci®c network design problem. 4. Test examples In order to demonstrate validity of the model and the solution method proposed here, two test examples are adopted. The ®rst example involves a small network presented in Fig. 2, which was originally used by Harker and Friesz (1984) for investigating the bounds for the CNDP and subsequently employed to test the dierent solution algorithms for the CNDP by Suwansirikul et al. (1987) and Friesz et al. (1992). The second example involves a middle size network shown
Fig. 2. Test network for Example 1.
96
Q. Meng et al. / Transportation Research Part B 35 (2001) 83±105
in Fig. 3, the network of Sioux Falls city, which was ®rst considered by Abdulaal and LeBlanc (1979). The network of Sioux Falls city may be the largest network that has ever been used extensively by researchers for the test and comparison of the various algorithms for the bilevel CNDP. The detailed data for the two test examples including OD demand and link travel cost function (BPR type) for capacity enhancement are presented in Suwansirikul et al. (1987).
Fig. 3. Test network for Example 2.
Q. Meng et al. / Transportation Research Part B 35 (2001) 83±105
97
Marcotte and Marquis (1992) provided some comparisons for the dierent heuristic methods used by Suwansirikul et al. (1987) with their own methods. Because of the similarity of these methods, here we only compare our algorithm with the methods cited and developed by Suwansirikul et al. (1987) and Friesz et al. (1992). For convenience, Table 1 lists the abbreviations for these methods. In the augmented Lagrangian algorithm for these examples, the following formula is used to update the penalty parameter. q
n1 b q
n :
43
The initial Lagrangian multiplier l is set to be 0 and the penalty parameter q is set to be 100 for the two examples. Parameter b is set to be 10 and 100 for Examples 1 and 2, respectively. The initial value of link capacity enhancement ya
0 ; a 2 A adopted in our algorithm is the same as Suwansirikul et al. (1987) used. The lower and upper bounds for the link capacity enhancement are shown in Tables 2±4 for the various cases. A predetermined step size ak 1=k is adopted in each step k for the convex combination method embedded in the inner loop of the augmented Lagrangian algorithm. Note that the outer iteration involves enlargement of the penalty term and should be terminated in a few limited iterations in actual applications (in fact, two outer iterations are performed in our test examples). Thus only the convergence of the inner loop of convex combination subroutine is a concern and examined here in the examples. The ®rst example contains only two OD pairs, 1 ! 6 and 6 ! 1. Table 2 summarizes two levels of demand for the two cases considered in this example, a higher level and a lower level of demand to look at the design eect of trac congestion. In the second example, the network contains 24 23 552 OD pairs and in total 10 links (link 16, 17, 19, 20, 25, 26, 29, 39, 48 and 74) are candidates selected for expansion. Tables 3 and 4 present comparisons of the augmented Lagrangian algorithm proposed here with the existing four algorithms for Example 1. The comparisons are made in terms of the objective function value, the resultant link capacity enhancements and the required number of solved Table 1 Abbreviations Abbreviation
Title of the algorithm
Source
MINOS HJ EDO SA AL
Modular in-core nonlinnear system Hooke±Jeeves algorithm Equilibrium decomposed optimization Simulated annealing algorithm Augmented Lagrangian algorithm
Suwansirikul et al. (1987) Abdulaal and LeBlanc (1979) Suwansirikul et al. (1987) Friesz et al. (1992) Presented in this paper
Table 2 The level of trac demand for Example 1 Case
Trac demand from 1 to 6
Trac demand from 6 to 1
Total trac demand
I II
5.0 10.0
10.0 20.0
15.0 30.0
98
Q. Meng et al. / Transportation Research Part B 35 (2001) 83±105
Table 3 Comparison of the algorithms for Case I of Example 1a Variable y1 y2 y3 y4 y5 y6 y7 y8 y9 y10 y11 y12 y13 y14 y15 y16 Value of objective function, SC Number of solved DUE, N a
MINOS
HJ
6.58
EDO
1.2
0.13
3.0
6.26
SA
AL
0.0062 3.1639
5.2631 0.0032
0.0064 7.01 0.22 211.25 ±
3.0 2.8 215.08 54
0.13 6.26 201.84 10
6.7240 198.10378 18,300
0.7171 6.7561 202.9913 2700
Note: Total trac demand is 15.0, the upper bound for all ya is 10.0 for the EDO, SA and AL algorithms.
Table 4 Comparison of the algorithms for Case II of Example 1a Variable y1 y2 y3 y4 y5 y6 y7 y8 y9 y10 y11 y12 y13 y14 y15 y16 Value of objective function, SC Number of solved DUE, N a
MINOS
HJ
EDO
SA
AL
4.61 9.86
5.40 8.18
4.88 8.59
10.1740
4.6153 9.8804
7.71
8.1
0.59
0.9
7.48 0.26 0.85
5.7769
7.5995 0.0016 0.6001 0.001 0.1130
1.32 19.14 0.85 557.14 ±
3.9 8.1 8.4 557.22 134
1.54 0.26 12.52 540.74 12
17.2786 528.497 24,300
1.3184 2.7265 17.5774 532.7100 4000
Note: Total trac demand is 30.0, the upper bound for all ya is 20.0 for the EDO, SA and AL algorithms.
user equilibrium problems. Table 5 presents such comparisons for Example 2. From these tables, we can observe that in both examples the solution obtained by our algorithm is quite close to the solution obtained by the simulated annealing approach, which may be considered to be a global
Q. Meng et al. / Transportation Research Part B 35 (2001) 83±105
99
Table 5 Comparison of the algorithms for Example 2a
a
Variable
HJ
HJ
EDO
SA
AL
Initial value of ya y16 y17 y19 y20 y25 y26 y29 y39 y48 y74 Value of objective function, SC Number of solved DUE, N
2.0 4.8 1.2 4.8 0.8 2.0 2.6 4.8 4.4 4.8 4.4 81.25 58
1.0 3.8 3.6 3.8 2.4 2.8 1.4 3.2 4.0 4.0 4.0 81.77 108
12.5 4.59 1.52 5.45 2.33 1.27 2.33 0.41 4.59 2.71 2.71 83.47 12
6.25 5.38 2.26 5.50 2.01 2.64 2.47 4.54 4.45 4.21 4.67 80.87 3900
12.5 5.5728 1.6343 5.6228 1.6443 3.1437 3.2837 7.6519 3.8035 7.3820 3.6935 81.752 2700
Note: The upper bound for all ya is 25.0 for the EDO, SA and AL algorithms.
optimum. On the other hand, the number of solved user equilibrium problems in our method is substantially less than that required by the simulated annealing algorithm. In comparison with the other three algorithms (MINOS, HJ and EDO) from these tables, our method requires a much larger number of user equilibrium problems to be solved, but gives an acceptable solution in terms of the resultant objective function value for these two examples. With the values of the Lagrangian multiplier and the penalty parameter obtained in the last step of the augmented Lagrangian algorithm, convergence of the inner loop method for Example 2 is shown in Figs. 4 and 5, in terms of the change of Zk
va ; ya shown in (27) and the value of objective
Fig. 4. Inner loop convergence of the augmented Lagrangian algorithm in terms of Zk
va ; ya for Example 2.
100
Q. Meng et al. / Transportation Research Part B 35 (2001) 83±105
Fig. 5. Change of the system cost objective function, SC
v; y, with the iteration of the inner loop of the augmented Lagrangian algorithm for Example 2.
Fig. 6. Changes of the gap function, h
v; y, with the iteration of the inner loop of the augmented Lagrangian algorithm for Example 2.
function SC, respectively. It can be seen that the convergence rate is fast initially and becomes slower as the iterative points move close to the local optimal solution. This is because the predetermined step size becomes smaller and smaller with the growing of the iteration number. Finally, Fig. 6 depicts a tendency of a sharp decrease of the gap function value for the DUE problem.
Q. Meng et al. / Transportation Research Part B 35 (2001) 83±105
101
5. Conclusions In this paper we have proposed an equivalent single level continuously dierentiable optimization model for the conventional bilevel CNDP. Although a complicated nonlinear, implicit constraint, which characterizes the DUE condition, is included in our formulation, we have shown that the value and the derivative of this constraint can be obtained readily by implementing a standard DUE trac assignment only. This property greatly facilitates the development and application of ecient existing nonlinear programming algorithms for our problem. We have employed an augmented Lagrangian algorithm to ®nd an exact local solution for the CNDP. In the inner loop of the algorithm, the descent direction is obtained by performing an all or nothing assignment. All these characteristics imply that our method has the potential for application to large-scale transportation networks. As shown from our comparisons with other existing algorithms using previous examples, the proposed method is at least not inferior to existing algorithms and proves to be an ecient alternative to the CNDP. Therefore, our study provides another theoretical scope to deal with the CNDP by identifying its dierentiable nature. Our methods can of course be applied to deal with other bilevel transportation optimization problems of similar nature, such as road pricing and signal control. Acknowledgements The authors thank Bingsheng He of Nanjing University, China, Terry Friesz of George Mason University, USA, M.J. Smith of University of York, UK, Sze Chun Wong of University of Hong Kong, China, and Hai-Jun Huang of Beijing University of Aeronautics and Astronautics, China, for their helpful comments. This study was supported by the Hong Kong Research Grants Council through a RGC-CERG Grant (HKUST719/96E). Appendix A Consider a real function F
x; y; x 2 X1 ; y 2 X2 , where X1 a bounded, closed convex subset of Rm and X2 is a closed subset of Rn . The following assumptions for this function are made. Assumption 1: F
x; y is a continuously dierentiable function with respect to y 2 X2 , for any ®xed x 2 X1 and its gradient is denoted by ry F
x; y. Assumption 2: F
x; y and each component of ry F
x; y are both continuous functions on the set X1 X2 . Assumption 3: There is a unique ®nite solution for the optimization problem minx2X1 F
x; y for any ®xed y 2 X2 and the solution is denoted by x
y. Assumption 4: For any y 2 X2 , 9d1 d1
y > 0 such that 8s 2 Rn and ksk 6 d1 , we have y s 2 X2 . De®ne the marginal function as follows: u
y min F
x; y: x2X1
Then we have the following properties for this function.
A1
102
Q. Meng et al. / Transportation Research Part B 35 (2001) 83±105
Theorem 1. 1. u
y is a continuous function on X2 . 2. Suppose that for each x 2 X1 function F
x; y is convex in y on a convex set X02 X2 , then u
y is convex on X2 . For a proof, see DemÕyanov et al. (1990, p. 187). Theorem 2. x
y is a continuous function with respect to y 2 X2 . Proof. Suppose the contrary holds. Then there exists a sequence of vectors yk ! y when k ! 1 and a number a > 0 such that kx
yk ÿ x
y k > a:
A2
Let xk x
yk . Since xk 2 X1 , we can assume without loss of generality that the sequence fxk g is convergent, and we let x1 denote its limit point. Since F
xk ; yk u
yk ;
A3
let k ! 1 for this equality, then the following statement is valid: F
x1 ; y u
y :
A4
By Assumption 3, we have x1 x
y , this contradicts (A2). The theorem is thus proved. Note that the proof is similar to Lemma 2.1 in DemÕyanov et al. (1990). Theorem 3. u
y is a continuously differentiable function in its domain X2 and its gradient is ru
y ry F
x
y; y:
A5
Proof. According to Assumptions 1 and 4, for 8e > 0; 9d2 min
d1 ; d
x
y; y > 0 such that 8s that satisfies ksk 6 d2 we have
F
x
y; y s ÿ F
x
y; y ÿ ry F
x
y; y; s 6 e:
A6 ksk For notational simplicity, let
u
y s ÿ u
y ÿ ry F
x
y; y; s : g
y; s ksk
A7
Based on the de®nition of the marginal function u
y, the following equalities can be easily derived:
F
x
y; y s ÿ F
x
y; y ÿ ry F
x
y; y; s 6 e:
A8 g
y; s 6 ksk
Q. Meng et al. / Transportation Research Part B 35 (2001) 83±105
103
On the other hand
F
x
y s; y s ÿ F
x
y s; y ÿ ry F
x
y; y; s : g
y; s P ksk
A9
According to the mean-value theorem, there exists a vector y~ y as; 0 < a < 1 such that D E F
x
y s; y s ÿ F
x
y s; y ry F x
y s; ~y ; s :
A10 Then the right hand side of inequality (A9) is equal to s ry F x
y s; ~ : y ÿ ry F
x
y; y; ksk
A11
From the continuity of ry F
x; y, 9d3 min
d1 ; d
x
y; y > 0, such that 8s that satis®es ksk 6 d3 , we have s ry F x
y s; ~ y ÿ ry F
x
y; y; P ÿ e:
A12 ksk Thus, for 8e > 0, 9d min
d2 ; d3 > 0 such that 8s that satis®es ksk 6 d2 , we have jg
y; sj 6 e:
A13
Namely, lims!0 g
y; s 0. Hence the gradient of u
y is ry F
x
y; y. Based on Theorem 2 and Assumption 3, each element of ry F
x
y; y is a continuous function with respect to y. Consequently, u
y is a continuously dierentiable function in y. References Abdulaal, M., LeBlanc, L.J., 1979. Continuous equilibrium network design models. Transportation Research B 13, 19±32. Bazaraa, M.S., Sherali, H.D., Shetty, C.M., 1993. Nonlinear Programming: Theory and Algorithms. Wiley, New York. Bertsekas, D.P., 1982. Constrained Optimization and Largrange Multiplier Methods. Academic Press, New York. Boyce, D.E., Janson, B.N., 1980. A discrete transportation network design problem with combined trip distribution and assignment. Transportation Research B 14, 147±154. Chen, M., Alfa, A.S., 1991. A network design algorithm using a stochastic incremental trac assignment approach. Transportation Science 25, 215±224. Chen, Y. 1992. Bilevel programming problems: analysis, algorithms and applications, Ph.D. dissertation. University de Montreal, Montreal, Canada. Cho, H.-J. 1988. Sensitivity analysis of equilibrium network ¯ows and its application to the development of solution methods for equilibrium network design problems, Ph.D Dissertation. University of Pennsylvania, Philadelphia. Dafermos, S., 1980. Trac equilibria and variational inequalities. Transportation Science 14, 42±54. Davis, G.A., 1994. Exact local solution of the continuous network design problem via stochastic user equilibrium assignment. Transportation Research B 28, 61±75. DemÕyanov, V.F., Malozemov, V.N., 1990. Introduction to minmax. Dover, New York. DemÕynaov, V.F., Malozemov, V.N., 1974. Introduction to Minmimax. Wiley, New York. Fisk, C.S., 1984a. Game theory and transportation systems modeling. Transportation Research B 18, 301±313.
104
Q. Meng et al. / Transportation Research Part B 35 (2001) 83±105
Fisk, C.S. 1984. Optimal signal controls on congested networks. In: Proceedings of the Ninth International Symposium on Transportation and Trac Theory. VNU Science Press, Delft, The Netherlands, pp. 197±216. Friesz, T.L., 1981. An equivalent optimization problem with combined multiclass distribution, assignment and modal split which obviates symmetry restriction. Transportation Research B 15, 361±369. Friesz, T.L., Anandalingham, A., Mehta, N.J., Nam, K., Shah, S.J., Tobin, R.L., 1993. The multiobjective equilibrium network design problem revisited: a simulated annealing approach. European Journal of Operational Research 65, 44±57. Friesz, T.L., Harker, P.T., 1985. Properties of the iterative optimization-equilibrium algorithm. Civil Engineering Systems 2, 142±154. Friesz, T.L., Cho, H.-J., Mehta, N.J., Tobin, R.L., Anandalingham, G., 1992. A simulated annealing approach to the network design problem with variational inequality constraints. Transportation Science 26, 18±26. Friesz, T.L., Tobin, R.L., Cho, H.-J., Mehta, N.J., 1990. Sensitivity analysis based heuristic algorithms for mathematical programs with variational inequality constraints. Mathematical Programming 48, 265±284. Gauvin, J., Dubeau, F., 1982. Dierentiable properties of the marginal function in mathematical programming. Mathematical Programming Study 19, 101±109. Harker, T.L., Friesz, T.L. 1984. Bounding the solution of the continuous equilibrium network design problem. In: Proceedings of the Ninth International Symposium on Transportation and Trac Theory. VNU Science Press, Delft, Netherlands. LeBlanc, L.J., 1975. An algorithm for the discrete network design problem. Transportation Science 9, 183±199. Luo, Z.-Q., Pang, J.-S., Ralph, D., 1996. Mathematical Programs with Equilibrium Constraints. Cambridge University Press, Cambridge, UK. Marcotte, P., 1983. Network optimization with continuous control parameters. Transportation Science 17, 181±197. Marcotte, P., 1986. Network design problem with congestion eects: a case of bi-level programming. Mathematical Programming 34, 142±162. Marcotte, P., Marquis, G., 1992. Ecient implementation of heuristic for the continuous network design problems. Annals of operation research 34, 163±176. Marcotte, P., Zhu, D.-L., 1996. Exact and inexact penality methods for the generalized bilevel programming problems. Mathematical Programming 74, 141±157. Murtagh, M., Saunders, M., 1982. A projected Lagrangian algorithm and its implementation for sparse nonlinear constraints. Mathematical Programming Study 16, 84±117. Patriksson, M. 1994. The trac assignment problem: models and methods. VSP BV, The Netherlands. Poorzahedy, P., Turnquist, M.A., 1982. Approximate algorithms for the discrete network design problem. Transportation Research B 16, 45±56. Powell, W.B., She, Y., 1982. The convergence of equilibrium algorithms with predetermined step size. Transportation Science 16, 45±55. Shimizu, K., Ishizhka, Y., Bard, J.F., 1997. Nondierentiable and two-level mathematical programming. Kluwer Academic Publishers, The Netherlands. Smith, M.J., 1979. Existence, uniqueness, and stability of trac equilibria. Transportation Research B 13, 295±304. Suwansirikul, C., Friesz, T.L., Tobin, R.L., 1987. Equilibrium decomposed optimization: a heuristic for the continuous equilibrium network design problem. Transportation Science 21, 254±263. Tan H.-N., Gershwin S.B., Athans, M. 1979. Hybrid optimization in urban trac networks. Report No. DOT-TSCRSPA-79-7, Laboratory for Information and Decision System, MIT, Cambridge, MA. Tobin, R.L., Friesz, T.L., 1988. Sensitivity analysis for equilibrium network ¯ows. Transportation Science 22, 242±250. Wong, S.C., Yang, H., 1997. Reserve capacity of a signal controlled road network. Transportation Research B 31, 397±402. Yang, H., 1995. Sensitivity analysis for queuing equilibrium network ¯ow and its application to trac control. Mathematical and Computer Modeling 22, 247±258. Yang, H., 1997. Sensitivity analysis for the elastic demand network equilibrium problem with applications. Transportation Research B 31, 55±70. Yang, H., Bell, M.G.H., 1998a. Models and algorithms for road network design: a review and some new developments. Transport Review 18, 257±278.
Q. Meng et al. / Transportation Research Part B 35 (2001) 83±105
105
Yang, H., Bell, M.G.H., 1998b. A capacity paradox in network design and how to avoid it. Transportation Research A 32, 539±545. Yang, H., Meng, Q. 2000. Highway pricing and capacity choice in a road network under a build-operate-transfer scheme. Transportation Research. Yang, H., Yagar, S., 1994. Trac assignment and trac control in general freeway-arterial corridor systems. Transportation Research B 28, 463±486. Yang, H., Yagar, S., Iida, Y., Asakura, Y., 1994. An algorithm for the in¯ow control problem on urban freeway networks with user-optimal ¯ows. Transportation Research B 28, 123±139.