A STRONGLY POLYNOMIAL CUT CANCELING ALGORITHM FOR MINIMUM COST SUBMODULAR FLOW ∗ SATORU IWATA
†,
S. THOMAS MCCORMICK
‡ , AND
MAIKO SHIGENO
§
Abstract. This paper presents a new strongly polynomial cut canceling algorithm for minimum cost submodular flow. The algorithm is a generalization of our similar cut canceling algorithm for ordinary min-cost flow. The algorithm scales a relaxed optimality parameter, and creates a second, inner relaxation that is a kind of submodular max flow problem. The outer relaxation uses a novel technique for relaxing the submodular constraints that allows our previous proof techniques to work. The algorithm uses the min cuts from the max flow subproblem as the relaxed most positive cuts it chooses to cancel. We show that this algorithm needs to cancel only O(n3 ) cuts per scaling phase, where n is the number of nodes. Furthermore, we also show how to slightly modify this algorithm to get a strongly polynomial running time. Finally, we briefly show how to extend this algorithm to the separable convex cost case, and that the same technique can be used to construct a polynomial time maximum mean cut canceling algorithm for submodular flow.
1. Introduction. A fundamental problem in combinatorial optimization is mincost network flow (MCF). It can be modeled as a linear program with guaranteed integral optimal solutions (with integral data), and many polynomial and strongly polynomial algorithms for it exist (see Ahuja, Magnanti, and Orlin [1] for background). Researchers in mathematical programming have developed a series of extensions of MCF having integral optimal solutions (see Schrijver [39]). Often the proofs of integrality are existential rather than algorithmic. We have long been interested in finding a generic (strongly) polynomial algorithm for such problems, i.e., one that is easily extended from the MCF case to more general cases. Finding such an algorithm would allow us to better understand what features of MCF algorithms depend on special structure, and which are more general. We hope that it would also allow researchers to be able to find algorithms for more general problems more quickly than in the past. A natural class of generic algorithm to consider is the class of (primal) cycle canceling algorithms, and (dual) cut canceling algorithms. These are natural in the sense that they take improving steps coming from the linear algebra of the constraints, and their sense of “improvement” depends only on the local effect on the objective value. It is (in theory) easy to figure out what these objects should look like for more general models than MCF. A natural first step in applying such a research program is to consider the case of submodular flow (SF), formally defined in Section 2. This problem very much resembles ordinary MCF, except that the usual conservation constraints have been relaxed into constraints on the total violation of conservation on any node subset. SF was shown to enjoy integral optimal solutions by Edmonds and Giles [4]. Early algorithmic contributions came from [3, 9, 10, 13]. Our algorithm uses many ideas from these papers. The first strongly polynomial algorithm for SF is due to Frank and Tardos [11], who generalized the strongly polynomial MCF algorithm by Tardos [42] ∗ An
extended abstract of a preliminary version of this paper appeared in [26]. of Mathematical Informatics, Graduate School of Information Science and Technology, University of Tokyo, Tokyo 113-8656, Japan (
[email protected]). ‡ The Sauder School of Business, University of British Columbia, Vancouver, BC V6T 1Z2 Canada. Research supported by an NSERC Operating Grant; part of this research was done during visits to RIMS at Kyoto University and SORIE at Cornell University (
[email protected]). § Institute of Policy and Planning Sciences, University of Tsukuba, Tsukuba, Ibaraki 305-8573, Japan (
[email protected]). † Department
1
to a fairly wide class of combinatorial optimization problems. A more direct generalization of the Tardos algorithm was presented by Fujishige, R¨ock, and Zimmermann [18] with the aid of the tree-projection method by Fujishige [15]. Unfortunately it has been surprisingly difficult to extend known MCF cycle and cut canceling algorithms to SF. One of the most attractive MCF canceling algorithms is the min mean cycle canceling algorithm of Goldberg and Tarjan [20] (and its dual, the max mean cut canceling algorithm of Ervolina and McCormick [5]). Cui and Fujishige [2] were able to show finiteness of this algorithm for SF, but no polynomial bound. McCormick, Ervolina, and Zhou [35] show that it is unlikely that a straightforward generalization of either min mean cycle or max mean cut canceling is polynomial for SF. Recently, Wallacher and Zimmermann [43] devised a weakly polynomial cycle canceling algorithm for SF based on the min ratio approach by Zimmermann [44]. This algorithm is provably not strongly polynomial, even for networks [36]. Another possibility, developed and analyzed in [41], is to cancel cycles or cuts which maximize the violation of complementary slackness with respect to a relaxed optimality parameter, the relaxed min/max canceling algorithms. These algorithms allow for a simpler computation of the cycle or cut to cancel, but otherwise enjoy the relatively simple analysis of min mean cycle/max mean cut canceling. The same authors found a way to generalize relaxed min cycle canceling to SF [28], which led to nearly the fastest weakly polynomial running time for SF. This algorithm leads to the same strongly polynomial bound as the fastest known one by Fujishige, R¨ock, and Zimmermann [18], and it can also be extended to SF with separable convex costs. The same ideas can be further extended to solve separable convex cost problems over any totally unimodular system [29]. However, in some circumstances cut canceling appears to be more generalizable than cycle canceling. One recent example of this is the problem of submodular function minimization (SFM; see [32] for a survey). One of the two basic SFM algorithms, by Iwata, Fleischer, and Fujishige [25], is a direct descendent of the cut canceling algorithm in the present paper, whereas no cycle canceling SFM algorithm is known. Hence it is important to also generalize cut canceling MCF algorithms to SF. This paper extends the relaxed most positive cut canceling algorithm for MCF of [41] in a non-trivial way to SF. Section 3 develops the concepts of dual approximate optimality we will need, and Section 4 develops the main weakly polynomial algorithm. The algorithm runs in O(n6 h log(nU )) time, where n is the number of vertices, h is the time for computing an exchange capacity, and U is the maximum absolute value of capacities. The next three sections generalize this algorithm in three directions: to a strongly polynomial algorithm (Section 5), to the case of separable convex costs (Section 6), and to a max mean cut variant (Section 7). The total running time bound of the strongly polynomial algorithm is O(n8 h log n). Finally, in Section 8, we describe how to apply our algorithms to the submodular flow problem with crossing submodular functions. The only other polynomial cut canceling algorithm we know for SF is the most helpful total cut canceling algorithm of [33]. However that algorithm is provably not strongly polynomial [6], and it needs an oracle to compute exchange capacities for a derived submodular function which appears to be difficult to derive from an oracle for the original submodular function. By contrast, the present algorithm can compute exchange capacities in derived networks using only an oracle for the original function, and is the first dual strongly polynomial algorithm for SF that we know of. Our algorithm also appears to be the first strongly polynomial algorithm (primal or dual) 2
for SF that avoids scaling or rounding the data. The original version of this paper [26] developed the first strongly polynomial cut canceling algorithm for SF. Subsequently, Fleischer, Iwata, and McCormick (FIM) [7] derived a new strongly polynomial SF algorithm, which improves the running time of the capacity scaling algorithm for SF in [24] by incorporating the cut canceling technique developed in the present paper. The FIM algorithm improves the cut canceling SF algorithm of this paper by updating dual prices using shortest path distances coming from a modified form of Dijkstra’s algorithm to effectively cancel many cuts at once. The FIM algorithm improves both the weakly and strongly polynomial bounds of this paper by a factor of n2 . The strongly polynomial bound matches the current best known bound, and the weakly polynomial bound is better than any previously known algorithm for the submodular flow problem when costs are large but capacities are not (again showing that cut canceling algorithms are sometimes superior). The FIM SF algorithm in turn then led to the IFF SFM algorithm [25], which solved one of the most important and long-standing open problems in combinatorial optimization. This development validates our claim that this research program can lead to quicker algorithm development, as the IFF algorithm followed the FIM algorithm relatively quickly, despite SFM seeming to be quite different from SF. 2. Submodular flow. An instance of submodular flow looks much like an instance of MCF. We are given a directed graph G = (N, A) with node set N of cardinality n and arc set A of cardinality m. We are also given lower and upper bounds `, u ∈ RA on the arcs, and costs c ∈ RA on the arcs. We need some notation to talk about relaxed conservation. If w ∈ RX and Y ⊆ X, P then we abbreviate y∈Y wy by w(Y ) as usual. For node subset S, define ∆+ S as {i → j ∈ A | i ∈ S, j ∈ / S}, ∆− S as {i → j ∈ A | i ∈ / S, j ∈ S}, and ∆S = ∆+ S ∪∆− S. We say that arc a crosses S if a ∈ ∆S. If ϕ is a flow on the arcs (i.e., ϕ ∈ RA ), for i ∈ N define the boundary ∂ϕ(i) of ϕ at i to be ϕ(∆+ {i}) − ϕ(∆− {i}), i.e., the net flow out of i. ThusP ∂ϕ is a vector in RN . We extend ∂ϕ to node subsets S ⊆ N as usual via ∂ϕ(S) = i∈S ∂ϕ(i) = ϕ(∆+ S) − ϕ(∆− S), which equals the net flow out of S. Later we will consider several auxiliary networks whose arc sets are supersets of A, so we will often subscript ∂ and ∆ by the set of arcs we want them to include at that point. Usual MCF conservation requires that ∂ϕ(i) = 0 for all i ∈ N , which implies that ∂ϕ(S) = 0 for all S ⊆ N . From this perspective it is natural to relax this constraint to ∂ϕ(S) ≤ f (S) for some set function f : D → R on some set family D ⊆ 2N with ∅, N ∈ D. Since ∂ϕ(∅) = ∂ϕ(N ) = 0 for any flow ϕ, we can and will assume that f (∅) = f (N ) = 0. For such a set function f , the base polyhedron B(f ) is defined by B(f ) = {x | x ∈ RN , x(N ) = f (N ) = 0, ∀S ∈ D : x(S) ≤ f (S)}. A vector in B(f ) is called a base. This is defined so that ∂ϕ(S) ≤ f (S) for all S ∈ D is equivalent to requiring that ∂ϕ ∈ B(f ). Then the (primal) optimization problem we would like to solve is: X Minimize ca ϕa a∈A
(1)
subject to `a ≤ ϕa ≤ ua ∂ϕ ∈ B(f ). 3
(a ∈ A)
This is a linear program with an exponential number of constraints. Note that for general D and f , B(f ) could be empty, and/or (1) could have fractional optimal solutions. In order to ensure that B(f ) 6= ∅ and that we have integral optimal solutions, it is necessary to require some structure on D and f . For S, T ⊆ N , we say that S and T intersect if S ∩ T 6= ∅, and cross if they intersect and also have S ∪ T 6= N . We say that D is a ring (resp. intersecting, crossing) family if S ∩ T, S ∪ T are also in D for all (resp. intersecting, crossing) pairs S, T ∈ D (note that a ring family is a distributive lattice). We say that f is ring (resp. intersecting, crossing) submodular on D if D is a ring (resp. intersecting, crossing) family and for all (resp. intersecting, crossing) pairs S, T ∈ D we have (2)
f (S) + f (T ) ≥ f (S ∪ T ) + f (S ∩ T ).
Note that every ring family is also an intersecting and crossing family, and so every ring submodular f is also intersecting and crossing submodular. Edmonds and Giles [4] originally formulated the submodular flow problem for crossing submodular functions, which can have empty base polyhedra. Fujishige [14] provided a necessary and sufficient condition for B(f ) to be nonempty. This condition can be checked efficiently by the bi-truncation algorithm of Frank and Tardos [12]. When f is crossing submodular on D and B(f ) 6= ∅, then (1) is a submodular flow problem. Edmonds and Giles showed that submodular flow problems are totally dual integral (TDI), and so always have integral optimal solutions with integer data. Furthermore, Fujishige [14] showed that a nonempty base polyhedron B(f ) of a crossing submodular function f on a crossing family D is identical to a base polyhedron e that contains D, which B(fe) of a ring submodular function fe on a ring family D was already implicit in Frank [9]. Therefore, theorems and algorithms concerning geometric properties of the base polyhedra of ring submodular functions carry over to crossing submodular functions. The cut canceling algorithms presented in this paper rely only on geometric properties such as the exchange capacity (defined below). This allows us to assume throughout most of the paper that D is a ring family and that f is ring submodular on D. We sketch out the minor modifications needed to our algorithms to make them work for the crossing submodular case in Section 8. The base polyhedron of a ring submodular function is always nonempty. Technically, the dual problem to (1) should have an exponential number of dual variables, one for each subset of N . However it is possible and much more convenient to simplify these to dual variables on nodes, or node potentials π ∈ RN . For arc a = i → j, we define the reduced cost on a by cπa = ca + πi − πj . Node potentials are widely used in MCF algorithms, where they arise naturally as LP dual variables. The first use of node potentials as a simpler substitute for the dual variables for the set constraints in submodular problems that we are aware of is made by Iri and Tomizawa [23] for the independent assignment problem. The weight splitting in the weighted matroid intersection algorithm of Frank [8] plays essentially the same role as the node potentials. These two equivalent matroid optimization problems are special cases of the 0-1 submodular flow problem, for which Frank [9] devised an efficient combinatorial algorithm using node potentials. We need some further concepts to characterize optimality for (1). If ϕ is a flow with boundary x = ∂ϕ such that x ∈ B(f ), then we say that subset S is ϕ-tight, or x-tight, or just tight if x(S) (= ∂ϕ(S)) = f (S). Ring submodularity implies that the union and intersection of tight sets are also tight. 4
If π is a set of node potentials with distinct values v1 > v2 > . . . > vh , then the kth level set of π is Lπk ≡ {i ∈ N | πi ≥ vk }, k = 1, 2, . . . , h. It will be convenient to let Lπ0 = ∅. Note that ∅ = Lπ0 ⊂ Lπ1 ⊂ Lπ2 ⊂ . . . ⊂ Lπh = N . Suppose Lπk ∈ D for k = 0, 1, . . . , h and define a function f π : D → R by f π (S) =
h X
{f ((S ∩ Lπk ) ∪ Lπk−1 ) − f (Lπk−1 )}
(S ∈ D).
k=1
Then f π is ring submodular, f π ≤ f , and if S nests with each Lπk (either S ⊆ Lπk or Lπk ⊆ S), then f π (S) = f (S). Also, x is in B(f π ) if and only if x ∈ B(f ) and every Lπk is x-tight for f , which is true if and only if x maximizes π T y over y ∈ B(f ), see [16, Theorem 3.15]. Cunningham and Frank [3] established an optimality criterion of the submodular flow problem in terms of node potentials, which can be reformulated as follows. See also [16, Theorem 5.3] for a direct proof. Theorem 2.1. A submodular flow ϕ is optimal if and only if there exists a node potential π : N → R such that cπa > 0 ⇒ ϕa = `a , cπa < 0 ⇒ ϕa = ua , and ∂ϕ ∈ B(f π ). Thus two equivalent ways to state the “submodular part” of the optimality condition are to require that an optimal ϕ and π satisfy that each Lπk is ϕ-tight, or that ∂ϕ maximizes the objective π T y over y ∈ B(f ). 3. Dual approximate optimality. We first cover the details of how to check node potentials π for optimality. We then show how to adapt checking exact optimality to checking a new kind of approximate optimality. This checking routine will form the core of our cut canceling algorithm, as it will produce the cuts that we cancel. For this paper a cut is just a non-empty, proper subset S ∈ D of N . We cancel a cut by increasing πi by some step length β for i ∈ S, which has the effect of increasing cπa on ∆+ S, decreasing cπa on ∆− S, and changing the level sets of π. Given a node potential π, we define modified bounds `π and uπ as follows. If cπa < 0 then `πa = uπa = `a ; if cπa = 0 then `πa = `a and uπa = ua ; if cπa > 0 then `πa = uπa = ua . Then Theorem 2.1 implies that π is optimal if and only if there is a feasible flow in the network Gπ with bounds `π , uπ , and f π . For the proof of correctness of our algorithm we will need some details of how feasibility of Gπ is checked. We do this using a variant of an algorithm of Frank [10] that we call the Feasibility Algorithm. It starts with any initial base x ∈ B(f π ) and any initial ϕ satisfying `π and uπ . We now define a residual network on N . If ϕij > `πij , then we have a backward residual arc j → i with residual capacity rij = ϕij − `πij ; if ϕij < uπij , then we have a forward residual arc i → j with rji = uπij − ϕij . To deal with relaxed conservation, for every i, j ∈ N with i 6= j and x ∈ B(f ), we define the exchange capacity w.r.t. x as re(x, i, j) = max{α | x + αχi − αχj ∈ B(f )}; thus re(x, i, j) > 0 means that there is no x-tight set containing i but not j. If we compute exchange capacity relative to f π instead of f , then we denote it by reπ . Since f π ≤ f , we have reπ (x, i, j) ≤ re(x, i, j). Make a jumping arc j → i with capacity 5
re(x, i, j) whenever re(x, i, j) > 0. Note that S is x-tight if and only if there are no jumping arcs w.r.t. x entering S. The residual network for π contains only the jumping arcs w.r.t. reπ . Also, [3, Theorem 7] shows that for any node potentials π ∈ RN and i, j ∈ N with πi = πj , the exchange capacity reπ (x, i, j) for f π is determined by a set S that nests with the Lπk , so that reπ (x, i, j) = f π (S) − x(S) = f (S) − x(S) ≥ re(x, i, j), and hence reπ (x, i, j) = re(x, i, j). The Feasibility Algorithm finds directed residual paths from N + = {i ∈ N | xi > ∂ϕ(i)} to N − = {i ∈ N | xi < ∂ϕ(i)} with a minimum number of arcs in the residual network. On each path it augments flow ϕ on the residual arcs, and modifies x as per the jumping arcs, which monotonically reduces the difference of x and ∂ϕ. By using a lexicographic selection rule, the algorithm terminates in finite time. At termination, either x coincides with ∂ϕ, which implies the optimality of π; or there is no directed path from N + to N − . In this last case, define T ⊆ N as the set of nodes from which N − is reachable by directed residual paths. No jumping arcs enter T , so it must belong to D, it must be tight for the final x, and it must contain all i with xi < ∂ϕ(i). Furthermore, we have ϕ(∆+ T ) = `π (∆+ T ) and ϕ(∆− T ) = uπ (∆− T ). Thus we get (3)
− π π V π (T ) ≡ `π (∆+ A T ) − u (∆A T ) − f (T ) = ∂ϕ(T ) − x(T ) > 0.
We call a node subset S with V π (S) > 0 a positive cut. Similar reasoning shows that for any other S ∈ D, we have V π (S) ≤ ∂ϕ(S) − x(S) ≤ ∂ϕ(T ) − x(T ) = V π (T ), proving that T is a most positive cut. Intuitively, V π (T ) measures how far away from optimality π is. We summarize as follows: Lemma 3.1. Node potentials π are optimal if and only if there are no positive cuts w.r.t. π. When π is not optimal, the output of the Feasibility Algorithm is a most positive cut. We denote the running time of the Feasibility Algorithm by FA. As usual, we assume that we have an oracle to compute exchange capacities, and denote its running time by h. Computing an exchange capacity is a submodular function minimization on a distributive lattice, which can be done via the ellipsoid method [21] or by combinatorial methods [25, 40] in strongly polynomial time. Cunningham and Frank [3, Theorem 7] show how to derive an oracle for computing exchange capacities for f π from an oracle for f with the same running time. Fujishige and Zhang [17] (see also [27]) show how to extend the push-relabel algorithm of Goldberg and Tarjan [19] to get FA = O(n3 h). Unfortunately, even for min-cost flow, an example of Hassin [22] shows that it may be necessary to cancel an exponential number of most positive cuts to achieve optimality. In [41] we get around this difficulty for MCF by relaxing the optimality conditions by a parameter δ > 0, and then applying scaling to δ. It then turns out that only a polynomial number of relaxed most positive cuts need to be canceled for a given value of the parameter, and only a polynomial number of scaled subproblems need to be solved. Until now it has been difficult to find a relaxation in the submodular flow case that would allow the same analysis to apply. Here we introduce a new relaxation that works. We relax the modified bounds on the arcs in A by widening them by δ just as π π,δ π in [5]. Define `π,δ a = `a − δ and ua = ua + δ for every arc a ∈ A. We also need to relax the submodular bounds on conservation by some function of δ. Define f π,δ (S) = f π (S) + δ|S| · |N − S|, which is closely related to a relaxation introduced in [24]. Since |S| · |N − S| is submodular, f π,δ is submodular. 6
We can now define the relaxed cut value of S as π,δ π,δ V π,δ (S) = `π,δ (∆+ (∆− (S). A S) − u A S) − f
We say that node potential π is δ-optimal if V π,δ (S) ≤ 0 holds for every cut S. Thus π is 0-optimal if and only if it is optimal. To check if π is δ-optimal, we need to see if there is a feasible flow in the network with bounds `π,δ , uπ,δ , and f π,δ . It turns out to be more convenient to move the δ b π = (N, A ∪ E) to be the directed relaxation off f π,δ onto a new set of arcs. Define G graph obtained from G by adding E = {i → j | i, j ∈ N, i 6= j} to the arc set. Extend the bounds `π , uπ , `π,δ and uπ,δ from A to E by setting `πe = uπe = 0, `π,δ =0 e and uπ,δ = δ for every e ∈ E. For convenience set I = A ∪ E and m0 = |A ∪ E| = e m + n(n − 1). Now checking `π,δ , uπ,δ and f π,δ for feasibility on G (with arc set A) is equivalent b π (with arc set I = A ∪ E). Also, to checking `π,δ , uπ,δ , and f π for feasibility on G the relaxed cut value on G can be re-expressed as π,δ π V π,δ (S) = `π,δ (∆+ (∆− I S) − u I S) − f (S).
b π , then π is δ-optimal. Otherwise, the If there is a feasible submodular flow in G Feasibility Algorithm gives us a node subset T that maximizes V π,δ (S), a relaxed most positive cut, or δ-MPC. When π is not δ-optimal, the Feasibility Algorithm also b π and a base x in B(f π ) such that xi ≤ ∂I ϕ(i) gives us an optimal “max flow” ϕ on G for i ∈ T and xi ≥ ∂I ϕ(i) for i ∈ / T. We also need to consider max mean cuts. The mean value of cut S is π
V (S) ≡
V π (S) . |∆A S| + |S| · |N − S| π
A max mean cut, or MMC, attains the maximum in δ(π) ≡ maxS V (S). By standard LP duality arguments δ(π) also equals the minimum δ such that there is a b π with bounds `π,δ , uπ,δ , and f π . Define U to be the maximum feasible flow in G absolute value of any `a , ua , or f ({i}). A max mean cut can be computed using log(nU ) O(min{m0 , 1+log log(nU )−log log n }) calls to the Feasibility Algorithm in the framework of discrete Newton’s Algorithm, see Ervolina and McCormick [34] or Radzik [37]. 4. Cut cancellation. We will start out with δ large and drive δ towards zero, since π is 0-optimal if and only if it is optimal. The next lemma says that δ need not start out too big, and need not end up too small. Its proof is similar to [5, Lemma 5.1]. Lemma 4.1. Suppose that `, u, and f are integral. Then any node potentials π are 2U -optimal. Moreover, when δ < 1/m0 , any δ-optimal node potentials are optimal. Our relaxed δ-MPC Canceling algorithm will start with δ = 2U and will execute scaling phases, where each phase first sets δ := δ/2. The input to a phase will be a 2δ-optimal set of node potentials from the previous phase, and its output will be a δ-optimal set of node potentials. Lemma 4.1 says that after O(log(nU )) scaling phases we have δ < 1/m0 , and we are optimal. Within a scaling phase we use the Feasibility Algorithm to find a δ-MPC T . We then cancel T by adding a constant step length β to πi for each node in T to get π 0 . In ordinary min-cost flow we choose the step length β based on the first arc in A whose reduced cost hits zero (as long as the associated flow is within the bounds, see 7
[41, Figure 2]). This bound on β is η ≡ min
π min{|cπa | | a ∈ ∆+ A T, ca < 0, ϕa ≥ `a } − π π min{ca | a ∈ ∆A T, ca > 0, ϕa ≤ ua }
,
(Note that optimality of T implies that every arc a with ϕa ≥ `a has negative reduced cost, and every arc a with ϕa ≤ ua has positive reduced cost.) Here we must also worry about the level set structure of π changing during the cancel, which is equivalent to the reduced cost of a jumping arc reaching zero. We 0 need to increase flows on jumping arcs leaving T so that the Lπk will be tight. We try to do this by decreasing flow on arcs of ∆− E T , all of which have flow δ since T is tight. If the exchange capacity of such an arc is at most δ we can do this. Otherwise, this E-arc will determine β, and the large exchange capacity will allow us to show that the algorithm makes sufficient progress. These considerations lead to the subroutine AdjustFlow below. It takes as input the optimal max flow ϕ from the Feasibility Algorithm, and its associated base x ∈ B(f π ). It computes an update ψ ∈ RE to ϕ and base x0 so that x0 will be the base associated with ϕ0 ≡ ϕ − ψ. Note that in AdjustFlow we always compute re(x, i, j) w.r.t. f , never w.r.t. f π , so that jumping arcs are w.r.t. f , not f π . We call an arc that determines β a determining arc. Algorithm AdjustFlow: begin ψe = 0 for e ∈ E; H := {j → i | i ∈ T, j ∈ N − T, 0 < πj − πi < η}; while H 6= ∅ do begin x0 := x − ∂E ψ; select e = j → i ∈ H with minimum πj − πi ; set H := H − {e}; if re(x0 , i, j) < δ then ψe := re(x0 , i, j); else return β := πj − πi [β is determined by jumping arc j → i]; end return β := η [β is determined by the A-arc determining η]; end. The full description of canceling a cut is now: Use the Feasibility Algorithm to find δ-MPC T , max flow ϕ, and base x. Run AdjustFlow to modify ϕ to ϕ0 = ϕ−ψ and x to x0 = x − ∂E ψ, and to choose β. Finally, increase πi by β for all i ∈ T . A scaling phase cancels δ-MPCs in this manner until π is δ-optimal. The next results show that canceling a δ-MPC cannot increase the δ-MPC value, and that V π,δ (S) will decrease by at least δ under some circumstances. Lemma 4.2. Suppose we cancel a δ-MPC T for π to get a node potential π 0 . Then π 0 ,δ V (S) ≤ V π,δ (T ) holds for any cut S. 0 0 Proof. It is easy to show similarly to [41, Lemma 5.3] that `πa ,δ ≤ ϕa ≤ uπa ,δ E holds for every a ∈ A. As a result of AdjustFlow, we obtain a flow ψ ∈ R with 0 ≤ ψe < δ for e ∈ E. Put ϕ0a = ϕa for a ∈ A and ϕ0e = ϕe − ψe for e ∈ E. To prove 0 that ϕ0 is feasible for B(f π ) we need: 0
Claim: x0 = x − ∂E ψ ∈ B(f π ). 8
Proof: If the initial H = ∅ in AdjustFlow (i.e., πj > πi implies that πj ≥ 0 πi + η), then ψ ≡ 0 and so x0 = x. In this case it can be shown that f π (S) = f π (S ∩ T ) + f π (S ∪ T ) − f π (T ). Since x(T ) = x0 (T ) = f π (T ), we have x0 (S) = 0 x0 (S ∩ T ) + x0 (S ∪ T ) − x0 (T ) ≤ f π (S ∩ T ) + f π (S ∪ T ) − f π (T ) = f π (S). Suppose instead that the initial H is non-empty. Here we know at least that x0 ∈ B(f ) (since ψji is never bigger than re(x0 , i, j)). Denote T ∩ (Lπk − Lπk−1 ) by Tk , and (Lπk − Lπk−1 ) − T by T k , so that the Tk partition T , and the T k partition N − T . Sq Sp Then a typical level set of π 0 looks like L0 = ( k=1 Tk ) ∪ ( k=1 T k ) = Lπp ∪ (T ∩ Lπq ) for some 0 ≤ p < q ≤ h. Now both T and Lπi are ϕ-tight for f π , so their union and intersection are both Si Si−1 ϕ-tight for f π . By the same reasoning, Tbi ≡ (T ∩ Lπi ) ∪ Lπi−1 = ( k=1 Tk ) ∪ ( k=1 T k ) is also ϕ-tight for f π . Since every arc of H starts in some T j and ends in some Ti with j < i, no arc of H crosses any Tbi , so each Tbi is also ϕ0 -tight for f π . Each Tbi nests with the Lπk , so we have f π (Tbi ) = f (Tbi ), so each Tbi is in fact ϕ0 -tight for f . This implies that the only possible jumping arcs entering level set L0 of π 0 are those from T j to Ti for p < j < i ≤ q. But these jumping arcs all belong to H and were removed 0 by AdjustFlow. Thus L0 is ϕ0 -tight for f , and so x0 ∈ B(f π ). Since ψe > 0 implies ϕe = δ, every e ∈ E satisfies 0 ≤ ϕ0e ≤ δ. Therefore we have 0
V π ,δ (S)
= ≤ = ≤ = =
0
0
0
π ,δ π `π ,δ (∆+ (∆− I S) − u I S) − f (S) 0 0 ∂I ϕ (S) − x (S) ∂I ϕ(S) − x(S) ∂I ϕ(T ) − x(T ) π,δ π `π,δ (∆+ (∆− I T) − u I T ) − f (T ) V π,δ (T )
0
(def’n of V π ,δ (S)) 0 (feas. of ϕ0 , x0 ∈ B(f π )) 0 0 (def’n of ϕ , x ) (T a δ-MPC) (T is ϕ-tight for f π ) (def’n of V π,δ (T )).
Corollary 4.3. With the same hypothesis and notation as Lemma 4.2, if the 0 determining arc a ∈ A ∪ E crosses S, then we have V π ,δ (S) ≤ V π,δ (T ) − δ. 0 π 0 ,δ Proof. If the determining arc a ∈ A, then `a + δ = `a ≤ ϕa ≤ ua = uπa ,δ − δ, and 0 + − + π 0 ,δ π 0 ,δ 0 hence ` (∆I S) − u (∆I S) ≤ ∂I ϕ (S) − δ. If a ∈ ∆E S, then ϕa = δ = `πa ,δ + δ, 0 0 − π ,δ 0 and hence `π ,δ (∆+ (∆− I S) − u I S) ≤ ∂I ϕ (S) − δ. If a = j → i ∈ ∆E S, then π0 0 π0 0 0 0 f (S) − x (S) ≥ re (x , i, j). Now a determining means that πi = πj , which implies 0 0 reπ (x0 , i, j) = re(x0 , i, j) ≥ δ, and hence f π (S) ≥ x0 (S) + δ. In all cases we have 0 0 V π ,δ (S) ≤ ∂I ϕ0 (S) − x0 (S) − δ, which implies V π ,δ (S) ≤ V π,δ (T ) − δ. Note that we have two different notions of “closeness to optimality” in this algorithm. At the outer level of scaling phases we drive δ towards zero (and so π towards optimality), and inside a phase we drive the δ-MPC value towards zero (and so π towards δ-optimality). Observe that the difference between ∂ϕ and x in the Feasibility Algorithm is an inner relaxation of the condition that π is δ-optimal, in that we keep a flow satisfying the bounds `π,δ and uπ,δ and a base x in B(f π ). As the phase progresses the difference is driven to zero, until the final flow proves δ-optimality of the final π. We can now prove the key lemma in the convergence proof of δ-MPC Canceling, which shows that we must be in a position to apply Corollary 4.3 reasonably often. Lemma 4.4. After at most n cut cancellations, the value of a δ-MPC decreases by at least δ. 9
Proof. Suppose that the first cancellation has initial node potentials π 0 and cancels δ-MPC T 1 to get π 1 , the next cancels T 2 to get π 2 , . . . , and the nth cancels T n to get π n . Each cancellation makes at least one arc into a determining arc. Consider the subgraph of these determining arcs. If a cut creates a determining arc whose ends are in the same connected component of this subgraph, then this cut must be crossed by a determining arc from an earlier cut. We can avoid this only if each new determining arc strictly decreases the number of connected components in the subgraph. This can happen at most n − 1 times, so it must happen at least once within n iterations. Let k be an iteration where T k shares a determining arc with an earlier cut T h . By h−1 h Corollary 4.3 applied to T = T h and S = T k , we have V π ,δ (T h ) ≥ V π ,δ (T k ) + δ. i−1 h k−1 Note that the δ-MPC value at iteration i is V π ,δ (T i ). If V π ,δ (T k ) ≥ V π ,δ (T k ), h−1 h k−1 0 Lemma 4.2 says that V π ,δ (T 1 ) ≥ V π ,δ (T h ) ≥ V π ,δ (T k ) + δ ≥ V π ,δ (T k ) + δ ≥ n−1 V π ,δ (T n ) + δ. h k−1 If instead V π ,δ (T k ) < V π ,δ (T k ), then let p be the latest iteration between h p−1 p and k with V π ,δ (T k ) < V π ,δ (T k ). The only way for the value of T k to increase like p−1 = 0 that crosses T k in the reverse orientation this is if there is an arc a ∈ A with caπ p−1 p p π ,δ k to its orientation in T , or f (T ) > f π ,δ (T k ) holds. Let ϕp be a flow proving − k − p + k p π p ,δ optimality of T p . If a ∈ ∆+ = A T ∩ ∆A T (a ∈ ∆A T ∩ ∆A T ) then we have ua p p−1 p ϕpa + 2δ ≥ ϕpa + δ (laπ ,δ = ϕpa − 2δ ≤ ϕpa − δ). When f π ,δ (T k ) > f π ,δ (T k ) holds, p there exist i, j ∈ N such that i ∈ T k \ T p , j ∈ T p \ T k . It follows from i → j ∈ ∆− I T p π p ,δ + p π p ,δ p and j → i ∈ ∆I T that ϕe = δ = le + δ for e = i → j and ϕe0 = 0 = ue0 − δ for p k π p ,δ k p k e0 = j → i. In either case, we have lπ ,δ (∆+ (∆− I T )−u I T ) ≤ ∂I ϕ (T ) − δ. Thus p π p−1 ,δ p equations in the proof of Lemma 4.2 apply, showing that V (T )−δ ≥ V π ,δ (T k ). 0 p−1 p Then Lemma 4.2 and the choice of p say that V π ,δ (T 1 ) ≥ V π ,δ (T p ) ≥ V π ,δ (T k )+ k−1 n−1 δ ≥ V π ,δ (T k ) + δ ≥ V π ,δ (T n ) + δ. Lemma 4.5. The value of every δ-MPC in a scaling phase is at most m0 δ. Proof. Let ϕ prove the 2δ-optimality of the initial π in a phase, so that ϕ violates the bounds lπ and uπ by at most 2δ on every arc of A ∪ E, and ϕsi = 0 for all i. To get an initial flow ϕ0 to begin the next phase that violates the bounds by at most δ, we need to change ϕ by at most δ on each arcPof A ∪ E. Then if we set ϕ0si = δ|∆+ I {i}| for all i, ϕ0 will satisfy f π also. The sum i∈N ϕ0si is at most m0 δ, and this is an upper bound on the value of the first δ-MPC in this phase. By Lemma 4.2, it is then also a bound on the value of every δ-MPC in the phase. Putting Lemmas 4.4 and 4.5 together yields our first bound on the running time of δ-MPC Canceling: Theorem 4.6. A scaling phase of δ-MPC Canceling cancels at most m0 n δ-MPCs. Thus the running time of δ-MPC Canceling is O(n3 log(nU )FA) = O(n6 h log(nU )). Proof. Lemma 4.5 shows that the δ-MPC value of the first cut in a phase is at most m0 δ. It takes at most n iterations to reduce this by δ, so there are at most m0 n = O(n3 ) iterations per phase. The time per iteration is dominated by computing a δ-MPC, which is O(FA). The number of phases is O(log(nU )) by Lemma 4.1. 5. A strongly polynomial bound. We first prove the following lemma. Such a result was first proved by Tardos [42] for MCF, and this version is an extension of dual versions by Fujishige [15] and Ervolina and McCormick [5]. Lemma 5.1. Suppose that flow ϕ0 proves that π 0 is δ 0 -optimal. If arc a ∈ A ∗ satisfies ϕ0a < ua − (m0 + 1)δ 0 , then all optimal π ∗ have cπa ≥ 0. If arc a ∈ A satisfies 10
∗
0
ϕ0a > `a + (m0 + 1)δ 0 , then all optimal π ∗ have cπa ≤ 0. If ∂I ϕ0 (T 0 ) > f π (T 0 ) + m0 nδ 0 for some π 0 and some T 0 ∈ D, then there is an E-arc i → j leaving some level set of π 0 such that all optimal π ∗ have πi∗ ≤ πj∗ . Proof. Note that ϕ0 proving δ 0 -optimality of π 0 implies that there is no jumping 0 0 arc entering any Lπk for any k. Since x0 ≡ ∂I ϕ0 ∈ B(f π ), this says that x0 is a 0 π -maximum base of B(f ). 0 0 0 0 Now change ϕ0 , a flow on A ∪ E feasible for `π ,δ and uπ ,δ , into flow ϕ b on just π0 π0 0 A, feasible for ` and u , by getting rid of ϕe for e ∈ E, and by changing ϕ0a by at most δ 0 on each a ∈ A. Note that this ϕ b will not in general satisfy ∂A ϕ b ∈ B(f ), but that ϕ b otherwise satisfies all optimality conditions: it satisfies the bounds and is b and N − = {i ∈ complementary slack with π 0 . Define N + = {i ∈ N | x0i > ∂A ϕ(i)} 0 b so that N | xi < ∂A ϕ(i)}, (4)
x0 (N + ) ≤ m0 δ 0 + ∂A ϕ(N b + ).
We now apply the successive shortest path algorithm for submodular flow of Fujishige [13] starting from ϕ. b (As originally stated, this algorithm is finite only for integer data, but the lexicographic techniques of Sch¨onsleben [38] and Lawler and Martel [31] show that it can be made finite for any data.) This algorithm looks for an augmenting path from a node i ∈ N + to a node j ∈ N − , where residual capacities on A-arcs come from ϕ, b and residual capacities on jumping arcs come from x0 . It chooses a shortest augmenting path (using the current reduced costs as lengths; such a path can be shown to always exist) and augments flow along the path, updating π 0 by the shortest path distances, and x0 as per the jumping arcs. This update maintains the properties that the current ϕ b satisfies the bounds and is complementary slack with the current π 0 , and the current x0 belongs to B(f ). The algorithm terminates with optimal flow ϕ∗ once the boundary of the current ϕ b coincides with the current x0 . By (4), the total amount of flow pushed by this algorithm is at most m0 δ 0 . This implies that for each a ∈ A, ϕ ba differs from ϕ∗a by at most m0 δ 0 , so ϕ0a differs 0 0 ∗ from ϕa by at most (m + 1)δ . In particular, if ϕ0a < ua − (m0 + 1)δ 0 , then ϕ∗a < ua , ∗ implying that cπa ≥ 0, and similarly for ϕ0a > `a + (m0 + 1)δ 0 . Suppose that P is an augmenting path chosen by the algorithm, and that flow is augmented by amount τP along P . Since P has at most n − 1 jumping P arcs, the boundary of any S ⊆ N changes by at most (n − 1)τP due to P . Since P τP ≤ m0 δ 0 by (4), the total change in ∂A ϕ(S) b during the algorithm is at most m0 (n − 1)δ 0 . Since 0 0 0 |∂I ϕ (S) − ∂A ϕ(S)| b ≤ m δ , the total change in ∂I ϕ0 (S) is at most m0 nδ 0 . Thus 0 0 0 π0 0 ∂I ϕ (T ) > f (T ) + m0 nδ 0 implies that ∂A ϕ∗ (T 0 ) > f π (T 0 ). This says that some 0 ∗ level set of π is not ϕ -tight. This implies that there is some E-arc i → j with πi0 > πj0 but πi∗ ≤ πj∗ . We now modify our algorithm a little bit. We divide our scaling phases into blocks of log2 (m0 + 1) phases each. At the beginning of each block we compute a max mean cut T 0 with mean value δ 0 = δ(π 0 ) and cancel T 0 (including calling AdjustFlow). This ensures that our current flow is δ 0 -optimal, so we set δ = δ 0 and start the block of scaling phases. It will turn out that only 2m0 blocks are sufficient to attain optimality. Theorem 5.2. This modified version of the algorithm takes O(n5 log nFA) = O(n8 h log n) time. Proof. Let T 0 be the max mean cut canceled at the beginning of the block, with associated node potential π 0 , flow ϕ0 , and mean value δ 0 . Complementary slackness 0 0 0 π0 0 for T 0 implies that ϕ0a = uπa + δ 0 for all a ∈ ∆− I T , that ϕa = `a − δ for all 11
0
0
+ 0 0 0 π 0 0 π 0 0 a ∈ ∆+ A T , that ϕa = `a for all a ∈ ∆I T , and that ∂I ϕ (T ) = f (T ). Let π , 0 0 ϕ , and δ be the similar values after the last phase in the block. Since each scaling phase cuts δ in half, we have δ 0 < δ 0 /(m0 + 1). − 0 0 0 0 0 π0 Subtracting ϕ0 (∆+ (T 0 ) yields I T ) − ϕ (∆I T ) − ∂I ϕ (T ) = 0 from V 0
(|∆A T 0 | + |T 0 | · |N − T 0 |)δ 0 = V π (T 0 ) (5)
0
0
− 0 0 0 π = (`π − ϕ0 )(∆+ I T ) + (ϕ − u )(∆I T ) 0
+ (∂I ϕ0 (T 0 ) − f π (T 0 )). 0
0 Now apply Lemma 5.1 to ϕ0 and π 0 . If the term for arc a of (`π − ϕ0 )(∆+ A T ) is at 0 0 0 π0 0 0 least δ > (m + 1)δ , then we must have that `a = ua and so ϕa < ua − (m + 1)δ 0 , ∗ 0 π0 and we can conclude that cπa ≥ 0. But each a in ∆+ A T had ca < 0, so this is a new 0 0 sign constraint on cπ . The case for terms of (ϕ0 − uπ )(∆− A T ) is similar. + 0 0 Suppose instead that all the terms in the ∆A T and ∆− A T sums of (5) are at 0 0 most δ . The total of all the E-arc terms is at most |T | · |N − T 0 |δ 0 . Therefore the only possibility left to achieve the large left-hand side of (5) is to have ∂I ϕ0 (T 0 ) > 0 f π (T 0 ) + (m0 n/2)δ 0 . Lemma 5.1 says that in this case there must be a jumping arc i → j leaving some level set of π 0 such that all optimal π ∗ have πi∗ ≤ πj∗ . Since πi0 was larger than πj0 , this is a new sign restriction on π. ∗ In either case each block imposes a new sign restriction on cπa for some I-arc a. At most 2m0 such sign restrictions can be imposed before π ∗ is completely determined, so after at most 2m0 = O(n2 ) blocks we must be optimal. Each block requires log(m0 + 1) = O(log n) scaling phases. The proof of Theorem 4.6 shows that each scaling phase costs O(n3 FA) time exclusive of computing the max mean cut. The time for computing a max mean cut is O(n2 FA), which is not a bottleneck.
6. Separable convex cost submodular flow. This section is devoted to a straightforward extension of our algorithm to the separable convex cost submodular flow problem. Of course, the linear case is a special case of the separable convex case, but this section is easier to understand after seeing the detailed proofs for the linear case of the previous section. Let ga : R → R be a convex cost function for an arc a ∈ A. Then the convex cost submodular flow problem is: X Minimize ga (ϕa ) a∈A
subject to ∂ϕ ∈ B(f ). (Here we take advantage of the ability to “hide” bounds in convex objective functions.) We denote by gaa (ξ) and ga` (ξ) respectively the left derivative and the right derivative of ga at ξ ∈ R. Then an optimality condition for this problem is: Theorem 6.1 ([16, Theorem 12.1]). A submodular flow ϕ is optimal if and only if there exists node potentials π : N → R such that gaa (ϕa ) ≤ πj − πi ≤ ga` (ϕa ) for every arc a = i → j, and ∂ϕ ∈ B(f π ). We define the modified bounds of a = i → j by uπa = sup{x | gaa (x) ≤ πj − πi },
`πa = inf{x | ga` (x) ≥ πj − πi }.
Then Theorem 6.1 implies that π is optimal if and only if there is a feasible flow in the network Gπ with bounds `π , uπ , and f π . From `π and uπ we get `π,δ = `π − δ 12
and uπ,δ = uπ + δ as before. Using `π,δ and uπ,δ we define δ-Most Positive Cuts and perform the δ-MPC Canceling algorithm just as in the linear cost case. The preliminary step length η is given by + min{−gaa (`π,δ a ) − πi + πj | a = i → j ∈ ∆ A T } η ≡ min , − min{ga` (uπ,δ a ) + πi − πj | a = i → j ∈ ∆ A T } and then β is obtained from AdjustFlow. We start with two technical lemmas showing that key properties of δ-MPC cancellation carry over to the separable convex case. Lemma 6.2. Suppose that T is a δ-MPC w.r.t. π, and we cancel T to get node 0 0 potential π 0 . If a ∈ ∆+ T (a ∈ ∆− T ), then `π,δ ≤ uπa (uπ,δ ≥ `πa ). When a is the a a 0 0 π π,δ π determining arc, `π,δ a ≥ `a (ua ≤ ua ). + Proof. If a ∈ ∆ T , we have πi < π 0 i , πj = π 0 j and a π,δ ) ≤ πj0 − πi0 . π 0 i − πi = β ≤ πj − πi − gaa (`π,δ a ) ⇒ g (` 0
0
It follows from the definition of uπa that gaa (uπa + ) + π 0 i − π 0 j > 0 for any > 0. 0 Hence it holds that gaa (uπa + ) > −π 0 i + π 0 j ≥ gaa (`π,δ a ), which, together with the 0 π0 π,δ convexity of ga implies uπa + > `π,δ , and thus u ≥ ` a a a . + 0 0 When a ∈ ∆ T is the determining arc, −π i + π j = gaa (`π,δ a ) holds. It follows 0 0 from the definition of `πa that ga` (`πa − ) + π 0 i − π 0 j < 0 for any > 0. Hence we have 0 0 ga` (`πa − ) ≥ gaa (`πa − ), which, together with the convexity of ga implies gaa (`π,δ a )> 0 − π0 π π,δ `π,δ a > `a − , and thus `a ≥ `a . The case a ∈ ∆ T is similar. Corollary 6.3. Suppose that T is a δ-MPC w.r.t. π, and we cancel T to get 0 node potential π 0 . For a max flow ϕ from the Feasibility Algorithm, we have `π ,δ ≤ 0 0 0 ϕ ≤ uπ ,δ . If a is the determining arc, `aπ ,δ + δ ≤ ϕa ≤ uπa ,δ − δ holds. 0 Proof. An arc a with `π,δ < `πa ,δ belongs to ∆− T , which implies ϕa = uπ,δ ≥ a a 0 0 + π0 δ , the arc a belongs to ∆ T , which > u `πa = `πa ,δ + δ. On the other hand, if uπ,δ a a 0 0 0 0 ≤ uπa = uaπ ,δ − δ. Hence, we have `π ,δ ≤ ϕ ≤ uπ ,δ . If a ∈ ∆+ T implies ϕa = `π,δ a 0 is the determining arc, ϕa ≤ uaπ ,δ − δ is shown above. From Lemma 6.2, we have 0 0 − `πa ,δ + δ = `πa ≤ `π,δ a ≤ ϕa . The case where a ∈ ∆ T is the determining arc is similar. It is easy to check that the above lemmas are enough to show that the analogues of Lemma 4.2, Corollary 4.3, and Lemma 4.4 hold for separable convex costs. When the ga are general, we must assume that we have an initial ϕ0 and bound B such that ϕ0 is B-optimal, and the best we can hope for is to compute a solution that is -optimal for some > 0 in time polynomial in log(B/). But if the ga are convex piecewise linear with integral breakpoints, then we can compute an exact optimal solution in polynomial time: Theorem 6.4. A scaling phase of δ-MPC Canceling cancels at most m0 n δ-MPCs. Thus δ-MPC Canceling takes O(n3 log(B/)FA) = O(n6 h log(B/)) time to find an -optimal solution. When `, u, and f are integer valued and each ga is a piecewise linear convex function with integral breakpoints, then δ-MPC Canceling finds an exact optimal solution in O(n6 log(nU )) time. Proof. Lemma 4.5 shows that the δ-MPC value of the first cut in a phase is at most m0 δ. It takes at most n iterations to reduce this by δ, so there are at most m0 n = O(n3 ) iterations per phase. The time per iteration is dominated by computing a δ-MPC, which is O(FA). When `, u and f are integer valued and each ga is a piecewise linear 13
function with integral breakpoints, then a proof similar to [5, Lemma 5.1] shows that a δ-optimal π with δ < 1/m0 is exactly optimal, so it suffices to choose = 1/m0 in this case. It is also possible to show similar to [30] that when each ga is piecewise linear and quadratic convex, that we can compute an exact optimal solution in polynomial time. Similarly, if each ga is piecewise linear, then we could compute an exact optimal solution in strongly polynomial time. 7. Max mean cut canceling. This section shows that machinery we developed for δ-MPC Canceling also leads to a polynomial max mean cut canceling algorithm for submodular flow. We are able to get this algorithm despite the negative result of [35] because we relax f by adding the arc set E, which is outside the class of algorithms considered in [35]. At each iteration, the max mean cut canceling algorithm finds a max mean cut (MMC) T , and calls AdjustFlow with T to get a step length β. We then cancel T by amount β, and repeat until no more positive cuts exist, at which point we are optimal, by Lemma 3.1. Lemma 7.1. Suppose we cancel an max mean cut T w.r.t. π with mean value δ(π) to get node potential π 0 . Then V
π0
π
(S) ≤ V (T ) holds for every cut S. b π with data `π,δ(π) , uπ,δ(π) , and f π . Since the Proof. Let ϕ be a feasible flow in G 0 π0 step length β is at most η, we have `a − δ(π) ≤ ϕa ≤ uπa + δ(π) for every a ∈ A. As a result of AdjustFlow, we obtain a flow ψ ∈ RE with 0 ≤ ψe < δ(π) for e ∈ E. Put ϕ0a = ϕa for a ∈ A and ϕ0e = ϕe − ψe for e ∈ E. As ∂I ϕ ∈ B(f π ), it 0 follows from the Claim in Lemma 4.2 that x0 ≡ ∂I ϕ0 ∈ B(f π ). Since ψe > 0 implies 0 0 ϕe = δ(π), every e ∈ E satisfies 0 = `πe ≤ ϕ0e ≤ δ(π) = uπe + δ(π). Therefore we have V
π0
0
(S)
0
0
=
− π π `π (∆+ A S) − u (∆A S) − f (S) |∆A S| + |S| · |N − S|
=
− π π `π (∆+ I S) − u (∆I S) − f (S) + − |∆A S| + |∆I S|
(`πe = uπe = 0)
≤
− 0 ∂I ϕ0 (S) + δ(π) · |∆+ A S| + δ(π) · |∆I S| − x (S) + − |∆A S| + |∆I S|
(feas. of ϕ0 )
0
0
(def’n of V
π0
(S))
0
0
0
(∂I ϕ0 (S) = x0 (S)) (def’n of δ(π)).
= δ(π) π = V (T )
Corollary 7.2. With the same hypothesis and notation as Lemma 7.1, if the π0
π
determining arc a ∈ A ∪ E crosses S, then we have V (S) ≤ (1 − 1/m0 )V (T ). 0 0 Proof. If the determining arc a ∈ A, then `πa = `a ≤ ϕa ≤ ua = uπa , and hence 0 0 − + + − − π 0 0 `π (∆+ A S) − u (∆A S) ≤ ϕ (∆A S) + δ(π) · |∆A S| − ϕ (∆A S) + δ(π) · |∆A S| − δ(π). + + + 0 π0 π0 0 If a ∈ ∆E S, then ϕa = δ(π) = `a + δ(π), and hence ` (∆E S) ≤ ϕ (∆E S) − δ(π). If π0 0 a ∈ ∆− E S, then as in the proof of Lemma 4.2 f (S) ≥ x (S) + δ(π). In all cases we 14
have V
π0
(S) ≤
+ − − 0 0 ϕ0 (∆+ I S) + δ(π) · |∆A S| − ϕ (∆I S) + δ(π) · |∆I S| − x (S) − δ(π) + − |∆A S| + |∆I S|
= δ(π) −
π δ(π) ≤ (1 − 1/m0 )V (T ). − + |∆I S|
|∆+ A S|
Now the analogue of Lemma 4.4 follows just as before: Lemma 7.3. After at most n cut cancellations, the value of an MMC decreases by a factor of at most (1 − 1/m0 ). This yields our first bound on the running time of MMC Canceling: Theorem 7.4. MMC canceling cancels at most O(m0 n log(nU )) = O(n3 log(nU )) MMCs. Proof. Let π 0 be the initial node potential. From Lemma 4.1, δ(π 0 ) ≤ 2U , and if δ(π) < 1/m0 then π is optimal. Hence there are O(m0 n log(nU )) iterations. The strongly polynomial argument in Section 5 applies to MMC Canceling as well: Theorem 7.5. MMC Canceling cancels at most O(n5 log n) MMCs and takes O(n10 h log n) time. Proof. Lemma 7.3 shows that after m0 ndlog(m0 + 1)e cancellations δ(π) decreases by a factor of at most 1/(m0 + 1). Applying Lemma 5.1 and the proof of Theo∗ rem 5.2, we can newly restrict the sign of a cπa for some I-arc a every O(m0 n log n) iterations. Therefore after O((m0 )2 n log n) = O(n5 log n) cancellations, we must be optimal. Since computing an MMC requires O(n2 FA) = O(n5 h) time, we obtain the claimed complexity. It seems likely that a faster version of MMC Canceling could be developed along the lines of the Dual Cancel and Tighten algorithm of [5], but this seems not so interesting in light of the faster algorithm of [7]. Also, it appears to be straightforward to combine the results of this section with the previous section to get a MMC Canceling algorithm for the separable convex costs case (as is done in [30]), but again this seems not so interesting since the algorithm of Section 6 is faster. 8. Crossing submodular functions. We return to the original Edmonds and Giles [4] assumption that f is only crossing submodular on D. Recall that in this case we can efficiently check whether B(f ) = ∅ by the bi-truncation algorithm of Frank and Tardos [12], and that Fujishige [14] showed that a nonempty base polyhedron B(f ) of a crossing submodular function f on a crossing family D is identical to the e that contains base polyhedron B(fe) of a submodular function fe on a ring family D D. We used this to partially justify that we could implement our algorithms using the e in place of the crossing submodular f and D. However, it ring submodular fe and D is not yet clear how to compute cut values, exchange capacities, etc. for fe when we are given only an oracle for f . All of our algorithms are based on applying the Feasibility Algorithm to various residual networks. The Feasibility Algorithm needs an initial base to work with. This can be computed using the bi-truncation algorithm. It also needs an oracle for computing exchange capacities re(x, i, j). Note that re(x, i, j) is the minimum value of f (S) − x(S) over S containing i but not j. The subfamily of members of D containing 15
i but not j is a ring family, on which f is ring submodular, and so it is reasonable to assume that we have an oracle for computing these exchange capacities. e This will have the effect of simulating the action of the algorithm on fe and D, π π e meaning that the cut values V in (3) will be computed relative to f instead of f π , and the Feasibility Algorithm provides a tight set for feπ , not for f π . However, since we are always computing exchange capacities w.r.t. f and not fe, in fact we never need to compute a value of fe during the algorithm (although the analysis uses fe). Acknowledgment. We heartily thank Lisa Fleischer for many valuable conversations about this paper and helpful comments on previous drafts of it, and an anonymous referee for suggesting several improvements to the paper. REFERENCES [1] R. K. Ahuja, T. L. Magnanti, and J. B. Orlin, Network Flows — Theory, Algorithms, and Applications, Prentice Hall, 1993. [2] W. Cui and S. Fujishige, A primal algorithm for the submodular flow problem with minimummean cycle selection, J. Oper. Res. Soc. Japan, 31 (1988), pp. 431–440. [3] W. H. Cunningham and A. Frank, A primal-dual algorithm for submodular flows, Math. Oper. Res., 10 (1985), pp. 251–262. [4] J. Edmonds and R. Giles, A min-max relation for submodular functions on graphs, Ann. Discrete Math., 1 (1977), pp. 185–204. [5] T. R. Ervolina and S. T. McCormick, Two strongly polynomial cut canceling algorithms for minimum cost network flow, Discrete Appl. Math., 46 (1993), pp. 133–165. [6] T. R. Ervolina and S. T. McCormick, Canceling most helpful total cuts for minimum cost network flow, Networks 23, (1993) pp. 41–52. [7] L. Fleischer, S. Iwata, and S. T. McCormick A faster capacity scaling algorithm for minimum cost submodular flow, Math. Programming, 92 (2002), pp. 119–139. [8] A. Frank, A weighted matroid intersection algorithm, J. Algorithms, 2 (1981), pp. 328–336. [9] A. Frank, An algorithm for submodular functions on graphs, Ann. Discrete Math., 16 (1982), pp. 97–120. [10] A. Frank, Finding feasible vectors of Edmonds-Giles polyhedra, J. Combinatorial Theory, B36 (1984), pp. 221–239. ´ Tardos, An application of simultaneous Diophantine approximation in [11] A. Frank and E. combinatorial optimization, Combinatorica, 7 (1987), pp. 49–65. ´ Tardos, Generalized polymatroids and submodular flows, Math. Program[12] A. Frank and E. ming, 42 (1988), pp. 489–563. [13] S. Fujishige, Algorithms for solving the independent-flow problems, J. Oper. Res. Soc. Japan, 21 (1978), pp. 189–204. [14] S. Fujishige, Structures of polyhedra determined by submodular functions on crossing families, Math. Programming, 29 (1984), pp. 125–141. [15] S. Fujishige, Capacity-rounding algorithm for the minimum-cost circulation problem: A dual framework of the Tardos algorithm, Math. Programming, 35 (1986), pp. 298–309. [16] S. Fujishige, Submodular Functions and Optimization, North-Holland, 1991. [17] S. Fujishige and X. Zhang, New algorithms for the intersection problem of submodular systems, Japan J. Indust. Appl. Math., 9 (1992), pp. 369–382. ¨ ck, and U. Zimmermann, A strongly polynomial algorithm for minimum [18] S. Fujishige, H. Ro cost submodular flow problems, Math. Oper. Res., 14 (1989), pp. 60–69. [19] A. V. Goldberg and R. E. Tarjan, A new approach to the maximum flow problem, J. ACM, 35 (1988), pp. 921–940. [20] A. V. Goldberg and R. E. Tarjan, Finding minimum-cost circulations by canceling negative cycles, J. ACM, 36 (1989), pp. 873–886. ¨ tschel, L. Lova ´sz, and A. Schrijver, Geometric Algorithms and Combinatorial [21] M. Gro Optimization, Springer-Verlag, 1988. [22] R. Hassin, Algorithm for the minimum cost circulation problem based on maximizing the mean improvement, Oper. Res. Lett., 12 (1992), pp. 227–233. [23] M. Iri and N. Tomizawa, An algorithm for finding an optimal “independent assignment”, J. Oper. Res. Soc. Japan, 19 (1976), pp. 32–57.
16
[24] S. Iwata, A capacity scaling algorithm for convex cost submodular flows, Math. Programming, 76 (1997), pp. 299–308. [25] S. Iwata, L. Fleischer, and S. Fujishige, A combinatorial strongly polynomial algorithm for minimizing submodular functions, J. ACM, 48 (2001), pp. 761–777. [26] S. Iwata, S. T. McCormick, and M. Shigeno, A strongly polynomial cut canceling algorithm for the submodular flow problem, Proceedings of the Seventh IPCO Conference (1999), G. Cornu´ ejols, R. E. Burkard, and G. J. Woeginger, eds., pp. 259–272. [27] S. Iwata, S. T. McCormick, and M. Shigeno, A fast cost scaling algorithm for submodular flow, Inform. Process. Lett., 74 (2000), pp. 123–128. [28] S. Iwata, S. T. McCormick, and M. Shigeno, Fast cycle canceling algorithms for minimum cost submodular flow, Combinatorica, 23 (2003), pp. 503–525. [29] S. Iwata, S. T. McCormick, and M. Shigeno, A relaxed cycle-canceling approach to separable convex optimization in unimodular linear space, manuscript. [30] A. V. Karzanov and S. T. McCormick, Polynomial methods for separable convex optimization in unimodular linear spaces with applications, SIAM J. Comput., 26 (1997), 1245– 1275. [31] E. L. Lawler and C. U. Martel, Computing maximal polymatroidal network flows, Math. Oper. Res., 7 (1982), pp. 334–347. [32] S. T. McCormick (2003). Submodular Function Minimization. A chapter in the forthcoming Handbook on Discrete Optimization, Elsevier, K. Aardal, G. Nemhauser, and R. Weismantel, eds. [33] S. T. McCormick and T. R. Ervolina, Cancelling most helpful total submodular cuts for submodular flow, Integer Programming and Combinatorial Optimization (Proceedings of the Third IPCO Conference), G. Rinaldi and L. A. Wolsey eds. (1993), pp. 343–353. [34] S. T. McCormick and T. R. Ervolina, Computing maximum mean cuts, Discrete Appl. Math., 52 (1994), pp. 53–70. [35] S. T. McCormick, T. R. Ervolina and B. Zhou, Mean canceling algorithms for general linear programs and why they (probably) don’t work for submodular flow, UBC Faculty of Commerce Working Paper 94-MSC-011 (1994). [36] S. T. McCormick and A. Shioura, Minimum ratio canceling is polynomial for linear programming, but not strongly polynomial, even for networks, Oper. Res. Lett., 27 (2000), pp. 199–207. [37] T. Radzik, Newton’s method for fractional combinatorial optimization, Proceedings of the 33rd IEEE Annual Symposium on Foundations of Computer Science (1992), pp. 659– 669; Parametric flows, weighted means of cuts, and fractional combinatorial optimization, Complexity in Numerical Optimization, P. Pardalos, ed. (World Scientific, 1993), pp. 351– 386. ¨ nsleben, Ganzzahlige Polymatroid-Intersektions Algorithmen, Dissertation, ETH [38] P. Scho Z¨ urich, 1980. [39] A. Schrijver, Total dual integrality from directed graphs, crossing families, and sub- and supermodular functions, Progress in Combinatorial Optimization, W. R. Pulleyblank, ed. (Academic Press, 1984), pp. 315–361. [40] A. Schrijver, A combinatorial algorithm minimizing submodular functions in strongly polynomial time, J. Combinatorial Theory, Ser. B, 80 (2000), pp. 346–355. [41] M. Shigeno, S. Iwata, and S. T. McCormick, Relaxed most negative cycle and most positive cut canceling algorithms for minimum cost flow, Math. Oper. Res., 25 (2000), pp. 76–104. ´ Tardos, A strongly polynomial minimum cost circulation algorithm, Combinatorica, 5 [42] E. (1985), pp. 247–255. [43] C. Wallacher and U. Zimmermann, A polynomial cycle canceling algorithm for submodular flows, Math. Programming, 86 (1999), pp. 1–15. [44] U. Zimmermann, Negative circuits for flows and submodular flows, Discrete Appl. Math., 36 (1992), pp. 179–189.
17