search

alteryx / AlteryxPrescriptive

R Package for Optimization in Alteryx
GNU General Public License v3.0
3 stars 6 forks source link

Add sensitivity analysis #2

Closed ramnathv closed 8 years ago

ramnathv commented 8 years ago

Here is one way of doing it using the glpkAPI package

library(Rglpk)
library(glpkAPI)
library(AlteryxPrescriptive)

mylp <- initProbGLPK()
readLPGLPK(mylp, getSampleData('cell_tower.lp'))
solFile <- tempfile(fileext = '.txt')
printSolGLPK(mylp, solFile)
file.edit(solFile)

@kuol Gurobi probably has a different format in the way it returns the sensitivity analysis report. I think we should write an S3 class that extracts the sensitivity analysis data for different solvers.

ramnathv commented 8 years ago 
Problem:    
Rows:       10
Columns:    15
Non-zeros:  36
Status:     UNDEFINED
Objective:  obj = 0 (MAXimum)

   No.   Row name   St   Activity     Lower bound   Upper bound    Marginal
------ ------------ -- ------------- ------------- ------------- -------------
     1 R0           B              0             0               
     2 R1           B              0             0               
     3 R2           B              0             0               
     4 R3           B              0             0               
     5 R4           B              0             0               
     6 R5           B              0             0               
     7 R6           B              0             0               
     8 R7           B              0             0               
     9 R8           B              0             0               
    10 R9           B              0                          20 

   No. Column name  St   Activity     Lower bound   Upper bound    Marginal
------ ------------ -- ------------- ------------- ------------- -------------
     1 r0           NL             0             0             1         < eps
     2 r1           NL             0             0             1         < eps
     3 r2           NL             0             0             1         < eps
     4 r3           NL             0             0             1         < eps
     5 r4           NL             0             0             1         < eps
     6 r5           NL             0             0             1         < eps
     7 r6           NL             0             0             1         < eps
     8 r7           NL             0             0             1         < eps
     9 r8           NL             0             0             1         < eps
    10 t0           NL             0             0             1         < eps
    11 t1           NL             0             0             1         < eps
    12 t4           NL             0             0             1         < eps
    13 t2           NL             0             0             1         < eps
    14 t3           NL             0             0             1         < eps
    15 t5           NL             0             0             1         < eps

Karush-Kuhn-Tucker optimality conditions:

KKT.PE: max.abs.err = 0.00e+00 on row 0
        max.rel.err = 0.00e+00 on row 0
        High quality

KKT.PB: max.abs.err = 0.00e+00 on row 0
        max.rel.err = 0.00e+00 on row 0
        High quality

KKT.DE: max.abs.err = 1.20e+03 on column 5
        max.rel.err = 9.99e-01 on column 5
        DUAL SOLUTION IS WRONG

KKT.DB: max.abs.err = 0.00e+00 on row 0
        max.rel.err = 0.00e+00 on row 0
        High quality

End of output
kuol commented 8 years ago

I have problems to get sensitivity reports on a very simple problem. The first part of code using glpkAPI gives me an pretty empty sensitivity analysis report. The second part of the code shows, Rglpk can give the correct answer. 

library(Rglpk)
library(glpkAPI)
library(AlteryxPrescriptive)

mylp <- initProbGLPK()
readLPGLPK(mylp, getSampleData('table_chair.lp'))
#solFile <- tempfile(fileext = '.txt')
printSolGLPK(mylp, "sol.txt")
file.edit("sol.txt")

## Useing Rglpk, it gives correct solution
file <- "table_chair.lp"
x <- Rglpk_read_file(file, type = "CPLEX_LP")
Rglpk_solve_LP(obj = x$objective, mat = x$constraints[[1]], dir = x$constraints[[2]],
               rhs = x$constraints[[3]], max = x$maximum)
kuol commented 8 years ago

Here's the file "table_chair.lp"

Maximize
  70 t + 50 c
Subject To
 c0: 4 t + 3 c <= 240
 c1: 2 t + c <= 100
End
kuol commented 8 years ago

I'm able to manually input the problem above and get the sensitivity report though. But man, I don't ever want to manually input -- it takes 40 lines of code to specify such a simple one!

## Table Chair with glpkAPI ----
prob <- initProbGLPK()
setProbNameGLPK(prob, "table_chair")
setObjDirGLPK(prob, GLP_MAX)
addRowsGLPK(prob, 2)
addColsGLPK(prob, 2)

setRowNameGLPK(prob, 1, "p")
setRowNameGLPK(prob, 2, "q")
setColNameGLPK(prob, 1, "x1")
setColNameGLPK(prob, 2, "x2")

setRowBndGLPK(prob, 1, GLP_UP, 0, 240)
setRowBndGLPK(prob, 2, GLP_UP, 0, 100)

lb <- c(0, 0)
ub <- c(240, 100)
type <- rep(GLP_UP, 2)
setRowsBndsGLPK(prob, 1:2, lb, ub, type)

#Set the type and bounds of the columns.
setColBndGLPK(prob, 1, GLP_LO, 0, 0)
setColBndGLPK(prob, 2, GLP_LO, 0, 0)
#Set the objective function.
setObjCoefGLPK(prob, 1, 70)
setObjCoefGLPK(prob, 2, 50)
#Set the type and bounds of columns and the objective function using a function which
#has the ability to work with vectors.
lb <- c(0, 0)
ub <- lb
type <- rep(GLP_LO, 2)
obj <- c(70, 50)
setColsBndsObjCoefsGLPK(prob, 1:2, lb, ub, obj, type)
#Load the constraint matrix.
ia <- c(1, 1, 2, 2)
ja <- c(1, 2, 1, 2)
ar <- c(4, 3, 2, 1)
loadMatrixGLPK(prob, 4, ia, ja, ar)
#Solve the problem using the simplex algorithm.
solveSimplexGLPK(prob)

printSolGLPK(prob, "table_sol.txt")
file.edit("table_sol.txt")
ramnathv commented 8 years ago 
mylp <- initProbGLPK()
readLPGLPK(mylp, getSampleData('table_chair.lp'))
#solFile <- tempfile(fileext = '.txt')
# you need this line here
solveSimplexGLPK(mylp)
printSolGLPK(mylp, "sol.txt")
file.edit("sol.txt")
kuol commented 8 years ago

Just to document it here, so we are aware of this: sensitivity analysis report doesn't always return Limiting variable. See the last line and the last column of the following. I'll need to modify my code for this :alien:

GLPK 4.47 - SENSITIVITY ANALYSIS REPORT                                                                         Page   1

Problem:    value
Objective:  296.2166065 (MINimum)

   No. Row name     St      Activity         Slack   Lower bound       Activity      Obj coef  Obj value at Limiting
                                          Marginal   Upper bound          range         range   break point variable
------ ------------ -- ------------- ------------- -------------  ------------- ------------- ------------- ------------
     1 yield        NS    2000.00000        .         2000.00000     1995.06864          -Inf     296.28365 bin3
                                           -.01360    2000.00000     2014.03479          +Inf     296.02579 cu

     2 fe           NU      60.00000        .               -Inf       55.89016          -Inf     306.77162 bin4
                                          -2.56823      60.00000       62.69978       2.56823     289.28294 bin3

     3 cu           BS      83.96751      16.03249          -Inf       93.88467       -.30613     270.51157 mn
                                            .          100.00000       79.98213        .21474     314.24798 bin5

     4 mn           NU      40.00000        .               -Inf       34.42336          -Inf     299.25255 bin4
                                           -.54440      40.00000       41.68691        .54440     295.29825 bin3

     5 mg           BS      19.96029      10.03971          -Inf       24.74427      -1.79618     260.36433 bin1
                                            .           30.00000        9.40292        .28757     301.95652 mn

     6 al           NL    1500.00000        .         1500.00000     1485.78425       -.25199     292.63444 cu
                                            .25199          +Inf     1504.92126          +Inf     297.45669 bin3

     7 si1          NL     250.00000        .          250.00000      235.32871       -.48520     289.09812 cu
                                            .48520          +Inf      255.06073          +Inf     298.67206 bin3

     8 si2          BS     250.00000      50.00000          -Inf      255.06073       -.48520     174.91697 si1
                                            .          300.00000      250.00000          +Inf          +Inf
ramnathv commented 8 years ago

@kuol I think we can close this, since you have already incorporated sensitivity analysis into the code.

kuol commented 8 years ago

yep definitely