Linear Programming

By the end of this lesson you will be able to read a linear program, put any of them into standard and slack form, solve a two-variable one by sliding a line across a picture, run pivot and simplex by hand on slack forms (with exact fractions, so nothing is rounded), start simplex from an infeasible-looking LP with initializeSimplex and its auxiliary LP, write down and read off the dual, and see how shortest paths and maximum flow become LPs — with every procedure also written in Dart and cross-checked against brute force.

What is a linear program? The 2-D picture

Picture a small workshop with limited hours of oven time and labour. You choose how many tables and chairs to build; each choice earns some profit, and the limited hours rule out some choices. A linear program (LP) is exactly that question in maths: choose numbers x₀, x₁, … to make a linear score as large (or small) as possible, while a list of linear rules (equalities and inequalities) stays true. "Linear" means variables are only multiplied by fixed numbers and added: 3x₀ + x₁ is linear, x₀x₁ or x₀² is not.

A point that obeys every rule is feasible; the score of a point is its objective value. There are exactly three outcomes (the fundamental theorem of linear programming): the LP is infeasible (no point obeys all rules), unbounded (feasible points exist but the score can grow forever), or has a finite optimal value. Real uses: airline crew scheduling, choosing oil-well sites, and the graph problems of the formulation section below. The algorithm of this lesson, simplex (Dantzig, 1947), is not polynomial in the worst case but is fast in practice.

With two variables you can see an LP. Every rule "ax₀ + bx₁ ≤ r" is a straight line with an allowed side; the points allowed by all rules form a convex polygon (the feasible region); and all points with the same score c₀x₀ + c₁x₁ = z lie on a straight line that slides across the plane as z grows. The best score is the last position at which that line still touches the polygon — always at a corner, or along an edge whose corners then tie. Everything simplex does later is a way of walking from corner to corner that works in any number of dimensions.

The geometric facts to keep: the feasible region is convex (it is an intersection of half-planes); an LP is infeasible when the region is empty and unbounded when it is open in a direction that improves the score; otherwise an optimum sits at a vertex. In n dimensions the "lines" become hyperplanes and a vertex is where n of them meet — hopeless to draw, but exactly what the algebra of slack forms tracks.

The fundamental theorem as code. fundamentalOutcome runs simplex and reports infeasible, unbounded, or "optimal at a vertex"; for an optimum it also checks that at most m of the n + m values (originals plus slacks) are nonzero, which is what "a vertex" means in slack form.

Input size → what's feasible. two variables (m ≤ 5 rows on this page): draw it. A few dozen variables and rows: the exact-fraction simplex below finishes in well under a second. Thousands of variables and rows: floating-point simplex or an interior-point solver, never vertex enumeration (C(n+m, m) vertices).

Why exact fractions? simplex compares numbers for equality all the time: "is c_j > 0?", "which delta[i] is the smallest, and is there a tie?", "is aux exactly 0?". With floating point, 1/3 + 1/3 + 1/3 may not equal 1, and a wrong tie or a wrong "is it zero?" changes the path — or the verdict "infeasible". So every animation on this page (and the Dart code) does its arithmetic with a fraction class: an integer numerator and a positive denominator, always reduced (Dart: BigInt so nothing overflows).

Standard and slack forms

LPs arrive in many shapes: "minimize", "at least", "exactly", numbers allowed to be negative. Like a kitchen with one standard pan size, algorithms want one shape. The standard form is: maximize cᵀx subject to Ax ≤ b and x ≥ 0 (A is an m × n table of numbers, b a list of m limits, c a list of n profits). Converting has exactly four possible obstacles, each with a fixed cure: (1) minimize → negate the objective; (2) a variable with no sign rule → write it as x′ − x″ with x′, x″ ≥ 0; (3) an equality f = b → two inequalities f ≤ b and f ≥ b; (4) a ≥ constraint → multiply by −1.

The second shape is the slack form, which is what simplex really manipulates. For each constraint Σⱼ aᵢⱼxⱼ ≤ bᵢ, add a new variable, the slack x_n+i = bᵢ − Σⱼ aᵢⱼxⱼ ≥ 0 (rows and variables are numbered from 0, so the slack of constraint i has the id n + i) — the amount of room still left in that rule. Now every rule is an equation. The variables on the left of "=" are basic (set B, m of them), those on the right and in the objective are nonbasic (set N, n of them). A slack form is the tuple (nonbasic, basic, a, b, c, v): the objective is z = v + Σ_j∈nonbasic c_jx_j and each basic variable satisfies x_i = b_i − Σ_j∈N a_ijx_j. Setting all nonbasic variables to 0 gives the basic solution: each basic x_i = b_i.

Sign convention: the stored a_ij has the minus sign already taken out, so the equation reads x_i = b_i − Σ a_ijx_j. The animations print the equations with ordinary signs (a coefficient a[i][j] = 2 shows as "− 2x_j"), but every caption that talks about a[i][j] means the stored value.

The four obstacles as code. toStandard applies the four cures mechanically (minimize → negate; free variable → two columns; equality → two rows; "≥" → multiply by −1) and returns the index maps that read the original x back from the new columns. Variables and constraints are numbered from 0, like every Dart list.

Why the four fixes are always enough. Each obstacle has a mechanical cure that never changes the answer: "minimize f" becomes "maximize −f" (negate the final optimum to translate back); a free variable x becomes x′ − x″ with x′, x″ ≥ 0 (every real number is a difference of two nonnegative ones); an equality becomes a "≤" plus a "≥" pair; and a "≥" is multiplied by −1. Every cure adds at most one variable or one constraint per original one, so the standard form is at most about twice the size, which is why "any LP can be solved by an algorithm for standard-form LPs" costs only a constant factor. Slack form is then just standard form with each "≤" turned into an equation by a new variable that measures the leftover room.

Standard form = maximize cᵀx, every rule "≤", every variable ≥ 0. Slack form = the same LP as equations x_i = b_i − Σ a_ijx_j for the basic variables B, expressed in the nonbasic variables N (set to 0 to read off a point). simplex only ever rewrites slack forms.

Variable numbering. The variables are x₀ … x_n+m−1 (originals first, then one slack per constraint) and the number of a variable is also its key in the Dart Maps: a[i][j] with i ∈ basic and j ∈ nonbasic is a map lookup, and the slack of constraint i (counting from 0) is variable n + i. The auxiliary variable of initializeSimplex is called aux (id −1).

Slack form as code. slackEquations prints the slack form of any LP the way the page does (the stored a_ij has the minus sign already taken out, so the text shows −a_ij); slackCoefficient shows where an input entry lives: the entry a[i][j] of the input matrix sits at s.a[n + i][j] of the slack form.

Input size → what's feasible. the conversion costs O(mn) (every entry is copied or negated once), so m, n up to a few thousand is instant; the work is in solving, not converting.

Formulating problems as linear programs

Being able to say "this is an LP" is worth a lot: it hands the problem to off-the-shelf solvers. The skill is translation — choose variables, write the goal as a linear score and every requirement as a linear rule.

Single-pair shortest path. One variable d_v per vertex; maximize d_t subject to d_v ≤ d_u + w(u,v) for every edge and d_s = 0 (|V| variables, |E| + 1 constraints). It maximizes, not minimizes, because the constraints are upper bounds: minimizing would just set every d_v = 0, while pushing d_t up as far as the bounds allow makes it equal to the length of the shortest chain of edges (this is Bellman–Ford's relaxation inequality). The Dart solution to interview question 14 solves it with the simplex below and checks it against Floyd–Warshall.
Maximum flow. One variable f_uv per edge; maximize the flow leaving s subject to f_uv ≤ c(u,v) and conservation (flow in = flow out) at every vertex except s and t. Interview question 17 implements it and checks it against Edmonds–Karp.
Minimum-cost flow. The same constraints, but demand exactly d units to leave s and minimize Σ a(u,v)f(u,v).
Multicommodity flow. One flow per commodity sharing the edge capacities; the only known polynomial algorithm is to write it as an LP (objective "minimize 0" just asks for feasibility).

Example 1 — shortest path. Watch the graph turn into inequalities, simplex solve them, and the tight roads spell out the shortest path (the last frame checks it against Dijkstra). The second player shows the edge case (target unreachable); the third takes your own graph.

The shortest-path LP as code (also interview question 14): n variables d_v, |E| + 1 rows.

Example 2 — maximum flow. Same recipe: one variable per pipe, capacity rows, conservation as a pair of inequalities, objective = net flow out of the source. The first player is a 6-vertex network with 9 pipes (max flow 15); the last frame shows the minimum cut that certifies it.

The max-flow LP as code (also interview question 17): one variable per edge, |E| capacity rows and 2(|V| − 2) conservation rows.

Minimum-cost flow as code. The same rows with "net outflow = +d at s, −d at t, 0 elsewhere" (two "≤" rows each) and the objective "maximize −cost". On our own 4-unit network (s→x cap 3 cost 2, s→y cap 3 cost 5, x→y cap 1 cost 1, x→t cap 2 cost 7, y→t cap 4 cost 1) the cheapest 4 units cost 22, 6 units cost 40, and 7 units do not fit; the verify program checks every number against successive shortest paths.

Multicommodity flow as code. One variable per (commodity, edge), shared capacity rows, per-commodity conservation, and the null objective "minimize 0": simplex only has to say optimal (feasible) or infeasible. Two commodities that each fit alone can still be infeasible together: commodities 0→3 and 1→3 of demand 2 and 2 over edges 0→2 (cap 2), 1→2 (cap 2), 2→3 (cap 3) need 4 units on the edge 2→3 of capacity 3.

Input size → what's feasible. the examples here have |V| ≤ 8, |E| ≤ 20 (a tableau of a few hundred entries). Real networks (10⁴ edges) are far beyond a hand-rolled exact tableau: use Dijkstra, Dinic or a network-simplex library for the single-commodity problems, and an industrial LP solver for multicommodity flow, where the LP is the only polynomial method known.

Two traps. (1) The shortest-path LP maximizes; minimizing d_t gives the useless answer 0. (2) The LP only fixes d_t: other vertices of the optimal solution may lie below their true distances, so read the path from the tight roads, not from arbitrary d values. Also, this simple version assumes lengths ≥ 0 because variables are ≥ 0; negative-length roads need d_v written as a difference of two nonnegative variables (see the conversion section).

Why these LPs are "nice". Shortest-path and max-flow LPs have constraint matrices that are totally unimodular (a network structure: each edge contributes a +1 at its head and a −1 at its tail), so their corners are already whole numbers when the data are whole numbers. That is why solving the LP gives an actual path or an actual integral flow instead of "half a road". The max-flow LP's dual is exactly the minimum-cut LP; the shortest-path LP's dual is a minimum-cost unit flow from s to t. Both are answered in the question bank.

If a problem becomes an LP with polynomially many variables and constraints, it can be solved in polynomial time by the ellipsoid or interior-point methods, even though simplex itself is exponential in the worst case. If the variables are additionally required to be integers the problem (integer linear programming) is NP-hard.

pivot(form, l, e)

Think of the slack form as a system of equations written "solved for" the basic variables. A pivot changes which variables you have solved for: you choose one basic variable x_l to give up its place (the leaving variable) and one nonbasic variable x_e to take it (the entering variable). Algebraically: take the equation of x_l, solve it for x_e, and substitute the result everywhere else. The LP itself does not change — it is written from a new point of view, whose basic solution is a neighbouring vertex.

Three easy mistakes: (1) the pair (l, e) must have l basic and e nonbasic, and a[l][e] ≠ 0, or the division on line 1 is meaningless (the custom player checks this); (2) pivot does no thinking about feasibility or optimality — choosing l by the smallest delta (done by simplex, not by pivot) is what keeps every b[i] ≥ 0; (3) the stored a[i][j] use the sign convention x_i = b[i] − Σ a[i][j]x_j, so a coefficient printed as "− 2x_j" in the animation is a[i][j] = +2.

Correctness. After pivot the nonbasic variables are again the ones set to 0, the entering variable has the value b2[e] = b[l]/a[l][e], and every other basic variable has b2[i] = b[i] − a[i][e]·b2[e]. The algebra only ever solves one equation for one variable and substitutes it, so the new form has exactly the same solutions as the old one.

Running time. The pivot row costs one division for b2[e] (line 2), n − 1 for lines 3–4 and one for line 5; each of the m − 1 other rows costs 1 + (n − 1) + 1 (lines 8–11); the objective costs 1 + (n − 1) + 1 (lines 13–16): Θ(mn) operations in total. Cost table for the main LP (m = n = 3): the pivot row needs 4 divisions (1 + 2 + 1), the 2 other rows need 4 updates each (8 in all), and the objective needs 4 updates: 16 arithmetic steps, which grows as m·n for a bigger LP.

pivot is Gaussian elimination on one column: it swaps the roles of one basic and one nonbasic variable, leaves the LP unchanged in meaning, and costs Θ(mn).

The cost table as code. pivotOperationCount walks the same loops as the numbered pseudocode and counts the arithmetic steps: (m + 1)(n + 1) = 16 for m = n = 3, between mn and 4mn for every m, n ≥ 1 (so Θ(mn)).

Input size → what's feasible. one pivot is Θ(mn): at m = n = 100 that is ≈ 10⁴ fraction operations (milliseconds), at m = n = 10⁴ it is 10⁸ (a second with doubles, hopeless with exact BigInt fractions).

Dart implementation (the data type, the conversion to slack form, and pivot itself):

simplex(a, b, c)

Imagine walking along the edges of a mountain-shaped polyhedron, always stepping to a neighbouring corner that is higher (or level) — and stopping when every neighbouring step goes down. Because the region is convex, a corner with no better neighbour is the best corner overall. In slack-form terms: while some nonbasic variable has a positive objective coefficient, raising it from 0 raises the score, so increase it until some basic variable would hit 0 — that is one pivot, one step to the next corner.

simplex starts from a slack form whose basic solution is feasible (initializeSimplex, described in its own section below, supplies one). Each iteration: choose an entering e with c_e > 0; compute for each basic row the limit delta[i] = b[i]/a[i][e] (only rows with a[i][e] > 0 limit x_e); the smallest delta gives the leaving row; if no row limits it, the LP is unbounded. We always use Bland's rule (smallest index for e, and for l among ties). Watch the "delta" chips in the animation: they are the decision at each step.

When the walk goes in circles. The next three players use an LP with several rows equal to 0 (Beale's example, built so that many constraints meet at the origin). The rule for choosing e decides everything: under the popular largest-coefficient rule the walk returns to its starting basis after 6 pivots and would loop forever; under Bland's rule (always the smallest index) the same LP finishes. Each frame names the rule's decision.

Bland's rule as code. blandEntering and blandLeaving are the two choices; largestEntering is the rival that cycles on Beale's LP; pivotCount runs either rule (and returns −1 if it is still running after the pivot limit). cyclePivots below is the cycle detector of interview question 29.

Degeneracy and cycling. If the leaving row has b[l] = 0, then delta = 0: the pivot changes B but not the point or the score (the custom example's redundant constraint does this — you will see a basic variable equal to 0 in the answer). A run of degenerate pivots can return to an earlier slack form and repeat forever. Since a slack form is determined by its basic set B a run longer than C(n+m, m) iterations must repeat; Bland's rule provably prevents that, so simplex stops within C(n+m, m) iterations.

Loop invariant. Before every iteration of the while loop: (1) the slack form is equivalent to the original LP; (2) b_i ≥ 0 for all i ∈ B; (3) the basic solution is feasible. It holds at the start because initializeSimplex returns a feasible form. The pivot chooses l by the smallest delta, so the new b2[i] = b[i] − a[i][e]·b2[e] stays ≥ 0 when a_ie > 0 (that is what "smallest ratio" guarantees) and when a_ie ≤ 0 it can only grow. On termination all c_j ≤ 0, so z = v + Σ c_jx_j ≤ v for every feasible point and the basic solution achieves v: optimal. If all delta = ∞, setting x_e = λ with the other nonbasic variables 0 gives feasible points with score → ∞: unbounded.

Running time. Each iteration costs Θ(mn) (one pivot plus the O(m) ratio scan). With Bland's rule the number of iterations is at most C(n+m, m), which is exponential; the Klee–Minty cube shows some inputs really need 2ⁿ − 1 pivots. In practice simplex needs a small multiple of m pivots. The ratio table for the main LP: 3 pivots × Θ(3·3) work each.

simplex = repeat { pick e with c_e > 0; pick l by the smallest ratio; pivot } until no c_j > 0 (optimal) or no ratio exists (unbounded). Each pivot never lowers z; only Bland's rule (or a similar tie-breaking rule) guarantees it stops.

Analysis visual — how many pivots could it take? A slack form is determined by its basis, so a run that never repeats visits at most C(n+m, m) different bases; Bland's rule guarantees no repeat. That bound is huge, and some inputs really do need exponentially many pivots. The table (computed by this page) compares the two numbers for a square LP with n = m:

The C(n+m, m) bound as code. choose computes the binomial exactly, countBases counts the m-subsets of the n + m variables directly (a slack form is determined by its basis), and kleeMinty builds the cube on which the largest-coefficient rule needs exactly 2ⁿ − 1 pivots (checked for n = 1..7 in the verify program, against C(2n, n) = 2, 6, 20, 70, 252, 924, 3432).

Input size → what's feasible. n, m ≤ 40 with exact fractions: about 100 pivots of 1,600 operations, well under a second. n, m ≈ 10³: floating-point simplex (tableau of 10⁶ numbers per pivot, a few thousand pivots). Larger or adversarial inputs: interior-point methods, because C(n+m, m) and 2ⁿ − 1 are exponential.

Dart implementation (the loop and the read-off of x̄ from lines 13–14; initializeSimplex is in its own section below):

Duality

Suppose you are the workshop owner, and an outside buyer offers to rent all your oven and labour hours. Each hour gets a price y_i. The buyer wants the total bill as small as possible, but you will only agree if, for every product you could have built, the hours it needs are priced at least as high as the profit you would have made. The cheapest such prices form the dual LP. The remarkable fact: the buyer's minimum bill equals your maximum profit.

Given the primal max{cᵀx : Ax ≤ b, x ≥ 0}, the dual is min{bᵀy : Aᵀy ≥ c, y ≥ 0}: one dual variable per primal constraint, one dual constraint per primal variable. Weak duality: every feasible x and feasible y satisfy c·x ≤ b·y — so if they are equal, both are optimal. Strong duality: if simplex finds an optimal x̄, then y[i] = −c[n+i] (if the slack x_n+i is nonbasic in the final form, else 0) is an optimal dual solution and c·x̄ = b·ȳ. The four possible outcome pairs are (optimal, optimal) with equal values; (unbounded, infeasible); (infeasible, unbounded or infeasible); never anything else.

Weak duality as code. primalFeasible and dualFeasible test the two kinds of point; weakDualityChain returns the four numbers of the proof of weak duality, c·x ≤ (Aᵀy)·x = y·(Ax) ≤ b·y, which can never decrease from left to right (4,000 random pairs are checked in the verify program).

Seeing strong duality happen. Run simplex on the primal and on the dual and plot their objective values after each pivot. Every primal value is a lower bound on the true optimum z*, every dual value an upper bound; the two curves approach each other and touch at z*. The edge-case player shows what a missing bound looks like (unbounded primal, no dual point at all).

Do not read the two curves as one algorithm moving in lockstep: they are two independent simplex runs (the dual run may need its own initializeSimplex, so it can have a different number of pivots), and the chart lines them up only by pivot number. Also, the dual value at an intermediate step is the objective of a feasible dual point, not the value the primal "should" have at that step.

Weak duality: c·x ≤ b·y for every feasible pair. Strong duality: at optimum they are equal. So any feasible dual point is a proof that the primal cannot exceed it, and a matching pair is a proof of optimality that anyone can check without rerunning simplex.

Strong duality as code. strongDuality solves the primal, reads y off the final slack form, solves the dual LP separately, and returns the three values, which are equal. dualityOutcome gives the (primal, dual) status pair: optimal/optimal, unbounded/infeasible, infeasible/unbounded or infeasible/infeasible.

Input size → what's feasible. the dual has the same size as the primal (n rows, m columns), so duality costs O(mn) to write down; checking a given primal/dual pair is O(mn) too, which is why a matching pair is a proof anyone can verify for m, n up to 10³ in milliseconds.

Correctness of the read-off. Write the final objective z = v′ + Σ c′_jx_j over all n + m variables (basic variables get coefficient 0). Since every slack form is equivalent to the original, this equals c·x for all x, so Σ_i a[i][j]·y[i] = c[j] − c′_j ≥ c[j] (each c′_j ≤ 0 at the optimum) — the dual constraints hold — and Σ b[i]·y[i] = v′, the primal optimum. Complementary slackness falls out: x_j > 0 forces the j-th dual constraint tight, and y[i] > 0 forces the i-th primal constraint tight. The Dart verification checks strong duality and complementary slackness on 300+ random LPs.

Dart implementation (the dual as a standard-form LP, and the read-off of y):

initializeSimplex(a, b, c)

simplex needs a corner to start walking from, and the obvious one — all variables 0 — is only legal when every b[i] ≥ 0. Otherwise you first need to find a legal point, or prove there is none. The trick: build a helper LP, the auxiliary LP, that is deliberately easy — it has one extra "fudge" variable aux that loosens every rule — and use simplex itself to squeeze the fudge back to 0. If the fudge can reach 0, the rules can be satisfied together; if it cannot, they contradict each other.

The auxiliary LP as code. auxiliaryLp adds the column for aux (the last column) with coefficient −1 in every row and the objective "maximize −aux"; auxStartPoint is the point x = 0, aux = −min b that makes it feasible; feasibleViaAux checks that the LP is feasible exactly when the optimum of the auxiliary LP is 0.

Input size → what's feasible. one extra column and one extra pivot, then an ordinary simplex run: n, m ≤ 40 exact (well under a second), or a few thousand with doubles. If every b_i ≥ 0 skip it entirely.

Correctness. The LP is feasible ⇔ the auxiliary LP has optimal value 0. If the LP is infeasible the procedure returns "infeasible"; if feasible, the single pivot on line 8 makes all b[i] ≥ 0 (aux = −b[k], and every other row becomes b[i] − b[k] ≥ 0), the simplex loop on the auxiliary LP finishes with aux = 0, and the returned slack form is equivalent to the original LP with a feasible basic solution. Together with the termination argument for Bland's rule and strong duality this proves the fundamental theorem of linear programming: every LP is infeasible, unbounded, or has a finite optimum, and simplex says which.

Running time. One extra variable and one extra pivot, then an ordinary simplex run on an LP with n + 1 variables and m constraints: the same Θ(mn) per pivot, and the same iteration bound.

Two subtleties. (1) If the optimum of the auxiliary LP is 0 but aux is still basic (sitting at value 0), one extra degenerate pivot must push it out before aux is deleted. (2) After deleting aux the objective of L is written in terms of the current nonbasic variables by substituting each basic original variable's equation; forgetting this substitution silently optimizes the wrong function.

initializeSimplex either returns a slack form with a feasible basic solution (so simplex can start) or proves the LP infeasible. Together with the loop, simplex decides all three outcomes: infeasible, unbounded, or an optimal vertex.

Dart implementation (with the wrapper that runs both stages):

Quiz

Interview questions

Cheat sheet

Procedure / idea	Time	Space	Key facts	When to use
Standard form	O(mn) to convert	O(mn)	max cᵀx, Ax ≤ b, x ≥ 0; four obstacles, four cures	First step for any LP
Graphical method (2 variables)	O(m³) by vertex enumeration	O(m)	Optimum at a vertex; empty region = infeasible; open region can be unbounded	Intuition, tiny LPs
pivot	Θ(mn)	Θ(mn)	Swaps one basic and one nonbasic variable; needs a_le ≠ 0; same LP, new view	Inner step of simplex
simplex	Θ(mn) per pivot; ≤ C(n+m, m) pivots with Bland; 2ⁿ − 1 worst case (Klee–Minty)	Θ(mn)	Returns optimal, or "unbounded"; degenerate pivots leave v unchanged; Bland's rule prevents cycling	General LPs, fast in practice
initializeSimplex	One pivot + one simplex run	Θ(mn)	auxiliary LP: maximize −aux; optimum 0 ⇔ LP feasible; returns a feasible slack form or "infeasible"	Whenever some b_i < 0
Duality	Free once simplex has run (O(m))	O(m)	Weak: c·x ≤ b·y. Strong: equal at optimum. y_i = −c′_n+i for nonbasic slacks	Optimality certificates, sensitivity ("shadow prices"), min-cut, Farkas
Other LP algorithms	Ellipsoid and interior-point: polynomial	—	Ellipsoid is slow in practice; interior-point can beat simplex on huge LPs	Very large LPs