Linear Programming
Introduction to AlgorithmsChapter 29: Linear Programming
Many problems take the form of optimizing a linear function subject to linear constraints on resources. If the objective and constraints are all linear, we have a linear-programming problem. Linear programs arise in operations research, scheduling, network flows, and many other practical domains.
Key definitions:
- A linear function on variables x₁, x₂, …, xₙ: f(x₁, …, xₙ) = a₁x₁ + a₂x₂ + ⋯ + aₙxₙ
- A linear equality: f(x₁, …, xₙ) = b
- A linear inequality: f(x₁, …, xₙ) ≤ b or f(x₁, …, xₙ) ≥ b
- Linear constraints: either linear equalities or linear inequalities
- Strict inequalities are not allowed in linear programming
A linear-programming problem is the problem of either minimizing or maximizing a linear function subject to a finite set of linear constraints.
Algorithms for LP:
- Simplex algorithm — classical method, not polynomial worst-case but fast in practice (invented by Dantzig, 1947)
- Ellipsoid algorithm — first polynomial-time algorithm (Khachian, 1979), slow in practice
- Interior-point methods — polynomial-time, competitive with simplex on large inputs
- Integer linear programming (requiring integer solutions) is NP-hard
29.1 Standard and Slack Forms
Standard Form
In standard form, we are given:
- n real numbers c₁, c₂, …, cₙ (objective coefficients)
- m real numbers b₁, b₂, …, bₘ (constraint bounds)
- mn real numbers aᵢⱼ (constraint coefficients)
We wish to find x₁, x₂, …, xₙ that:
maximize Σⱼ cⱼxⱼ (objective function)
subject to Σⱼ aᵢⱼxⱼ ≤ bᵢ for i = 1, 2, …, m (constraints)
xⱼ ≥ 0 for j = 1, 2, …, n (nonnegativity constraints)
Compact (matrix) notation using matrix A = (aᵢⱼ), vectors b = (bᵢ), c = (cⱼ), x = (xⱼ):
maximize cᵀx
subject to Ax ≤ b
x ≥ 0
A standard-form LP is specified by the tuple (A, b, c).
Terminology:
- Feasible solution: a setting x̄ that satisfies all constraints
- Infeasible solution: a setting that violates at least one constraint
- Objective value: cᵀx̄ for a given solution x̄
- Optimal solution: a feasible solution with maximum objective value
- Optimal objective value: the objective value of an optimal solution
- Infeasible LP: no feasible solutions exist
- Unbounded LP: feasible solutions exist but no finite optimal objective value
Converting to Standard Form
Any LP can be converted to standard form. Four possible issues and their fixes:
-
Minimization instead of maximization: Negate the objective function coefficients.
- minimize c₁x₁ + c₂x₂ → maximize −c₁x₁ − c₂x₂
-
Variables without nonnegativity constraints: Replace each unconstrained variable xⱼ with two new nonneg variables: xⱼ = x'ⱼ − x''ⱼ where x'ⱼ, x''ⱼ ≥ 0.
-
Equality constraints: Replace f(x) = b with two inequalities: f(x) ≤ b and f(x) ≥ b.
-
Greater-than-or-equal-to constraints: Multiply through by −1 to convert ≥ to ≤.
Two LPs are equivalent if their feasible solutions correspond with equal (or negated, for min↔max) objective values.
Slack Form
To use the simplex algorithm, we convert to slack form where all constraints (except nonnegativity) are equalities.
For each inequality constraint Σⱼ aᵢⱼxⱼ ≤ bᵢ, introduce a slack variable s:
s = bᵢ − Σⱼ aᵢⱼxⱼ
s ≥ 0
We use xₙ₊ᵢ as the slack variable for the i-th constraint:
xₙ₊ᵢ = bᵢ − Σⱼ aᵢⱼxⱼ (along with xₙ₊ᵢ ≥ 0)
Key concepts:
- Basic variables: variables on the left-hand side of the equalities (set B, |B| = m)
- Nonbasic variables: variables on the right-hand side (set N, |N| = n)
- N ∪ B = {1, 2, …, n + m}
A slack form is concisely defined by the tuple (N, B, A, b, c, ν):
z = ν + Σⱼ∈N cⱼxⱼ (objective function)
xᵢ = bᵢ − Σⱼ∈N aᵢⱼxⱼ for i ∈ B (constraints)
All variables are constrained to be nonneg. Note: the values aᵢⱼ stored are the negatives of the coefficients as they appear in the slack form (because we subtract the sum).
Basic solution: set all nonbasic variables to 0, then xᵢ = bᵢ for each i ∈ B, and z = ν.
29.2 Formulating Problems as Linear Programs
Once a problem is cast as a polynomial-sized LP, it can be solved in polynomial time via the ellipsoid algorithm or interior-point methods.
Shortest Paths
The single-pair shortest-path problem (source s, destination t, weight function w) can be formulated as:
maximize d_t
subject to d_v ≤ d_u + w(u,v) for each edge (u,v) ∈ E
d_s = 0
- Variables: d_v for each vertex v ∈ V (|V| variables)
- Constraints: |E| + 1
We maximize d_t because the shortest-path constraints force each d_v to be the largest value ≤ min{d_u + w(u,v)} over incoming edges. Minimizing would trivially yield d_v = 0 for all v.
Maximum Flow
Given directed graph G = (V, E) with capacities c(u,v), source s, sink t:
maximize Σ_v f_{sv} − Σ_v f_{vs}
subject to f_{uv} ≤ c(u,v) for each u,v ∈ V (capacity)
Σ_v f_{vu} = Σ_v f_{uv} for each u ∈ V−{s,t} (conservation)
f_{uv} ≥ 0 for each u,v ∈ V (nonnegativity)
- |V|² variables, 2|V|² + |V| − 2 constraints (can be reduced to O(V + E))
Minimum-Cost Flow
Given additionally a cost a(u,v) per unit of flow on each edge and a flow demand d:
minimize Σ_{(u,v)∈E} a(u,v) · f_{uv}
subject to f_{uv} ≤ c(u,v) for each u,v ∈ V
Σ_v f_{vu} − Σ_v f_{uv} = 0 for each u ∈ V−{s,t}
Σ_v f_{sv} − Σ_v f_{vs} = d
f_{uv} ≥ 0 for each u,v ∈ V
Multicommodity Flow
Given k commodities Kᵢ = (sᵢ, tᵢ, dᵢ), each with its own source, sink, and demand. All commodities share edge capacities. The aggregate flow on each edge must not exceed capacity.
minimize 0 (feasibility problem only)
subject to Σᵢ f_{i,uv} ≤ c(u,v) for each u,v ∈ V
Σ_v f_{i,vu} − Σ_v f_{i,uv} = 0 for each commodity i, each u ∈ V−{sᵢ,tᵢ}
Σ_v f_{i,sᵢ,v} − Σ_v f_{i,v,sᵢ} = dᵢ for each commodity i
f_{i,uv} ≥ 0 for each u,v ∈ V, each i
The only known polynomial-time algorithm for multicommodity flow is via LP.
29.3 The Simplex Algorithm
The simplex algorithm is the classical method for solving LPs. It operates on the slack form, iteratively improving the objective value by moving between vertices of the feasible region.
Core Idea
Each iteration is associated with a basic solution (set nonbasic variables to 0, compute basic variables from equality constraints):
- Choose a nonbasic variable xₑ with positive coefficient in the objective (increasing it improves the objective)
- Increase xₑ until some basic variable hits 0 (the tightest constraint)
- Pivot: exchange the roles of the entering variable xₑ and the leaving variable xₗ
- Repeat until no coefficient in the objective is positive → optimum found
If no constraint limits the increase of xₑ → the LP is unbounded.
Worked Example
Starting LP in standard form:
maximize 3x₁ + x₂ + 2x₃
subject to x₁ + x₂ + 3x₃ ≤ 30
2x₁ + 2x₂ + 5x₃ ≤ 24
4x₁ + x₂ + 2x₃ ≤ 36
x₁, x₂, x₃ ≥ 0
Slack form (introducing x₄, x₅, x₆):
z = 3x₁ + x₂ + 2x₃
x₄ = 30 − x₁ − x₂ − 3x₃
x₅ = 24 − 2x₁ − 2x₂ − 5x₃
x₆ = 36 − 4x₁ − x₂ − 2x₃
Basic solution: (0, 0, 0, 30, 24, 36), z = 0.
Iteration 1: Enter x₁ (positive coeff 3). Tightest constraint: x₆ limits x₁ to 9. Leave x₆.
z = 27 + x₂/4 + x₃/2 − 3x₆/4
x₁ = 9 − x₂/4 − x₃/2 − x₆/4
x₄ = 21 − 3x₂/4 − 5x₃/2 + x₆/4
x₅ = 6 − 3x₂/2 − 4x₃ + x₆/2
Basic solution: (9, 0, 0, 21, 6, 0), z = 27.
Iteration 2: Enter x₃ (positive coeff 1/2). Tightest: x₅ limits x₃ to 3/2. Leave x₅.
z = 111/4 + x₂/16 − x₅/8 − 11x₆/16
x₁ = 33/4 − x₂/16 + x₅/8 − 5x₆/16
x₃ = 3/2 − 3x₂/8 − x₅/4 + x₆/8
x₄ = 69/4 + 3x₂/16 + 5x₅/8 − x₆/16
Basic solution: (33/4, 0, 3/2, 69/4, 0, 0), z = 111/4.
Iteration 3: Enter x₂ (positive coeff 1/16). Tightest: x₃ limits x₂ to 4. Leave x₃.
z = 28 − x₃/6 − x₅/6 − 2x₆/3
x₁ = 8 + x₃/6 + x₅/6 − x₆/3
x₂ = 4 − 8x₃/3 − 2x₅/3 + x₆/3
x₄ = 18 − x₃/2 + x₅/2
All objective coefficients are now negative → optimal.
Optimal solution: x₁ = 8, x₂ = 4, x₃ = 0, objective value = 28.
Pivoting
The PIVOT procedure takes a slack form (N, B, A, b, c, ν), a leaving variable xₗ, and an entering variable xₑ, and returns the new slack form.
PIVOT(N, B, A, b, c, ν, l, e)
// Compute coefficients for new basic variable xₑ
b̂ₑ = bₗ / aₗₑ
for each j ∈ N − {e}
âₑⱼ = aₗⱼ / aₗₑ
âₑₗ = 1 / aₗₑ
// Compute coefficients of remaining constraints
for each i ∈ B − {l}
b̂ᵢ = bᵢ − aᵢₑ · b̂ₑ
for each j ∈ N − {e}
âᵢⱼ = aᵢⱼ − aᵢₑ · âₑⱼ
âᵢₗ = −aᵢₑ · âₑₗ
// Compute new objective function
ν̂ = ν + cₑ · b̂ₑ
for each j ∈ N − {e}
ĉⱼ = cⱼ − cₑ · âₑⱼ
ĉₗ = −cₑ · âₑₗ
// Update basic and nonbasic sets
N̂ = N − {e} ∪ {l}
B̂ = B − {l} ∪ {e}
return (N̂, B̂, Â, b̂, ĉ, ν̂)
Lemma 29.1 (Effect of PIVOT on basic solution):
- x̄ⱼ = 0 for each j ∈ N̂
- x̄ₑ = bₗ / aₗₑ
- x̄ᵢ = bᵢ − aᵢₑ · b̂ₑ for each i ∈ B̂ − {e}
The Formal Simplex Algorithm
SIMPLEX(A, b, c)
(N, B, A, b, c, ν) = INITIALIZE-SIMPLEX(A, b, c)
let δ be a new vector of length m
while some index j ∈ N has cⱼ > 0
choose an index e ∈ N for which cₑ > 0
for each index i ∈ B
if aᵢₑ > 0
δᵢ = bᵢ / aᵢₑ
else δᵢ = ∞
choose an index l ∈ B that minimizes δₗ
if δₗ == ∞
return "unbounded"
else (N, B, A, b, c, ν) = PIVOT(N, B, A, b, c, ν, l, e)
for i = 1 to n
if i ∈ B
x̄ᵢ = bᵢ
else x̄ᵢ = 0
return (x̄₁, x̄₂, …, x̄ₙ)
Key properties:
- Entering variable xₑ: a nonbasic variable with cₑ > 0
- Leaving variable xₗ: the basic variable that most tightly limits the increase of xₑ
- If no constraint limits xₑ (all δᵢ = ∞), the LP is unbounded
- When all cⱼ ≤ 0 for j ∈ N, the current basic solution is optimal
Correctness (Lemma 29.2)
If INITIALIZE-SIMPLEX returns a feasible slack form, then:
- If SIMPLEX returns a solution in line 17, it is a feasible solution
- If SIMPLEX returns "unbounded," the LP is indeed unbounded
Loop invariant (maintained each iteration):
- The slack form is equivalent to the initial one
- bᵢ ≥ 0 for each i ∈ B
- The basic solution is feasible
Termination and Degeneracy
- Degeneracy: a pivot that leaves the objective value unchanged (occurs when bₗ = 0)
- Cycling: degeneracy can cause the simplex algorithm to revisit the same slack form, looping forever
- Cycling is the only reason SIMPLEX might not terminate
Lemma 29.4: The slack form is uniquely determined by the set B of basic variables.
Lemma 29.5: If SIMPLEX fails to terminate in at most C(n+m, m) iterations, it cycles.
Bland's rule: Always choose the variable with the smallest index when breaking ties for entering and leaving variables. This prevents cycling (Lemma 29.6).
Lemma 29.7: With Bland's rule, SIMPLEX either reports "unbounded" or terminates with a feasible solution in at most C(n+m, m) iterations.
29.4 Duality
Duality provides a way to prove that a solution is optimal. Given a maximization problem (the primal), we define a related minimization problem (the dual) whose optimal value equals that of the primal.
Constructing the Dual
Given a primal LP in standard form:
Primal:
maximize Σⱼ cⱼxⱼ
subject to Σⱼ aᵢⱼxⱼ ≤ bᵢ for i = 1, …, m
xⱼ ≥ 0 for j = 1, …, n
The dual LP is:
Dual:
minimize Σᵢ bᵢyᵢ
subject to Σᵢ aᵢⱼyᵢ ≥ cⱼ for j = 1, …, n
yᵢ ≥ 0 for i = 1, …, m
Transformation rules:
- Maximization → minimization
- Objective coefficients cⱼ ↔ right-hand sides bᵢ
- ≤ constraints → ≥ constraints
- Each primal constraint i has a dual variable yᵢ
- Each dual constraint j has a primal variable xⱼ
Weak Duality (Lemma 29.8)
For any feasible primal solution x̄ and any feasible dual solution ȳ:
Σⱼ cⱼx̄ⱼ ≤ Σᵢ bᵢȳᵢ
The primal objective value never exceeds the dual objective value.
Corollary 29.9: If a feasible primal solution x̄ and a feasible dual solution ȳ have equal objective values, then both are optimal for their respective LPs.
Strong Duality (Theorem 29.10 — Linear-Programming Duality)
If SIMPLEX returns a solution x̄ for the primal LP (A, b, c), then there exists a dual solution ȳ such that:
Σⱼ cⱼx̄ⱼ = Σᵢ bᵢȳᵢ
Both x̄ and ȳ are optimal for their respective LPs.
Reading the dual solution from the final slack form:
Let the final slack form have nonbasic set N and coefficients c'ⱼ. Then:
ȳᵢ = −c'_{n+i} if (n + i) ∈ N
ȳᵢ = 0 otherwise
The negatives of the slack variable coefficients in the final objective function give the optimal dual values.
Example: For the LP with optimal z = 28 − x₃/6 − x₅/6 − 2x₆/3:
- ȳ₁ = 0 (since n+1 = 4 ∈ B)
- ȳ₂ = −c'₅ = 1/6
- ȳ₃ = −c'₆ = 2/3
- Dual objective: 30(0) + 24(1/6) + 36(2/3) = 28 ✓
29.5 The Initial Basic Feasible Solution
The Problem
The initial basic solution (set all original variables to 0) may not be feasible if some bᵢ < 0. We need INITIALIZE-SIMPLEX to either:
- Determine the LP is infeasible, or
- Produce a slack form with a feasible basic solution
The Auxiliary Linear Program
Given primal LP L in standard form, introduce a new variable x₀ and form Lₐᵤₓ:
maximize −x₀
subject to Σⱼ aᵢⱼxⱼ − x₀ ≤ bᵢ for i = 1, …, m
xⱼ ≥ 0 for j = 0, 1, …, n
Lemma 29.11: L is feasible if and only if the optimal objective value of Lₐᵤₓ is 0.
INITIALIZE-SIMPLEX Procedure
INITIALIZE-SIMPLEX(A, b, c)
let k be the index of the minimum bᵢ
if b_k ≥ 0 // initial basic solution is feasible
return ({1,…,n}, {n+1,…,n+m}, A, b, c, 0)
// Form auxiliary LP: add −x₀ to each constraint, objective = −x₀
form L_aux by adding −x₀ to the left-hand side of each constraint
and setting the objective function to −x₀
let (N, B, A, b, c, ν) be the resulting slack form for L_aux
l = n + k
// Pivot to make basic solution feasible for L_aux
(N, B, A, b, c, ν) = PIVOT(N, B, A, b, c, ν, l, 0)
// Solve L_aux using simplex
iterate the while loop of SIMPLEX until optimal solution to L_aux is found
if optimal solution sets x̄₀ = 0
if x₀ is basic
perform one degenerate pivot to make it nonbasic
remove x₀ from constraints
restore original objective function of L,
replacing each basic variable by the right-hand side of its constraint
return the modified final slack form
else
return "infeasible"
Key steps:
- If all bᵢ ≥ 0, the initial basic solution is already feasible — return immediately
- Otherwise, form Lₐᵤₓ and pivot x₀ in for the most negative basic variable
- After the pivot, all bᵢ become nonneg → feasible basic solution for Lₐᵤₓ
- Solve Lₐᵤₓ with simplex
- If optimal value is 0 → L is feasible; remove x₀ and restore original objective
- If optimal value is negative → L is infeasible
Correctness (Lemma 29.12)
- If L has no feasible solution → INITIALIZE-SIMPLEX returns "infeasible"
- If L has a feasible solution → returns a valid slack form with a feasible basic solution
The single pivot in step 3 works because:
- The leaving variable xₗ has the most negative bₗ
- aₗₑ = −1 (coefficient of x₀), so b̂ₑ = bₗ/(−1) > 0
- For all other i: b̂ᵢ = bᵢ − aᵢₑ(bₗ/aₗₑ) = bᵢ − bₗ ≥ 0 (since bₗ is the minimum)
Fundamental Theorem of Linear Programming
Theorem 29.13: Any LP L in standard form either:
- Has an optimal solution with a finite objective value
- Is infeasible
- Is unbounded
SIMPLEX correctly handles all three cases:
- If infeasible → returns "infeasible" (via Lemma 29.12)
- If unbounded → returns "unbounded" (via Lemma 29.2)
- Otherwise → returns an optimal solution (via Theorem 29.10 and Lemma 29.7)
Previous chapter
Matrix OperationsNext chapter
Polynomials and the FFT