GitHub

Note: This is early days for this module. Anything before 0.1.0 is subject to massive changes.

SimplexTableaux

This module can be used to show how to solve linear programming problems using the simplex method by pivoting tableaux.

This is an illusration project for solving feasible optimization problems of the form $\min c^T x$ subject to $Ax ≥ b$ and $x \ge 0$ (canonical form) and of the form $\min c^T x$ subject to $Ax = b$ and $x \ge 0$ (standard form).

Caveats

As a demonstration project, this is not suitable for solving actual linear programming (LP) problems. At present it will fail if:

  • The LP is infeasible.
  • The LP is unbounded.
  • Other unidentified reasons. (In other words, still buggy.)

This module is set up for minimization problems only.

This module solves [not yet!!] LPs using the simplex algorithm which is not the most performant method. Further, all data is stored using arbitrary precision integers (that is, Rational{BigInt}) which gives exact results, but is much slower than floating point arithmetic. These issues are negligible for small problems.

Setting up a Tableau

Canonical form LPs

A linear program in canonical form is $\min c^T x$ s.t. $Ax ≥ b$, $x ≥ 0$.

For example, let A, b, and c be as follows:

julia> A = [3 10; 5 6; 10 2];

julia> b = [100, 100, 100];

julia> c = [25, 10];

julia> Tableau(A, b, c)
┌──────────┬───┬─────┬─────┬─────┬─────┬─────┬─────┐
│          │ z │ x_1 │ x_2 │ x_3 │ x_4 │ x_5 │ RHS │
│ Obj Func │ 1 │ -25 │ -10 │   0 │   0 │   0 │   0 │
├──────────┼───┼─────┼─────┼─────┼─────┼─────┼─────┤
│   Cons 1 │ 0 │   3 │  10 │  -1 │   0 │   0 │ 100 │
│   Cons 2 │ 0 │   5 │   6 │   0 │  -1 │   0 │ 100 │
│   Cons 3 │ 0 │  10 │   2 │   0 │   0 │  -1 │ 100 │
└──────────┴───┴─────┴─────┴─────┴─────┴─────┴─────┘

Notice that extra variables $x_3$, $x_4$, and $x_5$ are added to the Tableau as slack variables to convert inequalities into equations. That is, canonical form LPs are automatically converted into standard form.

Standard form LPs

A linear program in standard form is $\min c^T x$ s.t. $Ax = b$, $x ≥ 0$. For example,

julia> A = [2 1 0 9 -1; 1 1 -1 5 1]
2×5 Matrix{Int64}:
 2  1   0  9  -1
 1  1  -1  5   1

julia> b = [9, 7]
2-element Vector{Int64}:
 9
 7

julia> c = [2, 4, 2, 1, -1]
5-element Vector{Int64}:
  2
  4
  2
  1
 -1

julia> T = Tableau(A, b, c, false)
┌──────────┬───┬─────┬─────┬─────┬─────┬─────┬─────┐
│          │ z │ x_1 │ x_2 │ x_3 │ x_4 │ x_5 │ RHS │
│ Obj Func │ 1 │  -2 │  -4 │  -2 │  -1 │   1 │   0 │
├──────────┼───┼─────┼─────┼─────┼─────┼─────┼─────┤
│   Cons 1 │ 0 │   2 │   1 │   0 │   9 │  -1 │   9 │
│   Cons 2 │ 0 │   1 │   1 │  -1 │   5 │   1 │   7 │
└──────────┴───┴─────┴─────┴─────┴─────┴─────┴─────┘

The fourth argument false means that the constraints are already equalities and slack variables should not be appended.

Operations

Operations on a Tableau may include row and/or column indices. A row index refers to the constraint number. That is, row 1 corresponds to the Cons 1 row. A column index refers to a decision variable. That is, column 1 corresonds to the $x_1$ column.

Pivot operations

Basic pivot

Use pivot!(T, i, j) to perform a pivot on row i and contraint j. The entry at that location must not be 0.

julia> T
┌──────────┬───┬─────┬─────┬─────┬─────┬─────┬─────┐
│          │ z │ x_1 │ x_2 │ x_3 │ x_4 │ x_5 │ RHS │
│ Obj Func │ 1 │  -2 │  -4 │  -2 │  -1 │   1 │   0 │
├──────────┼───┼─────┼─────┼─────┼─────┼─────┼─────┤
│   Cons 1 │ 0 │   2 │   1 │   0 │   9 │  -1 │   9 │
│   Cons 2 │ 0 │   1 │   1 │  -1 │   5 │   1 │   7 │
└──────────┴───┴─────┴─────┴─────┴─────┴─────┴─────┘

julia> pivot!(T,1,4)
┌──────────┬───┬───────┬───────┬─────┬─────┬──────┬─────┐
│          │ z │   x_1 │   x_2 │ x_3 │ x_4 │  x_5 │ RHS │
│ Obj Func │ 1 │ -16/9 │ -35/9 │  -2 │   0 │  8/9 │   1 │
├──────────┼───┼───────┼───────┼─────┼─────┼──────┼─────┤
│   Cons 1 │ 0 │   2/9 │   1/9 │   0 │   1 │ -1/9 │   1 │
│   Cons 2 │ 0 │  -1/9 │   4/9 │  -1 │   0 │ 14/9 │   2 │
└──────────┴───┴───────┴───────┴─────┴─────┴──────┴─────┘

Basis pivot

Let B be a list of columns (with as many columns specified as there are constraints). basis_pivot!(T, B) makes that selection of columns into basic variables.

julia> T
┌──────────┬───┬─────┬─────┬─────┬─────┬─────┬─────┐
│          │ z │ x_1 │ x_2 │ x_3 │ x_4 │ x_5 │ RHS │
│ Obj Func │ 1 │  -2 │  -4 │  -2 │  -1 │   1 │   0 │
├──────────┼───┼─────┼─────┼─────┼─────┼─────┼─────┤
│   Cons 1 │ 0 │   2 │   1 │   0 │   9 │  -1 │   9 │
│   Cons 2 │ 0 │   1 │   1 │  -1 │   5 │   1 │   7 │
└──────────┴───┴─────┴─────┴─────┴─────┴─────┴─────┘

julia> B = [3,5];

julia> basis_pivot!(T, B)
┌──────────┬───┬─────┬─────┬─────┬─────┬─────┬─────┐
│          │ z │ x_1 │ x_2 │ x_3 │ x_4 │ x_5 │ RHS │
│ Obj Func │ 1 │  -6 │  -7 │   0 │ -20 │   0 │ -23 │
├──────────┼───┼─────┼─────┼─────┼─────┼─────┼─────┤
│   Cons 1 │ 0 │  -3 │  -2 │   1 │ -14 │   0 │ -16 │
│   Cons 2 │ 0 │  -2 │  -1 │   0 │  -9 │   1 │  -9 │
└──────────┴───┴─────┴─────┴─────┴─────┴─────┴─────┘

Other operations

Swap rows

Use swap!(T, i, j) to swap rows i and j of the Tableau.

julia> T
┌──────────┬───┬─────┬─────┬─────┬─────┬─────┬─────┐
│          │ z │ x_1 │ x_2 │ x_3 │ x_4 │ x_5 │ RHS │
│ Obj Func │ 1 │  -2 │  -4 │  -2 │  -1 │   1 │   0 │
├──────────┼───┼─────┼─────┼─────┼─────┼─────┼─────┤
│   Cons 1 │ 0 │   2 │   1 │   0 │   9 │  -1 │   9 │
│   Cons 2 │ 0 │   1 │   1 │  -1 │   5 │   1 │   7 │
└──────────┴───┴─────┴─────┴─────┴─────┴─────┴─────┘

julia> swap!(T, 1, 2)
┌──────────┬───┬─────┬─────┬─────┬─────┬─────┬─────┐
│          │ z │ x_1 │ x_2 │ x_3 │ x_4 │ x_5 │ RHS │
│ Obj Func │ 1 │  -2 │  -4 │  -2 │  -1 │   1 │   0 │
├──────────┼───┼─────┼─────┼─────┼─────┼─────┼─────┤
│   Cons 1 │ 0 │   1 │   1 │  -1 │   5 │   1 │   7 │
│   Cons 2 │ 0 │   2 │   1 │   0 │   9 │  -1 │   9 │
└──────────┴───┴─────┴─────┴─────┴─────┴─────┴─────┘

Scaling a row

Use scale!(T,r,s) to multiply all elements of row r by s.

julia> T
┌──────────┬───┬─────┬─────┬─────┬─────┬─────┬─────┐
│          │ z │ x_1 │ x_2 │ x_3 │ x_4 │ x_5 │ RHS │
│ Obj Func │ 1 │  -2 │  -4 │  -2 │  -1 │   1 │   0 │
├──────────┼───┼─────┼─────┼─────┼─────┼─────┼─────┤
│   Cons 1 │ 0 │   2 │   1 │   0 │   9 │  -1 │   9 │
│   Cons 2 │ 0 │   1 │   1 │  -1 │   5 │   1 │   7 │
└──────────┴───┴─────┴─────┴─────┴─────┴─────┴─────┘

julia> scale!(T,2,5)
┌──────────┬───┬─────┬─────┬─────┬─────┬─────┬─────┐
│          │ z │ x_1 │ x_2 │ x_3 │ x_4 │ x_5 │ RHS │
│ Obj Func │ 1 │  -2 │  -4 │  -2 │  -1 │   1 │   0 │
├──────────┼───┼─────┼─────┼─────┼─────┼─────┼─────┤
│   Cons 1 │ 0 │   2 │   1 │   0 │   9 │  -1 │   9 │
│   Cons 2 │ 0 │   5 │   5 │  -5 │  25 │   5 │  35 │
└──────────┴───┴─────┴─────┴─────┴─────┴─────┴─────┘

Return to start

Use restore!(T) to return T to its original values. 

julia> T
┌──────────┬───┬─────┬─────┬─────┬─────┬─────┬─────┐
│          │ z │ x_1 │ x_2 │ x_3 │ x_4 │ x_5 │ RHS │
│ Obj Func │ 1 │  -2 │  -4 │  -2 │  -1 │   1 │   0 │
├──────────┼───┼─────┼─────┼─────┼─────┼─────┼─────┤
│   Cons 1 │ 0 │   2 │   1 │   0 │   9 │  -1 │   9 │
│   Cons 2 │ 0 │   1 │   1 │  -1 │   5 │   1 │   7 │
└──────────┴───┴─────┴─────┴─────┴─────┴─────┴─────┘

julia> pivot!(T,2,4)
┌──────────┬───┬──────┬───────┬───────┬─────┬───────┬───────┐
│          │ z │  x_1 │   x_2 │   x_3 │ x_4 │   x_5 │   RHS │
│ Obj Func │ 1 │ -9/5 │ -19/5 │ -11/5 │   0 │   6/5 │   7/5 │
├──────────┼───┼──────┼───────┼───────┼─────┼───────┼───────┤
│   Cons 1 │ 0 │  1/5 │  -4/5 │   9/5 │   0 │ -14/5 │ -18/5 │
│   Cons 2 │ 0 │  1/5 │   1/5 │  -1/5 │   1 │   1/5 │   7/5 │
└──────────┴───┴──────┴───────┴───────┴─────┴───────┴───────┘

julia> restore!(T)
┌──────────┬───┬─────┬─────┬─────┬─────┬─────┬─────┐
│          │ z │ x_1 │ x_2 │ x_3 │ x_4 │ x_5 │ RHS │
│ Obj Func │ 1 │  -2 │  -4 │  -2 │  -1 │   1 │   0 │
├──────────┼───┼─────┼─────┼─────┼─────┼─────┼─────┤
│   Cons 1 │ 0 │   2 │   1 │   0 │   9 │  -1 │   9 │
│   Cons 2 │ 0 │   1 │   1 │  -1 │   5 │   1 │   7 │
└──────────┴───┴─────┴─────┴─────┴─────┴─────┴─────┘

Copies of the original matrix A and the vectors b and c used to create a Tableau can be extracted from T by calling get_Abc(T).

Non-modifying versions of functions

The functions basis_pivot!, pivot!, restore!, scale!, and swap! all modify the tableau to which they apply. There are versions of all of these that do not modify the tableau; instead, they return a copy of the tableau that is the result of the operation, leaving the tableau unchanged. The non-modifying versions have the same names, but with the ! removed.

Solving Linear Programs

Manually

User can use pivot! to emulate the Simplex Algorithm. We plan to automate this, but that will take a bit of time.

Feasible basic vector

Assuming the Tableau has been pivoted to a basic vector, check that the current state is feasible (i.e., the RHS for all the constraints is non-negative).

julia> T
┌──────────┬───┬─────┬─────┬─────┬─────┬─────┬─────┐
│          │ z │ x_1 │ x_2 │ x_3 │ x_4 │ x_5 │ RHS │
│ Obj Func │ 1 │  -2 │  -4 │  -2 │  -1 │   1 │   0 │
├──────────┼───┼─────┼─────┼─────┼─────┼─────┼─────┤
│   Cons 1 │ 0 │   2 │   1 │   0 │   9 │  -1 │   9 │
│   Cons 2 │ 0 │   1 │   1 │  -1 │   5 │   1 │   7 │
└──────────┴───┴─────┴─────┴─────┴─────┴─────┴─────┘

julia> basis_pivot!(T,[2,4])
┌──────────┬───┬───────┬─────┬───────┬─────┬──────┬──────┐
│          │ z │   x_1 │ x_2 │   x_3 │ x_4 │  x_5 │  RHS │
│ Obj Func │ 1 │ -11/4 │   0 │ -43/4 │   0 │ 29/2 │ 37/2 │
├──────────┼───┼───────┼─────┼───────┼─────┼──────┼──────┤
│   Cons 1 │ 0 │  -1/4 │   1 │  -9/4 │   0 │  7/2 │  9/2 │
│   Cons 2 │ 0 │   1/4 │   0 │   1/4 │   1 │ -1/2 │  1/2 │
└──────────┴───┴───────┴─────┴───────┴─────┴──────┴──────┘

julia> is_feasible(T)
true

julia> basis_pivot!(T,[1,3])
┌──────────┬───┬─────┬──────┬─────┬──────┬──────┬──────┐
│          │ z │ x_1 │  x_2 │ x_3 │  x_4 │  x_5 │  RHS │
│ Obj Func │ 1 │   0 │   -4 │   0 │    7 │   -3 │    4 │
├──────────┼───┼─────┼──────┼─────┼──────┼──────┼──────┤
│   Cons 1 │ 0 │   1 │  1/2 │   0 │  9/2 │ -1/2 │  9/2 │
│   Cons 2 │ 0 │   0 │ -1/2 │   1 │ -1/2 │ -3/2 │ -5/2 │
└──────────┴───┴─────┴──────┴─────┴──────┴──────┴──────┘

julia> is_feasible(T)
false

Using a solver

The function lp_solve finds a numerical solution to the linear program using a Julia solver (default: HiGHS).

julia> T
┌──────────┬───┬─────┬─────┬─────┬─────┬─────┬─────┐
│          │ z │ x_1 │ x_2 │ x_3 │ x_4 │ x_5 │ RHS │
│ Obj Func │ 1 │  -2 │  -4 │  -2 │  -1 │   1 │   0 │
├──────────┼───┼─────┼─────┼─────┼─────┼─────┼─────┤
│   Cons 1 │ 0 │   2 │   1 │   0 │   9 │  -1 │   9 │
│   Cons 2 │ 0 │   1 │   1 │  -1 │   5 │   1 │   7 │
└──────────┴───┴─────┴─────┴─────┴─────┴─────┴─────┘

julia> lp_solve(T)
Optimum value = -0.14285714285714324

5-element Vector{Float64}:
 0.0
 0.0
 0.0
 1.1428571428571428
 1.285714285714286

LaTeX output

Using LatexPrint users can get the code for pasting into a LaTeX document.

julia> using LatexPrint

julia> T
┌──────────┬───┬─────┬─────┬─────┬─────┬──────┬─────┐
│          │ z │ x_1 │ x_2 │ x_3 │ x_4 │  x_5 │ RHS │
│ Obj Func │ 1 │   0 │  -3 │  -2 │   8 │    0 │   9 │
├──────────┼───┼─────┼─────┼─────┼─────┼──────┼─────┤
│   Cons 1 │ 0 │   1 │ 1/2 │   0 │ 9/2 │ -1/2 │ 9/2 │
│   Cons 2 │ 0 │   0 │ 1/2 │  -1 │ 1/2 │  3/2 │ 5/2 │
└──────────┴───┴─────┴─────┴─────┴─────┴──────┴─────┘


julia> lap(T)
\begin{tabular}{|c|ccccc|c|}\hline 
{\Large\strut}$z$ &$x_{1}$ & $x_{2}$ & $x_{3}$ & $x_{4}$ & $x_{5}$ & RHS \\
{\Large\strut}$1$ & $0$ & $-3$ & $-2$ & $8$ & $0$ & $9$\\
\hline 
{\Large\strut}$0$ & $1$ & $\frac{1}{2}$ & $0$ & $\frac{9}{2}$ & $\frac{-1}{2}$ & $\frac{9}{2}$\\
{\Large\strut}$0$ & $0$ & $\frac{1}{2}$ & $-1$ & $\frac{1}{2}$ & $\frac{3}{2}$ & $\frac{5}{2}$\\
\hline 
\end{tabular}

Here is the LaTeX output: