Mathematical Programming 71 (1995) 153-177
An efficient cost scaling algorithm for the
assignment problem
Andrew V. Goldberg
,,1,
Robert Kennedy
2
Computer Science Department, Stanford University, Stanford, CA 94305-2140, United States
Received 18 October 1993; revised manuscript received 19 January 1995
Abstract
The cost scaling push-relabel method has been shown to be efficient for solving minimum-cost
flow problems. In this paper we apply the method to the assignment problem and investigate
implementations of the method that take advantage of assignment's special structure. The results
show that the method is very promising for practical use.
Keywords: Network optimization; Assignment problem; Algorithms; Experimental evaluation; Cost scaling
1. Introduction
Significant progress has been made in the last decade on the theory of algorithms for
network flow problems. Some of the algorithms that came out of this research have been
shown to have practical impact as well. In particular, the push-relabel method [ 11, 16]
is the best currently known way for solving the maximum flow problem [2, 8, 23]. This
method extends to the minimum-cost flow problem using cost scaling [ 11, 17]. Earlier
implementations of this method [5, 14] performed well on some problems but relatively
poorly on others. A later implementation [ 12] has been shown very competitive on a
wide class of problems. In this paper we study efficient implementations of the cost
scaling push-relabel method for the assignment problem.
* Con'esponding author. Present address: NEC Research Institute, 4 Independence Way, Princeton, NJ 08540,
United States, e-mail: avg@research.nj.nec.com.
1 This author's research was supported in part by ONR Young Investigator Award N00014-91-J-1855, NSF
Presidential Young Investigator Grant CCR-8858097 with matching funds from AT&T, DEC and 3M, and a
grant from the Powell Foundation.
2 This author's research was supported by the above-mentioned ONR and NSF grants.
0025-5610 (~) 1995--The Mathematical Programming Society, Inc. All rights reserved
SSD1 0025-5610(95)00008-9
154
A. V. Goldberg, R. Kennedy~Mathematical Programming 71 (1995) 153-177
The most relevant theoretical results on the assignment problem are as follows. The
best currently known strongly polynomial time bound of
O(n(m + n
log n)) is achieved
by the classical Hungarian method of Kuhn [21]. Here n denotes the number of nodes
in the input network and rn denotes the number of edges. (Implementations of the
Hungarian method and related algorithms are discussed in [7].) Under the assumption
that the input costs are integers in the range [-C ..... C], Gabow and Tarjan [10]
use cost scaling and blocking flow techniques to obtain an
O(v/nmlog(nC))
time
algorithm. Algorithms with the same running time bound based on the push-relabel
method appeared in [ 15, 24].
In this paper we study implementations of the scaling push-relabel method in the
context of the assignment problem. We use the ideas behind the minimum-cost flow
codes [5, 12, 14], the assignment codes [3,4,6], and the ideas of theoretical work on
the push-relabel method for the assignment problem [ 15], as well as new techniques
aimed at improving practical performance of the method. We develop several CSA
(Cost Scaling Assignment) codes based on different heuristics which improve the code
performance on many problem classes. The "basic" code CSA-B does not use any
heuristics, the CSA-Q code uses a "quick-minima" heuristic, and the CSA-S code
uses a "speculative arc fixing" heuristic. Another outcome of our research is a better
understanding of cost scaling algorithm implementations, which may lead in turn to
improved cost scaling codes for the minimum-cost flow problem.
We compare the performance of the CSA codes to four other codes: SFR10, an
implementation of the auction method for the assignment problem [6]; SJV and DJV,
implementations of Jonker and Volgenant's shortest augmenting path method [ 19] tuned
for sparse and dense graphs respectively; and ADP/A, an implementation of the interior-
point method specialized for the assignment problem [25]. We make the comparison
over classes of problems from generators developed for the First DIMACS Implemen-
tation Challenge [18] 3 and on problems obtained from digital images as suggested
by Knuth [20]. Of our codes, CSA-Q is best overall. This code outperforms ADP/A
on all problem instances in our tests, outperforms SFR10 on all except one class, and
outperforms SJV and DJV on large instances in every class. Although our second-best
code, CSA-S, is somewhat slower than CSA-Q on many problem classes, it is usually
not much slower and it outperforms CSA-Q on three problem classes, always outper-
forms ADP/A, is worse than SFR10 by only a slight margin on one problem class and
by a noticeable margin on only one problem class, and loses to the Jonker-Volgenant
codes only on one class and on small instances from two other classes. While we use
the CSA-B code primarily to gauge the effect of heuristics on performance, it is worth
noting that it outperforms ADP/A in all our tests, the Jonker-Volgenant codes on large
instances from all but one class, and SFR10 on all but one class of problems we tested.
This paper is organized as follows. Section 2 gives the relevant definitions. Sec-
tion 3 outlines the scaling push-relabel method for the assignment problem. Section 4
3 The DIMACS benchmark codes, problem generators, and other information we refer to are available by
anonymous ftp from diraacs, rutgers, odu.
A.V. Goldberg, R. Kennedy/Matheraatical Programming 71 (1995) 153-177
155
discusses our implementation, in particular the techniques used to improve our code's
practical performance. Section 5 describes the experimental setup. Section 6 gives the
experimental results. In Section 7, we give concluding remarks.
2. Background
Let G = (V = X U Y,E) be an undirected bipartite graph and let n = I V] and m = IEI .
A matching
in G is a subset of edges M C ~ that have no node in common. The
cardinality
of the matching is IMI . Given a weight function ~ : E -+ R, we define the
weight of M to be the sum of weights of edges in M. The
assignment problem
is to
find a maximum cardinality matching of maximum weight. We assume that the weights
are integers in the range [-C ..... C]. To simplify the presentation, we assume that
IXI = ]YI, G has a perfect matching (i.e., a matching of cardinality IX1), and every node
degree in G is at least two. We can dispense with these last assumptions without any
significant decrease in performance by using a slightly more complicated reduction to
the transportation problem than the one described below.
Our implementation reduces the assignment problem to the
transportation problem
defined as follows. Let G = (V,E) be a digraph with a real-valued
capacity u(a)
and
a real-valued
cost c(a)
associated with each arc 4 a and a real-valued
supply d(v)
associated with each node v. We assume that
~v d(v) = O. A pseudoflow
is a function
f • E -+ R+ satisfying the
capacity
constraints: for each
a G E, f(a) <~ u(a).
For
a pseudoflow f and a node v, the
excess flow into
v, ef(v), is defined by el(v) =
d(u) + ~(~,,~,)cE
f(u,
v) -
~-~4,,,w)cE
f(v, w). A node v with e/(v) > 0 is called
active.
Note that
~vcv ef(v) = O.
A flow
is a pseudoflow f such that, for each node
v, ef(v)
= 0. Observe that a
pseudoflow f is a flow if and only if there are no active nodes. The
cost
of a pseudoflow
f is given by cost(f)
= ~a~E c(a)f(a).
The transportation problem is to find a flow
of minimum cost.
We use a slight variation of the standard reduction from the assignment problem to
the minimum-cost flow problem (see, e.g., [22] ). Given an instance of the assignment
problem (G,~), we construct a transportation problem instance (G ---
(V,E),c,u)
as
follows. We define V = V = X U Y. For every edge {v, w} E E such that v C X and
w E Y, we add the arc
(v,w)
to E and define
c(v,w) = --d(v,w)
and
u(v,w) = 1.
Finally we define
d(v)
= 1 for all v E X and
d(w)
= -1 for all w E Y. Note that the
graph G is bipartite.
For a given pseudoflow f, the
residual capacity
of an arc a ~ E is uf(a) = u(a) --
f(a).
The set of
residual arcs Ef
contains the arcs a C E with
f(a) < u(a)
and the
reverse arcs a R, for every a E E with
f(a)
> 0. The
residual graph Gf = (V, EI)
is the
graph induced by the residual arcs. For a E E, we define
c(a R) = -c(a).
Note that if
4 Sometimes we refer to an arc a by its end points, e.g., (v, w). This is ambiguous if there are multiple arcs
fiom v to w. An alternative is to refer to v as the tail of a and to w as the head of a, which is precise but
inconvenient.
156
A.V. Goldberg, R. Kennedy~Mathematical Programming 71 (1995) 153-177
procedure
MIN-COST( V,E, u, c) ;
[initialization]
e ~-- C; Vv, p(v) *'-
O;
[loop]
while
e /> 1/n do
(e, f,p) *-
REFINE(E,p);
return (f);
end.
Fig. 1. The cost scaling algorithm.
G is obtained by the above reduction, then for any integral pseudoflow f and for any
arc
a ~ E, u(a),f(a) E
{0, 1}.
A price funetion
is a function p : V ---+ R. For a given price function p, the
reduced
cost
of an arc (v, w) is
C p ( V , w) = c( v , w) + p ( v ) - p ( w )
and the
partial reduced cost
is
c~,(v,w) = c(v,w) - p(w).
For a constant • ~> 0, a pseudoflow f is said to be
•-optimal with respect to a price
function p
if, for every residual
arc a C Ef,
we have
a E E ~ cp(a) >/ 0, aR E E ~ Cp(a) >/--•.
A pseudoflow f is
•-optimal
if f is •-optimal with respect to some price function p. If
the arc costs and capacities are integers and • <
l/n,
any •-optimal flow is optimal.
For a given f and p, an arc a E
Ef
is
admissible
iff
aEEandep(a) < ½• or a REEandcp(a)
<-½e.
The
admissible graph
GA = (V, EA)
is the graph induced by the admissible arcs.
3. The method
First we give a high-level description of the successive approximation algorithm (see
Fig. 1 ). For a detailed presentation of the successive approximation framework and the
associated proofs, see [17]. The algorithm starts with • = C and
p(v)
= 0 for all
v C V. At the beginning of every iteration, the algorithm divides • by a constant factor
o~ and sets f to the zero pseudoflow. The iteration modifies f and p so that f is an
(•/a)-optimal flow with respect to p. When • <
I/n, f
is optimal and the algorithm
terminates. The number of iterations of the algorithm is 1 + Llog,(nC)J.
Reducing • is the task of the subroutine
refine.
The input to
refine
is • and p such
that there exists a flow f that is e-optimal with respect to p. The output from
refine
is
• ~ = •/ce, a flow f, and a price function p such that f is •t-optimal with respect to p.
The generic
refine
subroutine (described in Fig. 2) begins by decreasing the value
of e, setting f to the zero pseudoflow (thus creating some excesses and making some
nodes active), and setting p(v) = -min(v,w)EE{Cp (v, w)} for every v C X. This converts
the f into an •-optimal pseudoflow (indeed, into a 0-optimal pseudoflow). Then the
subroutine converts f into an •-optimal flow by applying a sequence
of push
and
relabel
A.V. Goldberg, R. Kennedy~Mathematical Programming 71 (1995) 153-177 157
procedure
REFINE( 6, p ) ;
[initialization]
e 6--- e/OL';
Va E E, f(a) +-- 0;
Vv E X, p(v) ~-- - min(v,w)~E c~p(v, w);
[loop]
while
f is not a flow
apply a push or a relabel operation;
return(e, f, p) ;
end.
Fig. 2. The generic refine subroutine.
operations, each of which preserves e-optimality. The generic algorithm does not specify
the order in which these operations are applied. Next, we describe the push and relabel
operations for the unit-capacity case (see Fig. 3).
A push operation applies to an admissible arc (v, w) whose tail node v is active. It
consists of pushing one unit of flow from v to w, thereby decreasing ef(u) by one,
increasing ef(w), and either increasing f(v,w) by one if (v,w) E E or decreasing
f(w,u) by one if (w,v) E E. A relabel operation applies to a node v. The opera-
tion sets p(u) to the smallest value allowed by the e-optimality constraints, namely
maxf.,w)EE:{p(w) --c(v,w)} if v E X, or max(..w)EEs{P(W) --c(v,w) --e} otherwise.
The analysis of cost scaling push-relabel algorithms is based on the following facts
[ 15, 17]. During a scaling iteration,
• the node prices monotonically decrease;
• for any v E V, p(v) decreases by O(ne).
4. Implementation and heuristics
In this section we discuss implementation issues and heuristics aimed at speeding up
the method.
The efficiency of a scaling implementation depends on the choice of scale factor ce.
Although an earlier study [6] suggests that the performance of scaling codes for the
assignment problem may be quite sensitive to the choice of scale factor, our observations
are to the contrary. Spot checks on instances from several problem classes indicated that
PUSH(U,W),
send a unit of flow from v to w.
end.
RELAI~EL(u).
if v E
X
then
replace p ( v ) by max(.,w)
E E.: {P ( w ) -- c ( v, w ) }
else
replace
p(v) by max(u,v)EF~: {p(u) + C(U,V) -- e}
end.
Fig. 3. The push and relabel operations.