The Simplex Method is a canonical algorithm for linear programming and optimization. Usually implemented with Bland's Finite Pivoting Method, this algorithm has been shown to be pseudo-finite terminating, by showing that it avoids cycles (see Bland's original paper for more details). This is an implementation of the algorithm in pure Haskell - to both show its feasability, and practicality as obviously correct.
We treat linear expressions particularly - specifically in
standard form, such that sets of variables and their coefficients are mappings
between the variable name and coefficient value: name => coeff
; thus a linear
expression could be designed as the pair (name => coeff, const)
. This is
exemplified under the LinExpr
data type.
Likewise, we imagine inequalities as the same structure, but with an additional tag representing it's constraint type:
ineq := equality -- =
| ltequal -- <=
| gtequal -- >=
Such that we imagine the variables on the left, and the constant on the
right. Thus, an inequality is the triple (ineq, name => coeff, const)
.
There is some machinery we need to account for internally in the algorithm; particularly in that the simplex method only operates on positive variables, which is rarely desires. So, we use a state monad to build a tableau which stores our constraints incrementally, while keeping track of the details.
myConstraints :: MonadState Tableau m => m ()
myConstraints = do
addConstraint ineq1
addConstraint ineq2
addConstraint ineq3
Then, to get the underlying Tableau, we need to run the state monad,
supplying an initial Tableau via newTableau
(where we also need an
objective function / cost function):
myTableau :: Tableau
myTableau = execState myConstraints (newTableau objectiveIneq)
Now we can get to optimization!
The tableau contains enough data to start optimizing. To kick the algorithm
into gear, just run primalSimplex
over the beast:
myOptimizedTableau :: Tableau
myOptimizedTableau = primalSimplex myTableau
Then you can peek at the supposedly optimal solution via solution
:
mySolution :: Map MainVarName Rational
mySolution = solution myOptimizedTableau
This should be optimal such that all other feasible solutions are less than or equal to this optimal solution.
stack test
- generate a convext polytope as the initial (solvable) constraint set