but fix with value

Author: scheinerman

src/Solver.jl

@@ -86,7 +86,7 @@ Solve `T` using the simplex method.
 """
 function simplex_solve!(T::Tableau, verbose::Bool=true)
     stat = status(T)
-    if stat == :infeasible 
+    if stat == :infeasible
         if verbose
             @info "Infeasible basis detected. Restoring tableau to its initial state."
         end

src/Status.jl

@@ -21,10 +21,11 @@ end
     status(T::Tableau)::Symbol
 
 Return an indicator for the status of the tableau `T` being one of:
-* `:no_basis` -- no basis has been established for this tableau.
-* `:optimal` -- the tableau has reached a global minimization point.
+* `:no_basis` -- no basis has been established for this tableau
 * `:feasible` -- the tableau is in a feasible state, but not optimal (rhs is nonnegative).
 * `:infeasible` -- the tableau is in an infeasible state (rhs contains negative values).
+* `:optimal` -- the tableau has reached a global minimization point.
+* `:unbounded` -- no pivots are possible; objective function can be arbitrarily negative.
 """
 function status(T::Tableau)::Symbol
     # check for a basis
@@ -40,6 +41,10 @@ function status(T::Tableau)::Symbol
 
     # if not optimal, but feasible
     if all(x .>= 0)
+        i, j = find_pivot(T)
+        if i==0
+            return :unbounded
+        end
         return :feasible
     end
 

src/Tableau.jl

@@ -120,6 +120,16 @@ Return the value of the LP in `T` at the point `x`.
 This can also be invoked as `T(x)`.
 """
 function value(T::Tableau, x::Vector)
+    n = T.n_vars
+
+    if length(x) > n
+        x = x[1:n]     # trim if too long
+    end
+
+    if length(x) < n   # pad if too short
+        x = vcat(x, zeros(Int, n-length(x)))
+    end
+
     return T.c' * x
 end
 
@@ -129,8 +139,7 @@ end
 Return the value of the LP at the current basic vector. 
 """
 function value(T::Tableau)
-    x = basic_vector(T)
-    return value(T, x)
+    return T.M[1, end]
 end
 
 function (T::Tableau)(x::Vector)

test/runtests.jl

@@ -18,6 +18,8 @@ using SimplexTableaux
 
     set_basis!(T, [3, 4])
     @test !in_feasible_state(T)
+
+    lp_solve(T)
 end
 
 @testset "Duality" begin
@@ -28,10 +30,25 @@ end
     T = Tableau(A, b, c)
     dT = dual(T)
 
-    simplex_solve!(T, false)
+    x = simplex_solve!(T, false)
     v = value(T)
+    @test T(x) == v
 
     simplex_solve!(dT, false)
     dv = value(dT)
     @test v == -dv
 end
+
+@testset "Status" begin
+    A, b, c = ([-1 -2 0 0; -3 -4 0 0; 0 0 1 2; 0 0 3 4], [2, 2, 4, 5], [2, 3, -4, -7])
+    T = Tableau(A, b, c)
+    @test status(T) == :no_basis
+    set_basis!(T, [2, 3, 4, 5])
+    @test status(T) == :infeasible
+
+    T = Tableau([2 1 0 9 -1; 1 1 -1 5 1], [9, 7], [2, 4, 2, 1, -1])
+    set_basis!(T, [1, 2])
+    @test status(T) == :feasible
+    simplex_solve!(T, false)
+    @test status(T) == :unbounded
+end