On this page you will find programs that I have written, including:
Linear programming:
These versions have been coded as Excel macros so they are easy to use.
The simplex problem is a series of Hilbert matrices. Hilbert matrices are known to be unstable. The code is good for degree 10.
The assignment or transportation problem uses the same simplex code internally:
Option Explicit
Sub Simplex()
Application.ScreenUpdating = False
Application.Calculation = xlCalculationManual
Application.EnableEvents = False
Const ZT = 1E-16
Dim Z As Double
Dim A As Double
Dim I As Integer
Dim J As Integer
Dim M As Integer
Dim N As Integer
Dim NM As Integer
Dim P As Integer
Dim Q As Integer
Dim Phase As Integer
Dim T() As Double
Dim R() As Double
Dim B() As Double
Dim X() As Double
' Read model size
N = ActiveSheet.Cells(1, 1).Value ' Variables (columns)
M = ActiveSheet.Cells(1, 2).Value ' Constraints (rows)
NM = N + M ' Total columns with slack variables
ReDim T(1 To M, 1 To NM) As Double
ReDim R(1 To NM) As Double
ReDim B(1 To M) As Double
ReDim S(1 To M) As Double
ReDim X(1 To N) As Double
' Read in model
' N=Variables, M=Constraints
' Maximise: R(i)
' Subject to: T(j,i) <= B(j)
' X(i) >= 0
For J = 1 To M
S(J) = 0 ' Used to record basis variables
For I = 1 To N
T(J, I) = ActiveSheet.Cells(J + 2, I).Value
Next I
B(J) = ActiveSheet.Cells(J + 2, N + 1).Value
Next J
For I = 1 To N
R(I) = ActiveSheet.Cells(2, I).Value
Next I
Z = 0 ' System value
' Add slack variables
For J = 1 To M
For I = 1 To M
T(J, N + I) = 0
Next I
R(N + J) = 0
T(J, N + J) = 1
Next J
Rem SIMPLEX - Two phase method
Phase = 1
While (Phase = 1) Or (Phase = 2)
If (Phase = 1) Then
Rem Test if feasible - Find smallest B
Q = 1
For J = 2 To M
If (B(Q) > B(J)) Then Q = J
Next J
If (B(Q) >= 0) Then
Rem Feasible but not necessarily optimal
Phase = Phase + 1
Else
Rem Not feasible - Select special pivot - Find smallest ratio
P = 0
For I = 1 To N + M
If (T(Q, I) < 0) Then
If (P = 0) Then P = I
If R(P) / T(Q, P) > (R(I) / T(Q, I)) Then P = I
End If
Next I
If (P = 0) Then
Rem Infeasible
Phase = -2
End If
End If
End If
If Phase = 2 Then
Rem Test if optimal - Find largest R
P = 1
For I = 2 To N + M
If (R(P) < R(I)) Then P = I
Next I
If (R(P) <= 0) Then
Rem Optimal
Phase = Phase + 1
Else
Rem Not optimal - Select pivot - Find smallest ratio
Q = 0
For J = 1 To M
If (T(J, P) > 0) Then
If (Q = 0) Then Q = J
If (B(Q) / T(Q, P) > B(J) / T(J, P)) Then Q = J
End If
Next J
If (Q = 0) Then
Rem Unbounded
Phase = -1
End If
End If
End If
Rem Row operations - Test for near zero here
If (Phase = 1) Or (Phase = 2) Then
S(Q) = P
A = T(Q, P)
For I = 1 To N + M
T(Q, I) = T(Q, I) / A
If (Abs(T(Q, I)) <= ZT) Then T(Q, I) = 0
Next I
T(Q, P) = 1
B(Q) = B(Q) / A
If (Abs(B(Q)) <= ZT) Then B(Q) = 0
For J = 1 To M
If (J <> Q) Then
For I = 1 To N + M
If (I <> P) Then
T(J, I) = T(J, I) - T(J, P) * T(Q, I)
If (Abs(T(J, I)) <= ZT) Then T(J, I) = 0
End If
Next I
B(J) = B(J) - T(J, P) * B(Q)
If (Abs(B(J)) <= ZT) Then B(J) = 0
End If
Next J
For J = 1 To M
T(J, P) = 0
Next J
T(Q, P) = 1
A = R(P)
For I = 1 To N + M
R(I) = R(I) - A * T(Q, I)
If (Abs(R(I)) <= ZT) Then R(I) = 0
Next I
R(P) = 0
Z = Z + A * B(Q)
End If
Wend
If (Phase = 3) Then
' Decode solution
For I = 1 To N
X(I) = 0
For J = 1 To M
If (S(J) = I) Then X(I) = B(J)
Next J
ActiveSheet.Cells(3 + M, I).Value = X(I)
Next I
ActiveSheet.Cells(3 + M, N + 1).Value = Z
Else
For I = 1 To N
ActiveSheet.Cells(3 + M, I).Value = ""
Next I
ActiveSheet.Cells(3 + M, N + 1).Value = "No Solution"
End If
Application.ScreenUpdating = True
Application.Calculation = xlCalculationAutomatic
Application.EnableEvents = True
End Sub
You can contact me at Master@AlanCooper.id.au if you feel the need.