Linear Programming

 

Cynthia G. Huyser, CS 5329

November 27, 2000

 

Abstract

 

Linear programming is a means by which a real-valued objective function expressible as a linear equation can be optimized subject to a set of constraints.  Applications include scheduling and materials use problems.  Linear programming is a special case of nonlinear programming, and the most efficient known algorithms for worst-case linear problems employ nonlinear techniques.  The classical Simplex algorithm solves many problems efficiently, but is known to have an exponential worst-case running time.  Since the 1970’s a number of polynomial worst-time algorithms have appeared.

 

 

1.     Introduction.

 

 

The concept of the linear program (see, e.g., Dantzig and Thapa, 1997) was introduced by George B. Dantzig in 1947.  Linear programming is a method of optimizing an objective function by solving a system of linear equations where the solution is subject to a set of constraints (see, e.g. Danzig 1949).  Prior to its introduction, optimization of problems involving such constraint matrices had not been undertaken.  Linear programming is restricted to problems where real-valued solutions are acceptable and where the equations and constraints are linear. Linear programming is a subset of mathematical programming, which includes integer programming, nonlinear programming, stochastic programming, combinatorial programming, and network flow maximization (Dantzig and Thapa, 1997).  Applications of linear programming include materials and schedule optimization (Dantzig and Thapa, 1997), and determination of fuzzy weighted averages (Yuh-Yuan, et al., 2001).

 

The objective function is posed in terms of economy in a matrix of activities, commodities, and the flow of commodities in a given time span. Dantzig (1949) defined the function as a set of activities A such that

 A0 + x1A1 + x2A2 + … + xnAn = 0                   xj>=0.                          [Equation 1]1

A feasible problem is given by a set of real values xj>=0 that satisfies this equation; if no such set of activities can be found, the problem is not feasible (Dantzig 1949).  If a system of linear inequalities is written so that all the inequalities are of the same type (either “greater than” or “less than”), and the equations can be manipulated to yield an equation of the form

0x1 + 0x2 + 0x3 + … + 0xn >= d, where d > 0,                         [Equation 2]

 

the set of inequalities cannot be solved, and is termed infeasible (Dantzig and Thapa, 1997).  An optimal solution for Equation 1 will minimize the number of activities at positive levels (xj>0)  and maximize those reduced to zero (xj=0)  (Dantzig, 1949).   While Dantzig’s original problem minimizes a linear function, a maximization constraint can also be satisfied.  Additionally, any linear minimization problem can be rewritten as a maximization problem (von Neumann, cited in Dantzig and Thapa, 1997).  The original problem is called the primal, and the rewriting its dual.  Systems of linear equalities and the question of their feasibility can be solved in polynomial time (see Section 3 of this paper).

 

2.     The Simplex Algorithm.

 

Dantzig’s original method for solving linear programming problems is known as the Simplex Method (see, e.g., Dantzig and Thapa, 1997).  This approach finds solutions to these problems by exploring the solution set along its boundaries, from vertex (extreme point) to vertex of the area bounded by the linear constraints. An optimal solution is determined (or is found not to exist) upon the application of a series of “pivot steps” to a system of equations in “standard form2.”  Since pivoting inserts and removes redundant equations, no net change occurs the system’s solution set (Dantzig and Thapa, 1997).

 

The worst-case complexity for the Simplex algorithm is exponential (Dantzig and Thapa, 1997). For example, Zadeh (1972) identifies a transportation problem containing (2n + 2) nodes requiring (2n + 2n-2 – 2) iterations of the Simplex algorithm.  In practice, the matrices encountered are usually much sparser than in this worst case, and the Simplex algorithm can efficiently solve large practical problems (Dantzig and Thapa, 1997). Gill, et al. (1986) state that the complexity of the Simplex algorithm is “roughly linear in the problem dimension,” i.e., the number of variables (n).  The Simplex algorithm also has a degenerative case, in which the value of a basic variable becomes 0 and the choice of pivot cycles (i.e., a pivot element previously chosen is revisited); randomizing the pivot chain can relieve this problem (Dantzig and Thapa, 1997).

 

3.     Polynomial-time algorithms.

 

Each of the polynomial-time algorithms presented in this paper is based on an interior point method.  While the Simplex algorithm searches for solutions along the vertices of the boundary created by the linear constraints, interior point methods begin at a point within the feasible region (i.e., the interior of the region bounded by the constraints), and proceed outward toward the minimum (maximum) constraint (Dantzig and Thapa, 1997).

 

The first polynomial-time algorithm for systems of linear equations was given by Khachiyan (1979).  In the worst case (i.e., where no feasible solution can be found), the time complexity of this algorithm is O(n6L2), where L is the length of the binary encoded linear system3, n is the number of elements in a row, and m is the number of rows in the system of equations (Dantzig and Thapa, 1997).  Gill, et al. (1986) note that though Khachiyan’s algorithm provides an improvement over the worst-case time complexity of the Simplex algorithm, it rarely demonstrates improved running time in practice.

 

Karmarkar’s (1984) linear programming algorithm establishes an interior point method with a worst-case time of O(n3.5L2ln(L)ln(ln(L))).  In this method, the solution space is explored by mapping a series of projective spheres inscribing and circumscribing the polytope defined by the linear constraints4.    The number of steps in the algorithm depends on the ratio of the inscribed to the circumscribed sphere (“r/R”).  Karmarkar’s algorithm improves the objective function by i/(n-1) = “r/R” each iteration, but Anstreicher (1989) shows that the function can be improved by 2/n with each iteration, and shows the worst–case step improvement in the transformation will be 0.72148. 

 

Gill, et al. (1986) note that linear programming is a special case of nonlinear programming, and show that a “projected Newton barrier” method5 from nonlinear programming is equivalent to Karmarkar’s algorithm; algorithms demonstrating better time complexity than Karmarkar’s use nonlinear programming approaches.  Gonzaga (1987) gives an algorithm with worst-case O(n3L) performance, while Anstreicher (1999) gives an O((n3/ln n)L) algorithm.  Anstreicher’s (1999) path-following algorithm uses a preconditioned conjugate gradient (PCG) method6 based on the work of Nesterov and Nemirovskii (1994) with a barrier method, and improves performance with the careful choice of a parameter r (with respect to the relationship between m and n).  Anstreicher describes the algorithm in terms of the solution of a linear dual, but notes that the method demonstrated here can be applied to improve the complexity of most O(n3) partial updating algorithms for linear programming.

 

4.     Conclusions.

 

While the Simplex algorithm provides a practical solution in most cases, its exponential time complexity makes it infeasible for dense problem matrices, i.e. ones where most rows are filled with non-zero elements.  More recent algorithms have polynomial worst-case complexity and are better suited for dense problems.  Because Karmarkar’s algorithm and subsequent methods can be seen as special cases of nonlinear programming, it appears that the use of nonlinear techniques specialized to linear programming problems may be the most efficient means of solving such problems.  To date, the best order of complexity for such algorithms is that obtained by Anstreicher, O((n3/ln(n))L).

References

 

Anstreicher, K. M. (1989). The worst-case step in Karmarkar’s algorithm. Mathematics of Operations Research. 14. 294 – 1989.

 

Anstreicher, K. M. (1999). Linear programming in O((n3/ln n)L) operations. SIAM Journal of Optimization. 9.  803-812.

 

Arbell, A. (1993).  Exploring Interior-Point Linear Programming.  Cambridge, Massachusetts.  The MIT Press.

 

Aspvall, B. and Stone, R. E. (1980). Khachiyan’s linear programming algorithm. Journal of Algorithms. 1. 1-13.

 

Dantzig, G. B. (1949). Programming of interdependent activities: II  mathematical model. Econometrica. 17. 200 – 211.

 

Dantzig, G. B. and Thapa, M. N. (1997). Linear Programming 1: Introduction. New York. Springer.

 

Gill, P. E., Murray, W., Saunders, M. A., Tomlin, J. A., and Wright, M. H. (1986). On projected Newton barrier methods for linear programming and an equivalence to Karmarkar’s projective method.  Mathematical Programming. 26.  183-209.

 

Gonzaga, C. C. (1987). An algorithm for solving linear programming problems in O(n3L) operations.  In N. Megiddo, (Ed.). Progress in Mathematical Programming. pp. 1 – 28. New York.  Springer-Verlag.

 

Guh Y., Hon C., and Lee E. S. (to appear in 2001). Fuzzy weighted average: the linear programming approach via Charnes and Cooper’s Rule.  Fuzzy Sets and Systems. 117 . 157-160.

 

Karmarkar, N. (1984). A new polynomial algorithm for linear programming. Combinatorica. 4. 373-395.

 

Khachiyan [Hajican], L. G. (1979). A polynomial algorithm for linear programming. Soviet Mathematics Doklady. 20. 191-194.

 

Nesterov, Y. and Nemirovskii, A. (1994). Interior-Point Polynomial Algorithms in Convex Programming. Philadelphia. Society for Industrial and Applied Mathematics.

 

Zadeh, N. (1972).  A bad network problem for the Simplex method and other minimum cost flow algorithms. Mathematical Programming. 5. 255 – 266.

 

Notes

 

1.     Equation 1 can also be written

minimize cTx

subject to                                                              

(see, e.g., Dantzig and Thapa, 1997). 

 

2.     Dantzig (1949) defines a minimization LP in standard form as:

 

c1x1 + c2x2  + … + cnxn = z (Min)                    where x1 >= 0, x2>=0, …, xn>=0

a11x1 + a12x2 + … + a1nxn = b1                                          

a21x1 + a22x2 + … + a2nxn = b2                                                                                   [Equation 3]

    .                 .                    .         .                                                   

    .                 .                    .         .

    .                 .                    .         .

am1x1 + am2x2 + … + amnxn = bm                  

                                                           

c1…cn in Equation 3 represents the solution for the equation, i.e. the quantities of each component x1…xn necessary to minimize z.  The linear equations that follow the first line represent the constraints on the solution.  The pivot element is a non-zero array element in a constraint equation that is used to derive equations that have the same solution set, but that do not contain the pivot element. All equations are divided by the amount necessary to make the coefficient of the pivot 1, and a multiple of the pivot equation is subtracted from each of the non-pivot equations so that the element corresponding to the pivot becomes 0. The pivot element, now unique in the equations, becomes a basic variable. When each non-solution row contains a basic variable, a solution is found for z by setting the non-basic variables to 0.  In a minimization problem, if all coefficients in the solution equation resulting from the last pivot are positive, the solution is optimal.  If not, an element in the solution equation having a negative coefficient is identified as the next pivot, and the pivoting process is repeated in such a way as to drive one of the previous pivots to 0; the zero pivot becomes non-basic, and the new pivot basic. 

 

3.    Khachiyan (1979) defines the length of a binary encoding of a system of m ³ 2 linear equalities for n ³ 2 real variables x1… xj … xn

 

ai1x1 + … + ainxn £ bi        i = 1, 2, …, m

 

with integer coefficients aij and bi to be

 

 

where lg(x) means log2(x).  So defined, L represents the sum of bits necessary to represent all constraint coefficients, plus the sum over each row of the bits required to write the integer coefficient bi plus the product of the row number and the number of variables in a row; the purpose of the final, additional bit is not clear from Khachiyan’s explication.

 

4.     The process is initialized with a point at the center of the polytope, which is then projected as the center of a unit sphere inscribed in the projected polytope.  A step is then taken in the direction of improving the objective function; when the point so derived is returned, the improvement it represents over the previous position is checked with the “potential function”

                                                [Equation 4]

 

If the improvement is less than half the size of the step taken to produce the new point, the solution is infeasible.  The optimality of feasible solutions is checked whenever the time since the last check exceeds the time to check it.

 

5.     The barrier method of nonlinear programming removes the inequality constraint (xi ³ 0) by adding a barrier term (Arbel, 1993), which decreases with each iteration (Gonzaga, 1987).  The barrier term maintains xi near the center of the polytope, and points generated with this method are said to follow a central trajectory (Arbel, 1993).  Algorithms that attempt to follow this central trajectory are known as path-following algorithms (Arbel, 1993).

 

6.     The preconditioned conjugate gradient method applies a conjugate gradient algorithm to preconditioned input (Anstreicher, 1999).  A conjugate gradient algorithm is an iterative method that involves guessing an initial solution, evaluating the error generated by it, and iteratively correcting it (Arbel, 1993).  The preconditioning suggested in Anstreicher (1999) is the application of a symmetric transformation to a matrix that stands for a system of equations.  Anstreicher (1999) describes a preconditioning matrix as the product

 

                        

 

where is an approximated diagonal solution matrix with values corresponding to the diagonal matrix S with values s1¼sn, and for r > 1,

 

,         i = 1, ¼, n, and r is not O(1).

 

The choice of r is critical to the improvement of complexity shown by this algorithm.  Anstreicher shows that if r is chosen between 0 and ½, the total number of operations is £ O, as opposed to O(n3L) when r is assumed to be 1.