Title Page
 Table of Contents
 Introduction and overview
 Problem formulation and literature...
 Transfer decomposition
 Applications of simplicial...
 Summary and future research
 Biographical sketch

Title: Decomposition techniques for the traffic assignment problem
Full Citation
Permanent Link: http://ufdc.ufl.edu/UF00089537/00001
 Material Information
Title: Decomposition techniques for the traffic assignment problem
Physical Description: v, 107 leaves : ill. ; 28 cm.
Language: English
Creator: Lawphongpanich, Siriphong, 1956-
Publication Date: 1983
Subject: Traffic assignment   ( lcsh )
Traffic engineering   ( lcsh )
Network analysis (Planning)   ( lcsh )
Industrial and Systems Engineering thesis Ph. D
Dissertations, Academic -- Industrial and Systems Engineering -- UF
Genre: bibliography   ( marcgt )
non-fiction   ( marcgt )
Thesis: Thesis (Ph. D.)--University of Florida, 1983.
Bibliography: Bibliography: leaves 102-106.
Statement of Responsibility: by Siriphong Lawphongpanich.
General Note: Typescript.
General Note: Vita.
 Record Information
Bibliographic ID: UF00089537
Volume ID: VID00001
Source Institution: University of Florida
Rights Management: All rights reserved by the source institution and holding location.
Resource Identifier: aleph - 000491142
oclc - 11941230
notis - ACQ9644

Table of Contents
    Title Page
        Page i
        Page ii
    Table of Contents
        Page iii
        Page iv
        Page v
    Introduction and overview
        Page 1
        Page 2
        Page 3
        Page 4
    Problem formulation and literature survey
        Page 5
        Page 6
        Page 7
        Page 8
        Page 9
        Page 10
        Page 11
        Page 12
        Page 13
        Page 14
        Page 15
        Page 16
        Page 17
    Transfer decomposition
        Page 18
        Page 19
        Page 20
        Page 21
        Page 22
        Page 23
        Page 24
        Page 25
        Page 26
        Page 27
        Page 28
        Page 29
        Page 30
        Page 31
        Page 32
        Page 33
        Page 34
        Page 35
        Page 36
        Page 37
        Page 38
        Page 39
        Page 40
        Page 41
        Page 42
        Page 43
        Page 44
        Page 45
        Page 46
        Page 47
        Page 48
        Page 49
        Page 50
        Page 51
        Page 52
        Page 53
        Page 54
        Page 55
        Page 56
        Page 57
    Applications of simplicial decomposition
        Page 58
        Page 59
        Page 60
        Page 61
        Page 62
        Page 63
        Page 64
        Page 65
        Page 66
        Page 67
        Page 68
        Page 69
        Page 70
        Page 71
        Page 72
        Page 73
        Page 74
        Page 75
        Page 76
        Page 77
        Page 78
        Page 79
        Page 80
        Page 81
        Page 82
        Page 83
        Page 84
        Page 85
        Page 86
        Page 87
        Page 88
        Page 89
        Page 90
        Page 91
        Page 92
        Page 93
        Page 94
    Summary and future research
        Page 95
        Page 96
        Page 97
        Page 98
        Page 99
        Page 100
        Page 101
        Page 102
        Page 103
        Page 104
        Page 105
        Page 106
    Biographical sketch
        Page 107
        Page 108
        Page 109
Full Text








The author wishes to express his gratitude and appreciation to the

members of his supervisory committee, Dr. H. P. Benson, Dr. D. J.

Elzinga, Dr. R. L. Francis, Dr. H. W. Hamacher, Dr. D. W. Hearn, Dr. R.

W. Swain, and Dr. S. Tufekci, for their assistance and guidance.

The author is especially grateful to his advisor, Dr. D. W. Hearn,

for suggesting the research topic and for his encouragement and

patience throughout the research.

This research was partially supported by the National Science

Foundation under NSF Grant ECS-79-25065.




ABSTRACT . . . . . . .



1.1 Introduction . . . .
1.2 Overview . . . .


2.1 Problem Formulation . . .
2.2 Literature Survey . . .

. 1 . 5
. . 11


3.1 Partitioning of the Feasible Region in
Transfer Decomposition . . .
3.2 Validation of Transfer Decomposition for the
3.3 Properties of the Mapping n(z) . .
3.4 Algorithms for the Master VIP . .
3.5 Applications: Link Abstraction and Sub-Area
Focusing . . . . .

. 18
VIP 22
. 28
. 38

. 48


4.1 The Restricted Simplicial Algorithm .
4.2 Simplicial Decomposition for the Asymmetric
Problem . . . . .

. 61

. 78


5.1 Summary of Results . . . .
5.2 Summary of Applications . . .
5.3 Future Research . . . .

REFERENCES . . . . . .

BIOGRAPHICAL SKETCH . . . . ... .. ... .107

. .102

. . 1
. . 2

. 5

Abstract of Dissertation Presented to the Graduate Council
in Partial Fulfillment of the Requirements
for the Degree of Doctor of Philosophy




August, 1983

Chairman: Dr. Donald W. Hearn
Co-Chairman: Dr. D. Jack Elzinga
Major Department: Industrial and Systems Engineering

This study is concerned with decomposition techniques for the

traffic assignment problem. In particular, it generalizes earlier

decomposition results for the standard traffic assignment problem in

which the travel cost function has a symmetric Jacobian to the case

where the Jacobian is asymmetric. The standard traffic assignment

problem is often formulated as an optimization problem, whereas the

asymmetric case can be formulated in several ways. However, the varia-

tional inequality formulation seems to be the most useful and is adopt-

ed for this study.

Decomposition ideas originally developed for the optimization

formulation are extended to the variational inequality problem. New

convergent algorithms for both forms of the traffic assignment problem

are given and explored computationally. Particular attention is given

to computer storage requirements and the application of the projection

technique to the master problem of simplicial decomposition. Transfer

decomposition is developed for problems which are only partially

asymmetric. This technique decomposes the problem into two variational

inequality problems, one of which is equivalent to the standard traffic

assignment problem. Furthermore, the impact of transfer decomposition

as a tool in network aggregation is also explored.



1.1 Introduction

The traffic assignment problem (TAP) is a multicommodity flow

problem which arises when users (or commodities) must share a common

resource, the street network. In this case, each user has a travel

plan; i.e., he must leave from one location, his origin, to another

location, his destination, in the network. It is generally assumed

that users always choose the best possible route under the prevailing

traffic conditions. However, each user's choice is dependent upon

other users' choices as well, because the travel time on each street

depends on the total number of cars on that street: more congestion

means longer travel times. The problem is to determine the route that

each user will choose to complete his trip, which, in turn, determines

the amount of traffic on each street. In transportation planning,

these numbers serve as a prediction of how the traffic as a whole

distributes itself over the entire network, and based on this predict-

ion decisions to expand the current street network are made. Present-

ly, the TAP is used as the core of many traffic planning systems.

Although reference to the TAP is used throughout this study,

there are other network models in which congestion plays an important

role: communication systems, railway systems, and military logistic

systems. Thus the discussion in this study also applies to other


nonlinear cost, multicommodity flow problems as well. One interesting

feature of TAP, not shared by other models, is that it permits several

equivalent formulations.

Since in practice even a moderately sized problem is too

large for many computers, we are interested in decomposition techniques

for the TAP. However, it should be emphasized that the decision to

decompose a problem should not depend solely on the size of the prob-

lem. Decomposition should also be considered when the problem exhibits

special structure and when it is desirable to use a special purpose

computer code to achieve some goal (e.g., high numerical accuracy) at

least for part of the network. Such special codes are often incapable

of solving large problems.

This study considers the extension of decomposition techni-

ques for optimization problems to the variational inequality formula-

tion of the TAP. As shall be seen, the variational inequality formu-

lation can be converted to an optimization problem whose objective

function is not necessarily convex or differentiable. However, this

establishes a relationship between a variational inequality and an opti-

mization problem which provides evidence for the possibility of extend-

ing decomposition techniques in nonlinear programming literature to the

variational inequality problem. Our main goal is to explore this pos-

sibility and to unify the ideas for decomposing the variational inequa-

lity formulation of the TAP.

1.2 Overview

Chapter Two discusses different formulations of the TAP and

surveys the past work. From the survey, it is evident that this disser-

tation is one of a few systematic studies of decomposition techniques

for variational inequality problems. Chapters Three and Four present

two complementary decomposition schemes: transfer and simplical decom-

position, respectively. Transfer decomposition is a resource-directed

decomposition which can be viewed as a special case of generalized Ben-

ders' decomposition [Geoffrion 1972]. It was first developed for the

convex programming formulation of the "standard" TAP. In Chapter 3,

this technique is extended to the variational inequality problem.

After the validation of this decomposition scheme, we investigate pro-

perties of the resulting master problem which later lead to a computa-

tional algorithm for its solution. In the final section of this chap-

ter, an application of transfer decomposition to a problem in transpor-

tation planning is discussed.

Chapter Four presents simplicial decomposition, a price-

directed scheme originally developed for convex programming problems.

Also, it is a generalization of Dantzig-Wolfe decomposition and has

been investigated by many authors in a nonlinear programming context.

First, we investigate further algorithmic enhancements for the existing

simplicial decomposition scheme with a special concern for large prob-

lems. In this restricted version, users may specify the number of grid

points to be retained from one iteration to the next. Based on our

computational results, this allows a trade-off between computer storage

and the progress of the algorithm. In addition, the modified version

also allows for approximation of the master problem's solution which in

effect provides theoretical justification for the ongoing implementa-

tions of simplicial decomposition. In the subsequent section, the

standard simplical method is extended to the variational inequality

problem. It is well known that the algorithm may cycle when column


dropping is permitted in every iteration. However, in view of the size

of problems normally considered by traffic planners, column dropping is

essential to any column generating scheme; therefore, we propose a col-

umn dropping rule. Under this rule, columns can be dropped whenever

the algorithm makes satisfactory progress toward the solution. As in

the restricted case, only an approximate solution to the master varia-

tional inequality is needed to establish convergence.

Finally, areas for possible future research are suggested in

Chapter Five.



2.1 Problem Formulation

The traffic network can be represented by a directed graph,

G(N, L), where N is the set of nodes with cardinality n and L is the

set of links or arcs with cardinality Nodes and arcs usually repre-

sent intersections and streets, respectively, and a subset of these

nodes serves as origins and destinations for travel demands. Associ-

ated with each arc in the network there is a travel cost function which

depends on the aggregate flow on the network. The object is to assign

the travel demands along paths connecting origins and destinations so

as to satisfy a certain criterion. Two widely used criteria first pro-

posed by Wardrop [1952] are the system and user optimized criteria.

The system optimized criterion assumes that there is a central authori-

ty which controls the traffic, e.g., in a railroad network; thus, it is

only logical that this authority distributes travel demands so that the

total travel cost or delay for the entire network is minimized. On the

other hand, the user optimized criterion assumes that there is no cen-

tral authority, and the users themselves decide the distribution of

flows. Naturally, users always choose the cheapest path to reach their

respective destinations. However, this does not mean that one user's

decision is totally independent of the others; in fact, one user's de-

cision may affect the decision of every other user since the travel


cost functions depend on the aggregate flow. Under this assumption, an

assignment is said to be a user optimized assignment if it prevents

users from wanting to change their assigned routes.

Since the two criteria do not generally yield the same solu-

tion, the selection of criteria depends on the underlying application

of the resulting model. For the traffic network, the solution obtained

from the model is usually used as a prediction of the traffic volume on

the streets. Decisions to expand existing streets or to construct new

streets often depend on this prediction. Since the transportation

planners have no control over users of the network (drivers), the user

optimized criterion is the most suitable choice.

Mathematically, one can formulate the traffic assignment

problem with user optimized criterion as a nonlinear complementarity

problem (NCP) or a variational inequality problem (VIP); they are equi-

valent to each other. Define
k k
(i) (s t ) = the k-th pair of origin and destination nodes

(OD pair) with nonzero demand, and there are K such

(ii) f. = the amount of travel demand for the k-th OD pair
which traverses arc i; in other words, it is the

flow on arc i due to the k-th OD pair. Let fk
k k k
denote the vector (fl f2 fk ).

(iii) F = the aggregate (link) flow vector (or flow

pattern) induced by fk, k = 1, ., K;

that is F = (FI, F2, ... F) and F. = f.
k k k
(iv) P = the set of paths from s to t which has
cardinalit k and P = M = k m
cardinality m Let P = U P and P = M m

(v) h





(viii) B (B)

= the path flow vector where the j-th element of h

represents the flow on the j-th path from the set P.
k k
= the (fixed) travel demand from s to t

= the travel cost vector and C(F) =

(C (F),C2(F),..., C (F)).

= the node-arc (arc-path) incidence matrix,


Bij =
13 0

I 1


if node i is incident to arc j


if arc i belongs to path j


Then, the vector F is said to be a feasible flow vector if it belongs

to the following set:
k k k k
(2.1) S = { F : F. = f. Bf = r f. > 0, for all i and k }
1 1
k k
where r the requirement vector for the k-th OD pair, has two nonzero

entries with values dk and -dk Also, the vector h is said to be a

feasible path flow vector if it satisfies the following equations:

(2.2) E h. = dk k = 1,..., K

(2.3) h. 2 0 V j P.


(2.4) F = B h,

and the travel cost function can be defined in terms of path flows as


(2.5) (h) = BT C(F)

The travel cost vector has, until recently, been assumed to

have the form C.(F) = C.(F.), i.e., the travel cost on link i depends

only the flow on link i, F.. Under this assumption, the Jacobian of

C(F) is diagonal, hence symmetric, and, as explained in the next

section, the TAP is equivalent to a convex optimization problem.

However, in a realistic situation, the travel cost on each link of the

network also depends upon the flow on other links as well. An example

of such a situation arises in uncontrolled traffic intersections where

the travel cost on a link depends on the traffic on all intersecting

links in an arbitrary manner, i.e., the influence of flow on link i on

the travel cost of link j is not necessarily the same as the influence

of flow on link j on the travel cost of link i. Therefore, the

Jacobian of C(F) should be asymmetric in order to have a more realistic

model. This model has been studied both as a nonlinear complementarity

and a variational inequality problem, and they are the subjects of the

next two subsections.

2.1.1 Nonlinear Complementarity Formulation

It can be argued (e.g., Aashtiani and Magnanti [1979]) that a

path flow vector satisfies the user optimized criterion if and only if

all utilized paths between a given OD pair have equal travel cost which

is lower than or equal to the lowest travel cost among the nonutilized

paths. Under this situation, no user has any incentive to take a dif-

ferent path, and the user optimized path flow can be viewed as the

steady state or equilibrium path flow vector. Formally, h is a user

optimized or an equilibrium path flow vector if it satisfies

(2.6) (C.(h) u ) = 0 if h. > 0 for j c P, k = 1, ..., K,
Sk 3

(2.7) (C.(h ) -uk) 0 if h = 0 for j e k = 1, ..., K,

(2.8) E h = d k = 1 ,..., K,

(2.9) h. > 0 V j C P,
k k
where uk is the equilibrium travel cost between s and t Let

x = (h,u) ,

(2.10) y.(x) = C(h) u j P and k = 1, ..., K,

(2.11) Zk(x) = Z h. d k = 1, ..., K,
jkpk 3


(2.12) G(x) = (Yl(x),...,yM(x),Zl(x),...,ZK(x)).

Then, the equations (2.6) to (2.9) are identical to

(2.13) G(x)x = 0

(2.14) G(x) > 0

(2.15) x 0 ,

which is the standard form of a NCP. (For convenience, we denote the

inner product of two vectors, a and b, by ab, i.e., ab = E a.b..)
i 1

2.1.2 Variational Inequality Formulation

Smith [1979] gave a dynamic interpretation to the user

optimized criterion as follows: Consider two time epochs today and

tomorrow. Users are not permitted to change their routes during today,

but may do so tomorrow provided that the total travel cost, based on

today's path time (cost), is reduced. Note that this interpretation

also permits each user to change his route if he individually would


reduce his travel cost, since any such change must reduce the total
travel cost. Thus, h is an equilibrium path flow vector if no user

can find a faster route based on today's travel time. Mathematically,
* '.
h must be an element of S = { h : h satisfies equations (2.2) and

(2.3) } and

(2.15) C(h )(h h ) V 0 V h c S

or, in terms of link flows, F must be an element of S and

(2.16) C(F )(F F ) 2 0 V F e S.

The problem of finding h (F ) satisfying the conditions stated above

is called a variational inequality problem and is denoted as

VIP [S,C] (VIP[S,C]).

The following theorem from Smith [1979] establishes the

equivalence of the NCP and VIP formulations:
Theorem 2.1 The vector h is a solution to VIP[S,C] if and only if

it satisfies equations (2.6) (2.9).

2.1.3 Equivalent Optimization Problems

It is well-known [Ortega and Rheinboldt 1970] that if C(F)

has a symmetric Jacobian, then C(F) is a gradient mapping for the

following functional:

(2.17) Z(F) = f1 (F F ) C(FO + t(F F )) dt

where FO is an arbitrary vector in the domain of C(F), and the NCP and

VIP formulations can be viewed as the Kuhn-Tucker conditions for the

following optimization problem [Dafermos 1971]:

min { Z(F) : F S } .



Since (h) = C'(F) B, where the matrices (h) and C'(F) denote the

Jacobian of C(h) and C(F), respectively, C(h) has a symmetric Jacobian

if and only if C(F) has one also. Therefore, an optimization problem

analogous to P-2.1 exists for the TAP formulated in path-flow space as


When the Jacobian of C(F) is asymmetric, Hearn and Nguyen
[1982] showed that F is an equilibrium flow if and only if F solves

P-2.2: min { G(F) : F S } ,

where G(F) = max { C(F)(F H) : H e S }. This function G(F) is

called the gap function in Hearn [1982]; generally, it is neither

convex nor differentiable. By further assuming that C(F) is monotone,


(2.18) (C(F1) C(F2))(F1 F2) Z 0 V Fl, F2 C S,

they also showed that the above result also applies to the following

optimization problem:

P-2.3: max { g(F) : F S } ,

where g(F) = min { C(H)(H F) : H E S }. Hearn et al. [1983] show

that G(F) is a dual function of g(F) and vice versa.

2.2 Literature Survey

The survey of decomposition techniques for the TAP is divided

into two subsections, one for the case of symmetric Jacobian and the

other for the asymmetric case. In all of the papers concerning decompo-

sition techniques for the symmetric case, the authors assumed that the

Jacobian of C(F) is diagonal, which implies that Z(F) is separable by

arcs. However, the results, in most cases, still hold when the


assumption is relaxed. For the asymmetric case, only the VIP formula-

tion seems to lend itself to decomposition, since it can be easily

transformed into optimization problems for which there exist many de-

composition techniques. Thus, the emphasis here will be on decomposi-

tion of the VIP formulation for TAP, a rather new area of research as

far as the transportation literature is concerned. It was only in 1979

that Smith first formulated the TAP as a VIP, and, as a result, there

are only three papers on the subject. In view of the rather small num-

ber of papers on this topic, the survey below will include papers on

algorithms for VIP in the transportation literature as well.

In this study, we think of decomposition as a technique which

divides the entire problem into two or more simpler or smaller prob-

lems, one of which is often designated as the master problem and the

rest as subproblems. Some methods which are referred to as decomposi-

tion in the TAP literature, such as Dafermos [1980a], are better clas-

sified as approximation methods, which are discussed below.

2.2.1 The Symmetric Case

Decomposition techniques for the symmetric TAP can be dichot-

omized into classes of resource and price directed techniques. The

latter class includes column generation and simplicial decomposition.

Leventhal et al. [1973] applied the standard column generating techni-

que to the path flow formulation and obtained a subproblem which gen-

erates paths violating the Kuhn-Tucker conditions. Later, Cantor and

Gerla [1974] applied simplicial decomposition to the arc flow formula-

tion, i.e., problem P-2.1, and obtained as a subproblem a linear cost,

uncapacitated TAP. Furthermore, simplicial decomposition can be viewed

as a generalization of the Frank-Wolfe algorithm or the Dantzig-Wolfe


decomposition, and it is also applicable to general nonlinear program-

ming problems with linear constraints. Table 2-1 summarizes the past

work on this technique.

All of the resource directed techniques are of the general-

ized Benders' type, and they partition the underlying network, which,

in this case, is the "resource" to be shared by users. These techni-

ques differ from one another in the manner by which the network is

divided. Barton and Hearn [1979] presented two methods for dividing

the network: one is to split it into two or more subnetworks at the

nodes, and the other is to partition the arc set into two or more sub-

sets which then define the subnetworks. They called these methods

subset and transfer decomposition, respectively. A special case of

subset decomposition is the geographical decomposition technique [Maier

and Robinson 1977] which first decomposes the problem into K TAPs, each

having only one OD pair with positive demand. Then subset decomposi-

tion is applied to each TAP problem with one OD pair, or to use a more

familiar name, a nonlinear minimum cost flow problem. To summarize,

all decomposition techniques in the resource directed class alter the

underlying network structure by decomposing it into two or more sub-

networks, and, for this reason, a more appropriate name for them would

be "structural decomposition."

2.2.2 The Asymmetric Case

Bertsekas and Gafni [1982] proposed a projection algorithm

for the path flow formulation of the TAP, and, noting that the number

of paths are too numerous to generate a priori, they devised a scheme

to generate them sequentially. The resulting algorithm is essentially

the simplicial scheme in the path flow space. Concurrent with our

Table 2-1: Summary of Simplicial

Obj. Function

Feasible Region

Master Problem


Column Dropping

# of major iter.

von Hohenbalken Holloway Sacher Cantor & Gerla Leventhal, et al.
[1975] & [1977] [1976) [1980] [19741 [1973)

pseudo-convex convex quadratic Klienrock's convex
(convex) delay function

polytope compact and polyhedron network network
convex set constraints constraints

nonlinear prog. nonlinear prog. Linear nonlinear prog. nonlinear prog.

van Moeseke's K-T cond. K-T cond. K-T cond. K-T cond.

follows from follows from Caratheodory Caratheodory Uniqueness of
monotonicity monotonicity Theorem Theorem Lagrange mult.

finite infinite finite infinite finite



study, Smith [in press a] independently proposed the simplicial

decomposition scheme in the link flow space and his own algorithm

[Smith in press b] to solve the master VIP. In these papers all gene-

rated extreme points are retained until the algorithm terminates,where-

as our algorithm emphasizes the elimination of superfluous extreme

points whenever possible. Subsequent to our study, Pang and Yu [1982]

viewed the convex hull of retained extreme points during each iteration

as a sequence of subsets of the feasible region and gave conditions for

global convergence of their method. Currently, the sequence of subsets

generated from our algorithm is the only one known to satisfy their


In the transportation literature, there are two major classes

of algorithms: one involves solving the VIP using an approximation of

the mapping C(F) and the other uses existing optimization algorithms to

solve the equivalent optimization problems formulated in Section 2.1.

The prototype algorithm for the first class of algorithm is as follows:

Given Ft e S, let Ft+l solve the following approximate

VIP [S, Ct]: find Ft+l S such that

(2.19) C(Ft+ )(F Ft+) 2 0 V F S ,

where Ct() is an approximation of C(.) at the point F.

Below is a list of possible approximations for the vector function,

C(F). Among them, the projection method, originally introduced by

Dafermos [1980b], is the only globally convergent algorithm; the others

only converge locally. All of these approximations are discussed in a

survey paper by Pang and Chan [1982] and a book by Ortega and Rhein-

boldt [1970]. Therefore, we only list below the references pertaining

to the TAP.


(i) Diagonalization [Abdulaal and Leblanc 1979; Dafermos

1980a, 1982, 1983; and Florian and Spiess 1982]

The function C (*) is defined to be a diagonal function

t t t t t
C. (F) = C. (F,..., F Fi, Fi+1--.. F ) i = 1 ..., m
1 1 l' m
and the approximate VIP becomes the standard TAP.

(ii) Linear Approximation [Aashtiani and Magnanti 1979]

The function Ct(*) has the form

Ct(F) = C(Ft) + M(Ft) (F Ft) ,

where M(Ft) is an m x m matrix. Aashtiani and Magnanti let M(Ft) be

the Jacobian of C(F) at F = F

Outside the transportation literature, the choices of M(Ft)

are more numerous, e.g.,

L(Ft) or

M(Ft) = D(Ft) or

(U(Ft) + D(Ft)) / w*

where D(Ft), L(Ft), and U(Ft) are respectively the diagonal, lower and
t *
upper triangular submatrices of the Jacobian of C(F ), and 0 < w < 2.

(iii) Projection Method [Dafermos 1980b, 1983; and Bertsekas and

Gafni 1982]

The projection method is a special case of the linear

approximation method where the matrix M(F ) is taken to be G, some

fixed symmetric and positive definite matrix. The reason for the term

projection is that the vector F solving the approximate VIP [S,Ct] is

precisely the projection of the point Ft G- (Ft) onto the set S,

where the projection is defined with respect to the G-norm, i.e.,

Ft+l = Pr (t G C(F))


where for a given vector Z, PrS (Z) is the (unique) vector solving the


min { Z Y G : Y S }

where lix |I= (xGx).

Chapter 8 of Auslender [1976] and a recent paper by Dupuis

and Nguyen [1982] discuss optimization algorithms to solve the dual

pair of problems presented in Section 2.1. The algorithms discussed

therein include steepest descent, penalty and barrier methods, cutting

plane and subgradient optimization. Also, Nguyen and Dupuis [1981]

improved the cutting plane algorithm originally developed by Zuhovicki

et al. [1969], and applied it to the TAP. Their paper contains the

most complete computational experiments, so to judge the efficiency of

some of our algorithms, we will use their results as a basis for the




Barton and Hearn [1979] developed the transfer decomposition

technique for the symmetric traffic assignment problem. Also, they

showed that it is a special case of generalized Benders' decomposition

[Geoffrion 1972]. In this chapter, we develop, extend, and generalize

the technique to VIPs over a polyhedral region. When transfer decompo-

sition is applied to a general convex programming problem with linear

contraints, the master problem is generally a nondifferentiable convex

programming problem even when the original problem is differentiable.

Analogously, the master VIP resulting from the extension is a general-

ized VIP, a harder class of VIPs. However, as in the symmetric case,

it can be shown that the master VIP resulting from the TAP reduces to

the standard VIP for which there exist efficient solution techniques.

This chapter is divided into four sections. The first sec-

tion discusses how transfer decomposition partitions the feasible re-

gion which is assumed to be in a certain form. The second section

validates transfer decomposition for VIPs. The third section discusses

properties concerning the master VIP. The last section discusses appli-

cations of transfer decomposition to two problems of traffic assignment

methodology, subarea focusing and network aggregation.

3.1 Partitioning of the Feasible Region in Transfer Decomposition

Consider the following VIP, VIP[S,(HI,..., H H ,., H )]:
1 p 1 q
^1 ^2 p, ^ 1 ^2 ^q
find (u u ,...,u v v ..., v) S such that




P q. q
P i i ^i 4 1 i i
E H. (u ) (u u ) + E H. (v) (v v ) a 0
i=l i=l
i=1 1=

V (U1,..., v1,..., vq) E S,

(1 1 V) which satisfies
where S is a set of (u ,...,u ,v v) which satisfies

i i
(i) A.u = E.z i = ,..., p
1 1

P i q i
(ii) Z E.z + E A.v = r
i=l i=l1

(iii) a.. 5 u < b.. ;

11 i a. v.
13 3 13

where A., A. and E. are
1 1 1
a.., a.., b.., b.., b..

The above VIP

Barton and Hearn. They


' i '
a.. < z. < b..
13 3 13

V i, j

V i, j ,

matrices, r is a vector of scalars, and a..,
are scalars.

includes the class of problems considered by

consider the following multicommodity flow

P-3.1 min E Z.(X.) + Z (Y.)
e jcL 2
jEL1 jL2
subject to
k k k
Bx + B2y = r k = 1, ..., K

X k k
k k
k k
x y k 0 k = 1, ..., K ,

where B = [B1, B2] is the node-arc incidence matrix and B1 and B2

contain columns corresponding to arcs in arc sets L1 and L2 which
k k
partition L. The vectors x and yk denote flows on the arcs in L and

L2, respectively. Problem P-3.1 can also be stated as a VIP: find a



feasible point (X Y ) such that

(3.2) E Z.(X *)(X. X.) + Z.(Y )(Y. Y.) 0 V (X,Y) feasible,
jeL13 3 3 jeL2 3 3 3
1 2

where Z.(*) (Z.(*)) denotes the derivative of Z. (*) (Z.(*)). In

terms of the earlier notation, f = (xk y ), F = (X, Y), C(F)
(C (X), C2(Y)), where C1(X) = (Z (X1), Z(X2),...) and C(Y) =

(Z1(Y1), Z2(Y2),...), and S is the set of (X,Y) feasible to problem


To decompose problem P-3.1 Barton and Hearn introduce artifi-
cial variables w k = 1,..., K, and a node-arc incidence matrix D and

rewrite the feasible region as

k k k k k k k
(3.3) S = ( (X,Y) : Dw + B2y = r ; B1x = Dw ; X = x Y = y,

W = E w and x w y 0

Intuitively, the equations in (3.3) describe the process of removing

arcs in L from the network G(N,L) and rejoining the severed network,

G(N N, L \ L ), by artificial arcs designated by w's. Here, N1

denotes the nodes adjacent to arcs in L \ L1. This process transforms

the original network into two networks, one (the subnetwork) consists

of arcs in.L1 and the other (the master network) consists of the remain-

ing arcs L \ L and the artificial arcs. Figure 3-1 depicts an example

of this transfer decomposition process.

In this example,

L1 = { (5,6),(5,7),(5,9),(6,5),(6,8),(6,9), (7,8),(8,7),(9,7),(9,8)}

L \ L1 = { (1,5),(1,6),(2,5),(2,6),(7,3),(7,4),(8,3),(8,4)}

N1 = {1, 2, 3, 4, 5, 6, 7, 8 }



(a) Original Network




Master Network






(b) Decomposed System

Figure 3-1: Transfer Decomposition

L 2





N2 = {5, 6, 7, 8, 9} = nodes in the subnetwork,

and the set of artificial arcs = {(5,7),(5,8),(6,7),(6,8)}.

Let A1, E and 1 be block diagonal matrices each having B1,
1 1 1
D, and B2, respectively, as its block, the vectors u z and v be
1 K 1 K 1 K
the vector (x ,..., x ), (y ,...,y ) and (w ,...,w ), and finally,

H (u1) and (v1) be the vectors (CI( E xk),..., C( x k)) and

(C2( E yk),...,C2( E yk)). Since

H (u) (u u) + (v v)

^k k ^k ^k k ^k
SE C( x) (x x ) + E C ( y)(y y
k k

= C1(X) (X X) + C2 (Y)(Y Y),

k k
where, as before, X = E x and Y = E y It is easy to establish the

equivalence between VIPs defined by equations (3.2) and (3.1) with

p = q = 1.

3.2 Validation of Transfer Decomposition for the VIP

It suffices to prove the validity of transfer decomposition

as applied to the VIP defined by equation (3.1) with p = q = 1, all

lower bounds equal zero, and all upper bounds equal m. In this case,

the VIP reduces to VIP[S, (H, H)] : find (u ,v ) e S such that

(3.4) H(u*) (u u*) + 1(v) (v v ) b 0 V (u, v) e S,

S = { (u,v) : Ez + Av = r, Au = Ez; 0 5 v., u., z. i m}.

Also, it is assumed that

(i) (H(u), f(v)) is continuous, and

(ii) the set

S(z) = { u : Au = Ez O 5 u. 5 m } ? p V (v,z) E S',

where S' = { (v,z) : Ez + Av = r, 0 < vi, z. < m } From (i), there


must exist at least one solution to VIP[S,(H,i)] (Theorem 4.1, Kinder-

lehrer and Stampacchia [1980]).

We decompose VIP [S,(H,i)] as follows:

The master VIP, MVIP [S',(nr)]:
T *
Find (z ,v ) e S' such that for some E T e C(z )
T V *
(E ) (z z) + H(v )(v v ) > 0 V (z,v) C S',
t T *
where n(z ) = {E : H(u ) A ~ X + a = 0, Xu = 0,
a.(u. m) = 0 V i, a X. a 0, and u c S (z ) } ,
11 i 1
and u solves the following VIP:

The Sub-VIP, SVIP[S(z ),H]:
Find u e S(z ) such that
H(u )(u u ) 2 0 V u e S(z ) .

If H(u) has a symmetric Jacobian, then n(z) is the subdiffer-

ential of the function

V(z) = min { h(u) : u e S(z ) } ,


h(u) = 1 (u u ) H(u + t(u u ))dt,

'i *
and u would solve the problem V(z ). That is, SVIP[S(z ),H] becomes

an optimization problem, and the above process provides a method which

extracts the optimization component from the original VIP. This could

lead to a significant computational improvement since optimization

problems can often be solved more efficiently than VIPs. The following

theorem (Theorem 2.1, Auslender [1976]) which relates Lagrange multi-

pliers (dual variables) to VIPs is essential to the validation of the

above decomposition process.


Theorem 3.1 Let S = { u: Au = r, u a 0 } and H(u) be continuous.

Then, u e S satisfies

H(u )(u u ) 0 V u e S,

if and only if there exists 5 and A such that

H(u ) A T X = 0

Au = 0

X 0.

Proof: Apply the Kuhn-Tucker Theorem to the following linear pro-

gramming problem:

min { H(u )(u u ) : u S } O

The following theorem verifies the above decomposition


Theorem 3.2 VIP[S,(H,H)] and the pair, MVIP[S',(n,H)] and

SVIP[S(z),H], have equivalent solutions.
Proof: First, assume that (u ,v ) solves VIP[S,(H,H)]. Then by
Theorem 3.1, the pair (u ,v ) must also solve the following linear

programming problem:

m *
P-3.2: min { H(u ) (u u ) + H(v )(v v ) : (u,v) S } .

Applying Benders' decomposition to problem P-3.2, we obtain the follow-

ing pair of problems:

P-3.3 min { V(z) + H(v )(v v ) : (z,v) e S' }


P-3.4 V(z) = min { H(u )(u u) : u e S(z) }


By Theorem 2.1 of Geoffrion [1972], there exists a z such that (z ,v )
solves the master problem, P-3.3, and u solves the subproblem, P-3.4.

Furthermore, the function V(z) is convex and by Theorem 3 in Section

9.3 of Lasdon [1970] and Theorem 23.9 of Rockafellar [1970], the

subdifferential, 3V(z), is

8V(z) = { E T : H(u) A T A + a = 0; Au = 0; a (u. m) = 0

V i, and a., X. 0 and u solves P-3.4 }.
1 1

Thus, both problems, P-3.3 and P-3.4, are convex programs and their

solutions must satisfy the following optimality conditions: for prob-
T *
lem P-3.3, there exists a vector E T E v(z ) such that

(t J *
(3.4) (E) )(z z ) + H(v )(v v ) 0 V (z,v) E S',

and for problems P-3.4,

(3.5) H(u )(u u ) Z 0 V u S(z ).

Since 8V(z ) = n(z ), the equations (3.4) and (3.5) define the master
and subvariational inequality problems, respectively. Hence, (z ,v )

and u solves MVIP[S',(n,H)] and SVIP[S(z ),H], respectively.
To show the converse, let (z ,v ) solve MVIP[S',(n,H)], i.e.,
T *
there exists a vector E T c n(z ) such that

T *
(E ) (z z ) + H(v )(v v ) a 0 V (z,v) E S'


(3.6) E(Ez Ez ) + H(v )(v v ) 2 0 V (z,v) E S'.

Since the set S(z) < for all (z,v) e S', there must exist at least
one u feasible to the set S(z), i.e., Au = Ez and 0 < u. < m. Let u


** *
solve SVIP[S(z ),H]. Then, substituting Au and Au for Ez and Ez,

respectively, in equation (3.6) gives the following:

S(Au Au ) + H(v )(v v ) L 0

(AT + ( *
(A ()(u u ) + H(v ) (v v ) 2 0

V (u,v) E S

V (u,v) E S.

T *
However, from the definition of n(z ) we have A T = H(u ) X + a

which can be substituted into equation (3.7) to obtain

[H(u ) A + a ](u u ) + H(v )(v v ) 2 0 V (u,v) E S


(3.8) H(u )(u u ) + H(v )(v v ) a (X a )(u u ) V (u,v) E S.

Furthermore, (X a )(u u ) = E (X. a.)(u. u.) and
2 1 1

(i) if u.


(ii) if u.


= 0, then a. = 0, X. 2 0 and

- a.)(u. u.) = Lu. 0.
2. 2. 2.2

= m, then a. >

- a.) (u. ui)
2. 2. 2.

0, X. = 0 and

-*(u )
= -a.(u. m) 0.

(iii) if 0 < u. < m, then a. = 0, X. = 0 and
1 1 1

( a. a.)(u. u.) = 0.

Thus, the right hand side of equation (3.11) is nonnegative and

(u ,v ) solves VIP[S,(H,H)]. O




3.2.1 An Example

To illustrate the result, consider the following VIP:

S = { (Ul, u2, V, V) : V + = 1, u + u2 = z, Vi, ui, z a-0 }

H(u) = (5u1 + 2u2; 3u1 + 6u2)T

H(v) = (3v1+ v2; 2v1 + v2)

The subvariational inequality problem becomes the following:

find u e S(z) = { u: uI + u2 = z; u. i 0 } such that

(5u1 + 2u2)(ul ul) + (3u1 + 6u2) (u2 u2) > 0 V u e S(z),

and the solution is u = 2z/3 and u2 = z/3. For any z e S', the

mapping n(z) is

n(z) = {A : 4u A + a = 0, 4u X + a2 = 0, au = 0, a. 2 0} = {4u},

and the master VIP can be written explicitly as follows:

find (z,v) e S' = { (z,v) : z + v + v2 = 1; z, vl, v2 0 } such that

4z (z z ) + (3v + v)(v 1) + (2v +.4v2)(v v2) a 0

V (z,v) E S

and the solution is z = 10/26, v1 = 12/26 and v2 = 4/26. This yields

u = 20/78 and u2 = 10/78. Note that

H(u*) = (40/26, 40/26)T and (v) = (40/20, 40/20)T;
H(u ) = (40/26, 40/26) and H(v ) = (40/20, 40/20)

thus, the obtained solution solves the VIP.


3.2.2 A Special Case H(u) and H(v) have Symmetric Jacobians

If both H(u) and H(v) are gradient mappings of convex

functions, h(u) and h(v), then n(z) is the subdifferential of the

following convex function:

P-3.5: V(z) = min h(u)

subject to

Au = Ez

0 5 u. 5 m,

and the master VIP reduces to

P-3.6: min V(z) + h(v)

subject to

Ez + Av = r

0 S vi., z. m.

This is exactly what Barton and Hearn [1979] obtain when generalized

Benders' decomposition is applied to the problem below:

P-3.7: min h(u) + h(v)

subject to

Ez + Av = r

Au = Ez

0 < u., v., z. < m.

3.3 Properties of the Mapping n(z)

The purpose of this section is to investigate properties of

the mapping n(z). As defined previously, n(z) is generally a point-to-

set mapping which further implies that the master VIP belongs to the


class of generalized VIPs [Fang and Peterson 1982]. The solution of a

generalized VIP is significantly more difficult than a standard VIP.

We know of only one proposed, but untested, method [Fang and Peterson

1982]. Therefore, it is advantageous to consider conditions under

which the master VIP reduces to an ordinary VIP, i.e., those for which

n(z) is point-to-point. An important result is that the transfer

decomposition of a TAP always results in a standard VIP for the master


In addition, we also investigate various monotonicity and

continuity properties of the mapping n(z). These properties are impor-

tant in the development of algorithms for the master VIP. The majority

of existing algorithms for VIPs require the mapping to be, at least,

monotone and continuous.

For the sake of completeness, some well known definitions

(e.g., Auslender [1976]) essential to the analysis are listed below.

Definition Let n( *) be a point-to-set mapping defined on S ( R then

(i) n(*) is monotone on S if

(S )(x y) 2 0 V x,y E S and 5 E (x),
x y x

C n(y).

(ii) n(') is strictly monotone on S if

(C E )(x y) > 0 V x,y E S and 5 e n(x),
x y x

e n(y).

(iii) n(*) is strongly monotone on S if there is a positive

scalar a such that


(5 ) (x y) 2 a(x y) (x y)
x y

V x, y e S and S e n (x), e (y).

If n (*) is a point-to-point mapping, E and E are replaced by n (x) and
x y
n (y), respectively, in the above definitions.

Theorem 3.3 If H(u) is monotone on the set of all feasible u, then

n(z) is also monotone on the set of all feasible z.

Proof: For all z, = z2 feasible and E TE E n(z ) and E TE e n(z2)

we have

(E TE E T2)(I z2)

= (1 2) (Ez1 Ez2)

= (5 E2) (Au Au2)

T T 1 2
= (AE A 2)(u u)

=(H(u ) 1 + a H(u2) + 12 a2) (u u2)

=(H(u ) H(u) (u u2)

1 1 2 2 2 1
+ (X a 2 a2) (u u1),

1 2
where u and u are solutions to the subvariational inequality prob-
1 1 2 2
lems, SVIP[S(z ),H] and SVIP[S(z ),H], and (X a ) and (X a ) are

the Lagrange multipliers associated with the solution vectors u and
u respectively. Since H(u) is monotone, the first term of the last

equation is nonnegative. It is shown below that the last term of the

last equation is also nonnegative, and thus n(z) is also monotone.


The last term can be written as

1 1 2 2 2 1
Z (X. a A. + a.)(u. u.),
i i i1 1 1 1

and we must show that each element of the summation is nonnegative.

There are nine cases to be considered:

1 2 1 1 2 2
(i) 0 < u. < m and 0 < u. < m: Here ai, a., Xi and a.

are zero, thus (1 a .2 + a2)(u u) = 0.
S i 1 1 1) 1

(ii) u. = 0 and 0

2 2 1
A., and a. are zero, thus (X.
1 1 1

< u.

- a1

Here 1

2 2
+ a.)(u.
1 1

is nonnegative and a.,

1 1 2
- u.) = X.u. > 0.
i i i

1 2 1 1
(iii) u. = m and 0 < u. < m: Here a. is nonnegative and .i,
1 1 1 i

2, and a2 are zero, thus (X1 a X + a2)(u2 u) = -a(u m) 0
. i i i i i i 1 11 1

(since u. < m).

Similarly, it can be shown that the term

1 1 2 2 2 1
(X. a. A2 + a.)(u. u.) is nonnegative for the following cases:

1 2
(iv) 0 < u. < m and u. = 0.
1 1

(v) u.

(vi) u.

(ii) u


= m and u. = 0.

= 0 and u. = 0.

= 0 and u. = m.


1 2
(viii) u. = m and u. = m.
1 1

1 2
(ix) 0 < u. < m and u. = m. 0
1 1

Corollary 3.4 Assume that S(z ) n S(z2) = ) n(z) is strictly mono-

tone if H(u) is strictly monotone.

Proof: From Theorem 3.3,

(ET1 E 2) (z1 z2) 2 (H(u ) H(u 2))(u u2

1 1 2 2 2 1
+ (1 a x + a )(u u )

a (H(u ) H(u 2)) (u u2) > 0.

The last strict inequality follows from the strict monotonicity of H(u)
1 2 1 2
and the fact that u 1 u (since u e S(z ) and u S(z2)).

Corollary 3.5 Assume that E has full column rank. The mapping n(z) is

strictly monotone if H(u) is strictly monotone.

Proof: We show that the sets S(z ) and S(z ) are disjoint for any

z1 f z2. Assume the contrary, i.e., there exist zI and z2 such that
1 2 1 2
z1 i z2 and u = u for some u e S(zl) and u e S(z2). Since
1 2
Au Ezl, and Au = Ez2,
1 2
0 = A(u u ) = E(z1 z2).

By the assumption that E has full column rank, z1 z2 = 0 or zI = z2.

This contradicts the assumption that z1 # z2. Therefore, the sets

S(z ) and S(z2) must be disjoint and the result follows directly from

the previous corollary. 0

Corollary 3.6 Assume that E and A are nonsingular. n(z) is strongly

monotone if H(u) is strongly monotone.


Proof: From Theorem 3.3,

(E T1 E T2) (zI z2) > (H(u ) H(u2)) (u1 u2),


1. 2
where u and u are solutions to SVIP[S(z ),H] and SVIP[S(z2),H],
1 2
respectively, and z1 7 z2. Since E has full column rank, u 1 u and

the strong monotonicity of H guarantees that the right hand side of

equation (3.9) is positive. Letting p = IIE-1A 112 > 0, we have

(3.10) (H(u ) H(U2))(u u2)

2 a ||u u |2|

2 2
= a lu u2

= 2




IIElA ul u2112

E1(Au1 Au2)112

1( z2

I(z. ) 2I2

1 2
where u E S(z ) and u2 S(z2).


(ET E T2 (zI -

Combine equations (3.9) and (3.10) to

z2~ 2 II z z2 I2'

3.3.1 A Special Case the TAP

Under the assumption that S is the feasible region of a TAP,

it is shown below that the mapping n(z) is point to point and Lipschitz



For the TAP we have the following relationships:

S1 K 1 K k k k k k k
S' = { (w ,..., w y ,., y ) : Dw + B y = r ; x y w k 0 ,

k 1,..., K },

S 1 K k k k
S(w) = (x... x) B = Dw ; x 0, k = 1,..., K } where
1 K
w = (w ,...,w ), and

n(w) = { (DT ,..., DTK): C(X*) Bk = 0; A = 0, xk 0,

k = ,..., K } .

where, as before, C(F) = (C1(X),C2(Y)) is a nonnegative travel cost


For the k-th OD pair, each element of the subvector D Tk is

the "time" to traverse the shortest path in the subnetwork which con-

nects a pair of split nodes, i.e., nodes which appear in both the

master and the subnetwork. The following lemma shows that these times

are the same for all OD pairs.

Lemma 3.6 D k = D V k,A = 1,..., K.

Proof: From Theorem 3.1, Sk is the vector of dual variables for the

following minimum cost flow problem

min C1 (X)x

subject to

k k
If a and b are any pair of nodes in the subnetwork, then b a is the
b a
length of the shortest path from a to b, where [C1(X )]i, the i-th

element of C1 (X ), is the cost associated with arc i. Also each


element of the vector DT has the form, ~ where a and b are

split nodes in the subnetwork.

Assume, without loss of generality, that for some particular

indices k and

[DT]i = a- < [D ]i = E- a.

Since [D Tk]i is the length of the shortest path from a to b in the


[DT]i k a = k[C(X )]j,


where P is the set of indices for arcs belonging to a shortest path
from a to b with respect to the k-th OD pair. From the definition of


[DTi] i a k([C(X*
[D i = = E ([C(X )]j A)
J ab

Sk Xj = E [C(X*)]j [D T ]i
j eP j eP
J ab J ab

= [D ]i [DT]i < 0.

This is a contradiction since X. > 0 V ij, 0

Theorem 3.7 If C (x) is strictly monotone, n(w) is a point-to-point


Proof: From Lemma 3.6, a vector in n(w) is completely determined by

any subvector D (T, for k e {1,...,K}. Furthermore, each component of

DT represents the shortest path time between pairs of connected split


nodes. Theorem 6.2 of Aashtiani and Magnanti [1981] guarantees that

these path times are unique. Thus, the vector D Tk is unique, i.e.,

the set n(w) is a singleton. O

Theorem 3.8 If C (X) is strongly monotone, then

x X* 11 II I I w w
X1 -X2 "2 I l 2 112 '

where X1 and X2 are aggregate equilibrium flows to the traffic assign-

ment problem for the subnetwork with the trip vectors, wl, and w2,

i.e., SVIP[S(wl),C1] and SVIP[S(w2),C], and B is a positive constant.

Proof: (Theorem 4.1, Dafermos and Nagurney [1982].) 0

Corollary 3.9 Theorem 3.8 holds for any standard norm.

Proof: The result follows from the comparability of any standard

norm. O

Theorem 3.10 If C1(X) is strongly monotone and Lipschitz continuous

with respect to the 1-norm, i.e.,

C1(X1) C (X2) I11| p 11 x- x2 Il p > 0

Then, n(w) is also Lipschitz continuous with respect to 1-norm.
Proof: By Theorem 3.7, n(w) is a vector function (a point-to-point

mapping). -Define

(i) X. to be the aggregate equilibrium flow vector with
respect to the subnetwork trip vector, w., and

(ii) P to be a collection of arcs belonging to a utilized
(shortest) path from p to q for the problem SVIP[S(wi),C ].



Iln(wl) n(w2 )1 = K

I [C (X1)]i -
i e

2[C1(X2)i I
i P

where the outer summation is over the set of artificial arcs, i.e.,

pairs of connected split nodes, in the master network. Let

E 1[C (X1)]i -


=max { I [C (X1)]i
(pq) i P 1

- 2[C (X2) I ,
i ep

and let 4 be the number of artificial arcs in the master network.


1In(w ) nl(w2) |1

S K E 1


- 2[C (X2)]i
iE ab

Without loss of generality assume that E 1[C (X1)]i
iEP .

. [C (X2 )]i,


Ih(w1) n (w,) 1i 15 K.? {

- [C (X2)]j

S1[C (X) ]i
i CP

< K.* Z {

5 K- E 2

[C1(X1)]j [Cl(X2)]j


- [C1(X2)]j

5 K I p Xl X2 1

5 K p 8|w w2 l1


In the above sequence of inequalities, the second inequality holds

because P2 is not necessarily the shortest path between a and b with
respect to the cost vector C (X ); the third inequality holds because
for all i c P

[C (X1)]i [C (X2)]i ; |[C1(X1)]i [C1(X2)]i ,

the fourth inequality holds because

E 2 [C(X 1 )]j [C1(X2)]j | [C (X1)]j [C1(X2)lj|
ePab jL1

= [C(x) Cl(X2) II|

1 II
I1 2

and the last inequality follows from Theorem 3.8. Since all scalars

K, 0, p and 8 are positive, the theorem is proved. O

Corollary 3.11 Theorem 3.10 holds for any standard norm.

Proof: The result follows from the comparability property of any

standard norm. O

3.4 Algorithms for the Master VIP

Most algorithms for VIPs require the mapping to be at least

a continuous point-to-point mapping. Thus, the discussion in this sec-

tion will be limited to cases where the mapping n(z) is a vector func-

tion, i.e., the traffic assignment problem.
1 K
It is clear that a feasible flow vector (x ,..., x ) with

respect to a subnetwork trip vector cannot be feasible to the subnet-

work with respect to another subnetwork trip vector. Thus, the sets
S(") and 2 b 2
S(w ) and S(w ) would be disjoint when w r w Also, assume that


C1(X) and C2(X) are strongly monotone, C1(X) is a Lipschitz continuous

vector function, and C (Y) is a continuous vector function. By Corol-

lary 3.4, Theorems 3.7 and 3.10, the mapping n(w) is a Lipschitz con-

tinuous and strictly monotone vector function, and the composed mapping

(n(w),C2(Y)) is a continuous and strictly monotone vector function. In

this case, the cutting plane algorithm [Zuhovicki et al. 1969 and

Nguyen and Dupuis 1981] can be applied to the master problem.


'kO 0
Step 0 Let (w ,Y ) E S' and set t = 0 .

Ai i
Step 1 Define Pt(w,Y) = min (n(w (w w) + C2(Y) (Y Y)

i = 1,..., t) and let

t+l t+l
P-3.8 Pt(w Y ) = max Pt(w,Y) : (w,Y) S' }

Step 2 If (w Y ) = (w Yj) for some j = 1,..., t, stop and
".t+l t+l
(w Y ) is the solution. Otherwise, set t = t+l and go

to Step 1. 0

In Step 1, the value of n(w ) is obtained by solving the sub-

variational inequality problem for each w and problem P-3.8 can be

reformulated as a linear programming problem [Nguyen and Dupuis 1981]

max Z

subject to

ZO < nT(w)( w w) + C2(Y)(Y1 Y) i = 1,..., t

(w,Y) E S',

where S' is the set of linear constraints described earlier.


Depending on the assumptions imposed on the original problem,

other algorithms mentioned in Chapter 2 can be applied to the master

problem as well. For example, if the subnetwork matrix, B1, is a node-

arc incidence matrix of a tree, then n(w) is strongly monotone whenever

C1(x) is strongly monotone. With the tree structure, the matrices B

and D have linearly independent columns, and one of the rows in each

matrix is redundant. By deleting the same row from both matrices, S(w)

is equivalent to

1 K k ^k k
{ (x..., x) : x= Dw x 0, k 1,..., K.

where B1 and D are matrices B1 and D with one row missing, and they are

both nonsingular. Thus, the assumptions of Corollary 3.6 are fulfilled

and the above claim is proven. In the case of a strongly monotone

travel cost vector, most algorithms in Chapter 2 would be applicable to

the master problem whenever the subnetwork is a tree.

3.4.1 The Standard Traffic Assignment Problem

When C(F) = (C1(X),C2(X)) has a diagonal and positive semidef-

inite Jacobian, the master and subvariational inequalities reduce to

the following pair of optimization problems of Barton and Hearn [1979]:

P-3.9: (Master problem)

min V(w) + E Z.(Y.)

subject to
k k k
Dw + B2y = r k = 1,..., K
k k
y ,w k 0
Y= Ey


P-3.10 (Subproblem)/

V(w) = min E Z.I(Y
jcL 3

subject to

k k
Bx = Dw k = 1,..., K

X = x k

x" 0,

where z.(X ) = / C (t)]. dt and Z.(Y) = Yj [C2(t)] dt. Barton and
I 1 0 '1 j J ) 0 2
Hearn solved both problems by the Frank-Wolfe algorithm and noted that

the time consuming step is the line search for the master problem.

Each evaluation of V(w) requires solving the subproblem, P-3.10, which

is another standard traffic assignment problem. To alleviate this

problem, they suggested a heuristic method that performs the line

search, not on the objective function of the master problem, but on a

function defined to be the maximum of all tangent planes generated in

prior iterations. Below, we describe a convergent algorithm based on

the suggestion given in Geoffrion [1972].

From duality theory,

Sk k k
V(w) = maK i f { E Z (X ) u (B x Dw ) } ,
u x ?0 jEL1 k

where u is the Lagrange multiplier associated with the flow conserva-

tion constraints for the k-th OD pair. Equivalently, the master problem

can be stated as follows [Geoffrion 1972]:

min Z

subject to

SE Z.(Y.) + if { E Z.(X) uk B x Dw) V uk
jeL2 x 0 j EL k


k k k

k k k
y ,w 0; Y = E y,

and in practice the following relaxed version is used:

min Z

subject to

Z E Z.(Y.) + if { E Z.(X ) uk(B1xk Dw }
j0 EL2 x3 jeL 3 k
t= 1,..., T

k k k
Dw + B2y = r k = ,..., K

kk k
y ,w k 0; Y = E y

where ukt is usually the optimal Lagrange multiplier with respect to

the t-th subproblem.

Algorithm [Geoffrion 1972]

Step 0 Let (w ,Y0) be an element of S'. Solve the subproblem V(w)

to obtain the Lagrange multiplier, uk'0 for all k. Set T = 0
0 '"0
and OBJ = E Z (Y.) + V(w ).
jEL2 Y 3

Step 1 Solve the following problem:

P-3.11: min Z

subject to

k,t k
Z 2 E Z.(Y.) + E u (Dw
jEL2 k

+ if { E Z (X.) ukt(B x) } t = 0,..., T
x 20 jEL


k k k
Dw + B2y = r
k k
y ,w 0.

rV* *
Let (ZO,w ,Y ) be an optimal solution. If OBJ Z0 < ,

4 4* k,T+1
Step 2 Solve the subproblem V(w ) and obtain u .Set T = T + 1

and go to Step 1. 0

The property below makes problem P-3.11 more appealing.

k,t 1 ""t k,t
Property 3.12 If xt, k = ,..., K solve V(w ), then x V k

also solve

P-3.12: inf { E Z.(X ) I uk (B x ) }
xK 0 jeL1 k
Proof: The vector uk must satisfy the following Kuhn-Tucker condi-
tions for the subproblem V(w ),

C1( E xt) + Btu xt = 0 k = 1,..., K

Bk,t = Dk,t k = 1,..., K

k,t k,t) =
(X )(x )=0

Xk't,xk't Z 0

However, the Kuhn-Tucker conditions for problem P-3.12 consist of only
the first and the last two sets of equations. Therefore, xkt must

also be a Kuhn-Tucker point for problem P-3.12, and by sufficiency of

the conditions xkt is also optimal. O


tk,t k k,t k,t
6 = if { Z Z.(X ) E u (Bx ) } = E Z(X ) u Blx ,
x 0 jeL1 3 k jeL1 k


kt k t k t" t k,t
and by Lemma 3.6, E uktDw = uD E w = u DW where uD = u D for
k k
some k. Thus, the first set of constraints in P-3.11 becomes

Z E Z.(Y.) + utDW + t t = 0,..., T
jeL2 3 3

and P-3.11 is equivalent to

P-3.13: min VT (W) + E Z.(Y.)
J L2

subject to
k k k
Dw + B2Y = r k = ,..., K
k k
.W = E w ; Y = E y
k k
k k
w ,y 0O,

where V (W) = max { utDW + 8t : t = 0,..., T } and P-3.13 could be

solved by any method for the standard TAP.

To illustrate, consider the symmetric traffic assignment

problem in Figure 3-2 [Barton and Hearn 1979]. Figure 3-3 gives two

possibilities for decomposing the original network. The resulting

subnetwork in Figure 3-3(a) is a tree, and since there is only one path

between any pair of nodes, the solution to the subproblem is trivial.

As for the subnetwork in Figure 3-3(b), the subproblem is a symmetric

traffic assignment problem, and in our implementation we "solve" it by

performing five Frank-Wolfe iterations. For both choices of the sub-

network, the approximate (or minmax) master problem, p-3.13, is

"solved" by performing ten Frank-Wolfe iterations. Table 3-1 summa-

rizes the results of the above computations. For comparison, it takes

the Frank-Wolfe algorithm fifty iterations or 1.84 CPU seconds (on the

IBM 3033N computer at the University of Florida) to achieve the

optimal flow


3 4
1 10 20

2 30 40

Trip Matrix

Zi(Fi) = F Ti(1+0.15(t/Ci) )dt.

Cost Function
Cost Function

Arc Ti Ci
1,5 5 10
1,6 6 16
2,5 3 35
2,6 9 .18
5,6 1 50
5,7 5 25
5,9 2 35
6,5 1 50
6,8 5 25
6,9 2 35
7,3 3 25
7,4 6 24
7,8 1 50
8,3 8 39
8,4 6 43
8,7 1 50
9,7 2 35
9,8 2 25

Volume Delay = T. ( +0.15(F./C.) )

Optimal objective value 1453.

Figure 3-2: An Example





Figure 3-3: Choices of the Subnetwork


Table 3-1: Computational Examples for the Minmax Master Problem

No. of Iter.

Figure 3-2(a)















Figure 3-2(b)

Obj .












* The times are in CPU seconds on the IBM 3033N.


objective function value of 1457.47. However, one would not employ

transfer decomposition with the objective of improving computational

times the methods of Chapter 4 are more efficient. The purpose of

this example is to indicate that the technique is implementable when


3.5 Applications: Link Abstraction and Sub-Area Focusing

Both link abstraction and sub-area focusing are concerned with

replacing the original network with a simplified (and smaller) network.

Link abstraction refers to a procedure of collapsing a set of detailed

links into another (smaller) set of aggregate links. The main goal is

to abstract from the physical street or road network a simplified, but

accurate, representation to be used for analysis. In sub-area focus-

ing, the transportation planner is interested in a set of proposed

changes to street network within a given area. The accurate, but expen-

sive, approach is to solve the entire network problem for each of the

proposed changes. The less expensive, inaccurate, and often used ap-

proach is to discard the arcs external to the given area, estimate a

subnetwork trip table, and solve the traffic assignment subproblem for

each proposed change. A compromise is to replace the external arcs

with a simplified set of aggregate arcs.

Thus, both problems deal with the same type of simplification of

the original network, the physical street network for link abstraction

and the external arcs for the sub-area focusing problem. The major

obstacle in both problems involves the process of incorporating the

original link attributes into the aggregate links. Many heuristic

methods [Chan 1976 and Hearn 1978] to construct the aggregate cost


functions for the aggregate links have been proposed, but none can

guarantee that the obtained solution reasonably approximates the true


Under the assumption that the travel cost vector is monotone

and diagonal, transfer decomposition provides a formal mechanism for

incorporating the original link attributes into the aggregate links.

In this case, the master and the subproblem are optimization problems,

i.e., problems P-3.9 and P-3.10. However, the network defining the two

problems are different under this application. Consider the network in

SFigure 3.4(a) as the entire network, and let the indicated part of the

network represent the area to be kept unaltered; this might represent

the network of principal arterials for the link abstraction process or

the target area for the sub-area focusing problem. We shall refer to

this unaltered area as the principal area of the network. It is clear

that the principal area must be kept in the master network since it

cannot be replaced by pseudo-arcs (the set of aggregate arcs) and the

rest of the network is left in the subnetwork(s). Thus, the resulting

problems are

P-3.12: (The master problem)

min Z.(X.) + V(w)
j EL3 3

subject to

k k k
B x + Dw = r k = 1,..., K

X = xk

k k
x ,w O0

k = ,..., K.


Principal Area

(a) The Entire Street Network




\ /

S/ /
\ /
\ /
/ \
/ \
/ \
i,---- -D

(b) The Simplified Problem

Figure 3-4: An Example of Network Simplification


P-3.13: (The subproblem)

V(w) = min E Z.(Y.)
iL2 1 1

subject to
k k
B2y = Dw

Y= E yk

y 0

k = 1,..., K

where B1 and L are the node-arc incidence matrix and the arc set, re-

spectively, for the principal area. Similarly, B2 and E2 are the node-

arc incidence matrix and the arc set, respectively, for the network
external to the principal area. As before, w and D are the artificial

vector and its node-arc incidence matrix. Furthermore, Barton and

Hearn [1979] showed that the objective function V(w) is also a func-

tion of the aggregate flow on artificial arcs; i.e.,

P-3.14: V(w) = V(W) = min Z.(Y.)

subject to

B2y(p,q) = Dw(pq) V artificial arcs (p,q)

Y = Z y(p'q)

y(P'q) r 0 V artificial arcs (p,q),

where W = E w Note that the OD pair index is different from before.

In the subproblem, the origin and destination nodes are the split nodes

and any of the original origin and destination nodes already included

in the subnetwork.

Since V(W) depends on the aggregate flow on the artificial

arcs, it must be the exact aggregate cost function for the aggregate


links in the above discussion. If V(W) is known explicitly, solving

P-3.12 with V(w) replaced by V(W) would produce the correct aggregate
flow in the principal area, i.e., X Unfortunately, V(W) is defined,

in terms of an optimization problem. There is, however, one very
special (and obvious) case where the function V(W) is explicitly known.

Specifically, assume that all external arcs have linear cost. Then, as

illustrated in Figure 3-5, linear path costs could easily be assigned

to the artificial arcs. These path costs are the length of the short-

est paths which connect the head and tail nodes of the artificial arcs

and consist of only arcs.external to the principal area.

Since the evaluation of the function V(W) generally requires

solving an optimization problem, any algorithm which requires many

evaluations of V(W) would be quite slow. One way to alleviate this

problem is to approximate V(W). One of the two easily implementable
methods is curve fitting and the other is using the function V (W).
As defined earlier, V (W) is the maximum of tangent planes generated by

the algorithm in the previous section. For curve fitting, one can

obtain different approximations of V(W) by varying the form and the

norm used in calculating the error of the approximate function. For

the example problem below, we consider the following four variations of

the curve fitting scheme:

(i) LSQ1: V(W) = a + bW + WAW.

Norm = Euclidean norm.

(ii) LSQ2: V(W) = b + WAW.

Norm = Euclidean Norm.

(iii) MINMAX1: V(W) = a + bW + WAW.

Norm = Chebyshev norm.


Cost function

Principal Area

(a) The Original Network

Cost function



I/ 10U

--- ->-(T)
^lO. v-

Principal Area

(b) An Equivalent Network

Figure 3-5: A Network with Linear External Arcs


(iv) MINMAX2: V(W) = bW + WAW.

Norm = Chebyshev norm.

where a is a scalar, b is a vector and A is a matrix.

Consider the symmetric traffic assignment problem in Figure

3-6. Assume that the original network is decomposed as in Figure 3-7,

then the subproblem consists of two independent symmetric traffic

assignment problems, denoted as V (W) and V (W). For the four curve

fitting schemes, we randomly generated fifty flow vectors for the

artificial arcs from the set { (W1, W2, W8) : 0 W1, W2 5 30,

0 W W3 W4 70, 0 ~ W5,-W7 S 40, 0 S W6, W8 : 60 } using a uniform

random number generator. For the second approximation method, the
functions, V T(W) and V T(W), used in the calculations below are

obtained by solving the decomposed problem using the method described

in the previous section where T is equal to twenty, i.e., V i(W) is

taken to be the maximum of the first twenty tangent planes generated by

the algorithm.

Table 3-2 displays the resulting flow vectors for arcs inside

the principal area and the errors due to each approximation. From

these results, it is encouraging that the errors in the solution are

reasonably small. Also, it is interesting to note that V (W) achieves

the smallest errors among the five approximations.

In spite of the fact that the correct aggregate function is

generally unavailable in a convenient form, the analysis in this chap-

ter provides the first rigorous basis for construction of an aggrega-

tion function in practice. The construction of suitable approximate

functions along the lines suggested above is a matter for extensive

computational experimentation and is left to future research.


3 4
1 10 20

2 30 40

Trip Matrix

F 4
Zi (F) = Ti( +0.15(T/Ci) )dt

Cost Function

Optimal Objective Value = 2750.86

Figure 3-6: Original Network

1 5
1 6
2 5
2 6
5 6
5 7
5 8
6 5
6 7
6 8
7 8
7 9
7 13
8 7
8 10
8 13
9 10
9 11
9 12
10 9
10 11
10 12
11 3
11 4
11 12
12 3
12 4
12 11
13 9
13 10


S(w -


- '


Figure 3-7: The Decomposed System

V (W)



Table 3-2:

Computational Results


7,8 5.84 0.0 0.0 0.0 9.16 0.0
7,9 23.55 22.85 22.33 25.97 20.42 24.94
7,13 47.99 43.13 44.22 47.32 40.83 48.61
8,7 0.0 0.0 0.0 0.0 0.0 0.0
8,10 3.25 0.0 0.0 0.0 9.16 0.0
8,13 25.21 34.02 33.46 26.71 29.58 26.45
9,10 0.0 0.0 0.0 0.0 0.0 0.0
10,9 0.0 0.0 0.0 3.61 0.0 0.0
13,9 40.25 45.33 43.45 45.68 44.18 43.62
13,10 32.94 31.82 34.22 28.35 26.24 31.44

Abs. error in the soln.

(1) Euclidean norm
(2) Chebyshev norm
% error in the soln.

(1) Euclidean norm
(2) Chebyshev norm
Error in obj. function

(1) Absolute error
(2) Percent error



16.67% 14.99%
18.36% 17.19%














Simplicial decomposition was originally developed for a

convex programming problem with linear constraints, i.e.,

P-4.1: min { Z(x) : Bx = b, x L 0 }

where Z(x) is a convex function, B is an m x n matrix, b is a constant

vector in Rm, and x is a.vector in Rn. To simplify the notation, the

feasible region S = { x : Bx = b, x a 0 } is assumed to be bounded, and

thus it is compact and representable as a convex hull of its extreme

points. Although the number of extreme points of S is finite, but

generally large, making it impractical to enumerate them a priori.

Therefore it is more efficient to obtain a small subset of extreme

points initially and generate others as required. This is exactly the

basic scheme of simplicial decomposition which can be precisely stated

as follows:

The Basic Algorithm
0 0 ^ 0
Step 0 Let x be an element of S. Set S = {x } and x = x .

Step 1 Let y solve the subproblem: min { VZ(x)y : y S and

set n = n U {y).

Step 2 Let x solve the master problem: minimize Z(x) over the con-

vex hull of 2. If x solves problem P-4.1, stop. Otherwise,

purge 1 of all extreme points with zero weight in the expres-

sion of x as a convex combination of extreme points in Q. Go

to Step 1.



The advantages of this basic algorithm are (i) all con-

straints are assigned to the subproblem; thus the iterates are always

primal feasible and no dual variables are necessary to link the master

to the subproblem, (ii) the number of extreme points is kept minimal by

dropping those with zero weight, (iii) the master problem is relatively

simple, and (iv) its convergence rate coincides with that of the algo-

rithm employed for the master problem. From a theoretical standpoint,

the disadvantages are the upper bound on the cardinality of 0 during

any iteration is n+l (Caratheodory's Theorem, e.g.; Bazaraa and Shetty

[1979]) and the master problem must be solved exactly. The latter does

not pose any real difficulty. In practice, the master problem can

never be solved exactly on a computer and the method still produces

good solutions. The requirement for the exact solution originates from

the argument often used to establish convergence. It rests on the fact

that the number of extreme points is finite and that to each simplex

one can associate with it the optimal objective function value of the

corresponding master problem. Since the optimal objective function

value decreases from one iteration to the next, the same simplex can

never be regenerated and the algorithm must terminate after a finite

number of iterations. However, there is an alternative argument

[Holloway 1974] which only requires a decrease in the objective func-

tion value (not necessarily optimal) to insure convergence. Using this

argument, the algorithm proposed below only requires an approximate

solution to the master problem.

As for the first disadvantage, it is impractical, if not im-

possible, in most applications to reserve enough computer memory to

store n + 1 extreme points. Often, one would allocate enough computer


memory for only a few extreme points and hope that it is sufficient.

Below, we propose a modified simplicial decomposition scheme, called

restricted simplicial decomposition, which allows the user to restrict

the number of extreme points to any desired number r Z 2. Computa-

tional experiments demonstrate an interesting and potentially very

useful trade-off between r and the progress of the algorithm.

After the treatment of restricted simplicial decomposition,

the same basic decomposition scheme is extended to the VIP. A straight-

forward application of the basic algorithm to the VIP does not work; a

counter example showing that the resulting algorithm cycles through

three extreme points is given in Magnanti [in press]. Most authors

[Bertsekas and Gafni 1982 and Smith in press a] circumvent the problem

by disallowing column (extreme point) dropping. In the algorithm for

VIP below, column dropping is allowed when the algorithm has made sig-

nificant progress toward the solution, and, as before, the exact solu-

tion to the master VIP problem is not required.

Before proceeding, it is important to note that the set of

feasible flow vectors for the TAP is an unbounded polyhedral set ex-

pressible in the form of the feasible region of problem P-4.1. For a

TAP with a nonnegative travel cost function C(*), it is easily seen

that the only flow vectors of interest are those that lie in the convex

hull of all extreme aggregate flows. An "extreme aggregate flow vector

(or pattern)" is so called because it is an extreme point of the set of

feasible flow vectors and is the sum of individual flow patterns which

correspond to minimum cost paths between OD pairs for some given value

of C(*). Since the convex hull of these extreme aggregate flows is

compact, the simplicial decomposition scheme is applicable to the TAP.


Henceforth, the feasible region S for the TAP is taken to be the above

convex hull.

4.1 The Restricted Simplicial Algorithm

This algorithm is a variation on the basic algorithm present-

ed earlier for problem P-4.1. It may be viewed as an application of

the projection technique developed by Goldstein [1964] and Levitin and

Polyak [1966] to the master problem, combined with a scheme for retain-

ing at most r grid points where 2 5 r n n + 1 is specified by the user.

Grid points are defined to be elements of S which define the simplices

in Step 2. In the basic-algorithm, grid points are taken to be extreme

points here they may be prior iterates (the x) as well. In fact, it

is crucial to the convergence of the method that the immediate prior

iterate be in the convex hull of the defining grid points in the next

iteration. Without this feature, monotonicity of the objective

function could not be guaranteed.

Although there are many methods for solving the master prob-

lem, the projection technique is a natural choice. It allows for the

approximation of the master problem solution, and the accuracy can be

controlled by varying the number of projections performed per itera-

tion. Also, the simplicity of the master problem's constraints permits

an almost closed form solution to the projection problem.

Bertsekas [1979] described the projection algorithm as

follows: assume that Z(x) in problem P-4.1 is a twice continuously dif-

ferentiable convex function. (For the TAP, the function Z(x) is as de-

fined in equation (2.17) and when C(*) is monotone and continuously

differentiable the function Z(x) is convex and twice continuously

differentiable.) Since the set S is assumed to be compact, there


exists a constant Y > 0 such that

(4.1) IlVZ(x) VZ(y) | ; y i x yl V x,y E S,

i.e., VZ(x) is also Lipschitz continuous [Bertsekas 1976]. In addi-

tion, assume that {Mk) is a sequence of matrices satisfying

(4.2) a llx 112 < x Mkx s5 lx 12 V x null space of B, and V k

where a and 8 are positive scalars and

(4.3) a > y/2.

Then an approximation of the objective function using Mk can be

written as

(4.4) Z(y I ) = Z(x) + VZ(x)(y x) + ((y )Mk(y x)

and the solution to

P-4.2: min { Z(y x) : y S

is the projection of the vector x Mk VZ(X) onto the set S with

respect to the norm lIx I= (x Mx) Under the given assumption, Bert-
sekas [1976] proved that x's generated successively by solving problem

P-4.2 converge to x the solution of problem P-4.1.

Since problem P-4.2 can be formulated in terms of extreme (or

grid) points of S, i.e.

P-4.3: min { Z(AX x) : 1. = 1, X. > 0 V i }
1 1

where columns of A represent the extreme (or grid) points of S, its

solution can be calculated quite easily when the matrix A MkA is


diagonal. The solution method is almost closed form because the objec

tive function of problem P-4.3 is a strictly convex quadratic in the A

variables with no crossproduct terms. The required calculations are

similar to those described in Section 5 of Held et al. [1974].

Let 0 denote the set of retained grid points and H(0) denote

the convex hull of f, i.e.,

H() = x : x = Z X.a E A. = 1, X. i 0 and a e }

Restricted Simplicial Decomposition Algorithm

Step 0 Let x be an element of S. Set 0 = {x0} and t = 0.

Step 1 Solve

P-4.4: VZ(xt)yt = min { VZ(xt )y : y S }.

If VZ(x )(y x ) 0, stop and x is a solution. Other-


(i) let St+l = t U {yt) if It I < r, or

(ii) let nt+l be St with two of its elements replaced by x

and yt if I t = r.


Step 2 Solve

P-4.5: Z(xt+lI t) = min { Z(x t) x E H(t+l) }

Purge t of all grid points with zero weight in the expres-
sion of x as a convex combination of the grid points. Set

t = t+l and go to Step 1. O

It is important to observe that x always belongs to H(t+l)

at the beginning of Step 2. In Step 1, if it I < r, then xt H(Qt)

implies that xt H( t U {yt} ); otherwise, x belongs to


H( t+) by construction. This observation is later used to prove that

the objective function decreases at each iteration.

Under the column dropping scheme in Step 1, the set of re-

tained grid points always consists of the newly generated extreme point

and the previous iterate for the case r = 2. This is exactly the

scheme in Frank and Wolfe [1959]. When r > 2, two of the retained grid

points are replaced by the same two points mentioned above if the total

number of retained grid points equals r. Otherwise, the newly generat-

ed extreme point is added to the current set of grid points.

We now prove the convergence of the restricted simplicial

algorithm. The argument is a synthesis of known results from Bertsekas

[1978] and Holloway [1974]. Lemmas 4.2 and 4.3 show that the approxi-

mate objective function decreases from iteration to iteration which, in

turn, implies that the same is true for the real objective function.

From this conclusion, the convergence proof follows.

Lemma 4.1 VZ(x ) (y x ) = 0 V y e S if and only if x is an opti-

mal solution to problem P-4.1.

Proof: This is a standard nonlinear programming result. O

Lemma 4.2 If x is not optimal to problem P-4.1, then

Z(xt+l x) < Z(x x = Z(xt).

t t+l t+l
Proof: By the above discussion, xt H( t). Since x minimizes
t t+l
Z(x I xt) over H(tl)

t+1 t t t
Z(x t xt (x xs

If Zt+l It t
If Z(x Ix ) = Z(x x ), then x is also a minimum and


VZ(xt [ x) (y x) = z(xt) (y xt) 0

which contradicts the fact that x is not optimal to problem P-4.1 and

the stopping criterion is not satisfied, i.e., VZ(xt) (yt x ) < 0. O

Lemma 4.3 Given that Mt, t = 1, 2,..., satisfy the assumptions above,

t+1 I t ^
Z(x I x) < Z(xt xt) implies Z(xt+) < Z(xt)

Proof: The argument is essentially from Bertsekas [1979]. Since
x solves problem P-4.5,

St+1+ t t+t t
VZ(xt+l xt) (x xt+l) 0 V x H(t)

(4.5) or [ VZ(xt) + Mt(xt+- x ) ] (x x t-l) a 0 V x e H() .

In particular, let x = x and (4.5) becomes

t t+1 t t+1 t
SVZ(x ) + Mt(xt x ) ] (x xt ) 0.

(4.6) -VZ(xt)(xtl x ) (x x )Mt(x x )

The second order Taylor expansion of the objective function at the
point x gives

Z(xt+) Z(x ) + VZ(x ) (xt+l x)

+ / (VZ(x + u(x x)) VZ(xt))(xt+l x)du

t t+1 t t+1 t
Z(xt) Z(xt+l) = -VZ(xt)(x+ xt)

t t+l t t t+1 t
f (VZ(xt + u(x x)) VZ(x ))(x tl x )du


t t+l t
S-VZ(xt) (x t+- x )

II Z(x + u(x -t+ x ) Z(xt) II Ix xt II du

S(x t )Mxt+1 x f u Y x +l t 2 du

Sa IIxt xt2 /2 -x xt 112 = (a /2)Ixt+l xt2 > 0.

The first inequality follows from the Schwartz inequality, the second

follows from the Lipschitz condition and (4.6), the third follows from
t+l t
the property of Mt and the last follows from the fact that x + x O

Having established the monotonicity in the objective function

values of the iterates, the convergence of the sequence {xt) to an op-

timal solution of problem P-4.1 follows from arguments similar to those

of Lemma 4.2 and Theorem 4.3 of Holloway [1974] (see also Zangwill

[1969]). Therefore, the proofs for the next two results are omitted.

Lemma 4.4 A sequence {xt}t generated by the algorithm cannot admit a

subsequence {xt}t T satisfying the following properties:

(i) x x t T

(ii) y + y" t e T

(iii) VZ(xO")(y x-) < 0.

Theorem 4.5 Under the stated assumptions, the algorithm converges to

an optimal solution of problem P-4.1.

The algorithm as stated only performs one projection per
t -1 t
iteration, and since the projection of the vector x M Vz(x ) onto

H(t+ ) is not necessarily optimal to the original master problem,



P-4.6 min { Z(x) : x e H(St+l) }

the same extreme point, y is often successively regenerated in Step

1. One way to insure that yt does not belong to t is to solve the

original problem exactly. However, the following theorem suggests that

only an approximate solution is required.

Theorem 4.6 Let x be a solution to problem P-4.6, but not a solution

to problem P-4.1, and 4 be the set of solutions to min { VZ(x)y : y S }.

Then there exist a 6 > 0 and a yE # such that y solves

min { VZ(x)y : y S I and Vf(x)(y x) < 0 for all x satisfying

lix x '5 6.

Proof: Since VZ(x) is continuous, we have that for all e > 0 there

exists a 61 > 0 such that IIVZ(x) VZ(x) I< e whenever ix x < 61.

Furthermore, the problem min { VZ(x)y : y c S } is a linear program and
there must exist an e > 0 such that there is an element y from the set

I which solves min { VZ(x)y : y E S I whenever IVZ(x) VZ(x) < e ,

or whenever x x I < 61'


g(x) = VZ(x)(y x),

then g(x) is continuous and g(x) < 0, since x is not a solution. Hence,
there exists a 62 > 0 such that g(x) < 0 for all x satisfying

Ix x I 6 To obtain the result, let 6 = min{ 1 6 .
2 1' 2
Therefore, preventing the regeneration of the same extreme

point requires a good approximate solution to the master problem which

can be obtained by performing many projections per iteration. That is,

Step 2 can be replaced by the following steps:


0 t
Step 2a Let w = x k = 0 .

Step 2b Solve
k+1 k w k t+1
Z(wk+ w ) = min { Z(w w) : H(l) }

Step 2c If w does not satisfy a stopping criterion, set k = k + 1
t+l k+l
and go to Step 2b. Otherwise, set x = w and purge

t+l of all grid points with zero weight in the expression of
x as a convex combination of grid points. Set t = t + 1

and go to Step 1. O

As long as w is not optimal to P-4.6, Lemmas 4.2 and 4.3

guarantee that
k+1 k
Z(wk+l) < Z(w ) V k = 1, 2,...,

which further implies that

Z(xt+1) < Z(x .

From this inequality, the results of Lemma 4.4 and Theorem 4.5 follow.

4.1.1 Computational Testing

Since the Frank-Wolfe direction finding problem, i.e., the

problem in Step 1, is generally nontrivial and the projection is rela-

tively easy to perform when A MkA is diagonal, the following stopping

rule of Step 2c is adopted:
%X t+l I t+1
Stop the projection process if G(w I| ) 5 e or

t+l > MAXP, where

SII ut+1 '
(i) G(x ) = VZ(x)(x y)/ VZ(x)y, and y solves

min { VZ(x)y : y e H(k+l) .


(ii) MAXP = maximum allowable number of projections.

The numerator of G(x I t+l) is the gap function [Hearn 1982] at x with
t+l I t+l
respect to the region H(n ), and G(x 0l ) is therefore the gap

function normalized by VZ(x)y. We refer to it as the relative gap

function. Instead of letting et be constant for all t, the sequence t

is defined as follows:

e = min { e ; 6* G(xt S) }
st tmin{ t_l;1

where 0 < 6 < 1 and EO is set a priori. This construction allows the

values of t to be relatively large in the early stages of the algo-

rithm since the main interest at this stage is to generate grid points

and not to solve the master problem. However, in the latter part of

the algorithmic process, the algorithm presumably would have already

obtained the simplex which contains a solution, and a good approximate

solution to the master problem would mean a good approximate solution

to problem P-4.1. So, by making the et small at the latter part of the

algorithm, more effort is put into solving the master problem. The

constant MAXP is only used to prevent an excessive number of projec-

tions during each iteration.

The test problems include three traffic assignment problems,

the symmetric problem (Figure 4-1) from Nguyen and Dupuis [1981], the

nine node problem (Figure 4-2) and the windowpane problem (Figure 4-3)

from Barton and Hearn [1979]. For the experiments the matrix AT MA is

taken to be the diagonal matrix whose diagonal elements all equal to
.75 times the gap value at w .

The program for restricted simplicial decomposition is writ-

ten in FORTRAN with variables r, 6 and MAXP as its parameters. Ini-

tially, a (free) PDP 11/34 computer was used to run the program for


3 4

1 400 800

2 600 200

Trip Table

Figure 4-1: Nguyen & Dupuis Symmetric Problem


S1 10 20

S/ 2 30 40
6 < 8 4
Trip Table

Figure 4-2: Nine Node Problem

8 4
1 2 3 4
1 10 20 10

2 10 20 10

3 10 10

4 10 10 20

Trip Table

Figure 4-3: Windowpane Problem


each problem and each r (the number of grid points) at various values

of the parameters 6 and MAXP in order to determine their "optimal" val-

ues, 6 and MAXP For each pair of 6 and MAXP values, the program was

run until it achieved an objective function value which is as good as

or better than the objective function value at the 100th iteration of

the Frank-Wolfe algorithm (written by Nguyen and James [1975]). For

each run, the following data were recorded: the number of iterations,

the number of projections performed, the CPU time and the final objec-

tive function value. As an example, Table 4-1 summarizes the results

for the nine-node problem with r = 4. The CPU times in this table are

estimated CPU times required to execute the same program on the IBM

3033N computer (the PDP 11/34 does not have a CPU clock). Using these

estimated times as the sole criterion, Table 4-2 lists the best values

of 6 and MAXP for each value of r, the number of retained grid points,

for all three problems.

For the final comparison, we ran the programs using the para-

meters listed in Table 4-2 on the IBM 3033N computer and plotted the

objective function values against the actual CPU times on the IBM

3033N. The results are summarized in Figures 4-4, 4-5 and 4-6.

From the results above it is clear that as r increases the

time required to achieve the objective function at the 100th iteration

of the Frank-Wolfe algorithm decreases. Furthermore, it is encouraging

that the value of r need not be too large to obtain fast convergence.

In all three problems, the largest value of r is 7, and this is quite

reasonable when compared to the Caratheodory bound 20 for the symme-

tric problem, 19 for the nine node problem and 25 for the windowpane



Table 4-1: Computational results for the nine node problem with r = 4.

MAXP 100








obj. fun.
(I iter., I proj.)
Est. CPU time (sec.)

The objective function value at the 100th iteration of Frank-Wolfe algorithm
is 1455.9588 (3.10 sec.) and the best lower bound is 1453.1054.

1455.8660 1455.7278 1455.5480
(13,1025) (11,1091) (11,1252)
.6966 .6710 .7297

1455.8140 1455.7813 1455.8455
(16.1087) (13,1092) (12,1068)
.7937 .7210 .6874--

1455.7771 1455.9149 1455.9149
(20,1124) (19,1047) (19,1047)
.9065 .8536 .8536

1455.8676 1455.9387 1455.9387
(31,1075) (30,1037) (30,1037)
1.1619 1.1232 1.1232

1455.9497 1455.9497 1455.9497
(42,996) (42,996) (42,996)
1.4053 1.4063 1.4063

1455.9374 1455.8856 1455.8856
(45,997) (43,1029) (43,1029)
1.4812 1.4432 1.4432

1455.9219 1455.9219 1455.9219
(67,1031) (67,1031) (67,1031)
2.0400 2.040 2.040


Table 4-2:

"Optimal" values of 6 and MAXP.

Nguyen & Dupuis Nine-node Problem Windowpane Prob.
Symmetric Problem
100th iter. of F-U 100th Iter. of F-W = 100th Iter. of F-W -
85100.07 (3.12 sec.) 1455.91 (3.10 sec.) 1945.99 (3.41 sec.)

85098.5 1455.85 1945.98
(0.1,200) (0.1, 200) (0.1, 100)
(40,3804, 1.94) (29, 1699, 12.9) (47, 3898, 2.77)

85099.68 1455.72 1945.86
(0.1, 100) (0.1, 200) (0.1, 200)
(22, 1963, 1.28) (11, 1091, 0.71) (6 5T,-07;'0.44

85097.5 1455.84 1945.7898
(0.2, 300) (0.1, 100) (0.4, 100)
(6, 223, 0.32) (9, 704, 0.55) (3, 44, 0.16)

85029.9375 1455.44
(0.1, 100) (0.1, 100) same as above
(5, 95, 0.22) (6, 171, 0.26)

same as above (0.5, 100) same as above
(6. 118, 0.24)

objective function
(6, MAXP)
(I iter., # proj., time in sec.)


obj. fun.



85,500 r-2 (Frank & Wolfe)




85,000 r6,7

0.5 1.0 1.5 sec.

Figure 4-4: Nguyen & Dupuis Symmetric Problem


-obj. fun.
1490. -




-2 (Frank & Wolfe)

1.5 sec.

Figure 4-5: Nine Node Problem


obj. fun.



1960- r=2 (Frank & Wolfe)



05 sec.

0.5 1.0 1.5 sec.

Figure 4-6: Windowpane Problem


In summary, the tests show empirically that there is a trade-

off between computer storage used and the progress of the algorithm.

The implication for practical applications is that choosing r > 2 can

offer an improvement over the Frank-Wolfe method (where r = 2).

4.2 Simplicial Decomposition for the Asymmetric Problem

A simple extension of the basic algorithm to the VIP[C,S]

would be to replace the optimization problem in Step 2 with the corres-

ponding VIP. However, as Magnanti [in press] has recently noted with

an example due to J. Hammond, the resulting algorithm cycles. To

prevent cycling, we only-allow column dropping when the algorithm has

made satisfactory progress toward the solution. For the symmetric

case, the objective function was used to monitor the progress of the

algorithm. To establish the convergence of the simplicial technique

for the VIP[C,S], some choice of monitoring function must be made since

no objective function exists. There are several choices: the gap func-

tion (G(F), Section 2.1), its dual function (g(F), Section 2.1), and

the user objective function [Smith 1983]. For simplicial decomposi-

tion, the gap function is the natural choice since it is the byproduct

of solving the subprogram in Step 1. The gap value is simply the quan-

tity C(x )x minus the optimal objective function of the subproblem.

In addition, it is as yet uncertain how the other two functions might

be evaluated computationally.

The algorithm of the VIP requires a small positive tolerance

6 and a convergent monotone sequence {t I where t > e > 0, and
St 0 as t t-
et + 0 as t + co


The SDVI Algorithm

Step 0: Let x0 be an element of S, and define 0 = {x0 DO = 0 and
G = Set t = 0.

Step 1: Solve

C(xt)yt = min { C(xt)y : y e S .

If G(x ) = C(x ) (x y ) = 0, stop and x is the solution.


(i) If G(xt) t -t 6, let tl = t U yt

(ii) If G(x ) < t 6, let nt+l = (st\D ) U y

-t+l t t
Let G = min { G G(x ) .

Step 2: Find xt+1 H(t+1) such that

t+l t+l t+l
'C(x )(x x ) 2 V x H(SI ) .

Let Dt+l be columns of 2t with zero weight in the expression
t t
of x as a convex combination of extreme points in t.

Increment t and go to Step 1. D

The parameters et allow for approximate solving of the master

variational inequality problem in Step 1. In fact, when C(x) is strong-

ly monotone, i.e., there exists an a > 0 such that

(C(x1) C(x2)) (x x2) 2 a 1 x211 V x1, x2 S

and if x solves VIP[H( t),C] exactly, then

a lxt x2 < (C(xt) C(x))(xt- x)

=C( t) (x x) C(x)(xt x)


t t ^.
5 C(x) (xt -x)

Smax C(x) (x w) et
weH( t)

.lxt ^112 c/

Therefore, depending on a, the sequence et can be chosen to
t ^
assure x as near x as desired.

4.2.1 Convergence Proof

For the convergence proof we say that the algorithm completes

one major iteration every time the gap, G(xt), decreases, and we refer

to the iterations between any two major iterations as minor iterations.

Define T to be
t t t -t
T = { t : C(x) (x y ) < G t = 1,2,3, ... } ,

then T is the index set of the major iterations.

It is clear that if the algorithm terminates at the k-th

iteration, x must be the solution of VIP[S,C]. Also k e T since

C(xk (x y) = 0 < min { C(xt (xt yt) : t < k } and C(xt) (xt y

> 0 V t < k. Therefore, in the case of finite termination, the algo-

rithm yields the solution at a major iteration. Otherwise, the algo-

rithm generates a sequence {x }, and we will show that every convergent

subsequence of {x t ^ converges to the solution of VIP[S,C]. It

should be noted that the algorithm does not drop columns at every major

iteration; instead, it only drops columns when there is a sufficient

decrease in the gap, i.e., by an amount 6. As we will show, the para-

meter 6 assures that ultimately no columns will be dropped, and this

leads to the proof of the following convergence theorem.

Theorem 4.7 Assume that the algorithm generates an infinite sequence


{x } and T is an infinite set, then any convergent subsequence of T

converges to the solution of VIP[S,C].

Proof: Since xt S for all t e T where S is compact, there must exist
t ^ *
a convergent subsequence {x } where T T. Let x be the limit
point of the subsequence {x t} Since the gap function, G(x), is a

nonnegative continuous function on S, there is an integer K such that

for all t e T and t > K

IG(xt) G(x*) < 6.

By monotonicity of the sequence {G(xt) tT' we have

t *
(i) G(x) G(x ) < 6 V t T and t > K

(ii) For any tl, t2 T such that t2 > t > K

t t2 t1 *
G(x ) -G(x )

Note that (ii) also implies that for any t such that t e T and tl < t < t2

t t t
(iii) G(x 1) G(x ) < Gx ) G(x 2) < 6.

From (ii), (iii), and the fact that there is no column dropping during

any iteration t T there must be no column dropping after iteration K,
and the sequence of sets n t > K is a monotonically increasing se-

quence. Since S has a finite number of extreme points, the sequence
t ^ *
{ }t>K must converge to some set 0, i.e., there exists an integer K

such that
t ^
0 = V t > K.

For all t,

C(xt) (x x) -t V x e H( t)


ort t
max C(xt) (x w) 5 e
t ^ ^ t
However, since t = Q for t > K, y must be an element of
A t
^ = t where

C(x) (x y ) = max C(x) (x y).

Thus; for all t e K

t t t t t t
G(xt) = C(xt) (x y ) = max C(x) (x w) S .
wEH(Qt )

Because t is a sequence converging to zero,
t *
lim G(x ) = 0 = G(x ) .

Therefore, x must be the solution of VIP[S,C]. Since the

subsequence defined by T is arbitrary, the proof is complete. 0

In the proof above we assume that T is an infinite index set

whenever the algorithm generates an infinite sequence. Now we show

that this assumption is in fact a true statement.

Lemma 4.8: If the algorithm generates an infinite sequence, then T is

an infinite set.

Proof: Assume otherwise; that is, there exists a K such that for all

t > K
t = C(xt t t -k
G(x) = C(x) (x y ) > G = 8 > 0.

But, this implies that there is no column dropping after iteration K

and using a similar argument as in the above proof it can be shown that

G(x ) < E V t > K.
Since converges to zero, we obtain a contradiction.
Since s converges to zero, we obtain a contradiction. O


4.2.2 Discussion of Algorithm

In principle, any convergent method for solving variational

inequalities could be applied to the master problem of Step 2. How-

ever, especially powerful methods, owing to the simplicity of the con-

straint structure, are projection methods. Furthermore, Bertsekas and

Gafni [1982] have recently developed a convergent projection method

which is applicable to our master problem. Specifically, they have

shown that their method is linearly convergent whenever C(x) is Lip-

schitz continuous and strongly monotone, even if there is a linear

transformation of the variables of the form x = AX. In other words,

the linear convergence result will hold if the problem is solved in the

space of X rather than of x. Their proof is important because the as-

sumed properties of C may not carry over to the variational inequality

problem in the X space (see the discussion in Bertsekas and Gafni

[1982]). In all of our experimentation we have solved the master prob-

lem of Step 2 in X space and defined Dt by those columns of ft for

which the associated X variable is zero.

Bertsekas and Gafni have also employed their method for sol-

ving the traffic assignment problem. However, there are important dif-

ferences between the strategy of their implementation and ours. Their

method works in the space of path flows, ours in the space of aggregate

extreme flows (as does the Frank-Wolfe method). Although, there are

generally more extreme aggregate flows than path flows, this does not

necessarily imply that more will be required in the expression of a

given flow vector F. To illustrate, consider a network with two OD

pairs and having three paths pi, p2' p3 joining the first pair and four

paths p 4 p', p6 and p7 joining the second pair, where the pi are of


dimension R. These give rise to 12 (3 times 4) aggregate extreme

flows. Now let F = .5p1 + .5P2 + Op3 + .3p4 + .2p5 + .5p6 + Op7 be a

flow pattern expressed in terms of five path flows. Since F = .3(pl +

p4) + .5(p2 + p6), the same flow pattern is expressible in terms of

just three aggregate extreme flow vectors (pl + 4 ), (p + p5), and

(P2 + P6). Thus working in the aggregate extreme flow space can save

computer storage, a paramount concern on large networks.

Another difference in the algorithms is that Bertsekas and

Gafni only do one projection iteration (in the path flow space) before

generation of new paths (via a subproblem which is equivalent to Step 1

of our algorithm). Our Step 2 does as many projection iterations as

necessary to achieve an et solution. This trade-off takes advantage of

the fact that on large networks the generation of shortest paths in

Step 1 is far more time consuming than the projection iterations.

The computational experiments summarized in the next section

show that even on small networks there are computational savings in

both storage and computer time.

4.2.3 Computational Testing

In this section, computational testing of the algorithm on

four small problems is summarized. In the proof of convergence, we

have assumed that et + 0 and that the stopping rule is G(xt) = 0. For

implementation purposes it is more practical to require an e-solution.

Accordingly, we chose et = 10 + 5 104/t and 6 = 104. Further-

more, the algorithm is terminated when the value of the relative gap,
t t t t t -4
defined to be C(x)(x yt)/C(x )y is less than 10-4. This is a

more appealing choice than G(xt) and is further justified by the re-

sults for problem ND1 (below). In all testing we solve the master


variational inequality problem in Step 2 by the projection method. The

matrix which defines the norm for the projection mapping does not

change from iteration to iteration, and it is taken to be a diagonal

matrix whose diagonal elements have the same positive value.

All programs are written in FORTRAN and were executed on the

IBM 3033N computer at the University of Florida. In order to make an

approximate comparison with the times reported by Nguyen and Dupuis

[1981], we must relate their times on a Cyber 173 machine to the IBM

3033N. The ratio of the machine cycle times for the IBM 3033N and the

Cyber 173 is 58/50 = 1.16. However, when we solved their symmetric

traffic equilibrium by the program TRAFFIC written by Nguyen and James

[1975] on the IBM 3033N, we found that the average ratio of the time on

the IBM 3033N and the Cyber 173 is approximately 3. We cannot make any

definite conclusion from these facts except that IBM 3033N is clearly a

slower machine.

Two algorithms are selected for the comparison, the algo-

rithm by Nguyen and Dupuis [1981] and a variation of the algorithm pre-

sented in Bertsekas and Gafni [1982J. Other techniques, which have as

their subproblem a symmetric traffic assignment problem, require a stand-

ard nonlinear programming algorithm, e.g., the Frank-Wolfe algorithm,

for its solution. Such techniques would construct many shortest path

trees per iteration and therefore are likely to be too time consuming

for a large network. The reader interested in these techniques is re-

ferred to the paper by Fisk and Nguyen [1982].

The algorithm by Nguyen and Dupuis [1981] is an interesting

modification of the algorithm for n-person games by Zuhovicki et al.

[1969]. It employs the cutting-plane algorithm to solve the optimi-

zation version of the variational inequality problem, VIP[C,S], i.e.,


max min C(x)(x u)
ucS xeS

A relaxed version of this max-min problem can be written as a linear

programming problem

max x

subject to

C(x (xi U) 2 x i = 1, 2, ..., n

U C S,

where x i = 1,..., n, are points belonging to S. During each iter-

ation, such a relaxed problem is solved by the simplex method, and the

resulting solution is used to generate another constraint (a cut),

n+l n+l
C(xn+ )(xn u) 2 x0,

to be included in the next relaxed linear programming problem. In
their computational testing x0, the value of the objective function at

iteration t, monitors the progress of the algorithm.

For comparison purposes we have also implemented a variation

of the Bertsekas and Gafni method which we refer to as the commodity

projection method. (As previously described, their method works in the

path flow space and does not involve column dropping.) In the commod-

ity projection method, a flow vector is retained for each "commodity"

which is defined to be the aggregate flow from an origin to all desti-

nations. This strategy only requires a projection for each origin

at each iteration, rather than for each origin-destination (OD) pair.

We only report here the results for the "all-at-once" version proposed

by Bertsekas and Gafni because it performed better than the "one-at-a-

time" version (see the discussion in Bertsekas and Gafni [1982]) on our


test problems. In implementing the commodity projection method the

choice of the projection matrix is as described before.

The first 3 problems, ND1, ND2, and ND3, are taken from

Nguyen and Dupuis [1981]; all of them have the same network topology

and trip matrix as shown in Figure 4-1. This network admits 8 paths

from node 1 to node 3, 6 paths from 1 to 4, 4 paths from 2 to 3 and 6

paths from 2 to 4. Thus, there are a total of = 768 possible

extreme flow patterns. The cost vectors for ND1, ND2 and ND3 corres-

pond to the symmetric, diagonally dominant and nondiagonally dominant

cost vectors in Nguyen and Dupuis [1981].

In Problem ND1 the Jacobian C(x) is diagonal; therefore an

equivalent optimization problem exists. Table 4-3 summarizes the com-

putation results for the simplicial algorithm as compared with the re-

sults reported by Nguyen and Dupuis. Judging from the objective val-

ues, the two methods produce equivalent solutions, with the-simplicial

algorithm time offering an improvement even though the IBM 3033N is the

slower machine. Of course this comparison is tentative because the

computers differ and because both algorithms can be improved by "fine

tuning." In particular, the Nguyen and Dupuis method could be modified

to solve the linear program approximately at each iteration and this

might improve its overall performance.

Note, also, that the relative gap provides a more accurate

measure of the solution than the gap itself.

Tables 4-4 and 4-5 summarize the results for Problem ND2 and

ND3, respectively. Note that out of 768 extreme flow patterns the

simplicial algorithm generates no more than 6 extreme flow patterns and

retains at most 4 patterns during any iteration. For problems ND2, the

Table 4-3: Computational Results for Problem ND1

Nguyen and Oupuis The Simplicial Algorithm'

t F(xt) (CPUmele.) t F(xt) Gap Rel. Gap p p (CTimec.)
1 0 1 CPUsec') 11 2 (CPU sec.)










































.1004 10"
.2528 105

.4844 104

.5183 104

.1835 104

.1849 104

.6565 103

.6744 102

.6750 101

.1136 101

.2823 100

.4951 10-1

.5343 10-1

.1855 10'1

.1873 10-1

.6558 10-2

.6711 10-3

.6713 10-4










Notes: 1: Time
2: Pl
3: p2

on Cyber 173
I of active extreme aggregate
# of projections

flow patterns

4: Time on IBM 3033
5: 150 is the maximum number of projections allowed per iteration

Table 4-4: Computational Results for Problem ND2

Nguyen and Oupuls The Simplicial Algorithm The Commodity Projection Method

t o Tie t Re. ap. P P Time2 t Rel. Gap. p5 Time2
0 Time t Rel. Gap. p1 p2 Time2 t IR. Gap. p3
.27 in1 .02












Notes: 1: Time on Cyber 173
2: Time on IBM 3033
3: pl a I of active extreme aggregate flow patterns

4: p2 -
5: p3 I

.2677 101

.2035 10-1

.5091 10'2

.1449 10.2

.4229 10-3

.1241 10-3

.9681 10-4

.9681 .7111

of projections
of active extreme commodity

flows for commodities 1 and 2







.2677 101

.3446 100

.1169 100

.1852 10-1

.6524 10-4










Table 4-5: Computational Results for Problem ND3

Nguyen and Dupuis The Simplicial Algorithm The Commodity Projection Method

t x Time1 t Rel. Gap. p p2 Time2 t Rel. Gap. p Time2
0 P3 t 3
















.1941 11n

.6743 10-1

.6768 10-2

.2409 10-2

.9164 10-3

.3479 10-3

.1326 10"3

.9026 10-4









on Cyber 173
on IBM 3033
I of active extreme aggregate flow patterns

: p I # of projections
: P3 of active extreme commodity
flows for commodities 1 and 2

Notes: 1: Time
2: Time
3: p, "


simplicial algorithm built 6 minimum path trees and performed 109 pro-

jections, whereas the commodity projection method built 26 minimum path

trees and performed 26 projections. Although the simplicial algorithm

performed many more projections, the CPU time indicates that the time

required to build the additional 20 minimum path trees more than off-

sets the time to do the extra projections. A similar conclusion can be

drawn from Table 4-5. The computation times of the simplicial algo-

rithm are the best on both ND2 and ND3, but, for the reasons already

given, comparison with the Nguyen and Dupuis method must be regarded as


The network topology and the trip matrix for Problem 4 are

shown in Figure 4-7. The network is from Fisk and Nguyen [1982], but

we constructed the cost function shown in Table 4-6. This network ad-

mits 4 paths between each OD pair with nonzero demand; in other words,

there are (4)12 = 16,777,216 possible extreme flow patterns. Table 4-7

summarizes the computational results of the simplicial algorithm versus

the commodity projection method. The conclusions from Problem ND2 and

ND3 also hold for this case. Considering the very large number of ex-

treme flow patterns it is encouraging that the simplicial algorithm

generated only 18 extreme flow patterns.

We feel the above results demonstrate the potential of the

simplicial algorithm with regard to computer time, and they show fur-

ther that the method obtains good solutions while generating only a

small number of extreme flow patterns.


Network Topology

5 6 7 8

1 0 35 30 20

2 40 0 35 25

3 30 20 0 35

4 35 25 40 0

Trip Matrix

Figure 4-7: Data for Problem 4

Table 4-6: Cost Vector for Problem 4

Arc 0 0 Ci(x) Arc 0 0 C (x)

1 1 9 f(x1 + 14) 19 14 11 f(1/2 x16 + 19 + x25 + 1/4 x26)
2 1 12 f(x2 + 1/2 x11 + 1/4 x12 + 1/4 x14) 20 14 18 f(x16 + x20)
3 2 10 f(1/4 x10 + x3 + x21 + 1/4 x22) 21 15 10 f(x3 x21)
4 2 16 f(xo0+ x4) 22 15 16 f(1/4 x10 + 1/4 x3 + x2)
5 3 19 f(x34 + 1/4 x33 + 1/4 x32 + x5) 23 16 14 f(x29 + 1/4 x30 + 1/4 x18 + x17 + x23)
6 3 22 f(x32 + x6) 24 16 20 f(x18 + x24)
7 4 17 f(x7 + x36) 25 17 11 f(x19 + x25)
8 4 21 f(x8 + 1/2 x27 + 1/4 x28 + 1/4 x36) 26 17 13 f(x15 + 1/4 x16 + 1/4 x19 + x20 + x26)
9 9 6 f(xg + x22) 27 18 8 f(1/4 x8 + x27 + x35 + 1/4 x36)
10 9 16 f(x10 + 1/2 x3 + 1/4 x4 + 1/4 x22) 28 18 21 f(xg + x28)
11 10 5 f(1/4 x2 + x11 + x13 + 1/4 x14) 29 19 14 f(x29 + x23)
12 10 12 f(x2 + x12) 30 19 15 f(x30 + 1/4 x18 + 1/2 x17 + 1/4 x23)
13 11 5 f(x11 + x13) 31 20 7 f(x33 + x31)
14 11 9 f(xI + 1/4 x2 + 1/4 xll + x14) 32 20 22 f(1/4 x33 + x32 + x6 + 1/4 x5)
15 12 13 f(x15 + x26) 33 21 7 f(x33 + 1/2 x32 + 1/4 x31 +1/4 xS)
16 12 18 f(x16 + 1/4 x19 + 1/2 x20 + 1/4 x26) 34 21 19 f(x34 + x5)
17 13 15 f(x30 + x17) 35 22 8 f(x27 + x35)
18 13 20 f(1/2 x30 + x18 + x24 + 1/4 x23) 36 22 17 f(x7 + 1/4 x8 + 1/4 x27 + x36)

f(Z) 3.0 (1 + 0.15(Z/35)4)

Table 4-7: Computational Results for Problem 4

The Simpllcial Algorithm The Commodity Projection Method

t Rel. Gap. p2 p2 Time t Rel. Gap. p3 Time1'




Notes: 1: Time on IBM 3033
2: pl I of active extreme aggregate flow patterns
3: P2 I of projections




4: p3 I of extreme commodity flows
5: 300 max. allowable number
of projections





5.1 Summary of Results

This study is concerned with decomposition techniques for the

variational inequality formulation of the traffic assignment problem.

This formulation includes the standard (optimization) formulation which

has been investigated previously by many authors. In this disserta-

tion, it is shown that two decomposition techniques, transfer and sim-

plicial decompositions, originally developed for the optimization

formulation, can be extended to the variational inequality formulation.

Chapter 3 discusses the extension of transfer decomposition

to variational inequality problems. First, the modified version of

transfer decomposition for the variational inequality problem is pre-

sented, and it is shown that transfer decomposition provides a tech-

nique to separate the symmetric component from the rest of the problem.

Then, several properties for the mapping which defines the master vari-

ational inequality problem are stated and proved. The main property is

that the mapping resulting from applying transfer decomposition to the

traffic assignment problem is always point-to-point. Using this prop-

erty, an existing algorithm can be applied to the master variational

inequality. In the last section of this chapter, it is shown that

transfer decomposition also provides a rigorous basis for constructing

aggregate functions used in network simplification.


University of Florida Home Page
© 2004 - 2010 University of Florida George A. Smathers Libraries.
All rights reserved.

Acceptable Use, Copyright, and Disclaimer Statement
Last updated October 10, 2010 - - mvs