Linear optimization

Linear programming using the Simplex method

Maximize M = 40 x1 + 60 x2 

subject to:
2 x1 + x2 ≤ 70
x1 + x2 ≤40
x1 + 3 x2 ≤90

x≥ 0

Here are the constraints

constraints=Plot[{70- 2 x,40-x,90/3-x/3},{x,0,40},PlotStyle->{Dashing[{.02,.02}],Dashing[{.03,.01}],Dashing[{1,0}]},AxesLabel->{"x1","x2"},AspectRatio->1]

We see that the region in which a valid solution exists is the inside region bounded by the constraints.

Here is a plot of the objective function

objective=Plot[{-2/3 x+500/60,-2/3 x+1000/60,-2/3 x+1500/60,
	-2/3 x+2100/60			},{x,0,40},PlotStyle->{ Dashing[{.04,.04}]},AspectRatio->1]

The objective function is increasing upward.  Looking at the constraints and the monotonically increasing objective function we note that the optimimum will occur on a boundary and will for all reasonably cases be a vertex.

Show[constraints,objective,AspectRatio->1]

We see where the optimimum is!, x1=15, x2=25. Of course we can't solve every problem graphically so we need a procedure.

Express the problem in an augmented matrix form with the vector being
{M, x1, x2, x3, x4, x5} and the augmented column {70, 40, 90, 0} the last entry (which will be M) being the current value of the objective function.

Start with the real variables set to zero (x1=x2=0) which is a known trivial solution.  The zero variables will be called nonbasic variables these will change as the procedure progresses.  To eliminate the need to deal with inequalities and because we note that if a constraint is limiting the solution, we will be at the limiting value (i.e., it will be an equality) we will introduce slack variables to make all of the inequalities into equations.  For this problem we need three slack variables, x3, x4 and x5.  If we start at x1 = x2 = 0, then the initial values of the slack variables are x3=70, x4=40 and x5=90.  This solution is found from inspection.  At any step the non zero variables will be called basic variables.

Here is the initial augmented matrix.  There is some convenience to entering it line by line.  (We can use a subscripted variable.)

z[1][1]= {0,2,1,1,0,0,70}

z[1][2]={0,1,1,0,1,0,40}

z[1][3]={0,1,3,0,0,1,90}

z[1][4]={1,-40,-60,0,0,0,0}

Here we write it as a matrix.

a1=Table[z[1][i],{i,1,4}]

We see, from the objective function, M = 40 x1 + 60 x2, that increasing x2 is the way to get the biggest impact.  The idea of the simplex method is to move only if the objective function will increase and then to always take the biggest increase.  If instead of looking that the original objective function we choose to look at the augmented matrix, we should choose to remove the largest negative in the last row and continue doing this until there are no more negatives.  

So, we will move to the next vertex.  To do this, which of x3, x4, x5 should be set to 0?  

ANS: The one that represents the strongest constraint on the value of x2!  

This is found by the getting the smallest (non positive) quotient between the x2 column (3rd column) and the b vector (7th column)

We make a label for the third column

kk=3;

Here we calculate the quotients which are the values of x2 at the intersection with the three constraints.

ztest=Table[z[1][i][[7]]/z[1][i][[kk]],{i,1,3}]

We see that we need to use the last constraint.  This means that x5 will be as small as possible = 0 (i.e., no slack!) and the third constraint becomes active.  Note that each vertex represents a point where two of the variables are 0 and the other three are nonzero. (Why is this?)  Thus any solution that we examine must have this form.

The new value of x2 = 90/3 = 30.

We need to identify this active constraint with an index value.  This command will find the position of the minimum of the smallest value in ztest.

j=Position[ztest, Min[ztest]][[1,1]]

Now we want a new set of equations that show the change of variables,
  i.e.,  showing that x2 has been increased in favor of x5 becoming 0.
  
 Keeping x1 = 0 but changing x2 requires that changes be made in x3, x4 and x5.
 x3 = 70 - x2 = 40,
 x4 = 40 - x2 = 10,
 x5 = 90 - 3 x2  = 0 
 
So that x1 = 0, x2 = 30, x3 = 40, x4 = 10, x5 =0 is a new basic feasible solution.

Before we continue on, we need to process the equations above to reflect the changes.

Since we know that each solution has 2 zero variables and 3 nonzero variables, we require that each basic (nonzero) variable has only 1 non zero entry in its column in the augmented matrix and this is not in the bottom row.  This way, each nonzero variable will be the only contribution to the value in its equation.  The current values of the nonzero variables will be immediately evident.  The same solution could exist, if say we added two rows and messed up our convention, but we would loose the convenience of always knowing the solution.  Thus we will not do this.   

To make x2 active (nonzero) with a values at the 3ed constraint, we do a Gaussian elimination using the x2 (i.e., kk=3)  from the third row, j=3.

Do[If[i!=j, z[2][i]=z[1][i]-z[1][j]  z[1][i][[kk]]/ z[1][j][[kk]],z[2][i]=z[1][i]/z[1][j][[kk]]],{i,1,4}]

Here is the new augmented matrix.

a2=Table[z[2][i],{i,1,4}]

We see immediately that x2=30, x3=40 and x4=10 (since x1=x5=0).  We are thus assured that we at a contrainst intersection as we wish.

Note that the objective function is 1800 as required and it is expressed by the non basic variables (which are 0)

Now we need to increase x1 as far as possible.  Check to see which constraint becomes active?

kk=2;

ztest=Table[z[2][i][[7]]/z[2][i][[kk]],{i,1,3}]

It will be the second one and we need to identify this active constraint with an index value.

j=Position[ztest, Min[ztest]][[1,1]]

The second constraint will become active first so that x4 will go to 0.

Do[If[i!=j, z[3][i]=z[2][i]-z[2][j]  z[2][i][[kk]]/ z[2][j][[kk]],z[3][i]=z[2][i]/z[2][j][[kk]]],{i,1,4}]

a3=Table[z[3][i],{i,1,4}]

So the objective function is 2100! There are no remaining "negatives" in the last row so that no improvement can be attained.

 
So we are done!  For x4=x5=0, the solution for {x1,x2,x3} = {15, 25, 15}.  We see that x2 
was reduced a little to get a big gain from x1.  Now x4 and x5=0 which correspond to the constrints,
x1 + x2 <=40
x1 + 3 x2 <=90 
being active.   

Here is a check of the objective function!

15*40+25*60

Of course, Mathematica did not really need our help

ConstrainedMax[40 x1 + 60 x2,{2 x1 + x2 <=70,
		x1 + x2 <=40, x1 + 3 x2<=90},{x1,x2}]

Maximize M = 20 x1 + 30 x2 

subject to:
.1x1 + .075 x2 <= 1500
.08 x1 + .10 x2 <=1200
.05 x1 + .15 x2 <=1500

x<= 0

Here are the constraints

constraints=Plot[{1500/.075- .10/.075 x,1200/.10-.08x/.10,1500/.15-.05 x/.15},{x,0,30000},PlotStyle->{Dashing[{.02,.02}],Dashing[{.03,.01}],Dashing[{1,0}]},AxesLabel->{"x1","x2"},AspectRatio->1]

We see that the region in which a valid solution exists is the inside region bounded by the constraints.

Here is a plot of the objective function

objective=Plot[{-2/3 x+1000,-2/3 x+10000,
	-2/3 x+20000			},{x,0,30000},PlotStyle->RGBColor[1,0,1],AspectRatio->1]

The objective function is increasing upward.  Looking at the constraints and the monotonically increasing objective function we note that the optimimum will occur on a boundary and will for all reasonably cases be a vertex.

Show[constraints,objective,AspectRatio->1]

We see where the optimimum is!, x1=15, x2=25. Of course we can't solve every problem graphically so we need a procedure.

Maximize M = 20 x1 + 30 x2 

subject to:
.1x1 + .075 x2 <= 1500
.08 x1 + .10 x2 <=1200
.05 x1 + .15 x2 <=1500

x<= 0

Here is the initial augmented matrix.

a={{0,.1,.075,1,0,0,1500},
{0,.08,.1,0,1,0,1200},{0,.05,.15,0,0,1,1500},
{1,-20,-30,0,0,0,0}}

We see, from the objective function, M = 40 x1 + 60 x2,  that increasing x2 is the way to get the biggest impact.  The idea of the simplex method is to only move if the objective function will increase and then to always take the biggest increase.

So, we will move to the next vertex.  To do this, which of x3, x4, x5 should be set to 0?  

ANS: The one that represents the strongest constraint on the value of x2!  

This is found by the getting the smallest (non positive) quotient between the x2 column (3rd column) and the b vector (7th column)

ztest=Table[a[[i,7]]/a[[i,3]],{i,1,3}]

We need to identify this active constraint with an index value.

j=Position[ztest, Min[ztest]][[1,1]]

We see that we need to use the last constraint.  This means that x5 will be as small as possible =0 (no slack!) and the third constraint becomes active.  Note that each vertex represents a point where two of the variables are 0 and the other three are nonzero. (Why is this?)  Thus any solution that we examine must have this form.

 The new value of x2 = 90/3 = 30.  
 
 Now we want a new set of equations that show the change of variables,
  i.e.,  showing that x2 has been increased in favor of x5 becoming 0.
  
 Keeping x1 =0, but changing x2 requires that changes be made in x3, x4 and x5.
 x3 = 70 - x2 = 40,
 x4 = 40 - x2 = 10,
 x5 = 90 - 3 x2  =0 
 
So that x1 = 0, x2 = 30, x3 = 40, x4 = 10, x5 =0 is a new basic feasible solution.

Before we continue on, we need to process the equations above to reflect the changes.

Since we know that each solution has 2 zero variables and 3 nonzero variables, we require that each basic (nonzero) variable has only 1 non zero entry in its column in the augmented matrix and this is not in the bottom row.  This way, eaxh nonzero variable will be the contribution to the value on its equation.  The current values of the nonzero variables will be immediately evident.   

Making these changes is the same as a Gaussian elimination using the x2 from the third row.

Here we do the elimination

Do[If[i!=j, zx[i]=a[[i]]-a[[j]] a[[i,j]]/a[[j,j]],zx[i]=a[[i]]/a[[j,j]]];
	Print[zx[i]],{i,1,4}]

a2=Table[zx[i],{i,1,4}]

We see immediately that x2=30, x3=40 and x4=10 (since x1=x5=0).  We are thus assured that we at a contrainst intersection as we wish.

Note that the objective function is 1800 as required and it is expressed by the non basic variables (which are 0)

Now we need to increase x1 as far as possible.  Check to see which constraint becomes active?

ztest=Table[a2[[i,7]]/a2[[i,2]],{i,1,3}]

FullForm[ztest]

Min[ztest]

The second constraint will become active first so that x4 will go to 0.

Now choose j=2

j=2;

Do[If[i!=j, z[s][i]=a2[[i]]-a2[[j]] a2[[i,j]]/a2[[j,j]],z[s][i]=a2[[i]]/a2[[j,j]]];
	Print[zx[i]],{i,1,4}]

a3=Table[z[5][i],{i,1,4}]

a3={a2[[1]]-a2[[2]] a2[[1,2]]/a2[[2,2]],a2[[2]]/a2[[2,2]],a2[[3]]-
a2[[2]] a2[[3,2]]/a2[[2,2]],a2[[4]]- a2[[2]] a2[[4,2]]/a2[[2,2]] }

There are no remaining "negatives" in the last row so that no improvement can be attained.

 
So we are done!!  For x4=x5=0, the solution for {x1,x2,x3} = {15,25,15}.  We see that x2 
was reduced a little to get a big gain from x1.  Now x4 and x5=0 which correspond to the constrints,
x1 + x2 <=40
x1 + 3 x2 <=90 
being active.   

Here is a check of the objective function!

15*40+25*60

Maximize M = 20 x1 + 30 x2 

subject to:
.1x1 + .075 x2 <= 1500
.08 x1 + .10 x2 <=1200
.05 x1 + .15 x2 <=1500

x<= 0

Here are the constraints

constraints=Plot[{1500/.075- .10/.075 x,1200/.10-.08x/.10,1500/.15-.05 x/.15},{x,0,30000},PlotStyle->{Dashing[{.02,.02}],Dashing[{.03,.01}],Dashing[{1,0}]},AxesLabel->{"x1","x2"},AspectRatio->1]

We see that the region in which a valid solution exists is the inside region bounded by the constraints.

Here is a plot of the objective function

objective=Plot[{-2/3 x+1000,-2/3 x+10000,
	-2/3 x+20000			},{x,0,30000},PlotStyle->RGBColor[1,0,1],AspectRatio->1]

The objective function is increasing upward.  Looking at the constraints and the monotonically increasing objective function we note that the optimimum will occur on a boundary and will for all reasonably cases be a vertex.

Show[constraints,objective,AspectRatio->1]

We see where the optimimum is!, x1=15, x2=25. Of course we can't solve every problem graphically so we need a procedure.

Maximize M = 20 x1 + 30 x2 

subject to:
.1x1 + .075 x2 <= 1500
.08 x1 + .10 x2 <=1200
.05 x1 + .15 x2 <=1500

x<= 0

Here is the initial augmented matrix.

a={{0,.1,.075,1,0,0,1500},
{0,.08,.1,0,1,0,1200},{0,.05,.15,0,0,1,1500},
{1,-20,-30,0,0,0,0}}

We see, from the objective function, M = 40 x1 + 60 x2,  that increasing x2 is the way to get the biggest impact.  The idea of the simplex method is to only move if the objective function will increase and then to always take the biggest increase.

So, we will move to the next vertex.  To do this, which of x3, x4, x5 should be set to 0?  

ANS: The one that represents the strongest constraint on the value of x2!  

This is found by the getting the smallest (non positive) quotient between the x2 column (3rd column) and the b vector (7th column)

ztest=Table[a[[i,7]]/a[[i,3]],{i,1,3}]

We need to identify this active constraint with an index value.

j=Position[ztest, Min[ztest]][[1,1]]

We see that we need to use the last constraint.  This means that x5 will be as small as possible =0 (no slack!) and the third constraint becomes active.  Note that each vertex represents a point where two of the variables are 0 and the other three are nonzero. (Why is this?)  Thus any solution that we examine must have this form.

 The new value of x2 = 90/3 = 30.  
 
 Now we want a new set of equations that show the change of variables,
  i.e.,  showing that x2 has been increased in favor of x5 becoming 0.
  
 Keeping x1 =0, but changing x2 requires that changes be made in x3, x4 and x5.
 x3 = 70 - x2 = 40,
 x4 = 40 - x2 = 10,
 x5 = 90 - 3 x2  =0 
 
So that x1 = 0, x2 = 30, x3 = 40, x4 = 10, x5 =0 is a new basic feasible solution.

Before we continue on, we need to process the equations above to reflect the changes.

Since we know that each solution has 2 zero variables and 3 nonzero variables, we require that each basic (nonzero) variable has only 1 non zero entry in its column in the augmented matrix and this is not in the bottom row.  This way, eaxh nonzero variable will be the contribution to the value on its equation.  The current values of the nonzero variables will be immediately evident.   

Making these changes is the same as a Gaussian elimination using the x2 from the third row.

Here we do the elimination

Do[If[i!=j, zx[i]=a[[i]]-a[[j]] a[[i,j]]/a[[j,j]],zx[i]=a[[i]]/a[[j,j]]];
	Print[zx[i]],{i,1,4}]

a2=Table[zx[i],{i,1,4}]

We see immediately that x2=30, x3=40 and x4=10 (since x1=x5=0).  We are thus assured that we at a contrainst intersection as we wish.

Note that the objective function is 1800 as required and it is expressed by the non basic variables (which are 0)

Now we need to increase x1 as far as possible.  Check to see which constraint becomes active?

ztest=Table[a2[[i,7]]/a2[[i,2]],{i,1,3}]

FullForm[ztest]

Min[ztest]

The second constraint will become active first so that x4 will go to 0.

Now choose j=2

j=2;

Do[If[i!=j, z[s][i]=a2[[i]]-a2[[j]] a2[[i,j]]/a2[[j,j]],z[s][i]=a2[[i]]/a2[[j,j]]];
	Print[zx[i]],{i,1,4}]

a3=Table[z[5][i],{i,1,4}]

a3={a2[[1]]-a2[[2]] a2[[1,2]]/a2[[2,2]],a2[[2]]/a2[[2,2]],a2[[3]]-
a2[[2]] a2[[3,2]]/a2[[2,2]],a2[[4]]- a2[[2]] a2[[4,2]]/a2[[2,2]] }

There are no remaining "negatives" in the last row so that no improvement can be attained.

 
So we are done!!  For x4=x5=0, the solution for {x1,x2,x3} = {15,25,15}.  We see that x2 
was reduced a little to get a big gain from x1.  Now x4 and x5=0 which correspond to the constrints,
x1 + x2 <=40
x1 + 3 x2 <=90 
being active.   

Here is a check of the objective function!

15*40+25*60

Maximize M = 40 x1 + 60 x2 

subject to:
2 x1 + x2 <= 70
x1 + x2 <=40
x1 + 3 x2 <=90

x<= 0

Here are the constraints

constraints=Plot[{1500/.07- .10/.07 x,1500/.10-.09x/.10,1900/.18-.06 x/.18},{x,0,30000},PlotStyle->{Dashing[{.02,.02}],Dashing[{.03,.01}],Dashing[{1,0}]},AxesLabel->{"x1","x2"},AspectRatio->1]

We see that the region in which a valid solution exists is the inside region bounded by the constraints.

Here is a plot of the objective function

objective=Plot[{-35/30 x+1000,-35/30 x+10000,
	-35/30 x+17900			},{x,0,30000},PlotStyle->RGBColor[1,0,0],AspectRatio->1]

The objective function is increasing upward.  Looking at the constraints and the monotonically increasing objective function we note that the optimimum will occur on a boundary and will for all reasonably cases be a vertex.

Show[constraints,objective,AspectRatio->1]

We see where the optimimum is!, x1=15, x2=25. Of course we can't solve every problem graphically so we need a procedure.

Express the problem in an augmented matrix form with the vector being
{M, x1, x2, x3, x4, x5} and the augmented column {70, 40, 90, 0} the last entry (which will be M) being the current value of the objective function.

Start with the real variables set to zero (x1=x2=0) which is a known trivial solution.  The zero variables will be called nonbasic variables these will change as the procedure progresses.  To eliminate the need to deal with inequalities and because we note that if a constraint is limiting the solution, we will be at the limiting value (i.e., it will be an equality) we will introduce slack variables to make all of the inequalities into equations.  For this problem we need three slack variables, x3, x4 and x5.  If we start at x1=x2=0, then the initial values of the slack variables are x3=70, x4=40 and x5=90.  This solution is found from inspection.  At any step the non zero variables will be called basic variables.

Here is the initial augmented matrix.  There is some convenience to entering it line by line.  (We can use a subscripted variable.)

z[1][1]= {0,.1,.07,1,0,0,1500}

z[1][2]={0,.09,.1,0,1,0,1500}

z[1][3]={0,.06,.18,0,0,1,1900}

z[1][4]={1,-35,-30,0,0,0,0}

a1=Table[z[1][i],{i,1,4}]

We see, from the objective function, M = 35 x1 + 30 x2, that increasing x1 is the way to get the biggest impact.  The idea of the simplex method is to move only if the objective function will increase and then to always take the biggest increase.  If instead of looking that the original objective function we choose to look at the augmented matrix, we should choose to remove the largest negative in the last row and continue doing this until there are no more negatives.  

So, we will move to the next vertex.  To do this, which of x3, x4, x5 should be set to 0?  

ANS: The one that represents the strongest constraint on the value of x1!  

This is found by the getting the smallest (non positive) quotient between the x1 column (2rd column) and the b vector (7th column)

We make a label for the third column

kk=2;

ztest=Table[z[1][i][[7]]/z[1][i][[kk]],{i,1,3}]

We see that we need to use the first constraint.  This means that x3 will be as small as possible = 0 (i.e., no slack!) and the third constraint becomes active.  Note that each vertex represents a point where two of the variables are 0 and the other three are nonzero. (Why is this?)  Thus any solution that we examine must have this form.

The new value of x1 = 1500/.01 = 15000.  
 
We need to identify this active constraint with an index value.

j=Position[ztest, Min[ztest]][[1,1]]

Now we want a new set of equations that show the change of variables,
  i.e.,  showing that x2 has been increased in favor of x5 becoming 0.
  
 Keeping x1 =0, but changing x2 requires that changes be made in x3, x4 and x5.
 x3 = 70 - x2 = 40,
 x4 = 40 - x2 = 10,
 x5 = 90 - 3 x2  =0 
 
So that x1 = 0, x2 = 30, x3 = 40, x4 = 10, x5 =0 is a new basic feasible solution.

Before we continue on, we need to process the equations above to reflect the changes.

Since we know that each solution has 2 zero variables and 3 nonzero variables, we require that each basic (nonzero) variable has only 1 non zero entry in its column in the augmented matrix and this is not in the bottom row.  This way, each nonzero variable will be the only contribution to the value in its equation.  The current values of the nonzero variables will be immediately evident.   

Making these changes is the same as a Gaussian elimination using the x1 from the first row.

Do[If[i!=j, z[2][i]=z[1][i]-z[1][j]  z[1][i][[kk]]/ z[1][j][[kk]],z[2][i]=z[1][i]/z[1][j][[kk]]],{i,1,4}]

a2=Table[z[2][i],{i,1,4}]

We see immediately that x1=15000, x4=150 and x5=1000 (since x2=x3=0).  We are thus assured that we at a contrainst intersection as we wish.

Note that the objective function is 525000 as required and it is expressed by the non basic variables (which are 0)

Now we need to increase x2 as far as possible.  Check to see which constraint becomes active?

kk=3;

ztest=Table[z[2][i][[7]]/z[2][i][[kk]],{i,1,3}]

We need to identify this active constraint with an index value.

j=Position[ztest, Min[ztest]][[1,1]]

The second constraint will become active first so that x4 will go to 0.

Do[If[i!=j, z[3][i]=z[2][i]-z[2][j]  z[2][i][[kk]]/ z[2][j][[kk]],z[3][i]=z[2][i]/z[2][j][[kk]]],{i,1,4}]

a3=Table[z[3][i],{i,1,4}]

So the objective function is 547297!!  There are no remaining "negatives" in the last row so that no improvement can be attained.
  

Here is a check of the objective function!

12162.2*35+4054.05*30

Of course, Mathematica did not really need our help

ConstrainedMax[35 x1 + 30 x2,{.1 x1 + .07x2 <=1500,
		.09x1 + .1x2 <=1500, .06x1 + .18 x2<=1900},{x1,x2}]

z[1][1]= {0,.1,.07,1,0,0,1500}

z[1][2]={0,.09,.1,0,1,0,1500}

z[1][3]={0,.06,.18,0,0,1,1900}

z[1][4]={1,-35,-30,0,0,0,0}

Maximize M = 35 x1 + 30 x2 

subject to:
.1x1 + .075 x2 <= 1500
.08 x1 + .10 x2 <=1500
.05 x1 + .15 x2 <=1900

x<= 0

Here are the constraints

constraints=Plot[{1500/.07- .10/.07 x,1500/.10-.09x/.10,1900/.18-.06 x/.18},{x,0,30000},PlotStyle->{Dashing[{.02,.02}],Dashing[{.03,.01}],Dashing[{1,0}]},AxesLabel->{"x1","x2"},AspectRatio->1]

We see that the region in which a valid solution exists is the inside region bounded by the constraints.

Here is a plot of the objective function

objective=Plot[{-35/30 x+1000,-35/30 x+10000,
	-35/30 x+17900			},{x,0,30000},PlotStyle->RGBColor[1,0,0],AspectRatio->1]

The objective function is increasing upward.  Looking at the constraints and the monotonically increasing objective function we note that the optimimum will occur on a boundary and will for all reasonably cases be a vertex.

Show[constraints,objective,AspectRatio->1]

Here is the initial augmented matrix.

a={{0,.1,.07,1,0,0,1500},
{0,.09,.1,0,1,0,1500},{0,.06,.18,0,0,1,1900},
{1,-35,-30,0,0,0,0}}

We see, from the objective function, M = 40 x1 + 60 x2,  that increasing x2 is the way to get the biggest impact.  The idea of the simplex method is to only move if the objective function will increase and then to always take the biggest increase.

So, we will move to the next vertex.  To do this, which of x3, x4, x5 should be set to 0?  

ANS: The one that represents the strongest constraint on the value of x2!  

This is found by the getting the smallest (non positive) quotient between the x2 column (3rd column) and the b vector (7th column)

ztest=Table[a[[i,7]]/a[[i,2]],{i,1,3}]

We set j=3, (for the contraint from the third row)

j=1;

We see that we need to use the last constraint.  This means that x5 will be as small as possible =0 (no slack!) and the third constraint becomes active.  Note that each vertex represents a point where two of the variables are 0 and the other three are nonzero. (Why is this?)  Thus any solution that we examine must have this form.

 The new value of x2 = 90/3 = 30.  
 
 Now we want a new set of equations that show the change of variables,
  i.e.,  showing that x2 has been increased in favor of x5 becoming 0.
  
 Keeping x1 =0, but changing x2 requires that changes be made in x3, x4 and x5.
 x3 = 70 - x2 = 40,
 x4 = 40 - x2 = 10,
 x5 = 90 - 3 x2  =0 
 
So that x1 = 0, x2 = 30, x3 = 40, x4 = 10, x5 =0 is a new basic feasible solution.

Before we continue on, we need to process the equations above to reflect the changes.

Since we know that each solution has 2 zero variables and 3 nonzero variables, we require that each basic (nonzero) variable has only 1 non zero entry in its column in the augmented matrix and this is not in the bottom row.  This way, eaxh nonzero variable will be the contribution to the value on its equation.  The current values of the nonzero variables will be immediately evident.   

Making these changes is the same as a Gaussian elimination using the x2 from the third row.

Here we do the elimination

a

j

Do[If[i!=j, zx[i]=a[[i]]-a[[j]] a[[i,j]]/a[[j,2]],zx[i]=a[[i]]/a[[j,2]]];
	Print[zx[i]],{i,1,4}]

a2=Table[zx[i],{i,1,4}]

We see immediately that x2=30, x3=40 and x4=10 (since x1=x5=0).  We are thus assured that we at a contrainst intersection as we wish.

Note that the objective function is 1800 as required and it is expressed by the non basic variables (which are 0)

Now we need to increase x1 as far as possible.  Check to see which constraint becomes active?

ztest=Table[a2[[i,7]]/a2[[i,2]],{i,1,3}]

Min[ztest]

The second constraint will become active first so that x4 will go to 0.

Now choose j=2

j=2;

Do[If[i!=j, z[s][i]=a2[[i]]-a2[[j]] a2[[i,j]]/a2[[j,j]],z[s][i]=a2[[i]]/a2[[j,j]]];
	Print[zx[i]],{i,1,4}]

a3=Table[z[5][i],{i,1,4}]

a3={a2[[1]]-a2[[2]] a2[[1,2]]/a2[[2,2]],a2[[2]]/a2[[2,2]],a2[[3]]-
a2[[2]] a2[[3,2]]/a2[[2,2]],a2[[4]]- a2[[2]] a2[[4,2]]/a2[[2,2]] }

There are no remaining "negatives" in the last row so that no improvement can be attained.

 
So we are done!!  For x4=x5=0, the solution for {x1,x2,x3} = {15,25,15}.  We see that x2 
was reduced a little to get a big gain from x1.  Now x4 and x5=0 which correspond to the constrints,
x1 + x2 <=40
x1 + 3 x2 <=90 
being active.

Duality

We had in the first example:

maximize M = cT x
s.t.
Ax <= b,  x>0

which we solved as x = [15, 25, 15]  It turns out that this problem has a dual

minimize m = bT y
s.t.
AT y >= c, y>= 0

y = [y1,y2,y3]

minimize m = 70 y1 + 40 y2 + 90 y3

subject to:
2 y1 + y2 + y3 >= 40
y1 + y2 + 3 y3 >= 60

y<= 0

Express the problem in an augmented matrix form with the vector being
{m,x1,x2,x3,x4,x5} and the augmented column {70,40,90,0} the last entry being the current value of the objective function.

a={{0,2,1,1,-1,0,40},
   {0,1,1,3,0,-1,60},
   {1,70,40,90,0,0,0}}

MatrixForm[a]

We need to find an initial allowed solution.

2 y1 + y2 + y3 == 40, how about y1=12, y2=0, y3=16.  Arrange the tableau

a2={a[[1]],a[[2]]-3 a[[1]],
a[[3]] -(90)a[[1]] }

MatrixForm[a2]

The second constraint will become active first so that x4 will go to 0.

a3={a2[[1]]+2 a2[[2]]/5,a2[[2]]/(-5),
a2[[3]] -(110)a2[[2]]/5 }

MatrixForm[a3]

So the objective function is 2280!!  There is a remaining negative

a4={a3[[1]]-  a3[[2]]/2,a3[[2]]/(2/5),
a3[[3]] +(6)a3[[2]]/(2/5) }

MatrixForm[a4]

We see that y1=0,y2=30,y3=10, is a better solution.

Now if the additional constraint is added, we would need a new original solution of x1=x2=x3=0 and then redoing the problem from scratch.  However with the dual, we can start from the final solution above, y2=30, y3=10.

15*40+25*60

a={{0,2,1,1,3/2, -1,0,40},
   {0,1,1,3,2,   0,-1,60},
   {1,70,40,90,70,0,0,0}}

MatrixForm[a]

We need to find an initial allowed solution.

2 y1 + y2 + y3 == 40, how about y1=12, y2=0, y3=16.  Arrange the tableau

a2={a[[1]]/2,a[[2]]-  a[[1]]/2,
a[[3]] -(70)a[[1]]/2 }

MatrixForm[a2]

a3={a2[[1]]-  a2[[2]]/5,a2[[2]]/(5/2),
a2[[3]] -(55)a2[[2]]/(5/2) }

MatrixForm[a3]

Now we have unfilled expectations,

a4={a3[[1]]/(2/5),a3[[2]]-a3[[1]]/2  ,
a3[[3]] +(6) a3[[1]]/(2/5) }

MatrixForm[a4]

So the objective function is 2100 which is the orginal solution  There is a remaining negative.  This means that a constraint is not satisfied or that the objective function is not optimized.  The first one will be active.

a5={a4[[1]]/(5/4),a4[[2]]-a4[[1]]/5,
a4[[3]] +(5/2)a4[[1]]/(5/4) }

MatrixForm[a5]

We see that y1=0,y2=4,y3=24, is a better solution.