Pagefy

Pagefy

Back to Data Structures and Algorithms

Linear Programming

Introduction to Algorithms

Chapter 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:

  1. Minimization instead of maximization: Negate the objective function coefficients.

    • minimize c₁x₁ + c₂x₂ → maximize −c₁x₁ − c₂x₂
  2. Variables without nonnegativity constraints: Replace each unconstrained variable xⱼ with two new nonneg variables: xⱼ = x'ⱼ − x''ⱼ where x'ⱼ, x''ⱼ ≥ 0.

  3. Equality constraints: Replace f(x) = b with two inequalities: f(x) ≤ b and f(x) ≥ b.

  4. 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):

  1. Choose a nonbasic variable xₑ with positive coefficient in the objective (increasing it improves the objective)
  2. Increase xₑ until some basic variable hits 0 (the tightest constraint)
  3. Pivot: exchange the roles of the entering variable xₑ and the leaving variable xₗ
  4. 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):

  1. x̄ⱼ = 0 for each j ∈ N̂
  2. x̄ₑ = bₗ / aₗₑ
  3. 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):

  1. The slack form is equivalent to the initial one
  2. bᵢ ≥ 0 for each i ∈ B
  3. 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:

  1. If all bᵢ ≥ 0, the initial basic solution is already feasible — return immediately
  2. Otherwise, form Lₐᵤₓ and pivot x₀ in for the most negative basic variable
  3. After the pivot, all bᵢ become nonneg → feasible basic solution for Lₐᵤₓ
  4. Solve Lₐᵤₓ with simplex
  5. If optimal value is 0 → L is feasible; remove x₀ and restore original objective
  6. 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:

  1. Has an optimal solution with a finite objective value
  2. Is infeasible
  3. 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)