Simplex Method for Linear programming

I am trying to perform simplex method of linear programming for the following:
image
The code does produces correct result on the first iteration. But it won't go beyond result4. Any help on how to perform this iteration would be helpful.
I tried applying while loop within the function to check if the first four columns of last row are negative using any() function. But it did not help. Thank you.

row1 <- c(-4,6,5,-4,1,0,0,0,20)
row2 <- c(3,-2,4,1,0,1,0,0,10)
row3 <- c(8,-3,3,2,0,0,1,0,20)
row4 <- c(-4,-1,-4,-5,0,0,0,1,0)
matrix <- rbind(row1,row2,row3,row4)
matrix
func_simplex <- function(matrix){
  
#finding the pivot column
  pvt_col <- which(matrix[4,1:4]==min(matrix[4,1:4]))

#after division by the coefficient 
  Coef_div <- matrix[,9]/matrix[,pvt_col]

#finding the minimum coef division and setting the pivot element 
  pvt_row <- which(Coef_div==min(Coef_div))

  pvt_ele <- matrix[pvt_row,pvt_col]


#Divide the row by Pivot Element value
  matrix[pvt_row, ] <- matrix[pvt_row, ]/pvt_ele

  for (i in 1:nrow(matrix)){
    if (i != pvt_row){
      matrix[i,] <- matrix[i,]-matrix[i,pvt_col]*matrix[pvt_row,]
    }
  }
  
  return(matrix)
}
result <- func_simplex(matrix)
result
result2 <- func_simplex(result1)
result2
result3 <- func_simplex(result2)
result3
result4 <- func_simplex(result3)
result4
result5 <- func_simplex(result4)
result5

gives

     [,1] [,2] [,3]   [,4]   [,5] [,6]   [,7]   [,8] [,9]
row1    0    0    1 -5.100  0.500    0  0.700  0.900   24
row2    0    0    0 20.775 -1.875    1 -3.175 -3.725  -91
row3    0    1    0  6.700 -0.500    0 -0.900 -1.300  -28
row4    1    0    0  4.675 -0.375    0 -0.475 -0.825  -17

p.s. the code you provided to the forum has a typographical error which means it wont work without correction.
that is you do not save the first result as result1 that you then will use towards result2. i..e you make result, and then dont use it going forwards. If this is issue is present in your true script, then it might explain your odd results as you may have edited your functoins and run things at different times, and gotten stuck with results that dont align with your changes.

1 Like

Right that was typograhical error. But i am still not getting the optimal value. The code keeps getting same result for result1 and result3, result2 and result4 and so on. So it is not going for the optimal value. Can you help me on this?

row1 <- c(-4,6,5,-4,1,0,0,0,20)
row2 <- c(3,-2,4,1,0,1,0,0,10)
row3 <- c(8,-3,3,2,0,0,1,0,20)
row4 <- c(-4,-1,-4,-5,0,0,0,1,0)
matrix <- rbind(row1,row2,row3,row4)
matrix
func_simplex <- function(matrix){

#finding the pivot column
  pvt_col <- which(matrix[4,1:8]==min(matrix[4,1:8]))
  pvt_col
#after division by the coefficient 
  Coef_div <- matrix[1:3,9]/matrix[1:3,pvt_col]
  Coef_div
#finding the minimum coef division and setting the pivot element 
  pvt_row <- which(Coef_div==min(Coef_div))
  pvt_row
  pvt_ele <- matrix[pvt_row[1],pvt_col]
  pvt_ele

#Divide the row by Pivot Element value
  matrix[pvt_row[1],] <- matrix[pvt_row[1],]/pvt_ele
  matrix
  for (i in 1:nrow(matrix)){
    if (i != pvt_row[1]){
      matrix[i,] <- matrix[i,]-matrix[i,pvt_col]*matrix[pvt_row[1],]
    }
  }

  return(matrix)

}

result <- func_simplex(matrix)
result
result1 <- func_simplex(result)
result1
result2 <- func_simplex(result1)
result2
result3 <- func_simplex(result2)
result3
result4 <- func_simplex(result3)
result4

Have you tried your func_simplex() function for a problem with a known solution? The lp() function from the lpSolve package reports an error with status 3, which I believe is for an unbounded solution. Changing the first constraint results in a solution.

library(lpSolve)

f.obj <- c(4, 1, 4, 5)
f.con <- matrix (c(-4, 6, 5, -4,
                   3, -2, 4, 1, 
                   8, -3, 3, 2), nrow=3, byrow=TRUE)
f.dir <- c("<=", "<=", "<=")
f.rhs <- c(20, 10, 20)
lp ("max", f.obj, f.con, f.dir, f.rhs)
#> Error: status 3

# changing the x4 coefficient in the first constraint from -4 to 4:
f.con <- matrix (c(-4, 6, 5, 4,
                   3, -2, 4, 1, 
                   8, -3, 3, 2), nrow=3, byrow=TRUE)

lp ("max", f.obj, f.con, f.dir, f.rhs)
#> Success: the objective function is 34
lp ("max", f.obj, f.con, f.dir, f.rhs)$solution
#> [1] 1 0 0 6

Created on 2023-06-20 with reprex v2.0.2

1 Like

Oh you are right! the original ones do not converge to the solution. I did change the coefficient in my code. But my code is still not producing the optimum result.

Classic trap set by R for programs experienced in languages like C and its progeny. The matrix that is being returned is that last matrix the pops out of the for() loop. Here is the pattern to escape from Las Vegas

v <- vector(length = length(object))
 for(i in seq_along(object)) v[i] = sqrt(i)

which will get each iteration through the loop save to a variable in the .Global environment.

1 Like

After a lot of hassle, I think I am there.

row1 <- c(-4, 6, 5, 4, 1, 0, 0, 20)
row2 <- c(3, -2, 4, 1, 0, 1, 0, 10)
row3 <- c(8, -3, 3, 2, 0, 0, 1, 20)
row4 <- c(-4, -1, -3, -5, 0, 0, 0, 0)
row0 <- c("C", "x1", "x2", "x3", "x4", "x5", "x6", "x7", "Z")
only_matrix <- rbind(row1, row2, row3, row4)
only_matrix
var <- c("x5", "x6", "x7", "Z")
matrix <- data.frame(var, only_matrix)
matrix

simplex_func <- function(matrix) {
  max_iterations <- 100
  counter <- 1
  
  while (any(matrix[nrow(matrix), -1] < 0) && counter <= max_iterations) {
    pvt_col <- which(matrix[nrow(matrix), -1] == min(matrix[nrow(matrix), -1])) + 1
    
    Coef_div <- rep(0, nrow(matrix) - 1)
    for (i in 1:(nrow(matrix) - 1)) {
      if (matrix[i, ncol(matrix)] != 0) {
        Coef_div[i] <- matrix[i, ncol(matrix)] / matrix[i, pvt_col]
      } else {
        if (matrix[i, pvt_col] < 0) {
          Coef_div[i] <- -1
        } else {
          Coef_div[i] <- 0
        }
      }
    }
    
    min_coef <- min(Coef_div[Coef_div >= 0])
    pvt_row <- which(Coef_div == min_coef)[1]
    pvt_ele <- matrix[pvt_row, pvt_col]
    
    matrix[pvt_row, -1] <- matrix[pvt_row, -1] / pvt_ele
    
    for (i in 1:nrow(matrix)) {
      if (i != pvt_row) {
        matrix[i, -1] <- matrix[i, -1] - matrix[i, pvt_col] * matrix[pvt_row, -1]
      }
    }
    
    matrix[pvt_row, 1] <- trimws(paste("x", pvt_col - 1, sep = ""))
    
    counter <- counter + 1
  }
  
  return(matrix)
}

result <- simplex_func(matrix)
result

This topic was automatically closed 42 days after the last reply. New replies are no longer allowed.

If you have a query related to it or one of the replies, start a new topic and refer back with a link.