A scaling algorithm for optimizing arbitrary functions over vertices of polytopes

In this paper we present a scaling algorithm for minimizing arbitrary functions over vertices of polytopes in an oracle model of computation which includes an augmentation oracle. For the binary case, when the vertices are 0–1 vectors, we show that the oracle time is polynomial. Also, this algorithm allows us to generalize some concepts of combinatorial optimization concerning performance bounds of greedy algorithms and leads to new bounds for the complexity of the simplex method.


Introduction and formulation of the main result
We consider the following optimization problem: min{ f (x) : x ∈ S}, (1.1) where f is a function and S is a finite subset of rational vectors in R n .We assume that f can be accessed via an oracle which uses a black-box data structure to compute f .For instance, in the linear case, this structure can e.g.contain the coefficient vector.We will say that f is polynomially computable if there is an algorithm (an evaluation oracle) such that both the number of arithmetic operations and binary sizes of the numbers in the course of that algorithm are bounded by a polynomial in the binary size of that data structure and in the binary size of the argument.In this case, the binary size of f (x) is also polynomially bounded in the mentioned sense for every x ∈ S.
B Sergei Chubanov sergei.chubanov@de.bosch.com Let O PT denote the optimal value.A feasible solution x is an ε-approximate solution to the above problem if We assume that f belongs to a class C, of functions defined on S, closed under additions of linear functions, i.e., if f ∈ C and h is linear over S, then f + h ∈ C.An algorithm computing g(x) for every g ∈ C and x ∈ S is further called an evaluation oracle.
Also, we assume that we have an augmentation oracle which, given x ∈ S and g ∈ C, either returns y ∈ S such that g(y) < g(x), or correctly decides that x minimizes g over S.
Further we consider a model of computation including arithmetic operations, additions of linear functions to functions of C, and calls to the augmentation oracle.By an oracle running time or oracle complexity of an algorithm we understand the number of calls to the augmentation oracle performed by this algorithm for a given input.The access to solutions in S is only possible via the augmentation oracle, i.e., it is only possible to check whether x ∈ S is optimal for a given function g ∈ C and, if this is not the case, to get another solution in S with a lower value of g.Also, there is no direct access to the evaluation oracle, which is used exclusively inside the augmentation oracle.That is, in this model of computation we are basically interested in the number of calls to the augmentation oracle and put aside any details related to the arithmetic model of computation.
In the subsequent sections we present a general scaling framework where the main idea is to perform a perturbation of the objective function by adding a linear function multiplied by a scaling parameter.The linear function which is to be added at each iteration depends on the current solution.The scaling parameter decreases by a factor of two whenever the current solution is optimal for the perturbed objective function, which is verified by means of the augmentation oracle.An advantage of this scaling technique is that it is applicable to arbitrary objective functions, in contrast to methods based on bit-scaling like that proposed by Schulz et al. [11], which are restricted to linear objective functions.
The choice of a linear function for perturbation of the original objective function depends on the problem in question.We will consider two cases; in the first case, S is a set of binary vectors and in the second case S is the vertex set of a polytope defined by a system of linear equations and nonnegativity constraints.
If S is a set of binary vectors, our algorithm runs in polynomial oracle time.The previous results on the oracle complexity of 0-1 problems focused on linear problems.Schulz et al. [11] proved that a linear 0-1 problem can be solved in a polynomial number of calls to an augmentation oracle.Chubanov [2] proposed a variant of the simplex method for integer linear problems given by verification oracles, i.e., by the oracles which only verify whether a given solution is optimal and which do not necessarily return an improvement direction if the given solution is not optimal.This variant of the simplex method visits O(n) vertices of the convex hull of the feasible set in the 0-1 case and O(nq) vertices in the more general case where the feasible set is contained in [0, q] n .Orlin et al. [10] considered the case when local optimality does not necessarily imply global optimality and proposed a fully polynomial time approximation scheme for finding an ε-local optimal solution, i.e., one whose objective value is within a factor of 1 + ε of the minimum objective value in the respective neighborhood of this solution.In the present paper we consider only the case when local optimality means global optimality; so the existence of an FPTAS in the sense of the concept of ε-local optimality introduced in [10] seems to be an open question in the nonlinear case.
For the well-known polynomially solvable cases like the assignment problem and the minimum-weight spanning tree problem, as well as the problem of finding a maximum-weight independent set of a matroid, for which simple augmentation oracles are known, our scaling algorithm represents an alternative polynomial-time algorithmic approach.The main advantage of our algorithm is that its iterations are easy to parallelize because any improvement of the current objective function (which can differ from the objective function of the problem in question) is sufficient at each iteration, in contrast to greedy algorithms, which require the best local improvement.
On the other hand, our algorithm implies a new complexity bound for the simplex method.For the polytopes of the form P = {x ∈ R n : Ax = b, x ≥ 0}, where A is a matrix and b is a vector of the respective dimension (0 denotes a zero vector), we show that there is a variant of the simplex method which solves a linear problem over P by visiting a number of vertices which is polynomially bounded in certain parameters related to the problem, in particular in the ratio between the largest and the smallest nonzero elements of a basic feasible solution.It should be noted that our bound depends on the dimension, i.e., on the number of variables, only logarithmically; the respective bound on the diameters of polytopes is presented in Corollary 5.1.Our bound can be viewed as an improvement of the bound obtained recently by Kitahara and Mizuno [7], which depends linearly on the dimension.
Another important consequence of our scaling algorithm is that any greedy algorithm for binary optimization, with the property that it is able to find an optimal solution for any function in the class C, must run in polynomial oracle time.Moreover, the greedy algorithm runs in strongly polynomial oracle time for the case of a linear objective function.To prove this, we use our scaling algorithm.
It is not clear if our algorithm can be generalized for the case when S is a finite subset of integer vectors which are not necessarily vertices of a polytope.In particular, the question is, whether there is an algorithm for (1.1) with pseudopolynomial oracle running time when S is a finite subset of integer vectors.More precisely, assuming that S is a set of nonnegative integer vectors such that x ≤ u for some integer vector u, the question is whether there exists an algorithm for finding an ε-approximate solution for (1.1) whose oracle complexity is polynomial exclusively in n, u ∞ and log 1 ε .A special case of (1.1) where the objective function is linear and S = {x ∈ Z n : Ax = b, 0 ≤ x ≤ u} was studied by De Loera, Hemmecke, and Lee [3]; they developed a pseudopolynomial oracle algorithm for this case.
In this paper, we do not consider any concrete class of nonlinear functions which would be closed under additions of linear functions and permit an efficient augmentation oracle.This question is very nontrivial and should be studied separately.

Scaling algorithm
Let us consider a family of functions parameterized by y ∈ S and δ ∈ R + , such that, for all y ∈ S, where c y ∈ R n is a vector depending on y such that Note that, according to (2.2), c y defines a linear function which is uniquely minimized by y.If all elements of S are vertices of the convex hull of S, then such c y exists for every y ∈ S. In the next sections we consider different special cases where a suitable c y is easy to find.Denote These values are positive according to our assumptions about c y .

Lemma 2.1
The following statements are true: (i) If y is optimal for g y,δ , then y is an (μ + δ)-approximate solution of (1.1).
(ii) If x = y and g y,δ (y) ≥ g y,δ (x), then f Proof Proof of (i): Let x * be an optimal solution of (1.1).Then, Proof of (ii): This implies (ii).
Lemma 2.1 naturally suggests the following algorithm, which we further call the scaling algorihtm: Input: A family of functions (2.1) and an initial solution x 0 ∈ S. Output: x * ∈ S optimal for (1.1).
1. Set k := 0 and δ := 1.While x 0 is not optimal for g x 0 ,δ , set δ := 2δ. 2. Call the augmentation oracle for f and x k ; if x k is optimal for f , then return x * := x k .
3. Call the augmentation oracle to find x k+1 ∈ S with g x k ,δ (x k+1 ) < g x k ,δ (x k ) or prove that x k is optimal for g x k ,δ .4. If x k is optimal for g x k ,δ , then set δ := δ/2 and x k+1 := x k . 5. Set k := k + 1 and go to Step 2.
Further, we refer to repetitions of steps 2-5 as to iterations of the scaling algorithm.
At each iteration, the scaling algorithm tries to improve the current value of function g x k ,δ associated with the current solution.This function belongs to C because it is obtained from f ∈ C by adding a linear function and C is closed under additions of linear functions.So it is important to keep in mind that the augmentation oracle should be applicable to every function in C.
Step 1 sets δ to 1 and then multiplies it by 2 in the while-loop until x 0 becomes optimal for g x 0 ,δ .The while-loop of Step 1 runs in oracle time iterations because if x 0 is not optimal for g x 0 ,δ then f (x 0 ) can be improved by at least μ − δ (which follows from (ii) of Lemma 2.1) and f Thus, δ is either 1 after Step 1, in which case x 0 is optimal for g x 0 ,1 , or δ > 1.In the latter case, after Step 1, for a solution x which is optimal for g x 0 ,δ/2 , we have Therefore, the initial value δ constructed in Step 1 satisfies Let γ > 0 such that γ ≤ μ + be a strict lower bound on the gap between two different objective values that can be taken at feasible solutions:

Theorem 2.1 The scaling algorithm finds an optimal solution in oracle time
Proof As already mentioned, Step 1 requires a number of oracle calls bounded by (2.4).The first division of δ occurs at exactly the first iteration because x 0 is optimal for g x 0 ,δ after Step 1.
The number of iterations between two consecutive divisions of δ is bounded by O(μ + /μ − ) due to Lemma 2.1.Indeed, let x k be optimal for g x k ,δ .Then (i) of Lemma 2.1 implies that x k is (μ + δ)-approximate for the current δ.Denote the current δ by δ .The new δ is equal to δ /2 according to Step 4. Statement (ii) of Lemma 2.1 implies that, in each iteration before the next division of δ, the objective value will be improved by at least μ − δ /2.The number of such improvements is bounded by If δ is divided by 2, the current solution x k is (μ + δ)-approximate for the current δ.If δ was less than γ /μ + , this would mean that the current solution would be γapproximate, i.e., optimal for f .In this case Step 2 would return an optimal solution and the algorithm would stop.It follows that δ is still not less than γ /μ + whenever a division of δ occurs.Therefore, (2.5), implies that the number of divisions of δ is bounded by the logarithm in (2.6).
Thus, we have proved that the scaling algorithm runs in oracle time bounded by (2.6).

Scaling algorithm for binary optimization
In this section, we assume that Let g y,δ be defined as where (−1) y is a componentwise power.That is, the ith component of (−1) y is 1 if y i = 0 and it is equal to −1 if y i = 1.Note that g y,δ (y) = f (y).When restricted to 0-1 vectors x and y, this function has the following interpretation: when computing g y,δ (x), the value δ is added to f (x) whenever x i = y i for an index i = 1, . . ., n.
That is, one can say that we additionally pay δ for each inversion of a binary digit when moving from y to x.So if x and y are 0-1 vectors and x = y, then Let m be an upper bound on the number of nonzeros in any feasible solution: Proof Note that μ + ≤ 2m and μ − ≥ 1 because c y = (−1) y and then apply Theorem 2.1.
It is not hard to see that the scaling algorithm is a polynomial-time method for such problems as minimum-cost spanning trees, the assignment problem and, more generally, for all combinatorial problems, with polynomially computable objective functions, for which there is a polynomial-time augmentation oracle for the class C. The advantage of the scaling algorithm for instance for the spanning tree problem is a possibility of parallelization.From the matroid theory it follows that a spanning tree is optimal if and only if a replacement of an edge of this tree by an edge not belonging to the tree does not lead to a better solution.An augmentation oracle in the course of the scaling algorithm does not need to find the best improvement, it only needs to find a better solution for the current objective function or to prove that the current solution is optimal for the current objective function.That is, we can decompose the neighborhood explored by the augmentation oracle according to the number of processors in a parallel model of computation, so that each processor would independently explore the respective subset of the neighborhood.

Greedy algorithm
Now we are going to prove that every greedy algorithm of a sufficiently general class must run in a polynomial oracle time.Let N be a mapping from S to 2 S , the set of all subsets of S, such that x ∈ N (x) for all x ∈ S. In other words, N defines a neighborhood for each x ∈ S. Assume that N is defined so that the following optimality condition holds for any g ∈ C : x * ∈ arg min  This condition says that x * is optimal for g if and only if x * is locally optimal for g.For instance, if S represents the set of spanning trees of a graph, then N (x) can be defined as the set consisting of x and all spanning trees which can be obtained by means of removing an edge of the current tree x and adding another edge of the graph.
(More generally, we can consider S representing the set of bases of a matroid.)A greedy algorithm based on such pairwise interchanges belongs to the class of greedy algorithms that we define below.Consider the following greedy algorithm: Input: Problem (1.1) and an initial solution x 0 ∈ S.

Find x k+1 in arg min
In the context of the greedy algorithm, by an augmentation oracle we mean one which is able to minimize a given function g ∈ C over the neighborhood N (x) of a given solution x ∈ S. To distinguish such oracles from a larger class of augmentation oracles we have already introduced, we call such oracles local optimization oracles.Whenever we are speaking about oracle running time of the greedy algorithm, we mean a model of computation with a local optimization oracle.Although the greedy algorithm does not use optimization of other functions than f in the neighborhood of the current solution at each iteration, condition (4.1) is important when estimating the running time: Proof First, let us prove the following statement: If then x k minimizes g x k ,δ .
The statement can be proved in the following way.Assume the contrary, i.e., taking into account (4.1), that there is x ∈ N (x k ) with Then, since x k+1 minimizes f over N (x k ), The above statement implies that if (4.2) takes place, then x k is a 2mδ-approximate solution.Indeed, using the fact that (4.2) implies that x k is optimal for g x k ,δ , we can write: where x * is optimal for f .Let δ be initialized exactly as in the scaling algorithm.Let δ be divided by 2 whenever the improvement of the objective function in the greedy algorithm is smaller than δ, i.e., whenever (4.2) holds.The number of iterations between two consecutive divisions of δ does not exceed 4m because if no division of δ occurs then the current value of f is improved by at least δ and a division of δ implies that the current solution is 2mδ-approximate for the current δ.Now, to estimate the number of divisions of δ, we observe that a γ -approximate solution is optimal, which implies the respective logarithmic bound on the number of divisions of δ.
If m is much smaller than n, say, when n is exponential in m, the above theorem improves the bound O(n) for the diameter of 0,1-polytopes obtained by Naddef [9]: Corollary 4.1 Let P be a where each vertex has no more than m nonzero components.Then, the diameter of P is bounded by O(m log n).
Proof Let w be a vertex of P and f be defined as f (x) = (1 − 2w) T x.This function is uniquely minimized by w.Let v be another vertex of P. Let S be the vertex set of P and, for each x ∈ S, N (x) consist of x and the vertices adjacent to x.The greedy algorithm starting at x 0 = v constructs a path from v to w in the 1-skeleton of P to minimize f over S. Now we apply Theorem 4.1 noting that γ can be chosen as γ = 1/2 and | f (x)| ≤ n for all x ∈ S in this case.
Actually, the greedy algorithm runs in strongly polynomial oracle time provided that the objective function is linear:

Corollary 4.2 Let f be a linear function of the form f (x) = c T x for some integer vector c. Then, if S is a set of 0-1 vectors, then the greedy algorithm runs in oracle time O(mn log m).
Proof We now formulate a new problem such that μ + , μ − , and the objective values of feasible solutions will be of strongly polynomial size and the greedy algorithm will construct exactly the same sequence of solutions as when solving the original problem.
Let D = be the subset of all solution pairs (x, y) ∈ S × S such that c T (x − y) = 0 and D < be the subset of all (x, y) ∈ S × S such that c T (x − y) < 0.
The integer vector c defining the objective function f is contained in the polyhedron which means that Q is not empty.The coefficient matrix of the system of linear inequalities defining Q is a −1, 0, 1-matrix.The absolute values of the subdeterminants of the coefficient matrix of the system of constraints defining Q are bounded by (2m) n because each row of this coefficient matrix contains no more than 2m nonzero entries.Consider a matrix obtained from the coefficient matrix by replacing one of the columns by the −1, 0-vector of the righthand sides of the constraints defining Q.The absolute values of the subdeterminants of this matix are bounded by (2m + 1) n .Standard linear algebra implies that a minimal face of Q contains a vector a such that a ∞ ≤ (2m + 1) n and a ∞ ≥ 1/(2m) n .
For every x ∈ S and y ∈ S, we have a T x < a T y if and only if c T x < c T y.It follows that ∀x ∈ S : arg min Let f be defined as f (x) = a T x.If f had a unique local optimum in N (x) for every x, then we could conclude immediately that the greedy algorithm would construct the same sequence of solutions when applied to f .In the general case, we need to construct a new local optimization oracle which, under the same circumstances, outputs exactly the same result for f as that for f .Thus, let ϕ : S −→ S be defined as follows: For each x ∈ S, let ϕ(x) := x if the original local optimization oracle decides that x is optimal for f and, otherwise, ϕ(x) := y, where y is a better solution returned by the original local optimization oracle with respect to f .Define a new local optimization oracle in the following way: Input: x ∈ S and g ∈ C; Output: y ∈ S such that g(y) < g(x) or a decision that x is optimal for g.
• if g = f then run the original oracle to minimize g over N (x); • else if ϕ(x) = x, then decide that x is optimal, else return ϕ(x).
Given x ∈ S and a linear function g ∈ C, the new oracle returns a solution in arg min y∈N (x) g(x).Indeed, this is the case if g = f .Otherwise, this is implied by ( Let the greedy algorithm use the new oracle and be applied to the new problem, where the objective is to minimize f over x ∈ S. In this case, the greedy algorithm constructs exactly the same sequence of solutions as when applied to the original problem, with the original oracle.To estimate the oracle running time for the new problem, set γ = 1/(2(2m) n ) and apply Theorem 4.1.Note that a T x 0 ≤ n(2m + 1) n and the optimal value of the new problem is not less than −n(2m + 1) n because |a T x| ≤ n a ∞ for all x ∈ S. Then Theorem 4.1 implies that the greedy algorithm runs in oracle time O(mn log n) for the new problem.Since the sequence of solutions considered in this case by the greedy algorithm is exactly the same as that for the original problem, it follows that the greedy algorithm solves the original problem in the same oracle time using the original oracle.
Remark.It should be emphasized that Corollary 4.2 does not assume any modification of the objective function, i.e., the greedy algorithm is applied to the original objective function, in contrast to the method of Frank and Tardos [5] which relies on simultaneous Diophantine approximation.Schulz et al. [11] also apply simultaneous Diophantine approximation at a preprocessing step.
If P is a nondegenerate 0,1-polytope being the set of feasible solutions of a linear program of the standard form, then the above complexity bounds hold true for the number of vertices visited by the simplex method equipped with Dantzig's best improvement rule.In this case, the bound of Corollary 4.2 follows from the analysis of the simplex method proposed by Kitahara and Mizuno [6].In the next section, we will consider a variant of the simplex method, for the general case, which can use any pivot rule modified so that the simplex method using it becomes our scaling algorithm.

Optimization of functions over vertex sets of polytopes
Let the following system define a polytope P : where A ∈ Z m×n and b ∈ Z m .In this section, S is assumed to be the vertex set of the polytope P.
Let Z (x) denote the index set of zero components of a vector x and 1 J denote the characteristic vector of a set J .
Let g y,δ be defined in the following way: Such a function has the following interpretation: when moving from y to x, we additionally pay δ • x i for each component x i such that y i = 0. Since x ≥ 0 for all x in P, the value added to f (x) is nonnegative: At this stage, in order to apply our scaling algorithm, we should still prove that (2.2) takes place for c y = 1 Z (y) , i.e., that the added value is positive if x = y.Let α > 0 be a lower bound and β > 0 be an upper bound on the nonzero components of the vertices of P, the polytope defined by (5.1).

Observation 5.1
The following statement is true: ∀x, y ∈ S, x = y : 1 T Z (y) (x − y) ≥ α. (5.4) Proof If y is a vertex of P, then the system Az = b, z i = 0, ∀i ∈ Z (y), (5.5) uniquely identifies y, i.e., the only solution of (5.5) is y.Then for each x ∈ S, x = y, there exists i ∈ Z (y) with x i ≥ α.Since x is nonnegative, (5.3) implies (5.4).
Thus, by the above observation we have the property (2.2) and can apply the scaling algorithm.

Theorem 5.1 If S is the vertex set of P, then the scaling algorithm solves the problem
where m is the number of rows of A.
Let be an upper bound on the absolute values of the subdeterminants of A. Proof For every vertex x ∈ S, let N (x) denote the set consisting of x and the vertices of P adjacent to x.For any linear objective function, a vertex x is optimal if and only if it is optimal in N (x).Let the augmentation oracle be an algorithm which either proves that x is optimal in N (x) or delivers a solution in N (x) which is better than x for a given linear objective function.Such an augmentation oracle can e.g.be based on any anticycling pivot rule according to which an entering variable has a negative reduced cost.Then, the augmentation oracle can take the following form: Given x ∈ S and g ∈ C, choose a basis defining x and use the pivot rule until a new basic feasible solution is constructed or it is proved that x is optimal for g.(The polytope P can be degenerate and the pivot rule may need more than one iteration to reach another vertex.)With such an oracle, the scaling algorithm becomes a variant of the simplex method.During each call to the oracle, the pivot rule can start with any basis defining the given vertex x or use the current basis obtained at the previous call to the oracle (in this case the scaling algorithm should store the respective basis after each call to the oracle).Let x and y, x = y, be adjacent vertices of P. By Observation 5.1, It is clear that c T x 0 ≤ β c 1 and O PT ≥ −β c 1 .Note that γ = 1 2 2 is a strict lower bound on |c T (x − y)| where x, y ∈ S, and c T x = c T y.Then Theorem 5.1 implies the required complexity bound.
The coefficient of the logarithmic term in the bound (5.7) does not depend on the dimension.That is, our approach can be useful when developing column generation algorithms, when the number of columns of A is much larger than the number of rows of A.
Obviously, the above complexity bound implies the following bound on diameters of polytopes: Proof Consider a vertex x 1 of P.This vertex uniquely minimizes the linear function 1 T Z (x 1 ) x over P. Now consider the problem of minimizing this function over the vertex set of P. The objective value at x 1 is zero while the objective values of the other vertices are not less than α, which follows from (5.4).The length of the path found by the variant of the simplex method described in the proof of Theorem 5.2, when starting from another vertex x 0 , is bounded by (5.8) according to (5.7).The other end of the path is x 1 because x 1 is the unique optimal solution for the linear function 1 T Z (x 1 ) x. Thus, we obtain the bound (5.8).
The bounds (5.7) and (5.8) depend on the dimension at most logarithmically (to see this for (5.7), observe that c 1 ≤ n c ∞ ).Recently, Kitahara and Mizuno [6] proposed a bound of the form O nm β α log β α , which depends linearly on the dimension.Bonifas et al. [1] proposed a bound of O( 2 n 3.5 log n ) for polytopes in an n-dimensional space, where is an upper bound on the absolute values of the subdeterminants of the coefficient matrix of a system of linear inequalities defining the respective polytope.This bound does not depend on the number of inequalities, but the dependence of this bound on is at least quadratic.For totally unimodular coefficient matrices, this bound is strongly polynomial and substantially improves the bound obtained by Dyer and Frieze [4].In our case, if A is totally unimodular, the estimate (5.8) can be competitive only if the ratio β/α is bounded by a polynomial of n whose degree is not so high.

Conclusions
We have proved that augmentation is polynomially equivalent to optimization for arbitrary functions over sets of binary vectors.Also, we have given complexity bounds for optimization of arbitrary functions over vertex sets of polytopes.As a consequence, we have obtained new bounds on diameters of polytopes and guaranteed bounds on the running time which can be achieved by the simplex method.For instance, the simple method introduced in the paper allows to estimate the quality of the current solution in the course of column generation methods, the respective estimates depending only logarithmically on the number of columns of the coefficient matrix.

Theorem 3 . 1
If S is a set of 0-1 vectors, then the scaling algorithm runs in oracle time bounded by O m log m( f (x 0 ) − O PT ) γ .

Theorem 4 . 1
Let (4.1) hold for all g ∈ C.Then, if S is a set of 0-1 vectors, the greedy algorithm runs in oracle time O m log m( f (x 0 ) − O PT ) γ .

Theorem 5 . 2
Let c ∈ Z n .There exists a variant of the simplex method which solves the linear program minimize c T x subject to(5.1)