A Simplex Command for Mathematica and a few ways it could possibly be utilized

© Copyright Joseph M. Herrmann 1993

Abstract:  
This notebook simplex.ma contains a simplex command which produces a simplex tableau for a linear programming problem.  Examples of its use to solve a standard maximization problem, find multiple optimal feasible solutions, solve linear programming problems by the Big M method, and do a sensitivity analysis are included.  The simplex command does not currently have error checking capabilities and assumes that the matrix and vector arguments are correct.

This notebook may be useful to those using the simplex method in their work and educators and students learning the simplex method to solve linear programming problems.

As you work through the notebook, please execute the Mathematica commands which appear in the smaller boldface type.  You should first execute the Simplex command which follows.

The Simplex Command

Simplex[A_,P_,V_]:=Module[{AT,ML,NL,VL,BT,B,BI,Aa,AaT,
	Bot,VT,Vc,VcL,Pa,var,vartop,S},
	AT=Transpose[A];
	ML=Length[A];
	NL=Length[AT];
	VL=Length[V];
	BT={};
	Do[BT=Append[BT,AT[[V[[i]]]]],{i,1,VL}];
	B=Transpose[BT];
	BI=Inverse[B];
	Aa=BI.A;
	AaT=Transpose[Aa];
	Bot={};
	Do[Bot=Append[Bot,P[[V[[i]]]]],{i,1,VL}];
	VT=Table[i,{i,1,NL-2}];
	Vc=Complement[VT,V];
	VcL=Length[Vc];
	Pa=Table[0,{i,1, NL}];
	Do[Pa[[Vc[[i]]]]=P[[Vc[[i]]]]-Bot.AaT[[Vc[[i]]]],
	  {i,1,VcL}];
	Pa[[NL]]=P[[NL]]-Bot.AaT[[NL]];
	Pa[[NL-1]]=1;
	var=Table["x"<>ToString[V[[i]]],{i,1, ML}];
	AaT=Prepend[AaT,var];
	vartop=Table["x"<>ToString[i],{i,1, NL-2}];
	vartop=Append[vartop,"P"];
	vartop=Append[vartop," "];
	vartop=Prepend[vartop," "];
	Aa=Transpose[AaT];
	Aa=Prepend[Aa,vartop];
	Pa=Prepend[Expand[Pa],"P"];
	S=Append[Aa,Pa];
	TableForm[S]]

This command assumes a standard maximization linear programming problem which consists of maximizing an objective function with inequalities of the form ≤ positive constants and all variables assumed nonnegative.  In Simplex[A_,P_,V_], the matrix A consists of the coefficients of the constraint equations  formed from the inequalities by introducing nonnegative slack variables with the last column of the matrix A the positive constants on the right hand side of the ≤ inequalties.  The vector P consists of the coefficients of the objective function written with all coefficients on one side of the equation and the variable to be maximized positive.   The vector V consists of the column number of the variables which are to be solved to be possibly positive.

For example, considering the standard maximization linear  programming problem

Objective function: 		Maximize  P=x+2y
subject to the constraints:    		  x+5y ≤ 85
         					2x+8y≤ 140
              				        x≤ 50	
					     x≥0, y≥0

is modified by the introduction of slack variables to 

Objective function: 		         Maximize  -x-2y+P=0
subject to the constraint equations    		  x+5y+s1= 85
						2x+8y+s2= 140
						         x+s3= 50	
					x≥0, y≥0, s1≥0, s2≥0, s3≥0

For the simplex, the variables are now considered as x1, x2, x3, x4, and x5 so that the constraint equations become:

                     x1+5x2+x3= 85
                   2x1+8x2+x4= 140
                             x1+x5= 50	

where all the variables are assumed nonnegative.  The objective function becomes:

                      -x1-2x2+P= 0.

We now let A={{1, 5, 1, 0, 0, 0, 85},{2, 8, 0, 1, 0, 0, 140},{1, 0, 0, 0, 1, 0, 50}} and P={-1, -2, 0, 0, 0, 1, 0} then Simplex[B,Q,{1,2,3}]
produces the following simplex tableau which solves the constraint equations for x1, x2, and x3 assuming  x4 and x5 are zero.

Note that if the simplex tableau begins to wrap around at the margin of this notebook, you may widen this notebook by pulling on the lower right-hand corner while holding down the left menu button on the mouse and dragging the mouse to the right.

For a brief introduction of how to the use of the simplex method to solve a standard maximization linear programming problem, see the following section.   If you are unfamiliar with the simplex method, you may want to consult a finite mathematics textbook such as Raymond A. Barnett and Michael R. Ziegler, Finite Mathematics for Management, Life, and Social Sciences, Fourth Edition, MacMillan Publishing Co, 1987, Chapter 2.  For further examples of the use of the Simplex command, see the sections below.

Calculation of a Simplex for a Standard Maximization
 Linear Programming Problem:  Maximizing an objective function with inequalites of the form ≤ a positive constant

A standard maximization linear programming problem consists of maximizing an objective function with inequalities of the form ≤ a positive constant and all variables assumed nonnegative.

We shall consider the following example:

Objective function: 		Maximize  P=x+2y
subject to the constraints:    		  x+5y ≤ 85
						2x+8y≤ 140
							x≤ 50	
						     x≥0, y≥0

If slack variables are introduced then we obtain

Objective function: 		Maximize  P=x+2y
subject to the constraints:    		  x+5y+s1= 85
						2x+8y+s2= 140
							x+s3= 50	
						x≥0, y≥0, s1≥0, s2≥0, s3≥0

To use the Simplex command we will rename x, y, s1, s2, and s3 as x1, x2, x3, x4, and x5 respectively.  Therefore, the constraint equations are now:

                     x1+5x2+x3= 85
                   2x1+8x2+x4= 140
                             x1+x5= 50	

with the objective function:

                      -x1-2x2+P= 0

Thus the matrix of coefficients for the constraint equations becomes:

B={{1, 5, 1, 0, 0, 0, 85},
   {2, 8, 0, 1, 0, 0, 140},
   {1, 0, 0, 0, 1, 0, 50}}

Note that if a variable doesn't appear in the constraint equation its coefficient is zero in the matrix.  In particular, note that since the variable P doesn't appear in the constraint equations that its column in the matrix of coefficients consists of zeros.

Similarly,  the vector of coefficients of the objective function is:

Q={-1, -2, 0, 0, 0, 1, 0}

Finally, if we wish to solve for x1, x2, and x3 while  x4 and x5 are zero then the following command will give us the resulting simplex:

Simplex[B,Q,{1,2,3}]

From this simplex, we see that x1=50, x2=5, x3=10, x4=0, x5=0, and P=60.  Furthermore, this is the optimal solution since the bottom row is entirely positive and all variables are nonnegative.

There are many possible ways this simplex command could be used.  The first is to use it to do each simplex in the simplex algorithm until a solution is determined.

For the example above this would be done as follows:

First we begin at the initial simplex when x=x1 and y=x2 are set to 0.  Thus we will determine the simplex where s1=x3, s2=x4, and s3=x5 are possibly nonzero.  This is called the initial tableau.

Simplex[B,Q,{3,4,5}]

We now consider the next simplex where x2 will become positive and x3 zero since the pivot element is in the 1st row, 2nd column position.  (If you are unfamiliar with how to find the pivot element or read the simplex tableau, please consult a textbook which discusses the simplex method such as mentioned above.)  The following command determines this next simplex.

Simplex[B,Q,{2,4,5}]

We now consider the next simplex where x1 will become positive and x4 zero since the pivot element is in the 2nd row, 1st column position.  The following command determines the next simplex.

Simplex[B,Q,{2,1,5}]

We now consider the next simplex where x3 will become positive and x5 zero since the pivot element is in the 3rd row, 3rd column position.  The following command determines the next simplex.

Simplex[B,Q,{2,1,3}]

We have now found the optimal solution we determined before.

How Multiple Optimal Feasible Solutions may be Identified Using the Simplex Command

Another useful way to use this command is to use it together with the Mathematica ConstrainedMax command as follows:

ConstrainedMax[x+2y,{x+5y<=85, 2x+8y<=140, x<=50},{x,y}]

Note the ConstrainedMax command determines the optimal solution as P=60, x=50, and y=5.  How is the simplex tableau useful?

The simplex tableau at an optimal feasible solution may alert us to multiple solutions which the ConstrainedMax command does not indicate.  Furthermore, the Simplex command allows the possibility of a sensitivity analysis of the amount of variation in the various elements of our matrices and vectors and still have the optimal feasible solution remain an optimal feasible solution.

As an example consider the following problem:

 		Maximize  P=4x+18y
such that                   4x+3y ≤ 42
                                 2x+9y ≤ 36
                                         x ≤ 10
                                     x≥0, y≥0

The ConstrainedMax command determines the following solution:

ConstrainedMax[4x+18y,{4x+3y<=42, 2x+9y<=36,x<=10},{x,y}]

However if we use the Simplex command (remember to introduce slack variables) to solve for x, y, and s3 then we find:

Simplex[{{4,3,1,0,0,0,42},{2,9,0,1,0,0,36},{1,0,0,0,1,0,10}},
	{-4,-18,0,0,0,1,0},{1,2,5}]

The extra zero occuring in the bottom row, x3 column alerts us to the fact that there are multiple solutions to this problem.  This extra zero shows that an optimal feasible solution is possible when x3 is positive.  If we pivot on the x3 column, we find the pivot element to be the 1st row, 3rd column position and so x1 will become zero as x3 becomes positive.  The solution for positive x2 and x3 can be found by:

Simplex[{{4,3,1,0,0,0,42},{2,9,0,1,0,0,36},{1,0,0,0,1,0,10}},
	{-4,-18,0,0,0,1,0},{3,2,5}]

Again this is an optimal feasible solution with P=72 but now x1=0 and x2=4.  So all points on the line segment between x1=0, x2=4 and x1=9, x2=2 will be optimal feasible solutions with P=72.

Sensitivity Analysis of a Linear Programming Problem Using the Simplex Command

A sensitivity analysis is also easy using the Simplex command.  Let's reconsider our original linear programming problem:

Objective function: 		Maximize  P=x+2y
subject to the constraints:    		  x+5y ≤ 85
						2x+8y≤ 140
							x≤ 50	
						     x≥0, y≥0

How much could the coefficients change in the objective function and still have the feasible solution x=50, y=5 be optimal?  Note this can be determined by letting the objective function be  P=ax+by and then solving the following simplex tableau.  Note the -a ,-b in the vector Q corresponding to the coefficients of this new objective function in the Simplex command.

B={{1, 5, 1, 0, 0, 0, 85},
   {2, 8, 0, 1, 0, 0, 140},
   {1, 0, 0, 0, 1, 0, 50}}
Q={-a,-b,0,0,0,1,0}
Simplex[B,Q,{1,2,3}]

If b/8 and a-b/4 are positive then this is the simplex tableau at the optimal feasible solution x=50 and y=5.  Thus if b is positive and a/b > 1/4, then x=50 and y=5 is the optimal solution.   So if each unit of x has a profit of $1, (i.e., a=1) then b<4, or equivalently,  the profit of each unit of y must be less than $4 if x=50 and y=5 are to remain the optimal solution to the above problem. 

Similarly, we can ask how the solutions will change if additional steel or glass is available or if our limit on the number of watches which can be produced is changed.  In this case we change the values 85,  140, and 50 at the end of the inequalities to a, b, and c with the corresponding change in the matrix B and execute the Simplex command again.

B={{1, 5, 1, 0, 0, 0, a},
   {2, 8, 0, 1, 0, 0, b},
   {1, 0, 0, 0, 1, 0, c}}
Q={-1,-2,0,0,0,1,0}
Simplex[B,Q,{1,2,3}]

In this case we see that as long as c,  b/8 - c/4, and a - 5b/8 + c/4 are positive then the simplex tableau is at an optimal feasible solution (since the values of x1, x2, and x3 are nonnegative and the bottom row is nonnegative).  Thus if c is positive, b/c ≥ 2 and a/c ≥ (5/8)(b/c) - 1/4, then x=c and y=b/8 - c/4 is the optimal solution.  If c=50 and b=140 and a≥75, then we have x=50 and y=5 is the optimal solution.

The Big M Method Using the Simplex Command

The simplex command will also work for the big M method.  Consider the linear programming problem:

 		Minimize  C=40x1+60x2+50x3
such that                  4x1+x2+x3 ≤ 30
                               2x1+x2+5x3 ≥ 20
                                          6x1-x2=10   
                                     x1≥0, x2≥0, x3≥0

First, we rewrite the problem in a form suitable for solution by the simplex method using the big M method.  After introduction of the slack variable  s1, the surplus variable  s2 and the artificial variables  a1 and  a2, the constraints become 

                                       4x1+x2+x3+s1= 30
                               2x1+x2+5x3-s2+a1 = 20
                                               6x1-x2+a2=10   
                             x1≥0, x2≥0, x3≥0, s1≥0, s2≥0, a1≥0, a2≥0
                             
with the original objective function C=40x1+60x2+50x3 which was to be minimized, changed to an objective function to be maximized by multiplying (-1)*(40x1+60x2+50x3).  Thus we could solve our original linear programming problem of minimizing the objective function C=40x1+60x2+50x3 by maximizing  P=-40x1-60x2-50x3.

The objective function P=-40x1-60x2-50x3 is modified by subtracting from the objective function each artificial variable multiplied by a large positive coefficient. 

		Maximize  P=-40x1-60x2-50x3-Ma1-Ma2

Note that the large positive coefficient was designated by the letter M.  If the linear programming problem has a solution and if M is chosen sufficiently large then maximizing the objective function will lead to the artificial variables becoming zero. 

If we begin at the corner point  s1= 30, a1 = 20, a2=10 we must have the objective function reflect how each of these variables depend on  x1, x2, x3, and s2.  However since  s1 doesn't appear in the objective function it has no effect.  So we need to solve only for  a1 and  a2 and we find 

                               a1 = 20 - 2x1- x2- 5x3 + s2
                               a2 = 10 - 6x1 + x2   

If we then substitute these expressions into the objective function, we find

	Maximize  P=-40x1-60x2-50x3-M(20 - 2x1- x2- 5x3 + s2)
	                              -M(10 - 6x1 + x2)
		                    =(8M-40)x1-60x2+(5M-50)x3-Ms2-30M

We are now ready to use the simplex command to solve for the variables s1, a1, and a2 in the fourth, sixth, and seventh columns.

B={{4, 1, 1,  1,  0, 0, 0, 0, 30},
   {2, 1, 5,  0, -1, 1, 0, 0, 20},
   {6,-1, 0,  0,  0, 0, 1, 0, 10}}
Q={-8M+40, 60, -5M+50, 0, M, 0, 0, 1, -30M}
Simplex[B,Q,{4,6,7}]

Since M is a very large positive number then 40-8M is the most negative entry in the bottom row and the pivot element is the 6 in the x1 column, x7 row.  Thus we shall solve for x1 and let x7=0.

Simplex[B,Q,{4,6,1}]

Since we are not at the optimal solution, we continue pivoting

Simplex[B,Q,{4,3,1}]

We now find that all numbers in the bottom row under the x1through the x7 columns are positive and so we have arrived at the optimal solution  x1=5/3,  x2=0, and  x3=10/3 with the optimal value of P=-700/3.  Since C=-P then the mimimum value of the  objective function of our original linear programming problem is 700/3 at this optimal solution.

Note that the ConstrainedMin command determines the same result.

ConstrainedMin[40x1+60x2+50x3,{4x1+x2+x3<=30,
	 2x1+x2+5x3>=20, 6x1-x2==10},{x1,x2,x3}]