Simplex Algorithm

[ftneek]

simplex2 incorporates Dual Simplex algorithm. Using (non dual) simplex, in the end if c row is all positive we assume we have reached the optimum. But if there is a negative in b column, the answer is considered infeasible as the assumption is all variables are non negative. Using dual simplex, if the "final tableau" would have a negative b iteration continues. If a c becomes negative again, regular simplex is used until all c are positive and there is a negative b.

So I think b>=0 requirement means only regular simplex will be used.

For whatever reason test 3 kept returning 0 on my virtual and physical calculator despite both sides of "==" being visually the same. I believe there to be some kind of issue, at least on my end. Why would "==" return 0 if the objects are the same? It seemed to work fine on the other 2 tests. Even if I were to reset my virtual and physical calcs, I hope to figure out the cause so that this issue doesn't randomly pop up again sometime in the future.

Using size(difference())==0 seems to work.. But I guess using this, it would be possible for a false pass if the values are all the same but ordered differently.

Code:

#cas

test_simplex():=
BEGIN
RETURN {test_case_1(),test_case_2(),test_case_3(),test_case_4(),test_case_5(),test​_case_6()};
END;


test_case_1():=
BEGIN
// Example 3.7.2
// artificial variables
// s.t.  x1,  x2,  x3,  x4 ≥ 0
//       x1 +2x3      + x4 = 20
//      2x1 + x2 + x3      = 10
//      -x1 +4x2 -2x3 +3x4 = 40
// min   x1 +4x2 +3x3 +2x4 = z
LOCAL expected,result;
expected:=[35,[5,0,0,15,0,0,0],[[1,1/2,1/2,0,3/4,0,-1/4,5],[0,0,0,0,-3/2,1,1/2,0],[0,3/2,-1/2,1,1/4,0,1/4,15],[0,1/2,7/2,0,-5/4,0,-1/4,-35]]];
result:=simplex2([[1,2,0,1,1,0,0,20],[2,1,1,0,0,1,0,10],[-1,4,-2,3,0,0,1,40],[1,4,3,2,0,0,0,0],[0,0,0,0,1,1,1,0]],{5,6,7},7,3,0);
RETURN size(difference(result,expected))==0;
END;

test_case_2():=
BEGIN
// -slack variables, artificial variables
// s.t.  x,  y ≥ 0
//      2x +4y ≥ 22
//      3x +2y ≥ 20
//      4x +5y ≥ 40
// min  6x +5y = z
LOCAL expected,result;
expected:=[320/7,[20/7,40/7,46/7,0,0,0,0,0],[[0,1,0,4/7,-3/7,0,-4/7,3/7,40/7],[1,0,0,-5/7,2/7,0,5/7,-2/7,20/7],[0,0,1,6/7,-8/7,-1,-6/7,8/7,46/7],[0,0,0,10/7,3/7,0,-10/7,-3/7,-320/7]]];
result:=simplex2([[2,4,-1,0,0,1,0,0,22],[3,2,0,-1,0,0,1,0,20],[4,5,0,0,-1,0,0,1,40],[6,5,0,0,0,0,0,0,0],[0,0,0,0,0,1,1,1,0]],{6,7,8},8,3,0);
RETURN size(difference(result,expected))==0;
END;

test_case_3():=
BEGIN
// slack varaibles
// s.t.  x,  y ≥ 0
//     -2x -4y ≤-22
//     -3x -2y ≤-20
//     -4x -5y ≤-40
// min  6x +5y = z
LOCAL expected,result;
expected:=[320/7,[20/7,40/7,46/7,0,0],[[0,0,1,6/7,-8/7,46/7],[1,0,0,-5/7,2/7,20/7],[0,1,0,4/7,-3/7,40/7],[0,0,0,10/7,3/7,-320/7]]];
result:=simplex2([[-2,-4,1,0,0,-22],[-3,-2,0,1,0,-20],[-4,-5,0,0,1,-40],[6,5,0,0,0,0]],{3,4,5},5,0,0);
RETURN size(difference(result,expected))==0;
END;

test_case_4():=
BEGIN
// s.t.  x,  y,  z ≥ 0
//      3x +2y + z = 10
//      2x +5y +3z = 15
// min -2x -3y -4z = Z
LOCAL expected,result;
expected:=[-130/7,[15/7,0,25/7,0,0],[[1,1/7,0,3/7,-1/7,15/7],[0,11/7,1,-2/7,3/7,25/7],[0,25/7,0,-2/7,10/7,130/7]]];
result:=simplex2([[3,2,1,1,0,10],[2,5,3,0,1,15],[-2,-3,-4,0,0,0],[0,0,0,1,1,0]],{4,5},5,2,0);
RETURN size(difference(result,expected))==0;
END;

test_case_5():=
BEGIN
// artificial variables, non zero offset
// s.t.  x1,  x2,  x3,  x4,  x5,  x6,  x7 ≥ 0
//      4x1 - x2       +x4 - x5 + x6          = 6
//     -7x1 +8x2 + x3                - x7     = 7
//       x1 + x2      +4x4 -4x5               = 12
// min -3x1 +2x2 + x3 - x4 + x5           +87 = z
LOCAL expected,result;
expected:=[1441/17,[131/85,189/85,0,35/17,0,0,0,0,0],[[1,0,1/17,0,0,32/85,-1/17,1/17,-8/85,131/85],[0,1,3/17,0,0,28/85,-3/17,3/17,-7/85,189/85],[0,0,-1/17,1,-1,-3/17,1/17,-1/17,5/17,35/17],[0,0,13/17,0,0,5/17,4/17,-4/17,3/17,-1441/17]]];
result:=simplex2([[4,-1,0,1,-1,1,0,0,0,6],[-7,8,1,0,0,0,-1,1,0,7],[1,1,0,4,-4,0,0,0,1,12],[-3,2,1,-1,1,0,0,0,0,-87],[0,0,0,0,0,0,0,1,1,0]],{6,8,9},9,2,0);
RETURN size(difference(result,expected))==0;
END;

test_case_6():=
BEGIN
// slack variables 'dual problem'
// s.t.  z1,z2,z3 ≥ 0
//      2z1 + 3z2 + 4z3 ≤ 6
//      4z1 + 2z2 + 5z3 ≤ 5
// max 22z1 +20z2 +40z3 = c
LOCAL expected,result;
expected:=[-320/7,[0,10/7,3/7,0,0],[[-6/7,1,0,5/7,-4/7,10/7],[8/7,0,1,-2/7,3/7,3/7],[46/7,0,0,20/7,40/7,320/7]]];
result:=simplex2([[2,3,4,1,0,6],[4,2,5,0,1,5],[-22,-20,-40,0,0,0]],{4,5},5,0,0);
RETURN size(difference(result,expected))==0;
END;

#end

edit: I added the rest of the examples shared in the threads as tests and added the problem statements to the comments. For any change to simplex2, test_simplex() should still return all 1's (at minimum). Feel free to add more examples as you encounter them. If you want to check an answer before adding it, lpsolve on xcas can be used as a reference. It might be nice to have a method that can read problems stored in a file, once we have a method to convert a system to standard form. That would probably make it easy to transfer problems found online to the calculator if needed. If anyone has a hint for not exposing all the test cases in the catalog, please let me know.