Lab3: Linear Programming Modeling

Case Study 1

Suppose a farmer has 75 acres on which to plant two crops: wheat and barley. To produce these crops, it costs the farmer (for seed, fertilizer, etc.) $120 per acre for the wheat and $210 per acre for the barley. The farmer has $15000 available for expenses. But after the harvest, the farmer must store the crops while awaiting avourable market conditions. The farmer has storage space for 4000 bushels. Each acre yields an average of 110 bushels of wheat or 30 bushels of barley. If the net profit per bushel of wheat (after all expenses have been subtracted) is $1.30 and for barley is $2.00, how should the farmer plant the 75 acres to maximize profit?

\(PROBLEM\) \(DEFINITION\)

We need to define the problem first:

• x as wheat
• y as barley

Then, we define the constrains:

• x + y \(\leq\) 75
• 120 x + 210 y \(\leq\) 15000
• 110 x + 30 y \(\leq\) 4000
• x \(\geq\) 0, y \(\geq\) 0

To get the maximize profit, we define the z

• \(z\) = \((110)(1.3)x\) + \((30)(2)y\)
• \(z\) = \(143x\) + \(60y\)

which means that \[\hat C = \begin{pmatrix} 143 \\ 60 \end{pmatrix}\].

Define the coefficient matrix

\[A = \begin{pmatrix} 1 & 1 & \\ 120 & 210 &\\ 110 & 30 & \end{pmatrix}\], and

\[B = \begin{pmatrix} 75 \\ 15000 \\ 4000 \end{pmatrix}\].

\(SOLUTION\) \(WITH\) \(LPSOLVE\)

Install the packages of lpsolve first and then use lpsolve library. lpsolve: is callable from R via an extension or module

library(lpSolve)

After that, we set the coefficients of the decision variables -> c

C <- c(143,60)
C

## [1] 143  60

Then, we create constraint matrix B

A <- matrix(c (1, 1 ,                            120, 210,  
                110, 30 ), nrow=3, byrow=TRUE)
A

##      [,1] [,2]
## [1,]    1    1
## [2,]  120  210
## [3,]  110   30

Right hand side of the constraints

B <- c(75,15000,4000)
B

## [1]    75 15000  4000

Direction of the constraints

constraints_direction <-c("<=", "<=","<=")
constraints_direction

## [1] "<=" "<=" "<="

\(THE\) \(OPTIMUM\) \(RESULT\)

we need to find the optimal solution

optimum <- lp(direction ="max",
              objective.in= C,
              const.mat = A,
              const.dir = constraints_direction,
              const.rhs = B,
              all.int = T)
str(optimum)

## List of 28
##  $ direction       : int 1
##  $ x.count         : int 2
##  $ objective       : num [1:2] 143 60
##  $ const.count     : int 3
##  $ constraints     : num [1:4, 1:3] 1 1 1 75 120 210 1 15000 110 30 ...
##   ..- attr(*, "dimnames")=List of 2
##   .. ..$ : chr [1:4] "" "" "const.dir.num" "const.rhs"
##   .. ..$ : NULL
##  $ int.count       : int 2
##  $ int.vec         : int [1:2] 1 2
##  $ bin.count       : int 0
##  $ binary.vec      : int 0
##  $ num.bin.solns   : int 1
##  $ objval          : num 6266
##  $ solution        : num [1:2] 22 52
##  $ presolve        : int 0
##  $ compute.sens    : int 0
##  $ sens.coef.from  : num 0
##  $ sens.coef.to    : num 0
##  $ duals           : num 0
##  $ duals.from      : num 0
##  $ duals.to        : num 0
##  $ scale           : int 196
##  $ use.dense       : int 0
##  $ dense.col       : int 0
##  $ dense.val       : num 0
##  $ dense.const.nrow: int 0
##  $ dense.ctr       : num 0
##  $ use.rw          : int 0
##  $ tmp             : chr "Nobody will ever look at this"
##  $ status          : int 0
##  - attr(*, "class")= chr "lp"

When we print the optimum$status and the output is 0 it is mean sucess, while 2 have a meaning of no feasible solution.

print(optimum$status)

## [1] 0

So, we got 0, which mean it is success

Now, we want to display the optimum values for x and y

best_solution <- optimum$solution
names(best_solution) <- c("x","y")
print(best_solution)

##  x  y 
## 22 52

So, we got the optimum values for x and y are 22 and 52 values

Now, we want to check value of objective function at optimal point

print(paste("Total profit: ", optimum$objval, sep=""))

## [1] "Total profit: 6266"

The total profit we got is 6266 dollars

Lastly, we need to disconnect from the model and the optimum solution

rm(optimum, constraints_direction,best_solution)

In conclusion, the farmer will get the maximum profit ($6266) by planting 22 acres of wheat and 52 acres of barley.

\(Solution\) \(with\) \(lpSolveAPI\)

The lpSoleveAPI R package is a second implementation of an interface of lpsolve to R.

First, we need to install the packages of lpSolveAPI

library(lpSolveAPI)

## Warning: package 'lpSolveAPI' was built under R version 3.6.3

Second, we need to set the constrain and the decision variables. The, we set the type of the problem we are trying to solve

lprec <- make.lp(nrow=3, ncol=2)
lp.control(lprec,sense="max")

## $anti.degen
## [1] "fixedvars" "stalling" 
## 
## $basis.crash
## [1] "none"
## 
## $bb.depthlimit
## [1] -50
## 
## $bb.floorfirst
## [1] "automatic"
## 
## $bb.rule
## [1] "pseudononint" "greedy"       "dynamic"      "rcostfixing" 
## 
## $break.at.first
## [1] FALSE
## 
## $break.at.value
## [1] 1e+30
## 
## $epsilon
##       epsb       epsd      epsel     epsint epsperturb   epspivot 
##      1e-10      1e-09      1e-12      1e-07      1e-05      2e-07 
## 
## $improve
## [1] "dualfeas" "thetagap"
## 
## $infinite
## [1] 1e+30
## 
## $maxpivot
## [1] 250
## 
## $mip.gap
## absolute relative 
##    1e-11    1e-11 
## 
## $negrange
## [1] -1e+06
## 
## $obj.in.basis
## [1] TRUE
## 
## $pivoting
## [1] "devex"    "adaptive"
## 
## $presolve
## [1] "none"
## 
## $scalelimit
## [1] 5
## 
## $scaling
## [1] "geometric"   "equilibrate" "integers"   
## 
## $sense
## [1] "maximize"
## 
## $simplextype
## [1] "dual"   "primal"
## 
## $timeout
## [1] 0
## 
## $verbose
## [1] "neutral"

Third, we need to set the type of decision variables

set.type(lprec,1:2, type= c("integer"))

Fourth, we set the objective function coefficients vector C

set.objfn(lprec,C)

Fifth, we add the constraints

add.constraint(lprec, A[1, ], "<=", B[1])
add.constraint(lprec, A[2, ], "<=", B[2])
add.constraint(lprec, A[3, ], "<=", B[3])

sixth, we display the lpsolve matrix

lprec

## Model name: 
##            C1   C2             
## Maximize  143   60             
## R1          0    0  free      0
## R2          0    0  free      0
## R3          0    0  free      0
## R4          1    1    <=     75
## R5        120  210    <=  15000
## R6        110   30    <=   4000
## Kind      Std  Std             
## Type      Int  Int             
## Upper     Inf  Inf             
## Lower       0    0

The, we solve the problem using solve function. If the result is 0 then it is success.

solve(lprec)

## [1] 0

After that, we get the decision variables values

get.variables(lprec)

## [1] 22 52

We get the decision variables values are 22 and 52

Thus, we get the value of the objective function

get.objective(lprec)

## [1] 6266

We get the value of the objective function is 6266

For your information, the default boundaries on the decision variable are c(0,0) and c(Inf, Inf)

get.bounds(lprec)

## $lower
## [1] 0 0
## 
## $upper
## [1] Inf Inf

\(In\) \(Conclusion\)

We can use lpSolve and lpSolveAPI to solve optimization.

The farmer will get the maximum profit ($6266) by planting 22 acres of wheat and 52 acres of barley.