Linear Programming: Method

Linear Programming: the Simplex Method, Duality, Sensitivity Analysis and Branch and Bound

The linear programming model

A linear program chooses the values of decision variables x₁, …, x_n that maximize or minimize a linear objective Z = c₁x₁ + … + c_nx_n subject to linear constraints, each of the form a_i1x₁ + … + a_inx_n ≤ b_i, ≥ b_i or = b_i, and usually to x_j ≥ 0. The coefficients c_j are unit profits or costs, the b_i are the resources available or the requirements to be met, and the a_ij are the amounts of resource i used per unit of activity j. The constraints describe a convex polyhedron, the feasible region; because the objective is linear, its optimum (if one exists) is attained at a corner point of that region, and the simplex method is a systematic walk from corner to corner, each step improving Z, until no neighbouring corner is better.

Standard (augmented) form

The simplex method works on equations. A ≤ constraint receives a slack variable s_i ≥ 0 that absorbs the unused resource; a ≥ constraint receives a surplus variable subtracted from its left-hand side; an equality receives none. A row whose right-hand side is negative is multiplied by −1 first, which reverses its inequality. A variable that may be negative (free) is written as the difference of two non-negative parts. Minimization is handled by maximizing −Z. With m equations in n + (slacks) unknowns, a basic solution fixes all but m variables at zero and solves for the m basic variables; when the basic variables are all non-negative it is a basic feasible solution, which is a corner point of the feasible region.

The simplex method

Each iteration starts from a basic feasible solution and its tableau, which expresses the basic variables and Z in terms of the nonbasic ones. The optimality test reads the Z row: in a maximization, a nonbasic variable whose entry z_j − c_j is negative would increase Z if raised from zero. If no entry is negative the current solution is optimal. Otherwise an entering variable is chosen, by convention the one with the most negative entry (Dantzig's largest-coefficient rule, the rule the textbook tableaus follow). The minimum ratio test decides how far it can be raised before a basic variable reaches zero: for every row with a positive coefficient in the entering column the ratio right-hand side / coefficient is formed, and the row with the smallest ratio gives the leaving variable. A Gauss-Jordan pivot on that element produces the next tableau. If the entering column has no positive coefficient, the entering variable can grow without bound and so can Z: the problem is unbounded, which in practice means a constraint is missing from the model.

Two-phase method

When the model has ≥ or = constraints the slack variables alone do not give a starting basic feasible solution. An artificial variable a_i ≥ 0 is added to each such row to play the role of a slack, and phase 1 minimizes their sum W (the tool maximizes −W) with the simplex method from this artificial basis. If the minimum of W is positive the original constraints cannot all be satisfied: the problem is infeasible. If it is zero a basic feasible solution of the real problem has been found; any artificial variable still in the basis at zero level is pivoted out (or its row is found redundant), and phase 2 continues with the original objective from that basis. The big-M method does the same with one combined objective; the two-phase form is numerically cleaner and is the one printed.

Degeneracy and Bland's rule

A basic feasible solution is degenerate when a basic variable is zero, which happens when more constraints pass through a corner than are needed to define it. A pivot can then leave the objective unchanged, and with an unlucky choice of entering and leaving variables the method can return to an earlier basis and cycle for ever. Bland (1977) proved that choosing, among the candidates, the entering variable with the smallest index and, in a tie of the ratio test, the leaving variable with the smallest index, makes cycling impossible. Bland's rule tends to need more iterations than the largest-coefficient rule, so the tool uses Dantzig's rule until a degenerate pivot occurs and Bland's rule from then on; the tableau remarks say which rule chose each entering variable. A degenerate optimum also means that the shadow prices may be different for an increase and for a decrease of a right-hand side, and the tool flags it.

Alternative optima

If at the optimum a nonbasic variable has a zero entry in the Z row, raising it changes the solution but not Z: every point on the edge between the current corner and the next one is optimal. The tool reports this as alternative optima. For a two-variable problem it shows on the drawing as the objective level line lying along an edge of the feasible region.

Duality

Every linear program (the primal) has a dual with one variable y_i per primal constraint and one constraint per primal variable. For a maximization with ≤ constraints and non-negative variables the dual is: minimize b₁y₁ + … + b_my_m subject to Σ_i a_ijy_i ≥ c_j for every j and y ≥ 0. The general rules (Winston, Table 6.x): in the dual of a max problem a ≤ constraint gives y_i ≥ 0, a ≥ constraint y_i ≤ 0, an equality a free y_i; a non-negative primal variable gives a ≥ dual constraint and a free primal variable an equality. For a min problem the dual is a max with the inequalities reversed. The weak duality theorem says that any feasible dual objective bounds any feasible primal objective, and the strong duality theorem that at the optimum the two are equal. The optimal dual values are read directly from the final primal tableau: they are the Z-row entries under the initial basic variables (the slack or artificial columns), which is also how the tool obtains them.

Shadow prices

The optimal dual value y_i is the shadow price of constraint i: the rate at which the optimal Z changes per unit increase of the right-hand side b_i, holding everything else fixed. It is the most the decision maker should pay for one more unit of that resource (or the saving from relaxing a requirement by one unit). A constraint with slack at the optimum has a zero shadow price, because more of a resource that is not fully used is worth nothing; a binding constraint usually has a positive price (negative for a constraint that works against the objective, such as a ≥ requirement in a maximization). The price holds only as long as the basis stays optimal, that is within the right-hand-side range of the sensitivity analysis.

Reduced costs

The reduced cost of a variable is c_j − z_j, with z_j = Σ_i a_ijy_i the value of the resources one unit of the activity would consume at the shadow prices. For a basic variable it is zero. For a nonbasic variable in a maximization it is negative (or zero with alternative optima): forcing one unit of the activity into the solution would lower Z by that amount, and the activity enters the optimal plan only if its unit profit rises by at least that amount. In a minimization the reduced cost is positive and is the amount the unit cost must fall. This is the sign convention of the Excel Solver sensitivity report; LINDO and some texts print the magnitude.

Sensitivity ranging

The sensitivity analysis asks how far one coefficient can change, the others fixed, before the present basis stops being optimal (Hillier and Lieberman ch. 7; Winston ch. 6).

Objective coefficient of a nonbasic variable. The variable enters the basis as soon as its reduced cost reaches zero: a maximization allows an increase of up to the magnitude of the reduced cost and any decrease.

Objective coefficient of a basic variable. A change Δ in c_j changes every Z-row entry of the final tableau by −Δ times the entry of that variable's row in the same column. The basis stays optimal while all entries keep their optimal sign, which gives a lower and an upper limit on Δ. Within the range the solution x does not change and Z changes by Δ x_j.

Right-hand side. The basic variables are x_B = B⁻¹b, with B⁻¹ the block of the final tableau under the initial basic columns. A change Δ in b_i moves them by Δ times column i of B⁻¹; a ratio test on those entries finds the largest increase and decrease that keep every basic variable non-negative. Within the range the basis and the shadow prices hold, the basic variables move linearly and Z changes by y_iΔ.

The ranges are stated as an allowable increase and an allowable decrease (as on the Solver report). The 100 percent rule covers simultaneous changes: if the sum over the changed coefficients of (change used)/(allowable change in that direction) does not exceed 1, the basis stays optimal; the rule is sufficient, not necessary, so a larger combined change may still be fine and must be checked by re-solving. For objective coefficients the rule guarantees the same solution; for right-hand sides it guarantees the same shadow prices.

Integer programming by branch and bound

When some variables must be integers (or binary, 0 or 1), the feasible set is no longer convex and the simplex method does not apply directly. The LP relaxation, the same problem with the integrality dropped, is solved first; its optimum bounds the integer optimum (from above in a maximization) and, if it happens to be integer, is the answer. Otherwise branch and bound (Land and Doig, 1960) picks a fractional integer variable with value v and creates two subproblems with the added constraint x ≤ ⌊v⌋ in one and x ≥ ⌈v⌉ in the other; between them they contain every integer point of the parent and exclude the fractional optimum. Each subproblem is an LP solved by the simplex method. A node is fathomed (not branched further) when its LP is infeasible, when its solution is integer (it becomes the incumbent if better than the current one), or when its LP value is no better than the incumbent, since no integer point below it can be better. The tool searches depth first, branching on the most fractional variable (the one nearest to a half), and prints the nodes explored. The tree can double at every level, so a problem with a few dozen integer variables can take thousands of LP solutions; the tool stops at 2,000 nodes or 20 seconds and reports the best solution found with a note that it is not proven optimal. Rounding the relaxation is not a substitute: the rounded point is often infeasible or far from the best integer point, as the California Manufacturing example of Hillier and Lieberman shows (relaxation (0.83, 1, 0, 1) with Z = 16.5; integer optimum (1, 1, 0, 0) with Z = 14).

The drawing for two variables

For two variables every constraint is a line and the feasible region the polygon on the admissible side of all of them. The objective Z = c₁x₁ + c₂x₂ is a family of parallel lines; the optimum is the last point of the region the level line touches as Z increases. The tool draws the constraint lines across a window slightly larger than the region, shades the polygon, marks the optimum and draws the level line through it; for an integer problem it adds the integer points of the region and the integer optimum.

Model syntax

The model is typed one statement per line. The first line is the objective: max or min (also maximize / minimize, optionally followed by Z =) and a linear expression. Terms may be written 3 x1, 3x1, 3*x1, -x2, + 2.5 y, 2/3 z; variable names are letters, digits and underscores starting with a letter. Each following line is a constraint: an optional label and colon (lumber: 8 x1 + 6 x2 <= 48), a linear expression, one of <=, >=, = (also =<, =>, ≤, ≥), and a linear expression; constants and variables may appear on both sides and are collected. Variables are non-negative by default, so x1, x2 >= 0 is accepted and ignored; x3 free or free x3 makes a variable unrestricted in sign; int x1, x2 makes variables integer and bin x3 binary (0 or 1, the bound x3 ≤ 1 is added as a constraint). A line with a single variable and a non-zero constant, such as x1 <= 10, is an ordinary constraint. Text after # is a comment, and lines reading subject to or s.t. are skipped. Up to 30 variables and 50 constraints; the tableaus are printed for up to six of each.

Limits

Dense floating-point arithmetic with a 10⁻⁹ tolerance: fine for the textbook and planning models the tool is meant for, not for badly scaled coefficients (10⁶ next to 10⁻⁶) or for hundreds of rows. Sensitivity ranges are those of the final basis found; at a degenerate optimum other bases give the same solution with different ranges and possibly different shadow prices. The ranges and shadow prices refer to the LP relaxation only when the model has integer variables. Branch and bound uses no cutting planes or presolve, so large integer models may hit the node limit.