Strongly polynomial algorithm for a class of minimum-cost flow problems with separable convex objectives
aa r X i v : . [ c s . D S ] J un A Strongly Polynomial Algorithm for a Class of Minimum-CostFlow Problems with Separable Convex Objectives
L´aszl´o A. V´egh ∗ May 24, 2018
Abstract
A well-studied nonlinear extension of the minimum-cost flow problem is to minimize theobjective P ij ∈ E C ij ( f ij ) over feasible flows f , where on every arc ij of the network, C ij is aconvex function. We give a strongly polynomial algorithm for the case when all C ij ’s are convexquadratic functions, settling an open problem raised e.g. by Hochbaum [16]. We also givestrongly polynomial algorithms for computing market equilibria in Fisher markets with linearutilities and with spending constraint utilities, that can be formulated in this framework (seeShmyrev [33], Devanur et al. [2]). For the latter class this resolves an open question raised byVazirani [37]. The running time is O ( m log m ) for quadratic costs, O ( n + n ( m + n log n ) log n )for Fisher’s markets with linear utilities and O ( mn + m ( m + n log n ) log m ) for spendingconstraint utilities.All these algorithms are presented in a common framework that addresses the general prob-lem setting. Whereas it is impossible to give a strongly polynomial algorithm for the generalproblem even in an approximate sense (see Hochbaum [16]), we show that assuming the ex-istence of certain black-box oracles, one can give an algorithm using a strongly polynomialnumber of arithmetic operations and oracle calls only. The particular algorithms can be derivedby implementing these oracles in the respective settings. Let us consider an optimization problem where the input is given by N numbers. An algorithmfor such a problem is called strongly polynomial (see [13]), if (i) it uses only elementary arithmeticoperations (addition, subtraction, multiplication, division, and comparison); (ii) the number ofthese operations is bounded by a polynomial of N (iii) if all numbers in the input are rational, thenall numbers occurring in the computations are rational numbers of size polynomially bounded in N and the maximum size of the input numbers. Here, the size of a rational number p/q is definedas ⌈ log ( p + 1) ⌉ + ⌈ log ( q + 1) ⌉ .The flow polyhedron is defined on a directed network G = ( V, E ) by arc capacity and nodedemand constraints; throughout the paper, n = | V | and m = | E | . We study the minimum costconvex separable flow problem : for feasible flows f , the objective is to minimize P ij ∈ E C ij ( f ij ),where on each arc ij ∈ E , C ij is a differentiable convex function. We give a strongly polynomialalgorithm for the case of convex quadratic functions, i.e. if C ij ( α ) = c ij α + d ij α with c ij ≥ ij ∈ E . We also give strongly polynomial algorithms for Fisher’s market with linear ∗ Department of Mathematics, London School of Economics. E-mail:
[email protected] . This work was done atthe College of Computing, Georgia Institute of Technology, supported by NSF Grant CCF-0914732. This is a revisionof the work published in STOC’12, May 19 - 22 2012, New York, NY, USA Copyright 2012 ACM 978-1-4503-1245-5/12/05 α (log α −
1) on certain arcs.These algorithms are obtained as special implementations of an algorithm that works for thegeneral problem setting under certain assumptions. We assume that the functions are representedby oracles (the specific details are provided later), and two further black-box oracles are provided.We give a strongly polynomial algorithm in the sense that it uses only basic arithmetic operationsand oracle calls; the total number of these operations is polynomial in n and m . We then verify ourassumptions for convex quadratic objectives and the Fisher markets, and show that we can obtainstrongly polynomial algorithms for these problems.Flows with separable convex objectives are natural convex extensions of minimum-cost flowswith several applications such as matrix balancing or traffic networks, see [1, Chapter 14] forfurther references. Polynomial-time combinatorial algorithms were given by Minoux [26] in 1986,by Hochbaum and Shantikumar [18] in 1990, and by Karzanov and McCormick [22] in 1997. Thelatter two approaches are able to solve even more general problems of minimizing a separable (notnecessarily differentiable) convex objective over a polytope given by a matrix with a bound on itslargest subdeterminant. Both approaches give polynomial, yet not strongly polynomial algorithms.In contrast, for the same problems with linear objectives, Tardos [35, 36] gave strongly polyno-mial algorithms. One might wonder whether this could also be extended to the convex setting. Thisseems impossible for arbitrary convex objectives by the very nature of the problem: the optimalsolution might be irrational, and thus the exact optimum cannot be achieved.Beyond irrationality, the result of Hochbaum [16] shows that it is impossible to find an ε -accurate solution in strongly polynomial time even for a network consisting of parallel arcs betweena source and a sink node and the C ij ’s being polynomials of degree at least three. This is basedon Renegar’s result [31] showing the impossibility of finding ε -approximate roots of polynomialsin strongly polynomial time. This is an unconditional impossibility result in a computation modelallowing basic arithmetic operations and comparisons; it does not rely on any complexity theoryassumptions.The remaining class of polynomial objectives with hope of strongly polynomial algorithms iswhere every cost function is convex quadratic. If all coefficients are rational, then the existence of arational optimal solution is guaranteed. Granot and Skorin-Kapov [12] extended Tardos’s method[36] to solving separable convex quadratic optimization problems with linear constraints, wherethe running time depends only on the entries of the constraint matrix and the coefficients of thequadratic terms in the objective. However, this algorithm is not strongly polynomial because ofthe dependence on the quadratic terms.The existence of a strongly polynomial algorithm for the quadratic flow problem thus remainedan important open question (mentioned e.g. in [16, 4, 17, 12, 34]). The survey paper [17] givesan overview of special cases solvable in strongly polynomial time. These include a fixed numberof suppliers (Cosares and Hochbaum, [4]), and series-parallel graphs (Tamir [34]). We resolve thisquestion affirmatively, providing a strongly polynomial algorithm for the general problem in time O ( m log m ).There is an analogous situation for convex closure sets: [16] shows that no strongly polynomialalgorithm may exist in general, but for quadratic cost functions, Hochbaum and Queyranne [15]gave a strongly polynomial algorithm.An entirely different motivation of our study comes from the study of market equilibriumalgorithms. Devanur et al. [5] developed a polynomial time combinatorial algorithm for a classical A solution x is called ε -accurate if there exists an optimal solution x ∗ with || x − x ∗ || ∞ ≤ ε . O ( n + n ( m + n log n ) log n ) for linear and O ( mn + m ( m + n log n ) log m ) forspending constraint utilities, with m being the number of segments in the latter problem. For thelinear case, Orlin [30] used the assumption m = O ( n ) and achieved running time O ( n log n ),the same as ours under this assumption. So far, no extensions of [30] are known for other marketsettings. For linear minimum-cost flows, the first polynomial time algorithm was the scaling method byEdmonds and Karp [6]. The current most efficient strongly polynomial algorithm, given by Orlin[29], is also based on this framework. On the other hand, Minoux extended [6] to the convexminimum-cost flow problem, first to convex quadratic flows [25], later to general convex objectives[26]. Our algorithm is an enhanced version of the latter algorithm, in the spirit of Orlin’s technique[29]. However, there are important differences that make the nonlinear setting significantly harder.Let us remark that Orlin’s strongly polynomial algorithm for linear Fisher market [30] is alsobased on the ideas of [29]. In what follows, we give an informal overview of the key ideas of thesealgorithms that motivated our result. For more detailed references and proofs, we refer the readerto [1].The algorithm of Edmonds and Karp consists of ∆-phases for a scaling parameter ∆. Initially,∆ is set to a large value, and decreases by at least a factor of two at the end of each phase. Anoptimal solution can be obtained for sufficently small ∆. The elementary step of the ∆-phasetransports ∆ units of flow from a node with excess at least ∆ to another node with demand atleast ∆. This is done on a shortest path in the ∆-residual network, the graph of residual arcswith capacity at least ∆. An invariant property maintained in the ∆-phase is that the ∆-residualnetwork does not contain any negative cost cycles. When moving to the next phase, the flow onthe arcs has to be slightly modified to restore the invariant property.Orlin’s algorithm ([29], see also [1, Chapter 10.6]) works on a problem instance with no uppercapacities on the arcs (every minimum-cost flow problem can be easily transformed to this form).The basic idea is that if the algorithm runs for infinite number of phases, then the solution convergesto an optimal solution; furthermore, the total change of the flow value in the ∆-phase and all3ubsequent phases is at most 4 n ∆ on every arc. Consequently, if an arc ij has flow > n ∆ inthe ∆-phase, then the flow on ij must be positive in some optimal solution. Using primal-dualslackness, this means that ij must be tight for an arbitrary dual optimal solution (that is, thecorresponding dual inequality must hold with equality). It is shown that within O (log n ) scalingphases, an arc ij with flow larger than 4 n ∆ appears.Based on this fact, [29] obtains the following simple algorithm. Let us start running theEdmonds-Karp algorithm on the input graph. Once there is an arc with flow larger than 4 n ∆,it is contracted and the Edmonds-Karp algorithm is restarted on the smaller graph. The method isiterated until the graph reduces to a single node. A dual optimal solution on the contracted graphcan be easily extended to a dual optimal solution in the original graph by reversing the contrac-tion operations. Provided a dual optimal solution, a primal optimal solution can be obtained bya single maximum flow computation. The paper [29] (see also [1, Chapter 10.7]) also contains asecond, more efficient algorithm. When an arc with “large” flow is found, instead of contractingand restarting, the arc is added to a special forest F . The scaling algorithm exploits properties ofthis forest and can thereby ensure that a new arc enters F in O (log n ) phases. The running timecan be bounded by O ( m log n ( m + n log n )), so far the most efficient minimum-cost flow algorithmknown.Let us now turn to the nonlinear setting. By the Karush-Kuhn-Tucker (KKT) conditions,a feasible solution is optimal if and only if the residual graph contains no negative cycles withrespect to the cost function C ′ ij ( f ij ). Minoux’s algorithm is a natural extension of the Edmonds-Karp scaling technique (see [25, 26], [1, Chapter 14.5]). In the ∆-phase it maintains the invariantthat the ∆-residual graph contains no negative cycle with respect to the relaxed cost function( C ij ( f ij + ∆) − C ij ( f ij )) / ∆. When transporting ∆-units of flow on a shortest path with respectto this cost function, this invariant is maintained. A key observation is that when moving to the∆ / /
2. Therole of the scaling factor ∆ is twofold: besides being the quantity of the transported flow, it alsoapproximates optimality in the following sense. As ∆ approaches 0, the cost of ij converges to thederivative C ′ ij ( f ij ). Consequently, the solution converges to a feasible optimal solution. A variantof this algorithm is outlined in Section 3. To formulate the exact assumptions needed for the general algorithm, several notions have to beintroduced. Therefore we postpone the formulation of our main result Theorem 4.5 to Section 4.2.Now we exhibit the main ideas on the example of convex quadratic functions. We only give an in-formal overview here without providing all technical details; the precise definitions and descriptionsare given in the later parts of the paper. Then in Section 6.1, we show how the general frameworkcan be adapted to convex quadratic functions.Let us assume that C ij ( α ) = c ij α + d ij α with c ij > ij ∈ E , and thereforeall cost functions are strictly convex. This guarantees that the optimal solution is unique. Thisassumption is made only for the sake of this overview, and not used in the formal presentationstarting in Section 2. However, it is useful as the uniqueness of the optimum enables certaintechnical simplifications. We discuss these simplifications at the end of the section. Our problem4an be formulated as follows. min X ij ∈ E c ij f ij + d ij f ij X j : ji ∈ E f ji − X j : ij ∈ E f ij = b i ∀ i ∈ Vf ≥ f ∗ be the optimal solution; it is unique by the strict convexity of the objective. Let F ∗ denote the support of f ∗ . An optimal solution can be characterized using the Karush-Kuhn-Tuckerconditions: for Lagrange multipliers π : V → R , we have π j − π i ≤ c ij f ∗ ij + d ij with equalitywhenever f ∗ ij >
0, that is, ij ∈ F ∗ . Consequently, if F ∗ is provided, then we can obtain f ∗ as theunique solution to the following system of linear equations (see Section 6.1 for details). π j − π i = 2 c ij f ∗ ij + d ij ∀ ij ∈ F ∗ X j : ji ∈ F ∗ f ∗ ji − X j : ij ∈ F ∗ f ∗ ij = b i ∀ i ∈ V (1) f ∗ ij = 0 ∀ ij ∈ E \ F ∗ We will assume the existence of the subroutine
Trial ( F, ˆ b ) (Oracle 2), where F ⊆ E is an arbitraryarc set and ˆ b : V → R such that the sum of the ˆ b i values is 0 in every undirected connectedcomponent of F . The subroutine solves the modification of (1) when F ∗ is substituted by F and b by ˆ b . The system is feasible under the above assumption on ˆ b , and a solution can be found in time O ( n . ) (see Lemma 6.1).Our starting point is a variant of Minoux’s nonlinear scaling scheme as described above, withthe only difference that the relaxed cost function is replaced by C ′ ij ( f ij + ∆) (see Section 3).Following Orlin [29], we can identify an arc carrying a “large” amount of flow in O (log n ) steps.The required amount, (2 n + m + 1)∆ at the end of the ∆-phase, is large enough that even if werun the algorithm forever and thereby converge to the optimal solution f ∗ , this arc must remainpositive. Consequently, it must be contained in F ∗ . However, we cannot simply contract suchan arc as in [29]. The reason is that the KKT-conditions give π j − π i = c ij f ∗ ij + d ij , a conditioncontaining both primal and dual (more precisely, Lagrangian) variables simultaneously.In every phase of the algorithm, we shall maintain a set F ⊆ F ∗ of arcs, called revealed arcs. F will be extended by a new arc in every O (log n ) phases; thus we find F ∗ in O ( m log n ) steps (seeTheorem 5.5). Given a set F ⊆ F ∗ , we introduce some technical notions; the precise definitionsand detailed discussions are given in Section 4.1. First, we waive the nonnegativity requirement onthe arcs in F : a vector E → R is called an F -pseudoflow , if f ij ≥ ij ∈ E \ F but the arcs in F are unconstrained.For an F -pseudoflow f and a scaling factor ∆ >
0, the (∆ , F )-residual graph E Ff (∆) containsall residual arcs where f can be increased by ∆ so that it remains an F -pseudoflow (that is, allarcs in E , and all arcs ji where ij ∈ F , or ij ∈ E \ F and f ij ≥ ∆.) We require that the flow f in this phase satisfies the (∆ , F ) -feasibility property: the graph E Ff (∆) contains no negative cycleswith respect to the cost function C ′ ij ( f ij + ∆).Let us now describe our algorithm. We start with F = ∅ and a sufficiently large ∆ value sothat the initial flow f ≡ , ∅ )-feasible. We run the Minoux-type scaling algorithm sendingflow on shortest paths in the (∆ , F )-residual graph from nodes with excess at least ∆ to nodes with5eficiency at least ∆. If there exist no more such paths, we move to the ∆ / / , F )-feasible one, on the cost of increasing thetotal excess by at most m ∆ / Adjust in Section 4.1). We include in F every edgewith f ij > (2 n + m + 1)∆ at the end of the ∆-phase.At the end of each phase when F is extended, we perform a special subroutine instead ofsimply moving to the ∆ / D F ( b ) defined as follows. Let D F ( b ) = max K | P i ∈ K b i | , where K ranges over the undirected connected components of F . (Notethat the requirement on ˆ b in the subroutine Trial ( F, ˆ b ) above was D F (ˆ b ) = 0.) If the discrepancy D F ( b ) is large, then it can be shown that F will be extended within O (log n ) phases as in Orlin’salgorithm (see the first part of the proof of Theorem 5.5).If the discrepancy is small, the procedure Trial-and-Error is performed, consisting of twosubroutines. First, we run the subroutine
Trial ( F, ˆ b ), where ˆ b is a small modification of b satisfying D F (ˆ b ) = 0. This returns an F -pseudoflow ˆ f , satisfying (1) with F in the place of F ∗ . (Thisstep be seen as “pretending” that F = F ∗ and trying to compute an optimal solution under thishypothesis.) The resulting ˆ f is optimal if and only if F = F ∗ . Otherwise, we use a second subroutine Error ( ˆ f , F ) (see Oracle 3), that returns the smallest value ˆ∆ > f is ( F, ˆ∆)-feasible.This subroutine can be reduced to a minimum cost-to-time ratio cycle problem (also known as thetramp streamer problem), see [1, Chapter 5.7]; a strongly polynomial time algorithm was given byMegiddo [23].If ˆ∆ < ∆ /
2, then we set ˆ∆ as our next scaling value and f = ˆ f as the next F -pseudoflow -we can proceed since ˆ f is ( F, ˆ∆)-feasible. Otherwise, the standard transition to phase ∆ / f . The analysis shows that a new arc shall be revealed in every O (log n )phases. The key Lemma 5.4 relies on the proximity of f and ˆ f , which implies that Trial-and-Error cannot return the same ˆ f if performed again after O (log n ) phases. Consequently, the set F cannot be the same, and has been therefore extended. Since | F | ≤ m , this shows that the totalnumber of scaling phases is O ( m log n ).Besides the impossibility of contraction, an important difference as compared ot Orlin’s al-gorithm is that F cannot be assumed to be a forest (in the undirected sense). There are simplequadratic instances with the support of an optimal solution containing cycles. In Orlin’s algorithm,progress is always made by connecting two components of F . This will also be an important eventin our algorithm, but sometimes F shall be extended with arcs inside a component.The main difference when applied to Fisher markets instead of quadratic costs is the imple-mentation of the black boxes Trial and
Error . These are implemented by a simple linear timealgorithm and the Floyd-Warshall algorithm, respectively. The description above made the simpli-fying assumption that c ij > ij ∈ E , that is, all cost functions are strictly convex, and thusthere is a unique optimal solution. This might not be true even for quadratic costs if c ij = 0 isallowed on certain arcs. An important difference between the description and the general algorithmis that in the general algorithm, the set F ∗ has to be more carefully defined; in particular, it willcontain the support of every optimal solution. We therefore have to introduce the additional notionof F -optimal solutions for F ⊆ F ∗ . The algorithm will find F -optimal solutions instead of optimalones; however, an F -optimal solution can be converted to an optimal solution via an additionalmaximum flow subroutine.The rest of the paper is organized as follows. Section 2 contains the basic definitions andnotations. Section 3 presents the simple adaptation of the Edmonds-Karp algorithm for convexcost functions, following Minoux [26]. Our algorithm in Section 4 is built on this algorithm withthe addition of the subroutine Trial-and-Error , that guarantees strongly polynomial runningtime. Analysis is given in Section 5. Section 6 adapts the general algorithm for quadratic utilities,6nd for Fisher’s market with linear and with spending constraint utilities. Section 7 contains afinal discussion of the results and some open questions. An Appendix contains the description ofthe shortest path subroutines used. A table summarizing notation and concepts can be found atthe end of the paper.
Let G = ( V, E ) be a directed graph, and let n = | V | , m = | E | . For notational convenience, weassume that the graph contains no parallel arcs and no pairs of oppositely directed arcs. Conse-quently, we can denote the arc from node i to node j by ij . All results straightforwardly extendto general graphs. We are given node demands b : V → R with P i ∈ V b i = 0. The flow is re-stricted to be nonnegative on every arc, but there are no upper capacities. On each arc ij ∈ E , C ij : R → R ∪ {∞} is a convex function. We allow two types of arcs ij : • Free arcs: C ij is differentiable everywhere on R . • Restricted arcs: C ij ( α ) = ∞ if α < C ij is differentiable on (0 , ∞ ) and has a left derivativein 0 that equals −∞ ; let C ′ ij (0) = −∞ denote this left derivative. Let us use the convention C ′ ij ( α ) = −∞ for α < C ′ ij is continuous on R for free and on [ ℓ ij , ∞ ) for restricted arcs. Restricted arcswill play a role in the Fisher market applications, where the function C ij ( α ) = α (log α −
1) will beused on certain arcs (with C ij (0) = 0 and C ij ( α ) = ∞ if α < X ij ∈ E C ij ( f ij ) X j : ji ∈ E f ji − X j : ij ∈ E f ij = b i ∀ i ∈ V (P) f ≥ ℓ ij ≤ f ij ≤ u ij . Onecan reduce general capacities to the above form via the following standard reduction (see e.g. [1,Sec 2.4]). For each arc ij ∈ E , let us add a new node k = k ij , and replace ij by two arcs ik and jk . Let us set b k = u ij − ℓ ij , C ik ( α ) = C ij ( α + ℓ ij ), C jk ≡
0. Furthermore, let us increase b i by ℓ ij and decrease b j by u ij . It is easy to see that this gives an equivalent optimization problem, and ifthe original graph had n ′ nodes and m ′ arcs, the transformed instance has n = n ′ + m ′ nodes and m = 2 m ′ arcs.Further, we may assume without loss of generality that G = ( V, E ) is strongly connected and(P) is always feasible. Indeed, we can add a new node t with edges vt , tv for any v ∈ V , withextremely high (possibly linear) cost functions on the edges. This guarantees that an optimalsolution shall not use such edges, whenever the problem is feasible. We will also assume n ≤ m .By a pseudoflow we mean a function f : E → R satisfying the capacity constraints. For theuncapacitated problem, it simply means f ≥
0. Let ρ f ( i ) := X j : ji ∈ E f ji − X j : ij ∈ E f ij , (2)7enote the flow balance at node i , and let Ex ( f ) = Ex b ( f ) := X i ∈ V max { ρ f ( i ) − b i , } denote the total positive excess. For an arc set F , let ←− F := { ji : ij ∈ F } denote the set of reversearcs, and let ←→ F = F ∪ ←− F . We shall use the vector norms || x || ∞ = max | x i | and || x || = P | x i | .Following [18] and [22], we do not require the functions C ij to be given explicitly, but assumeoracle access only. Oracle 1.
We are given an oracle, that we will refer to as the differential oracle , satisfying eitherof the following properties.(a) For every arc ij ∈ E , the oracle returns the value C ′ ij ( α ) in O (1) time for every α ∈ R . If α isrational then C ′ ij ( α ) is also rational.(b) For every arc ij ∈ E , the oracle returns the value e C ′ ij ( α ) in O (1) time for every α ∈ R . If α isrational then e C ′ ij ( α ) is also rational. These two options are tailored to main the applications. The more natural Oracle 1(a) holds forquadratic objectives, where C ′ ij ( α ) = 2 c ij α + d ij for the cost function C ij ( α ) = c ij α + d ij . Option(b) is needed for Fisher markets, where C ′ ij ( α ) = log α for cost functions of the form C ij ( α ) = α (log α − C ′ ij ( α ) = − log U ij for the other type of cost function, C ij ( α ) = − α log U ij , for arational U ij . Note that we do not assume an evaluation oracle returning C ij ( α ) or e C ij ( α ) - thesevalues are not needed for the algorithm.The next assumption slightly restricts the class of functions C ij for technical reasons.Each cost function C ij ( α ) is either linear or strictly convex, that is, C ′ ij ( α ) is either constant or strictly monotone increasing. ( ⋆ )Arcs with C ij ( α ) linear are called linear arcs , the rest is called nonlinear arcs . Let m L and m N denote their numbers, respectively. We use the terms linear and nonlinear for the correspondingreverse arcs as well. If C ij does not satisfy this restriction, R can be decomposed into intervalssuch that C ′ ij is either constant or strictly monotone increasing on each interval. We can replace ij by a set of paths of length two (to avoid adding parallel arcs) with appropriately chosen capacitiesand cost functions all of which satisfy the assumption. Indeed, the piecewise linear utility functionsin Fisher markets with spending constraint utilities will be handled in a similar way. If the costfunctions are explicitly given, for example, the slope of every linear segment is part of the input,then the size of the resulting network still only depends on the input size (that includes all numbersin the input). Hence a strongly polynomial algorithm on this instance will be strongly polynomialwith respect to the original instance as well. This does not hold however if the functions C ij isgiven in some different, implicit way. ∆ -feasibility Given a pseudoflow f , let us define the residual graph E f as E f := E ∪ { ij : ji ∈ E, f ij > } . Arcs in E are called forward arcs , and those in the second set backward arcs . Recall our assumptionthat the graph contains no pairs of oppositely directed arcs, hence the backward arcs are not8ontained in E . We will use the convention that on a backward arc ji , f ji = − f ij , and C ji ( α ) = C ij ( − α ), also convex and differentiable. The residual capacity is ∞ on forward arcs and f ij on thebackward arc ji .The Karush-Kuhn-Tucker conditions assert that the feasible solution f to (P) is optimal if andonly if there exists a potential vector π : V → R such that π j − π i ≤ C ′ ij ( f ij ) ∀ ij ∈ E f . (3)This is equivalent to asserting that the residual graph contains no negative cost directed cycleswith respect to the cost function C ′ ij ( f ij ).For a value ∆ >
0, let E f (∆) = E ∪ { ij : ji ∈ E, f ij ≥ ∆ } denote the subset of arcs in E f that have residual capacity at least ∆. We say that the pseudoflow f is ∆ -feasible , if there exists a potential vector π : V → R such that π j − π i ≤ C ′ ij ( f ij + ∆) ∀ ij ∈ E f (∆) . (4)Equivalently, f is ∆-feasible if and only if E f (∆) contains no negative cycles with respect to thecost function C ′ ij ( f ij + ∆). If ji is a reverse arc, then (4) gives C ′ ij ( f ij − ∆) ≤ π j − π i .We note that our notion is different (and weaker) than the analogous conditions in [26] and in[18], where ( C ij ( f ij + ∆) − C ij ( f ij )) / ∆ is used in the place of C ′ ij ( f ij + ∆). Algorithm 1Subroutine
Adjust (∆ , ¯ f ) INPUT
A 2∆-feasible pseudoflow ¯ f and a potential vector π satisfying (4) with ¯ f and 2∆. OUTPUT
A ∆-feasible pseudoflow f such that π satisfies (4) with f and ∆. for all ij ∈ E doif C ′ ij ( ¯ f ij + ∆) < π j − π i then f ij ← ¯ f ij + ∆. elseif ¯ f ji ≥ ∆ and π j − π i < C ′ ij ( ¯ f ij − ∆) then f ij ← ¯ f ij − ∆. else f ij ← ¯ f ij . return f .The subroutine Adjust (∆ , f ) (see Algorithm 1) transforms a 2∆-feasible pseudoflow to a ∆-feasible pseudoflow by possibly changing the value of every arc by ± ∆. Lemma 2.1.
The subroutine
Adjust ( ∆ , f ) is well-defined and correct: it returns a ∆ -feasiblepseudoflow with ( f, π ) satisfying (4). Further, Ex ( f ) ≤ Ex ( ¯ f ) + m N ∆ (recall that m N is thenumber of nonlinear arcs).Proof. First we observe that the “if” and “elseif” conditions cannot hold simultaneously: C ′ ij ( ¯ f ij +∆) < π j − π < C ′ ij ( ¯ f ij − ∆) would contradict the convexity of C ij . Consider the potential vector π satisfying (4) with ¯ f and 2∆. We prove that π satisfies (4) with f and ∆ as well.First, take a forward arc ij ∈ E with C ′ ij ( ¯ f ij + ∆) < π j − π i . By 2∆-feasibility we know π j − π i ≤ C ′ ij ( ¯ f ij + 2∆). These show that setting f ij = ¯ f ij + ∆ satisfies (4) for both ij and ji , using C ′ ij ( f ij − ∆) ≤ C ′ ij ( f ij ) = C ′ ij ( ¯ f ij + ∆) < π j − π i ≤ C ′ ij ( ¯ f ij + 2∆) = C ′ ij ( f ij + ∆) . Next, assume ¯ f ji ≥ ∆ and π j − π i < C ′ ij ( ¯ f ij − ∆). Note that f ij satisfies (4) by π j − π i
0. That is, for output f , thereis an optimal solution f ∗ such that k f − f ∗ k ∞ < ε . Algorithm 2Algorithm
Basic f ← ; ∆ ← ∆ ; do // ∆ -phase do //main part S (∆) ← { i ∈ V : ρ f ( i ) − b i ≥ ∆ } ; T (∆) ← { i ∈ V : ρ f ( i ) − b i ≤ − ∆ } ; P ← shortest s − t path in E f (∆) for the cost C ′ ij ( f ij + ∆) with s ∈ S (∆), t ∈ T (∆);send ∆ units of flow on P from s to t ; while S (∆) , T (∆) = ∅ ; Adjust (∆ / , f );∆ ← ∆ / while ∆ > ε/ (2 n + m N + 1); Return f .We start with the pseudoflow f ≡ and an initial value ∆ = ∆ . We assume that the value∆ is provided in the input so that is a ∆ -feasible and Ex ( ) ≤ (2 n + m )∆ ; in the enhancedalgorithm we shall specify how such a ∆ value can be determined. The algorithm consists of∆-phases, with ∆ decreasing by a factor of two between two phases; the algorithm terminates once∆ < ε/ (2 n + m N + 1).In the main part of phase ∆, let S (∆) = { i ∈ V : ρ f ( i ) − b i ≥ ∆ } and T (∆) = { i ∈ V : ρ f ( i ) − b i ≤ − ∆ } , the set of nodes with excess and deficiency at least ∆. As long as S (∆) = ∅ , T (∆) = ∅ , send ∆ units of flow from a node s ∈ S (∆) to a node t ∈ T (∆) on a shortest path in E f (∆) with respect to the cost function C ′ ij ( f ij + ∆). (Note that there must be a path connectingnodes in S (∆) and T (∆), due to our assumption that the graph G = ( V, E ) is strongly connected,and E ⊆ E f (∆).)The main part finishes once S (∆) = ∅ or T (∆) = ∅ . The ∆-phase terminates by performing Adjust (∆ / , f ) and proceeding to the next phase with scaling factor ∆ / E f (∆) for the cost function10 ′ ij ( f ij + ∆). This can be done only if there is no negative cost cycle. ∆-feasibility is exactlythis property and is maintained throughout (see Lemma 3.2 below). Details of the shortest pathcomputation will be given in Section 5.1(ii), for the enhanced algorithm. Theorem 3.1.
The Basic algorithm delivers an ε -accurate solution in O (log((2 n + m N + 1)∆ /ε ) phases, and every phase comprises at most O (2 n + m N ) flow augmentations. An appropriate ∆ can be chosen to be polynomial in the input size, hence this gives a weaklypolynomial running time bound. We now state the two simple lemmas needed to prove this theorem.The first lemma verifies the correctness and efficiency of the algorithm, showing that ∆-feasibilityis maintained throughout and the number of flow augmentations is linear in every ∆-phase. Weomit the proof; its analogous counterpart for the enhanced algorithm will be proved in Lemma 5.1. Lemma 3.2. (i) In the main part of the ∆ -phase, the pseudoflow is an integer multiple of ∆ oneach arc, and consequently, E f (∆) = E f .(ii) ∆ -feasibility is maintained when augmenting on a shortest path.(iii) At the beginning of the main part, Ex ( f ) ≤ (2 n + m N )∆ , and at the end, Ex ( f ) ≤ n ∆ .(iv) The main part consists of at most n + m N flow augmentation steps. Our second lemma asserts the proximity of a current flow to all later flows during the algorithm.If we let the algorithm run without ever terminating, it will converge to an optimal solution. Hencethe lemma justifies that the algorithm obtains an ε -accurate solution as claimed in Theorem 3.1.Moreover, it also helps to identify edges which must be contained in the support of an optimalsolution. The proof is also omitted; see Lemma 5.2 and the first part of the proof of Theorem 5.5.This is essentially the same argument that was used by Orlin (e.g. [1, Lemma 10.21]). Lemma 3.3.
Let f be the pseudoflow at the end of the main part of the ∆ -phase and f ′ in anarbitrary later phase. Then || f − f ′ || ∞ ≤ (2 n + m + 1)∆ . If f ij > (2 n + m + 1)∆ at the end of the ∆ -phase, then this property is maintained in all later phases, and there exists an optimal solution f ∗ with f ∗ ij > . For all such arcs, we can conclude π j − π i = C ′ ij ( f ∗ ij ) for an optimal solution f ∗ . It will belong tothe set of revealed arcs, defined in the next section. The overall aim of the algorithm is to identify alarge enough set of revealed arcs containing the support of an optimal solution. The above lemmaguarantees that the first such arc can be identified in a strongly polynomial number of steps in theBasic algorithm. We will however need to modify the algorithm in order to guarantee that the setof revealed arcs is always extended in a strongly polynomial number of steps. Let F ∗ denote the set of arcs that are tight in every optimal solution (note that in general, we donot assume the uniqueness of the optimal solution). This arc set plays a key role in our algorithm.Formally, F ∗ := { ij ∈ E : π j − π i = C ′ ij ( f ij ) holds ∀ f optimal to (P), ∀ π : V → R , s.t. ( f, π ) satisfies the inequalities (3) } . (5)11he next lemma shows that F ∗ contains the support of every optimal solution. Lemma 4.1.
Let f be an arbitrary optimal solution to (P), and f ij > for some ij ∈ E . Then ij ∈ F ∗ . The proof needs the following notion, also used later. Let x, y : E → R be two vectors. Let usdefine the difference graph D x,y = ( V, E x,y ) with ij ∈ E x,y if ij ∈ E and x ij > y ij or if ji ∈ E and x ji < y ji . Using the convention x ji = − x ij , y ji = − y ij it follows that x ij > y ij for every ij ∈ E x,y .We will need the following simple claim. Claim 4.2.
Assume that for two vectors x, y : E → R , ρ x = ρ y holds (recall the definition of ρ in(2)). Then every arc in the difference graph E x,y must be contained in a cycle in E x,y .Proof. For ij ∈ E x,y , let us set z ij = x ij − y ij if x ij > y ij . The assumption ρ x = ρ y implies that z ij is a circulation in E x,y with positive value on every arc. As such, it can be written as a nonnegativecombination of incidence vectors of cycles. Therefore every ij ∈ E x,y must be contained in acycle. Proof of Lemma 4.1.
Let f ∗ be another arbitrary optimal solution, and consider potentials π and π ∗ with both ( f, π ) and ( f ∗ , π ∗ ) satisfying (3). We shall prove that π ∗ j − π ∗ i = C ′ ij ( f ∗ ij ). Since ( f ∗ , π ∗ )is chosen arbitrarily, this will imply ij ∈ F ∗ . If f ∗ ij >
0, then ji ∈ E f ∗ and thus π ∗ j − π ∗ i = C ′ ij ( f ∗ ij )must hold.Assume now f ∗ ij = 0. Consider the difference graph D f,f ∗ . Since f ij > f ∗ ij , it follows that ij ∈ E f,f ∗ . Because of ρ f ∗ ≡ ρ f ≡ b , Claim 4.2 is applicable and provides a cycle C in E f,f ∗ containing ij . For every arc ab ∈ C , f ab > f ∗ ab and thus ab ∈ E f ∗ and ba ∈ E f . By (3),0 = X ab ∈ C π ∗ b − π ∗ a ≤ X ab ∈ C C ′ ab ( f ∗ ab ) and0 = X ab ∈ C π a − π b ≤ X ab ∈ C C ′ ba ( f ba ) = − X ab ∈ C C ′ ab ( f ab ) . The convexity of C ab and f ab > f ∗ ab give C ′ ab ( f ab ) ≥ C ′ ab ( f ∗ ab ). In the above inequalities, equalitymust hold everywhere, implying π ∗ j − π ∗ i = C ′ ij ( f ∗ ij ) as desired.We shall see that using Oracle 2 (to be described later), finding the set F ∗ enables us tocompute an optimal solution in strongly polynomial time. In the Basic algorithm, F = { ij ∈ E : f ij > (2 n + m + 1)∆ } is always a subset of F ∗ according to Lemmas 3.3 and 4.1. Furthermore, oncean edge enters F , it stays there in all later phases. The Enhanced algorithm provides a modificationof the basic algorithm with the guarantee that within every O (log n ) phases, a new arc enters F .In each step of the enhanced algorithm, there will be an arc set F , called the revealed arc set ,which is guaranteed to be a subset of F ∗ . We remove the lower capacity 0 from arcs in F and allowalso negative values here.Formally, for an edge set F ⊆ E , a vector f : E → R is an F -pseudoflow , if f ij ≥ ij ∈ E \ F (but it is allowed to be negative on F ). For such an f , let us define E Ff := E f ∪ ←− F = E ∪ ←− F ∪ { ji : ij ∈ E \ F, f ij > } . (6)If ij ∈ F , then the residual capacity of ji is ∞ . In every phase of the algorithm, we maintain an F -pseudoflow f for a revealed arc set F ⊆ F ∗ .Provided the revealed arc set F ⊆ F ∗ , we will aim for F -optimal solutions as defined below; weprove that finding an F -optimal solution is essentially equivalent to finding an optimal one. We12ay that f : E → R is F -optimal , if it is an F -pseudoflow with ρ f ≡ b and there exists a potentialvector π : V → R with π j − π i ≤ C ′ ij ( f ij ) ∀ ij ∈ E Ff . (7)This is stronger than the optimality condition (3) in that it also requires the inequality on arcs in ←− F . On the other hand, it does not imply optimality as it allows f ij < ij ∈ F . Nevertheless,it is easy to see that every optimal solution f ∗ is also F -optimal for every F ⊆ F ∗ . This is dueto the definition of F ∗ as the set of arcs satisfying π j − π i = C ′ ij ( f ij ) whenever ( f, π ) satisfies(3). Conversely, we shall prove that provided an F -optimal solution, we can easily find an optimalsolution by a single feasible circulation algorithm, a problem equivalent to maximum flows (see [1,Chapters 6.2, 7]). Lemma 4.3.
Assume that for a subset F ⊆ F ∗ , an F -optimal solution f is provided. Then anoptimal solution to (P) can be found by a feasible circulation algorithm. Further, ij ∈ F ∗ whenever f ij > .Proof. Assume f and ¯ f are both F -optimal solutions, that is, for some vectors π and ¯ π , the pairs( f, π ) and ( ¯ f , ¯ π ) both satisfy (7). We prove that (i) f ij = ¯ f ij whenever ij is a nonlinear arc; and (ii) if ij is a linear arc with f ij = ¯ f ij , then π j − π i = C ′ ij ( f ij ) = C ′ ij ( ¯ f ij ) = ¯ π j − ¯ π i .Note that (i) and (ii) immediately imply the second half of the claim as it can be applied for f and an arbitrary optimal (and consequently, F -optimal) solution ¯ f .The proof uses the same argument as for Lemma 4.1. W.l.o.g. assume f ij > ¯ f ij for an arc ij ,and consider the difference graph D f, ¯ f . Since ρ f ≡ ρ ¯ f ≡ b and f ij > ¯ f ij , Claim 4.2 is applicableand shows that ij must be contained on a cycle C ⊆ E f, ¯ f . For every arc ab ∈ C , ab ∈ E F ¯ f and ba ∈ E Ff follows (using ←→ F ⊆ E F ¯ f ∩ E Ff ). By (7),0 = X ab ∈ C ¯ π b − ¯ π a ≤ X ab ∈ C C ′ ab ( ¯ f ab ) and0 = X ab ∈ C π a − π b ≤ X ab ∈ C C ′ ba ( f ba ) = − X ab ∈ C C ′ ab ( f ab ) . Now convexity yields C ′ ab ( f ab ) = C ′ ab ( ¯ f ab ) for all ab ∈ C . The condition ( ⋆ ) implies that all arcs in C are linear, in particular, ij is linear. This immediately proves (i) . To verify (ii) , observe that allabove inequalities must hold with equality.This suggests the following simple method to transform an F -optimal solution f to an optimal f ∗ of (P). For every nonlinear arc ij , we must have f ∗ ij = f ij . Let H ⊆ E be the set of linear arcssatisfying π j − π i = C ′ ij ( f ij ). Consider the solutions h of the following feasible circulation problem: h ij = f ij ∀ ij ∈ E \ H X j : ji ∈ E h ji − X j : ij ∈ E h ij = b i ∀ i ∈ Vh ≥ f ∗ is an optimal solution, then (i) and (ii) imply that f ∗ ij = f ij for all ij ∈ E \ H and ij ∈ H for every arc with f ij = f ∗ ij . The degree conditions are satisfied because of ρ f ∗ ≡ ρ f ≡ b .Conversely, every feasible circulation h is an optimal solution to (P), since ( h, π ) satisfies (3).13n every step of our algorithm we will have a scaling parameter ∆ ≥ F ⊆ F ∗ . The Basic algorithm used the notion of ∆-feasibility; it has to be modified according to F . Let E Ff (∆) denote the set of arcs in E Ff with residual capacity at least ∆. That is, E Ff (∆) := E f (∆) ∪ ←− F = E ∪ ←− F ∪ { ji : ij ∈ E \ F, f ij ≥ ∆ } . (8)We say that the F -pseudoflow f is (∆ , F ) -feasible , if there exists a potential vector π : V → R sothat π j − π i ≤ C ′ ij ( f ij + ∆) ∀ ij ∈ E Ff (∆) . (9)This is equivalent to the property that E Ff (∆) contains no negative cycle with respect to the costfunction C ′ ij ( f ij + ∆).In accordance with (∆ , F )-feasibility, we have to modify the subroutine Adjust . The modifiedsubroutine, denoted by
Adjust (∆ , f, F ) is shown in Algorithm 3. The only difference from Algo-rithm 1 is that the condition (4) is replaced by (9), and that in the second condition, “ ¯ f ji ≥ ∆”is replaced by “ ¯ f ji ≥ ∆ or ij ∈ F ”. The following lemma can be proved by the same argument asLemma 2.1. Algorithm 3Subroutine
Adjust (∆ , ¯ f , F ) INPUT
A (2∆ , F )-feasible pseudoflow ¯ f and a potential vector π satisfying (9) with ¯ f and 2∆. OUTPUT
A (∆ , F )-feasible pseudoflow f such that π satisfies (9) with f and ∆. for all ij ∈ E doif C ′ ij ( ¯ f ij + ∆) < π j − π i then f ij ← ¯ f ij + ∆. elseif ( ¯ f ji ≥ ∆ or ij ∈ F ) and π j − π i < C ′ ij ( ¯ f ij − ∆) then f ij ← ¯ f ij − ∆. else f ij ← ¯ f ij . return f . Lemma 4.4.
The subroutine
Adjust (∆ , f, F ) is well-defined and correct: it returns a (∆ , F ) -feasible pseudoflow with ( f, π ) satisfying (9). Further, Ex ( f ) ≤ Ex ( ¯ f ) + m N ∆ . Finally, we say that a set F ⊆ E is linear acyclic , if F does not contain any undirected cyclesof linear arcs (that is, no cycle in F may consist of linear arcs and their reverse arcs). We shallmaintain that the set of revealed arcs, F is linear acyclic.This notion is motivated by the following: assume there exists a cycle consisting of linear arcsand their reverses. Given an F -pseudoflow, we could modify it by sending an arbitrary amount offlow around this cycle. Hence we would not be able to derive our proximity result Lemma 5.6 andLemma 5.4 that relies on it. On the other hand, we can pick an arbitrary arc on a cycle of lineararcs, remove it from F , an reroute its entire flow on the rest of the cycle. Given the set F ⊆ F ∗ of revealed arcs, we will try to find out whether F already contains thesupport of an optimal solution. This motivates the following definition. We say that the (notnecessarily nonnegative) vector x : E → R is F -tight , if x ij = 0 whenever ij / ∈ F and there existsa potential vector π : V → R with π j − π i = C ′ ij ( x ij ) ∀ ij ∈ F. (10)14or example, any optimal solution is F ∗ -tight by Lemma 4.1. Notice that an F -tight vector f isnot necessarily F -optimal as (7) might be violated for edges in E Ff \ ←→ F and also since Ex b ( f ) > ρ f ≡ b is equivalent to Ex b ( f ) = 0). Conversely, an F -optimal vector is notnecessarily F -tight as it can be nonzero on E \ F .Given F and some node demands ˆ b : V → R , we would like to find an F -tight x with Ex ˆ b ( x ) = 0.This is equivalent to finding a feasible solution ( x, π ) to the following system: π j − π i = C ′ ij ( x ij ) ∀ ij ∈ F X j : ji ∈ F x ji − X j : ij ∈ E x ij = ˆ b i ∀ i ∈ V (11) x ij = 0 ∀ ij ∈ E \ F Let us define the discrepancy D ˆ b ( F ) of F as the maximum of | P i ∈ K ˆ b i | over undirected connectedcomponents K of F . A trivial necessary condition for solvability is D ˆ b ( F ) = 0: indeed, summingup the second set of equalities for a component K , we obtain 0 = P i ∈ K ˆ b i . Oracle 2.
Assume we have a subroutine
Trial ( F, ˆ b ) so that for any linear acyclic F ⊆ E and anyvector ˆ b : V → R satisfying D ˆ b ( F ) = 0 , it delivers an F -tight solution x to (11) with ρ x ≡ ˆ b instrongly polynomial running time ρ T ( n, m ) . For quadratic cost functions and also for Fisher markets, this subroutine can be implementedby solving simple systems of equations (for quadratic, this was already outlined in Section 1.2).Consider now an F -tight vector f , and let err F ( f ) := inf { ∆ : f is (∆ , F )-feasible } . (12)Recall the definition (8) of the edge set E Ff (∆). As f is assumed to be F -tight and therefore f ij > ij ∈ F , we get that E Ff (∆) = E ∪ ←− F . Consequently, E Ff (∆) is independent of the value of∆. Because of continuity, this infimum is actually a minimum whenever the set is nonempty. If f is not (∆ , F )-feasible for any ∆, then let err F ( f ) = ∞ . f is F -optimal if and only if f is a feasibleflow (that is, Ex b ( f ) = 0) and err F ( f ) = 0. Oracle 3.
Assume a subroutine
Error ( f, F ) is provided, that returns err F ( f ) for any F -tightvector f in strongly polynomial running time ρ E ( n, m ) . Further, if err ∅ ( ) = ∞ , then (P) isunbounded. This subroutine seems significantly harder to implement for the applications: we need to solvea minimum cost-to-time ratio cycle problem for quadratic costs and all pairs shortest paths for theFisher markets.Having formulated all necessary assumptions, we are finally in the position to formulate themain result of the paper.
Theorem 4.5.
Assume Oracles 1-3 are provided and ( ⋆ ) holds for the problem (P) in a network on n nodes and m arcs, m N among them having nonlinear cost functions. Let ρ T ( n, m ) and ρ E ( n, m ) denote the running time of Oracle 2 and Oracle 3, and let ρ S ( n, m ) be the running time needed fora single shortest path computation for nonnegative arc lengths. Then an exact optimal solution canbe found in O (( n + m N )( ρ T ( n, m ) + ρ E ( n, m )) + ( n + m N ) ρ S ( n, m ) log m ) time.This gives an O ( m log m ) algorithm for quadratic convex objectives. For Fisher markets, weobtain O ( n + n ( m + n log n ) log n ) running time for linear and O ( mn + m ( m + n log n ) log m ) for spending constraint utilities. .3 Description of the enhanced algorithm Algorithm 4Algorithm
Enhanced Convex FlowError ( , ∅ ); f ← ; ∆ ← max { err ∅ ( ) , Ex b ( ) / (2 n + m N ) } ; F ← ∅ ; repeat // ∆ -phase do //main part S (∆) ← { i ∈ V : ρ f ( i ) − b i ≥ ∆ } ; T (∆) ← { i ∈ V : ρ f ( i ) − b i ≤ − ∆ } ; P ← shortest s − t path in E Ff (∆) for the cost C ′ ij ( f ij + ∆) with s ∈ S (∆), t ∈ T (∆);send ∆ units of flow on P from s to t ; while S (∆) , T (∆) = ∅ ; Extend (∆ , f, F ); if ( F was extended) and ( D b ( F ) ≤ ∆) then Trial-and-Error ( F ) else Adjust (∆ / , f, F );∆ ← ∆ / Subroutine
Extend (∆ , f, F ) for all ij ∈ E \ F , f ij > (2 n + m + 1)∆ doif F ∪ { ij } is linear acyclic then F ← F ∪ { ij } else P ← path of linear arcs in ←→ F between i and j ;send f ij units of flow on P from i to j ; f ij ← f = , ∆ = max { err ∅ ( ) , Ex b ( ) / (2 n + m N ) } and F = ∅ . Thealgorithm consists of ∆-phases. In the ∆-phase, we shall maintain a linear acyclic revealed arc set F ⊆ F ∗ , and a (∆ , F )-feasible F -pseudoflow f . The algorithm will always terminate during thesubroutine Trial-and-Error .The main part of the ∆-phase is the same as in the Basic algorithm. Let S (∆) = { i ∈ V : ρ f ( i ) − b i ≥ ∆ } and T (∆) = { i ∈ V : ρ f ( i ) − b i ≤ − ∆ } . As long as S (∆) = ∅ , T (∆) = ∅ , send ∆units of flow from a node s ∈ S (∆) to a node t ∈ T (∆) on a shortest path in E Ff (∆) with respectto the cost function C ′ ij ( f ij + ∆). (The existence of such a path P is guaranteed by our assumptionthat the graph G = ( V, E ) is strongly connected.)After the main part (the sequence of path augmentations) is finished, the subroutine
Ex-tend (∆ , f, F ) adds new arcs ij ∈ E \ F with f ij > (2 n + m + 1)∆ to F maintaining the linearacyclic property. This is achieved as follows: we first add all nonlinear such arcs to F . We add alinear arc to F if it does not create any (undirected) cycles in F . If adding the linear arc ij wouldcreate a cycle, we do not include it in F , but reroute the entire flow from ij using the (undirected)path in F between i and j .If no new arc enters F , then we perform Adjust (∆ / , f, F ) and move to the next scaling phasewith the same f and set the scaling factor to ∆ /
2. This is done also if F is extended, but it has ahigh discrepancy: D b ( F ) > ∆.Otherwise, the subroutine Trial-and-Error ( F ) determines the next f and ∆. Based on thearc set F , we find a new F -pseudoflow f and scaling factor at most ∆ /
2. The subroutine may also16erminate with an F -optimal solution, which enables us to find an optimal solution to (P) by amaximum flow computation due to Lemma 4.3. Theorem 5.5 will show that this is guaranteed tohappen within a strongly polynomial number of steps. The Trial-and-Error subroutine
The subroutine assumes that the discrepancy of F is small: D b ( F ) ≤ ∆. Step 1.
First, modify b to ˆ b : in each (undirected) component K of F , pick a node j ∈ K andchange b j by − P i ∈ K b i ; leave all other b i values unchanged. Thus we get a ˆ b with D ˆ b ( F ) = 0. Trial ( F, ˆ b ) returns an F -tight vector ˆ f . Step 2.
Call the subroutine
Error ( ˆ f , F ). If b = ˆ b and err F ( ˆ f ) = 0, then ˆ f is F -optimal.An optimal solution to (P) can be found by a single maximum flow computation, as described inthe proof of Lemma 4.3. In this case, the algorithm terminates. If err F ( ˆ f ) ≥ ∆ /
2, then keep theoriginal f , perform Adjust (∆ / , f, F ) and go to the next scaling phase with scaling factor ∆ / f = ˆ f and define the next scaling factor as∆ next = max { err F ( ˆ f ) , Ex b ( ˆ f ) / (2 n + m N ) } . The details how the shortest path computations are performed will be discussed in Section 5.1; inthe following analysis, we assume it can be efficiently implemented. At the initialization, err ∅ ( )must be finite or the problem is unbounded as assumed in Oracle 3. Trial-and-Error replaces f by ˆ f if err F ( ˆ f ) ≤ ∆ / f otherwise. Thefirst case will be called a successful trial, the latter is unsuccessful . The following is (an almostidentical) counterpart of Lemma 3.2. Lemma 5.1. (i) In the main part of the ∆ -phase, the F -pseudoflow f is an integer multiple of ∆ on each arc ij ∈ E \ F , and consequently, E Ff (∆) = E Ff .(ii) (∆ , F ) -feasibility is maintained in the main part and in subroutine Extend (∆ , f, F ) .(iii) At the beginning of the main part, Ex ( f ) ≤ (2 n + m N )∆ , and at the end, Ex ( f ) ≤ n ∆ .(iv) The main part consists of at most n + m N flow augmentation steps.(v) The scaling factor ∆ decreases by at least a factor of 2 between two ∆ -phases.Proof. For (i) , f is zero on every arc in E \ F at the beginning of the algorithm and after everysuccessful trial. In every other case, the previous phase had scaling factor 2∆, and thus by induction,the flow is an integer multiple of 2∆ at the end of the main part of the 2∆-phase, a propertyalso maintained by Extend (2∆ , f, F ). The 2∆-phase finishes with
Adjust (∆ , f, F ), possiblymodifying the flow on every arc by ± ∆. In the main part of the ∆-phase, the shortest pathaugmentations also change the flow by ± ∆. This implies E Ff (∆) = E Ff .For (ii) , P is a shortest path if there exists a potential π satisfying (9) with π j − π i = C ′ ij ( f ij +∆)on each arc ij ∈ P (see also Section 5.1). We show that when augmenting on the shortest path P , (9) is maintained with the same π . If ij, ji / ∈ P , then it is trivial as the flow is left unchangedon ij . Consider now an arc ij ∈ P ; the next argument applies both if ij is a forward or a reversearc. The new flow value will be f ij + ∆, hence we need π j − π i ≤ C ′ ij ( f ij + 2∆), obvious as C ′ ij is monotonely increasing. We next verify (9) for the backward arc ji ∈ E Ff (∆). This gives17 i − π j ≤ C ′ ji (( f ji − ∆) + ∆), that is equivalent to C ′ ij ( f ij ) ≤ π j − π i , again a consequence ofmonotonicity.In subroutine Extend , we reroute the flow f ij from a linear arc ij if ←→ F contains a directedpath P from i to j . This cannot affect feasibility since the C ′ ij ’s are constant on linear arcs. Alsonote that arcs in ←→ F have infinite residual capacities.For (iii) , Ex ( f ) ≤ n ∆ as the main part terminates with either S (∆) = ∅ or T (∆) = ∅ .Lemma 4.4 shows that Adjust (∆ / , f, F ) increases the excess by at most m N ∆ /
2. Consequently, Ex ( f ) ≤ (2 n + m N )(∆ /
2) at the beginning of the ∆ / next . By definition, the newexcess is at most (2 n + m N )∆ next .Further, (iii) implies (iv) , as each flow augmentation decreases Ex ( f ) by ∆. Finally (v) isstraightforward if the next value of the scaling factor is set as ∆ /
2. This is always the case, exceptif
Trial-and-Error is called and err F ( ˆ f ) ≤ ∆ /
2, when the next scaling factor is set as themaximum of err F ( ˆ f ) and Ex b ( ˆ f ) / (2 n + m N ). We show that this second term is also at most ∆ / f was obtained by Trial ( F, ˆ b ), and therefore ρ ˆ f ( i ) − b i = ˆ b i − b i ≤ ∆ due to the definitionof ˆ b and D b ( F ) ≤ ∆. It follows that Ex b ( ˆ f ) ≤ n ∆, and thus Ex b ( ˆ f ) / (2 n + m N ) < ∆ / Lemma 5.2. F ⊆ F ∗ holds in each step of the algorithm.Proof. The proof is by induction. A new arc ij may enter F if f ij > (2 n + m + 1)∆ after themain part of the ∆-phase. We shall prove that f ∗ ij > F -optimal solution f ∗ , and thusLemma 4.3 gives ij ∈ F ∗ .After the phase when ij entered, let us continue with the following modified algorithm: do notextend F and do not perform Trial-and-Error anymore, but always choose the next scalingfactor as ∆ /
2, and keep the algorithm running forever. (This is almost the same as the Basicalgorithm, with the difference that we have a revealed arc set F .)Let ∆ = ∆ and ∆ t = ∆ / t denote the scaling factor in the t ’th phase of this algorithm (withphase 0 corresponding to the ∆-phase). Consider any ∆ t -phase ( t ≥ n + m N )∆ t during the main part by Lemma 5.1(iv) and by ∆ t / Adjust (∆ t / , f, F ),amounting to a total modification ≤ (2 n + m N + )∆ t . Consequently, the total modification in the∆ t phase and all later phases is bounded by (2 n + m N + ) P ∞ k = t ∆ k ≤ n + m + )∆ t .We may conclude that when running forever, the flow f converges to an F -optimal solution f ∗ .Indeed, let f ( t ) denote the F -pseudoflow at the end of the t ’th phase. By the above observation, || f ( t ) − f ( t ′ ) || ∞ ≤ n + m + )∆ t for any t ′ ≥ t ≥
0. Consequently, on every arc ij ∈ E , thesequence f ( t ) ij converges; let f ∗ denote the limit. We claim the f ∗ is F -optimal.Firstly, f ∗ is clearly an F -pseudoflow. Property (7) is equivalent to the property that E Ff doesnot contain any negative cycle w.r.t. C ′ ij ( f ij ). This follows from the fact that E Ff (∆ t ) does notcontain any negative cycle w.r.t. C ′ ij ( f ( t ) ij ) due to the (∆ t , F )-feasibility of f ( t ) . Finally, Ex b ( f ∗ ) =lim t →∞ Ex b ( f ( t ) ) ≤ lim t →∞ n ∆ t = 0, and therefore Ex b ( f ∗ ) = 0.To finish the proof, we observe that f ∗ ij >
0. Indeed, f ij > (2 n + m + 1)∆ after the mainpart of the ∆-phase, and hence f ij > (2 n + m + )∆ at the end of the ∆-phase (after performing Adjust (∆ / , f, F )). By the above argument, the total change in all later phases is ≤ n + m + )∆ = (2 n + m + )∆, yielding the desired conclusion.Recall the characterization of arcs to free and restricted. Free arcs are differentiable on theentire R , whereas for a restricted arc ij , we have C ′ ij ( α ) = −∞ for α <
0. Therefore we have toavoid the flow value becoming negative even if ij ∈ F for a restricted arc.18 laim 5.3. f ij ≥ holds for every restricted arc ij during the entire algorithm, even if ij ∈ F .Proof. f ij ≥ f ij < Adjust subroutine (
Extend may not modify f ij as ij is a nonlinear arc). In case of a path augmentation, ji is contained on the shortest path P , and therefore π j − π i = C ′ ij ( f ij − ∆) must hold for a potential π (see the proof of Lemma 5.1).This is a contradiction as f ij − ∆ < C ′ ij ( f ij − ∆) = −∞ . A similar argument works for Adjust . Lemma 5.4.
When
Trial-and-Error ( F ) is performed in the ∆ -phase, err F ( ˆ f ) ≤ m + 1) ∆ holds. This lemma is of key importance. Before proving it, we show how it provides the stronglypolynomial bound. The main idea is the following: in
Trial-and-Error ( F ), we replace f by ˆ f and ∆ by a new value instead of ∆ / err F ( ˆ f ) < ∆ /
2; otherwise, we ignore ˆ f and proceedto the next phase as usual. Whereas err F ( ˆ f ) ≥ ∆ / Trial-and-Error ( F ) depends only onthe revealed arc set F . Consequently, if we had err F ( ˆ f ) ≥ ∆ /
2, then by the time the scaling factorreduces to a smaller value ∆ ′ such that 6( m + 1) ∆ ′ < ∆ /
2, the set F must have been extended. Theorem 5.5.
The enhanced algorithm terminates in at most O (( n + m N ) log m ) scaling phases.Proof. The set of revealed arcs can be extended at most m N + n − n −
1) linear arcs because of the linear acyclic property. We shall show that after any∆-phase, a new arc is revealed within 2 ⌈ log T ⌉ phases, for T = 24( m + 1) .As ∆ decreases by at least a factor of two between two phases, after ⌈ log T ⌉ steps we have∆ T ≤ ∆ /T . Assume that in the ∆ T phase, we still have the same revealed arc set F as in the∆-phase. Case I. D b ( F ) > ∆ . At the end of the main part of the ∆ T -phase, D b ( F ) > m + 1) ∆ T . Thusthere is an undirected connected component K of F with | P i ∈ K b i | > m + 1) ∆ T . Let ρ f ( K )denote the total f value on arcs entering K minus the value on arcs leaving K , that is, ρ f ( K ) := X ij ∈ E : i/ ∈ K,j ∈ K f ij − X ij ∈ E : i ∈ K,j / ∈ K f ij . We have | ρ f ( K ) | = (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)X i ∈ K ρ f ( i ) (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) = (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)X i ∈ K ( ρ f ( i ) − b i + b i ) (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) ≥ (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)X i ∈ K b i (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) − Ex b ( f ) . The last part is derived from the simple inequality | β + α + + α − | ≥ | β |− γ , whenever α + , α − , β, γ ∈ R with − γ ≤ α − ≤ ≤ α + ≤ γ . In our setting, β = P i ∈ K b i , α + = P i ∈ K max { ρ f ( i ) − b i , } , α − = P i ∈ K min { ρ f ( i ) − b i , } , and γ = Ex b ( f ). The conditions hold since γ = Ex b ( f ) = X i ∈ V max { ρ f ( i ) − b i , } = − X i ∈ V min { ρ f ( i ) − b i , } . For the second equality, note that P i ∈ V b i = P i ∈ V ρ f ( i ) = 0. Now we may conclude | ρ f ( K ) | ≥ (cid:12)(cid:12)(cid:12)(cid:12)(cid:12)X i ∈ K b i (cid:12)(cid:12)(cid:12)(cid:12)(cid:12) − Ex b ( f ) > m + 1) ∆ T − n ∆ T > (2 n + m + 1) m ∆ T . ij entering or leaving K with f ij > (2 n + m + 1)∆ T , acontradiction as at least one such arc must have been added to F in Extend (∆ T , f, F ). Note thatthe first such arc examined during Extend (∆ T , f, F ) does keep the linear acyclic property as itconnects two separate connected components of F . Case II. D b ( F ) ≤ ∆ . We may assume that either we are at the very beginning of the algorithmwith F = ∅ , or in a phase when F just has been extended; otherwise, we could consider an earlierphase with this property. We can interpret the initial solution and initial ∆ as the output of Trial-and-Error ( ∅ ). Case IIa. D b ( F ) > ∆ T . The argument of Case I, applied for ∆ T instead of ∆, shows thatwithin ⌈ log T ⌉ phases after the ∆ T phase, F shall be extended, showing that a new arc was revealedwithin 2 ⌈ log T ⌉ phases after the ∆-phase. Case IIb. D b ( F ) ≤ ∆ T . Recall the assumption that F has not changed between phases ∆and ∆ T , and thus D b ( F ) has not changed its value either. Let us apply the analysis of the Trial-and-Error subroutine for the ∆ T -phase. (Even if the subroutine is not actually performed, itsanalysis is valid provided that D b ( F ) ≤ ∆ T .)Let ˆ f be the arc set found by Trial ( F, ˆ b ). Let us assume that b is modified to ˆ b always the sameway for the same F ; with this assumption, the output of the subroutine is the same whether calledin the ∆ or in the ∆ T -phase. In the event of an unsuccessful trial in the ∆-phase, ∆ / ≤ err F ( ˆ f ).Using Lemma 5.4 for the ∆ T -phase, err F ( ˆ f ) ≤ m + 1) ∆ T ≤ ∆ / ≤ err F ( ˆ f ) / , a contradiction. On the other hand, if we had a successful trial in the ∆-phase, then ∆ T ≤ next /T , as ∆ T is the scaling factor T − next -phase. Lemma 5.4 and Ex b ( ˆ f ) ≤ nD b ( F ) ≤ n ∆ T together yield∆ next = max { err F ( ˆ f ) , Ex b ( ˆ f ) / (2 n + m N ) } ≤ m + 1) ∆ T ≤ ∆ next / , a contradiction again.Some preparation is needed to prove Lemma 5.4. We note that the linear acyclic property isimportant due to the following lemma; if F may contains undirected cycles of linear arcs, the claimis not true. Lemma 5.6.
For a linear acylic arc set F ⊆ E , let x and y be two F -tight vectors. Then || x − y || ∞ ≤|| ρ x − ρ y || holds.Proof. First, we claim that the difference graph D x,y = ( V, E x,y ) is acyclic. Indeed, if there existeda cycle C ⊆ E x,y , then we get 0 = P ab ∈ C C ′ ab ( x ab ) = P ab ∈ C C ′ ab ( y ab ) as in the proof of Lemma 4.1.Since x ab > y ab for every ab ∈ C , this is only possible if all arcs of C are linear ( ⋆ ), contradictingthe linear acyclic property of F . (Note that E x,y ⊆ ←→ F , since by definition, every F -tight vector issupported on F ).Define the function z by z ij = x ij − y ij > ij ∈ E x,y (again with the convention x ji = − x ij , y ji = − y ij if ij ∈ E ). ρ z ≡ ρ x − ρ y , therefore we have to prove z ij ≤ || ρ z || for ij ∈ E x,y . Thisproperty indeed holds for every positive z with acyclic support.Consider a reverse topological ordering v , . . . , v n of V , where v a v b ∈ E x,y implies a > b . Forthe arc ij ∈ E x,y , let i = v t ′ and j = v t ( t ′ > t ). Let V t = { v , . . . , v t } . V t is a directed cut in E x,y ,thus X p>t ≥ q z v p v q = X p ≤ t ρ z ( v p ) . z is positive on all arcs, this implies z v a v b ≤ P p ≤ t ρ z ( v p ) ≤ || ρ z || for all such arcs, in particular,for ij . Claim 5.7. If f and ˆ f are F -pseudoflows with ˆ f ij = 0 for ij ∈ E \ F , and f is (∆ , F ) -feasible,then ˆ f is (∆ + || f − ˆ f || ∞ , F ) -feasible.Proof. There is a potential π so that f and π satisfy (9), that is, π j − π i ≤ C ′ ij ( f ij +∆) if ij ∈ E Ff (∆).For α = || f − ˆ f || ∞ , we have f ij + ∆ ≤ ˆ f ij + ∆ + α . Consequently, (9) is satisfied for ( ˆ f ij , π ) and∆ + α for every arc in E Ff (∆).By the assumption that ˆ f is zero outside F , we have E F ˆ f (∆ + α ) = E ∪ ←− F ⊆ E Ff (∆) and thusthe claim follows. Proof of Lemma 5.4.
When
Trial-and-Error is applied, f is (∆ , F )-feasible with some potential π and Ex b ( f ) ≤ n ∆. We claim that there is an F -tight ¯ f so that | ¯ f ij − f ij | ≤ ∆ for every ij ∈ F ,and Ex b ( ¯ f ) ≤ (2 n + m + 2) m ∆.Indeed, (∆ , F )-feasibility gives C ′ ij ( f ij − ∆) ≤ π j − π i ≤ C ′ ij ( f ij + ∆) ∀ ij ∈ F. If ij is a free arc (that is, differentiable on the entire R ), then C ′ ij is continuous, so there must bea value f ij − ∆ ≤ β ≤ f ij + ∆ with C ′ ij ( β ) = π j − π i . This also holds if ij is a restricted arc, sinceby Claim 5.3, f ij ≥ C ′ ij is continuous on (max { , f ij − ∆ } , f ij + ∆), and C ′ ij (0) = −∞ . Letus set ¯ f ij = β . This increases Ex b ( f ) by at most | F | ∆.Let us set ¯ f ij = 0 for ij ∈ E \ F . Note that f ij ≤ (2 n + m + 1)∆ if ij / ∈ F (every arc with f ij > (2 n + m + 1)∆ is either added to F or is modified to f ij = 0 in the subroutine Extend ).Further, Ex b ( f ) ≤ n ∆, and thus we obtain an F -tight ¯ f with Ex b ( ¯ f ) ≤ n ∆ + | F | ∆ + (2 n + m + 1)( m − | F | )∆ ≤ (2 n + m + 2) m ∆ . On the other hand, Ex b ( ˆ f ) ≤ nD b ( F ) ≤ n ∆, since Ex ˆ b ( ˆ f ) = 0 and ˆ b is obtained from b bymodifying certain values by ≤ D b ( F ). Consequently, || ρ ¯ f − ρ ˆ f || ≤ || ρ ¯ f − b || + || ρ ˆ f − b || = 2 Ex b ( ¯ f ) + 2 Ex b ( ˆ f ) ≤ n + m + 3) m ∆ ≤ m ( m + 1)∆ . Applying Lemma 5.6 for x = ¯ f and y = ˆ f gives || ˆ f − ¯ f || ∞ ≤ m ( m + 1)∆. We also have || f − ¯ f || ∞ ≤ (2 n + m + 1)∆ ≤ (3 m + 1)∆ by the construction, and therefore || f − ˆ f || ∞ ≤ || f − ¯ f || ∞ + || ¯ f − ˆ f || ∞ < m + 1) ∆ − ∆Applying Claim 5.7 for f and ˆ f we conclude that ˆ f is 6( m + 1) ∆-feasible; recall that f was (∆ , F )feasible when we applied Trial-and-Error . Theorem 5.8.
Let ρ S ( n, m ) be the running time needed for one shortest path computation fornonnegative lengths. Then the running time of the algorithm is bounded by O (( n + m N )( ρ T ( n, m ) + ρ E ( n, m )) + ( n + m N ) ρ S ( n, m ) log m ) . Proof.
By Theorem 5.5, there are at most ( n + m N ) log m scaling phases, each dominated by O ( n + m N ) shortest path computations. The subroutine Trial-and-Error is performed onlywhen F is extended, that is, at most n + m N times, and comprises the subroutines Trial and
Error . 21 .1 Shortest path computations
For the sake of efficiency, we shall maintain a potential vector π during the entire algorithm suchthat ( f, π ) satisfies the condition (9) on (∆ , F )-feasibility.For the initial ∆ value, ∆ ≥ err ∅ ( ), and the latter value is computed by Error ( , ∅ ). Thismeans that f = is (∆ , ∅ )-feasible. Similarly, after every successful trial we have a new flow ˆ f computed by Error ( f, F ) and new scaling factor value ∆ next ≥ err F ( ˆ f ). In the applications, thissubroutine will also return a potential vector π such that ( f, π ) satisfies (9).Alternatively, such a potential vector may be obtained by the standard label correcting algo-rithm (see [1, Chapter 5.5]), since it is a dual proof of the fact that the graph E Ff (∆) contains nonegative cycles with respect to the cost function C ′ ij ( f ij + ∆); we have access to these values viathe value oracle (Oracle 1).In the main part of the ∆-phase, we may apply a variant of Dijkstra’s algorithm (see [1, Chapter4.5]) to compute shortest paths. This needs a nonnegative cost function, but instead of the original C ′ ij ( f ij + ∆) that may take negative values, we shall use C ′ ij ( f ij + ∆) − π j + π i , a nonnegativefunction by (9); the set of shortest paths is identical for the two costs. This subroutine can beimplemented by updating the potentials π , so that (∆ , F )-feasibility is maintained, and we obtain C ′ ij ( f ij + ∆) = π j − π i on every arc of every shortest path. For the sake of completeness, we describethis subroutine in the Appendix.As shown in the proof of Lemma 5.1(ii), once we have a potential π such that C ′ ij ( f ij + ∆) = π j − π i on every arc of a shortest path P , then sending ∆-units of flow on P maintains (9) for( f, π ). It is also maintained in Extend (∆ , f, F ) since flow values are modified only on arcs with C ′ ij constant. Finally, Adjust (∆ / , f, F ) modifies the flow so that (9) is maintained for the same π and ∆ / Trial-and-error returns a rational flow vector f and a rational value ∆. Since flowwill always be modified in units of ∆ in all other parts of the algorithm, we may conclude thata rational f will be maintained in all other parts. Under Oracle 1(a) (i.e., quadratic objectives),we shall maintain a rational potential vector π , while under Oracle 1(b) (i.e., Fisher markets),we shall maintain the rationality of the e π i values; during the computations, we shall use therepresentation of these values instead of the original π . For this aim, we will use a multiplicativevariant of Dijkstra’s algorithm, also described in the Appendix. We shall also verify that in thecorresponding applications, the subroutine Error ( f, F ) returns a potential vector π so that ( f, π )satisfies (9), with the π i or the e π i values being rational, respectively.Finally, it is easy to verify that whereas we are working on a transformed uncapacitated instance,we may use the complexity bound of the original instance, as summarized in the following remark. Remark 5.9.
A shortest path computation can be performed in time ρ S ( n, m ) = O ( m + n log n )using Fibonacci heaps, see [9]. Recall that the original problem instance was on n ′ nodes and m ′ arcs, and it was transformed to an uncapacitated instance on n = n ′ + m ′ nodes and m = 2 m ′ arcs. However, as in Orlin’s algorithm [29], we can use the bound O ( m ′ + n ′ log n ′ ) instead of O ( m ′ + m ′ log n ′ ) because shortest path computations can be essentially performed on the originalnetwork. 22 Applications
Assume that C ij ( α ) = c ij α + d ij α for each ij ∈ E , with c ij ≥
0. This clearly satisfies theassumption in Oracle 1(i) since C ′ ij ( α ) = 2 c ij α + d ij . Also, ( ⋆ ) is satisfied: ij is linear if c ij = 0.The subroutine Trial ( F, b ) can be implemented by solving a system of linear equations. π j − π i = 2 c ij x ij + d ij ∀ ij ∈ F X j : ji ∈ F x ji − X j : ij ∈ F x ij = b i ∀ i ∈ V (13) x ij = 0 ∀ ij ∈ E \ F The conditions in Oracle 2 is verified by the next claim.
Lemma 6.1.
Let F be linear acyclic (that is, there is no undirected cycle of arcs with c ij = 0 ) with D b ( F ) = 0 . Then (13) is feasible and a solution can be found in ρ T ( n, m ) = O ( n . + m ) time.Proof. Clearly, we can solve the system separately on different undirected connected components of F . In the sequel, let us focus on a single connected component; for simplicity of notation, assumethis component is the entire V .Consider first the case when all arcs are linear. Then we can solve the equalities correspondingto edges and nodes separately. As F is assumed to be linear acyclic, it forms a tree. If we fix one π j value arbitrarily, it determines all other π i values by moving along the edges in the tree. The x ij ’s can be found by solving a flow problem on the same tree with the demands b i . This is clearlyfeasible by the assumption D b ( F ) = 0, that is, P i ∈ V b i = 0 (note that we do not have nonnegativityconstraints on the arcs). Both tasks can be performed in linear time.Assume next both linear and nonlinear arcs are present, and let T be an undirected connectedcomponent of linear arcs. As above, all π j − π i values for i, j ∈ T are uniquely determined. Ifthere is a nonlinear arc ij ∈ F with i, j ∈ T , then x ij = ( π j − π i − d ij ) / (2 c ij ) = α is also uniquelydetermined. We can remove this edge by replacing b i by b i + α and b j by b j − α . Hence we mayassume that the components of linear arcs span no nonlinear arcs.Next, we can contract each such component T to a single node t by setting b t = P i ∈ T b i andmodifying the d ij values on incident arcs as follows. Let t correspond to a fixed node in T , andconsider an arc with i ∈ T , j / ∈ T . Let α denote the sum of d ab values on the t − i path in T ; letus add α to d ij . Similarly for an arc ij entering T we must subtract the sum of the costs on the t − j path from d ij . A solution to the contracted problem can be easily extended to the originalinstance.For the rest, we can assume all arcs are nonlinear, that is, c ij > ij ∈ F . Let A be thenode-arc incidence matrix of F : A i,ij = − A i,ji = 1 for all ij ∈ F , and all other entries are 0. Let C be the | F | × | F | diagonal matrix with C ij,ij = − c ij . (13) can be written in the form (cid:18) A T C A (cid:19) ( π, x ) T = (cid:18) db (cid:19) . This can be transformed into (cid:18) A T CL (cid:19) ( π, x ) T = (cid:18) db ′ (cid:19) , L is the weighted Laplacian matrix with L ii = P j : ij ∈←→ F c ij , L ij = L ji = − c ij if ij ∈ F and L ij = 0 otherwise, and b ′ is an appropriate vector with P i ∈ V b ′ i = 0.The main task is to solve the system Lπ = b ′ . It is well-know (recall that V is assumed to bea single connected component) that L has rank | V | − P i ∈ V b ′ i = 0. A solution can be found in O ( n . ) time [3]. All previously described operations(eliminating nonlinear arcs spanned in components of linear arcs, contracting components of lineararcs) can be done in O ( m ) time, hence the bound ρ T ( n, m ) = O ( n . + m ).To implement Error ( f, F ), we have an F -tight vector f , and we need to find the minimum∆-value such that there exists a π potential with π j − π i ≤ (2 c ij f ij + d ij ) + 2 c ij ∆ ∀ ij ∈ E ∪ ←− F . (14)We show that this can be reduced to the minimum-cost-to-time ratio cycle problem, defined asfollows (see [1, Chapter 5.7]). In a directed graph, there is a cost function p ij and a time τ ij ≥ C minimizing ( P ij ∈ C p ij ) / ( P ij ∈ C τ ij ). A stronglypolynomial algorithm was given by Megiddo [23, 24] that solves the problem in min { O ( n log n ) ,O ( n log n ( n + m log log n )) } time. The problem can be equivalently formulated asmin µ s. t. there are no negative cyclesfor the cost function p ij + µτ ij . (15)Our problem fits into this framework with p ij = 2 c ij f ij + d ij and τ ij = 2 c ij . In (15), the optimal µ value is − ∆. However, [23] defines the minimum ratio cycle problem with τ ij > ij ∈ E .This property is not essential for Megiddo’s algorithm, which uses a parametric search method for µ to solve (15) under the only (implicit) restriction that the problem is feasible.In our setting τ ij > τ ij = 0 for linear arcs. Also, there can becycles C with P ij ∈ C τ ij = 0. (This can happen even if F is linear acyclic, as C can be any cycle in E ∪ ←− F .) If we have such a cycle C with P ij ∈ C p ij <
0, then (15) is infeasible. In every other case,the problem is feasible and thus Megiddo’s algorithm can be applied.For this reason, we first check whether there is a negative cycle with respect to the p ij ’s in theset of linear arcs in E ∪ ←− F . This can be done via the label correcting algorithm in O ( nm ) time([1, Chapter 5.5]). If there exists one, then (14) is infeasible, thus err F ( f ) = ∆ = ∞ , and (P) isunbounded as we can send arbitrary flow around this cycle. Otherwise, we have P ij ∈ C τ ij > P ij ∈ C p ij <
0, and consequently, there exists a finite ∆ satisfying (14).Consequently, ρ T ( n, m ) = min { O ( n log n ) , O ( n log n ( n + m log log n )) } . Theorem 5.8 givesthe following running time bound. Theorem 6.2.
For convex quadratic objectives on an uncapacitated instance on n nodes and m arcs, the algorithm finds an optimal solution in O ( m ( n log n + m log m ( m + n log n ))) time. Fora capacitated instance, the running time can be bounded by O ( m log m ) . The bottleneck is clearly the m minimum-cost-to-time computations. As in Remark 5.9, it islikely that one can get the same running time O ( m ( n log n + m log m ( m + n log n ))) for capacitatedinstances via a deeper analysis of Megiddo’s algorithm.Let us verify that the algorithm is strongly polynomial. It uses elementary arithmetic operationsonly, and the running time is polynomial in n and m , according to the above theorem. It is left toverify requirement (iii) on strongly polynomial algorithms (see the Introduction): if all numbers in24he input are rational, then every number occurring in the computations is rational and is of sizepolynomially bounded in the size of the input.At the initialization and in every successful trial, we compute a new flow f by solving (13) asdescribed in Lemma 6.1, and compute the new ∆ and π values by Megiddo’s algorithm. Theseare strongly polynomial subroutines and return rational values of size polynomially bounded inthe input. Namely, solving (13) requires first contracting components of linear arcs and modifyingcosts and demands by additive terms. In the contracted instance, we need to solve a system oflinear equations by exact arithmetics. This can be done by maintaining that the sizes of numbersin the output are polynomially bounded in the input size, see e.g. [32, Chapter 3]. The new ∆ and π are obtained using Megiddo’s strongly polynomial parametric search algorithm. It is immediatethat ∆ will be of polynomial encoding size, since it equals the cost-to-time ratio of a certain cycle,with both costs and times of polynomial encoding size.Consider now the phases between any two successful trials (or between the initialization andthe first successful trial); the bound on the number of such phases is O (log m ). The value of ∆decreases by a factor of 2 at the end of each phases, and the value of f is modified by ± ∆ in pathaugmentations and by ± ∆ / Adjust subroutine. Consequently, the flow remains an integermultiple of ∆ on the arcs ij ∈ E \ F up to the Adjust subroutine(see also Lemma 5.1(i)). On arcs ij ∈ F , it will be the sum of the value returned by Trial-and-Error , plus an integer multiple of∆. The bound O ( n + m N ) on the number of path augmentations, and the bound O (log m ) on thenumber of phases guarantees that the numerators also remain polynomially bounded. In the linear Fisher market model , we are given a set B of buyers and a set G of goods. Buyer i has a budget m i , and there is one divisible unit of each good to be sold. For each buyer i ∈ B andgood j ∈ G , U ij ≥ i for one unit of good j . Let n = | B | + | G | ; let E be the set of pairs ( i, j ) with U ij > m = | E | . We assume that there is at least one edgein E incident to every buyer and to every good.An equilibrium solution consist of prices p j of the goods and allocations x ij , so that (i) all goodsare sold, (ii) all money of the buyers is spent, and (iii) each buyer i buys a best bundle of goods,that is, goods j maximizing U ij /p j .The classical convex programming formulation of this problem was given by Eisenberg and Gale[7]. Recently, Shmyrev [33] gave the following alternative formulation. The variable f ij representsthe money paid by buyer i for product j .min X j ∈ G p j (log p j − − X ij ∈ E f ij log U ij X j ∈ G f ij = m i ∀ i ∈ B X i ∈ B f ij = p j ∀ j ∈ Gf ij ≥ ∀ ij ∈ E Let us construct a network on node set B ∪ G ∪ { t } as follows. Add an arc ij for every ij ∈ E , andan arc jt for every j ∈ G . Set b i = − m i for i ∈ B , b j = 0 for j ∈ G and b t = P i ∈ B m i . Let all lowerarc capacities be 0 and upper arc capacities ∞ . With p j representing the flow on arc jt , and f ij the flow on arc ij , the above formulation is a minimum-cost flow problem with separable convexobjective. (The arc jt is restricted, with extending the functions p j (log p j −
1) to take value 0 in 025nd ∞ on ( −∞ , p j = f jt shall be used for some pseudoflow f in the above problem.Let us justify that an optimal solution gives a market equilibrium. Let f be an optimal solutionthat satisfies (3) with π : B ∪ G ∪ { t } → R . We may assume π t = 0. C ′ jt ( α ) = log α implies π j = − log p j . On each ij ∈ E we have π j − π i ≤ − log U ij with equality if f ij >
0. With β i = e π i ,this is equivalent to U ij /p j ≤ β i , verifying that every buyer receives a best bundle of goods.Oracle 1(b) is a valid assumption, since the derivatives on arcs ij between buyers and goodsare − log U ij , while on an arc jt it is log f jt . The property ( ⋆ ) is straightforward.Let us turn to Oracle 2. When the subroutine Trial is called, we transform b to ˆ b by changingthe value at one node of each component K of F . For simplicity, let us always modify b t if t ∈ K ,and on an arbitrary node for the other components. We shall verify the assumptions in Oracle 2only for such ˆ b ’s; the argument can easily be extended to arbitrary ˆ b (although it is not necessaryfor the algorithm). Let us call the component K containing t the large component .In Trial ( F ), we want to find a potential π : B ∪ G ∪ { t } → R ∪ {∞} , money allocations f ij for ij ∈ F , i ∈ B , j ∈ G , and prices p j = f jt for jt ∈ F such that π j − π i = − log U ij ∀ ij ∈ F, i ∈ B, j ∈ Gπ t − π j = log p j ∀ jt ∈ F X j ∈ G,ij ∈ F f ij = ˆ b i ∀ i ∈ B X i ∈ B,ij ∈ F f ij = p j ∀ jt ∈ F X i ∈ B,ij ∈ F f ij = ˆ b j ∀ jt ∈ E \ F We may again assume π t = 0. Let P j = e − π j for j ∈ G and β i = e π i for i ∈ B . With this notation, U ij /P j = β i for ij ∈ F . If jt ∈ F , then P j = p j .Finding f and π can be done independently on the different components of F . For any com-ponent different from the large one, all edges are linear. Therefore we only need to find a feasibleflow on a tree, and independently, P j and β i values satisfying U ij /P j = β i on arcs ij in this com-ponent. Both of these can be performed in linear time in the number of edges in the tree. Notethat multiplying each P j by a constant α > β i by the same α yields anotherfeasible solution.Let T , . . . , T k be the components of the large component after deleting t . If T ℓ contains a singlegood j , then we set p j = P j = 0 ( π j = ∞ ). If T ℓ is nonsingular, then F restricted to T ℓ forms aspanning tree. The equalities U ij /P j = β i uniquely define the ratio P j /P j ′ for any j, j ′ ∈ G ∩ T ℓ .We have that p j = P j and P i ∈ B ∩ T ℓ m i = P j ∈ G ∩ T ℓ p j by the constraints on the buyers in B ∩ T ℓ and goods in G ∩ T ℓ ; note that ˆ b i = − m i for all buyers in B ∩ T ℓ . Hence the prices in T ℓ are uniquelydetermined. Then the edges in F simply provide the allocations f ij . All these computations canbe performed in ρ T ( n, m ) = O ( m ) time.For Oracle 3, we show that Error ( f, F ) can be implemented based on the Floyd-Warshall algo-rithm (see [1, Chapter 5.6]). Let π be the potential witnessing that f is (∆ , F )-feasible. Assuming π t = 0, and using again the notation P j = e − π j for j ∈ G and β i = e π i for i ∈ B , we get U ij /P j ≤ β i if i ∈ B, j ∈ G, ij ∈ E, with equality if ji ∈ E Ff . (16)Furthermore, we have p j − ∆ ≤ P j ≤ p j + ∆ if p j > P j ≤ ∆ if p j = 0.26et us now define γ : G × G → R as γ jj ′ = max (cid:26) U ij ′ U ij : i ∈ B, ji, ij ′ ∈ E Ff (cid:27) . If no such i exists, define γ jj ′ = 0; let γ jj = 1 for every j ∈ G . Claim 6.3.
Assume we are given some P j values, j ∈ G . There exists β i values ( i ∈ B ) satisfying(16) if and only if P j ′ ≥ P j γ jj ′ holds for every j, j ′ ∈ G .Proof. The condition is clearly necessary by the definition of γ jj ′ . Conversely, if this conditionholds, setting β i = max j ∈ G U ij /P j does satisfy (16).If there is a directed cycle C with Π ab ∈ C γ ab >
1, then f cannot be (∆ , F )-feasible for any ∆.Otherwise, we may compute ˜ γ jj ′ as the maximum of Π ab ∈ P γ ab over all directed paths P in E Ff from j to j ′ (setting the value 0 again if no such path exists). This can be done by the multiplicativeversion of the Floyd-Warshall algorithm in O ( n ) time (note that this is equivalent to findingall-pair shortest paths for − log γ ab ).For (∆ , F )-feasibility, we clearly need to satisfy( p j − ∆)˜ γ jj ′ ≤ P j ˜ γ jj ′ ≤ P j ′ ≤ p j ′ + ∆ . Let us define ∆ as the smallest value satisfying all these inequalities, that is,∆ = max (cid:26) , max j,j ′ ∈ G p j ˜ γ jj ′ − p j ′ ˜ γ jj ′ + 1 (cid:27) . (17)We claim that f is (∆ , F )-feasible with the above choice. For each j ∈ G , let P j = max h ∈ G ˜ γ hj ( p h − ∆). It is easy to verify that these P values satisfy P j ′ ≥ P j γ jj ′ , and p j − ∆ ≤ P j ≤ p j + ∆. Thecondition (16) follows by Claim 6.3.The complexity of Error ( f, F ) is dominated by the Floyd-Warshall algorithm, O ( n ) [8]. Theproblem is defined on an uncapacitated network, with the number of nonlinear arcs m N = | G | < n .Thus Theorem 5.8 gives the following. Theorem 6.4.
For Fisher’s market with linear utilities, the algorithm finds an optimal solution in O ( n + n ( m + n log n ) log n ) . The algorithm of Orlin [30] runs in O ( n log n ) time, assuming m = O ( n ). Under this assump-tion, we get the same running time bound.To prove that the algorithm is strongly polynomial, let us verify the nontrivial requirement (iii)(see the Introduction). As discussed in Section 5.1, if the input is rational, we shall maintain that f ,∆ and the e π i values are rational; the latter are used in the computations instead of the π i ’s. At theinitialization and in every successful trial, the subroutines described above are strongly polynomialand therefore return rational f , ∆ and e π i values, of size polynomially bounded in the input (notethat the e π i values above are denoted by P i for i ∈ G and β i for i ∈ B , and e π t = 1). Betweentwo successful trials, we can use the same argument as in Section 6.1 for quadratic costs: there are O (log m ) such iterations, ∆ is divided by two at the end of every phase, the path augmentationschange f by ± ∆ and Adjust by ± ∆ /
2. The multiplicative Dijkstra algorithm described in theAppendix also maintains rational e π i values of polynomial encoding length.27 .3 Fisher’s market with spending constraint utilities The spending constraint utility extension of linear Fisher markets was defined by Vazirani [37].In this model, the utility of a buyer decreases as the function of the money spent on the good.Formally, for each pair i and j there is a sequence U ij > U ij > . . . > U ℓ ij ij > L ij , . . . , L ℓ j ij >
0. Buyer i accrues utility U ij for every unit of j he purchased by spendingthe first L ij dollars on good j , U ij for spending the next L ij dollars, etc. These ℓ ij intervalscorresponding to the pair ij are called segments . ℓ ij = 0 is allowed, but we assume P j ∈ G ℓ ij > i ∈ B and P i ∈ B ℓ ij > j ∈ G . Let n = | B | + | G | denote the total number of buyersand goods, and m denote the total number of segments. Note that m > n is also possible.No extension of the Eisenberg-Gale convex program is known to capture this problem. Theexistence of a convex programming formulation is left as an open question in [37]. This was settledby Devanur et al. [2], giving a convex program based on Shmyrev’s formulation. Let f kij representthe money paid by buyer i for the k ’th segment of product j , 1 ≤ k ≤ ℓ ij .min X i ∈ G p j (log p j − − X i ∈ B,j ∈ G, ≤ k ≤ ℓ ij f kij log U kij X j ∈ G, ≤ k ≤ ℓ ij f kij = m i ∀ i ∈ B X i ∈ B, ≤ k ≤ ℓ ij f kij = p j ∀ j ∈ G ≤ f kij ≤ L kij ∀ ij ∈ E. This gives a convex cost flow problem again on the node set B ∪ G ∪ { t } , by adding ℓ ij parallelarcs from i ∈ B to j ∈ G , and arcs jt for each j ∈ G . The upper capacity on the k ’th segmentfor the pair ij is L kij . To apply our method, we first need to transform it to an equivalent problemwithout upper capacities. This is done by replacing the arc representing the k ’th segment of ij by anew node ( ij, k ) and two arcs i ( ij, k ) and j ( ij, k ). The node demand on the new node is set to L kij ,while on the good j , we replace the demand 0 by − P i,k L kij , the negative of the sum of capacitiesof all incident segments. The cost function on i ( ij, k ) is − log U kij α , while the cost of j ( ij, k ) is 0.Let S denote the set of the new ( ij, k ) nodes. This modified graph has n ′ = n + m + 1 nodes and m ′ = 2 m + | G | arcs.Assumption ( ⋆ ) is clearly valid. Oracle 1(b) is satisfied the same way as for linear Fishermarkets, using an oracle for the e C ′ ij ( α ) values.In Trial ( F ), we want to find an F -tight flow f ′ on the extended network, witnessed by thepotential π : B ∪ S ∪ G ∪ { t } → R . We may assume π t = 0. Let P j = e − π j for j ∈ G and β i = e π i for i ∈ B and S kij = e − π ( ij,k ) . For the k ’th segment of ij , U kij /S kij = β i if i ( ij, k ) ∈ F and S kij = P j if j ( ij, k ) ∈ F .As for linear Fisher markets, if a component of F does not contain t , we can simply computeall potentials and flows as F is a spanning tree of linear edges in this component.For the component K with t ∈ K , let T ℓ be a component of K − t . F is a spanning tree of linearedges in T ℓ as well, therefore the ratio P j /P j ′ is uniquely defined for any j, j ′ ∈ G ∩ T ℓ . On the otherhand, we must have P j = p j , and we know that P j ∈ G ∩ T ℓ p j = − P v ∈ T ℓ b v by flow conservation.These determine the P j = p j values, and thus all other β i and S kij values in the component as well.The support of the flow f ij is a tree and hence it can also easily computed. The running time of Trial is again linear, ρ T ( n ′ , m ′ ) = O ( m ′ ) = O ( m ).28 rror ( f, F ) can be implemented the same way as for the linear Fisher market. We shall definethe values γ : G × G → R so that P j ′ ≥ P j γ jj ′ must hold, and conversely, given P j prices satisfyingthese conditions, we can define the β i and S kij values feasibly. Let γ jj ′ = max n U k ′ ij ′ U kij : i ∈ B,j ( ij, k ) , ( ij, k ) i, i ( ij ′ , k ′ ) , ( ij ′ , k ′ ) j ′ ∈ E Ff o . Given these γ jj ′ values, the ˜ γ jj ′ values can be computed by the Floyd-Warshall algorithm and theoptimal ∆ obtained by (17) as for the linear case.Finding the γ jj ′ values can be done in O ( m ′ ) time, and the Floyd-Warshall algorithm runs in O ( | G | ). This gives ρ E ( n ′ , m ′ ) = O ( m ′ + | G | ) = O ( m + n ). From Theorem 5.8, together withRemark 5.9, we obtain: Theorem 6.5.
For an instance of Fisher’s market with spending constraint utilities with n = | B | + | G | and m segments, the running time can be bounded by O ( mn + m ( m + n log n ) log m ) . It can be verified that the algorithm is strongly polynomial the same way as for the linear case.
We have given strongly polynomial algorithms for a class of minimum-cost flow problems withseparable convex objectives. This gives the first strongly polynomial algorithms for quadraticconvex cost functions and for Fisher’s market with spending constraint utilities. For Fisher’smarket with linear utilities, we get the same complexity as in [30].The bottleneck in complexity of all applications is the subroutine
Trial . However, the exactvalue of err F ( f ) is not needed: a constant approximation would also yield the same complexitybounds. Unfortunately, no such algorithm is known for the minimum cost-to-time ratio cycleproblem that would have significantly better, strongly polynomial running time. Finding such analgorithm would immediately improve the running time for quadratic costs.A natural future direction could be to develop strongly polynomial algorithms for quadraticobjectives and constraint matrices with bounded subdeterminants. This would be a counterpartof Tardos’s result [36] for linear programs. Such an extension could be possible by extending ourtechniques to the setting of Hochbaum and Shantikumar [18].The recent paper [39] shows that linear Fisher market, along with several extension, can becaptured by a concave extension of the generalized flow model. A natural question is if there is anydirect connection between the concave generalized flow model and the convex minimum cost flowmodel studied in this paper. Despite certain similarities, no reduction is known in any direction.Indeed, no such reduction is known even between the linear special cases, that is, generalizedflows and minimum-cost flows. The perfect price discrimination model [11], and the Arrow-DebreuNash-bargaining problem [38], are instances of the concave generalized flow model, but they arenot known to be reducible to convex cost flows. On the other hand, the spending constraint utilitymodel investigated in this paper is not known to be reducible to concave generalized flows.The algorithm in [39] is not strongly polynomial. Even for linear generalized flows, the firststrongly polynomial algorithm was only given very recently [40]. One could try to extend this to aclass of concave generalized flows in a similar manner as in the current paper, i.e. assuming certainoracles. This could lead to strongly polynomial algorithms for the market problems that fit intothis model. 29 related problem is finding a strongly polynomial algorithm for minimizing a separable convexobjective over a submodular polyhedron. Fujishige [10] showed that for separable convex quadraticcosts, this is essentially equivalent to submodular function minimization. Submodular utility al-location markets by Jain and Vazirani [21] also fall into this class, and are solvable in stronglypolynomial time; see also Nagano [27]. Other strongly polynomially solvable special cases are givenby Hochbaum and Hong [14].A common generalization of this problem and ours is minimizing a separable convex objectiveover a submodular flow polyhedron. Weakly polynomial algorithms were given by Iwata [19] and byIwata, McCormick and Shigeno [20]. One might try to develop strongly polynomial algorithms forsome class of separable convex objectives; in particular, for separable convex quadratic functions. Acknowledgment
The author is grateful to an anonymous referee for several suggestions that helped to improve thepresentation.
References [1] R. K. Ahuja, T. L. Magnanti, and J. B. Orlin.
Network Flows: Theory, Algorithms, andApplications . Prentice-Hall, Inc., Upper Saddle River, New Jersey, feb 1993.[2] B. Birnbaum, N. R. Devanur, and L. Xiao. Distributed algorithms via gradient descent forFisher markets. In
Proceedings of ACM EC , pages 127–136, 2011.[3] D. Coppersmith and S. Winograd. Matrix multiplication via arithmetic progressions.
Journalof Symbolic Computation , 9(3):251–280, 1990.[4] S. Cosares and D. S. Hochbaum. Strongly polynomial algorithms for the quadratic transporta-tion problem with a fixed number of sources.
Mathematics of Operations Research , 19(1):94–111, 1994.[5] N. R. Devanur, C. H. Papadimitriou, A. Saberi, and V. V. Vazirani. Market equilibrium via aprimal–dual algorithm for a convex program.
Journal of the ACM (JACM) , 55(5):22, 2008.[6] J. Edmonds and R. M. Karp. Theoretical improvements in algorithmic efficiency for networkflow problems.
Journal of the ACM (JACM) , 19(2):248–264, 1972.[7] E. Eisenberg and D. Gale. Consensus of subjective probabilities: The pari-mutuel method.
The Annals of Mathematical Statistics , 30(1):165–168, 1959.[8] R. Floyd. Algorithm 97: shortest path.
Communications of the ACM , 5(6):345, 1962.[9] M. L. Fredman and R. E. Tarjan. Fibonacci heaps and their uses in improved network opti-mization algorithms.
Journal of the ACM (JACM) , 34(3):596–615, 1987.[10] S. Fujishige. Submodular systems and related topics.
Mathematical Programming at Oberwol-fach II , (22):113–131, 1984.[11] G. Goel and V. V. Vazirani. A perfect price discrimination market model with production, anda (rational) convex program for it.
Mathematics of Operations Research , 36:762–782, 2011.3012] F. Granot and J. Skorin-Kapov. Towards a strongly polynomial algorithm for strictly con-vex quadratic programs: An extension of Tardos’ algorithm.
Mathematical Programming ,46(1):225–236, 1990.[13] M. Gr¨otschel, L. Lov´asz, and A. Schrijver.
Geometric Algorithms and Combinatorial Opti-mizations . Springer-Verlag, 1993.[14] D. Hochbaum and S. Hong. About strongly polynomial time algorithms for quadratic opti-mization over submodular constraints.
Mathematical Programming , 69(1):269–309, 1995.[15] D. Hochbaum and M. Queyranne. Minimizing a convex cost closure set.
SIAM Journal onDiscrete Mathematics , 16:192, 2003.[16] D. S. Hochbaum. Lower and upper bounds for the allocation problem and other nonlinearoptimization problems.
Mathematics of Operations Research , 19(2):390–409, 1994.[17] D. S. Hochbaum. Complexity and algorithms for nonlinear optimization problems.
Annals ofOperations Research , 153(1):257–296, 2007.[18] D. S. Hochbaum and J. G. Shanthikumar. Convex separable optimization is not much harderthan linear optimization.
Journal of the ACM (JACM) , 37(4):843–862, 1990.[19] S. Iwata. A capacity scaling algorithm for convex cost submodular flows.
Mathematical Pro-gramming , 76(2):299–308, 1997.[20] S. Iwata, S. McCormick, and M. Shigeno. Fast cycle canceling algorithms for minimum costsubmodular flow.
Combinatorica , 23(3):503–525, 2003.[21] K. Jain and V. V. Vazirani. Eisenberg-Gale markets: Algorithms and game-theoretic proper-ties.
Games and Economic Behavior , 70(1):84–106, 2010.[22] A. V. Karzanov and S. T. McCormick. Polynomial methods for separable convex optimizationin unimodular linear spaces with applications.
SIAM J. Comput. , 26(4):1245–1275, 1997.[23] N. Megiddo. Combinatorial optimization with rational objective functions.
Mathematics ofOperations Research , 4(4):414–424, 1979.[24] N. Megiddo. Applying parallel computation algorithms in the design of serial algorithms.
Journal of the ACM (JACM) , 30(4):852–865, 1983.[25] M. Minoux. A polynomial algorithm for minimum quadratic cost flow problems.
EuropeanJournal of Operational Research , 18(3):377–387, 1984.[26] M. Minoux. Solving integer minimum cost flows with separable convex cost objective polyno-mially.
Mathematical Programming Study , 25:237, 1985.[27] K. Nagano. On convex minimization over base polytopes.
Integer Programming and Combi-natorial Optimization , pages 252–266, 2007.[28] N. Nisan, T. Roughgarden, E. Tardos, and V. Vazirani.
Algorithmic Game Theory . CambridgeUniversity Press New York, NY, USA, 2007.[29] J. B. Orlin. A faster strongly polynomial minimum cost flow algorithm.
Operations Research ,41(2):338–350, 1993. 3130] J. B. Orlin. Improved algorithms for computing Fisher’s market clearing prices. In
Proceedingsof the 42nd ACM Symposium on Theory of Computing (STOC) , pages 291–300. ACM, 2010.[31] J. Renegar. On the worst-case arithmetic complexity of approximating zeros of polynomials.
Journal of Complexity , 3(2):90–113, 1987.[32] A. Schrijver.
Theory of Linear and Integer Programming . John Wiley & Sons, 1998.[33] V. I. Shmyrev. An algorithm for finding equilibrium in the linear exchange model with fixedbudgets.
Journal of Applied and Industrial Mathematics , 3(4):505–518, 2009.[34] A. Tamir. A strongly polynomial algorithm for minimum convex separable quadratic cost flowproblems on series-parallel networks.
Mathematical Programming , 59:117–132, 1993.[35] ´E. Tardos. A strongly polynomial minimum cost circulation algorithm.
Combinatorica ,5(3):247–255, 1985.[36] ´E. Tardos. A strongly polynomial algorithm to solve combinatorial linear programs.
OperationsResearch , 34(2):250–256, 1986.[37] V. V. Vazirani. Spending constraint utilities with applications to the adwords market.
Math-ematics of Operations Research , 35(2):458–478, 2010.[38] V. V. Vazirani. The notion of a rational convex program, and an algorithm for the Arrow-Debreu Nash bargaining game.
Journal of ACM (JACM) , 59(2), 2012.[39] L. A. V´egh. Concave generalized flows with applications to market equilibria.
Mathematics ofOperations Research , 39(2):573–596, 2014.[40] L. A. V´egh. Strongly polynomial algorithm for generalized flow maximization. In
Proceedingsof the 46th ACM Symposium on Theory of Computing (STOC) , 2014.
Appendix
In this Appendix we describe two variants of Dijkstra’s algorithm that are used for the shortestpath computations in our algorithm. This is an equivalent description of the well-known algorithm,see e.g. [1, Chapter 4.5]. The first, standard version is shown in Algorithm 5. We start from acost function c on a digraph D = ( V, A ) and a potential vector π with c ij − π j + π i ≥ S and T . The set R is initialized as R = S , and denotes in everyiteration the set of nodes that can be reached from S on a tight path, that is, all arcs of the pathsatisfying c ij − π j + π i = 0. Every iteration increases the potential on V \ R until some new tightarcs enter. We terminate once R contains a node in T ; a shortest path between S and T can berecovered using the pointers pred ( i ).In our algorithm, this subroutine will be applied if Oracle 1(a) holds. In the ∆-phase, weapply it for the digraph E Ff (∆) and the cost function c ij = C ′ ij ( f ij + ∆), and the potential π as inthe algorithm. Note that if the initial π is rational, and all c ij values are rational, the algorithmterminates with a π that is also rational. Oracle 1(a) guarantees that if f ij and ∆ are rationalnumbers, then so is c ij .Algorithm 6 shows a multiplicative variant of the previous algorithm; they are identical aftersubstituting c ij = log γ ij and π i = log µ i . This variant shall be applied under Oracle 1(b). We shall32ssume that every e π i value is rational, and set µ i = e π i , and γ ij = e C ′ ij ( f ij +∆) . The assumptionguarantees that if f ij and ∆ are rational numbers, then so is γ ij . Consequently, the rationality ofthe e π i values is maintained during the computations.33 lgorithm 5Subroutine Shortest Paths
INPUT
A digraph D = ( V, A ), disjoint subsets
S, T ⊆ V , a cost function c : A → R and a potential vector π : V → R with c ij − π j + π i ≥ ij ∈ A . OUTPUT
A shortest path P between a node in S and a node in T and a π ′ : V → R with c ij − π ′ j + π ′ i ≥ ij ∈ A , and equality on every arc of P . R ← S ; for i ∈ S do pred ( i ) ← N U LL ; while R ∩ T = ∅ do α ← min { c ij − π j + π i : ij ∈ A, i ∈ R, j ∈ V \ R } ; for j ∈ V \ R do π j ← π j + α ; Z ← { j ∈ V \ R : ∃ ij ∈ A, i ∈ R such that c ij − π j + π i = 0 } ; for j ∈ Z do pred ( j ) ← i ∈ R such that ∃ ij ∈ A : c ij − π j + π i = 0; R ← R ∪ Z ; π ′ ← π ; Algorithm 6Subroutine
Multiplicative Shortest Paths
INPUT
A digraph D = ( V, A ), disjoint subsets
S, T ⊆ V , a cost function γ : A → R and a potential vector µ : V → R with γ ij µ i µ j ≥ ij ∈ A . OUTPUT