CODEWITHSUNDEEP
rmatrixoptimization
Andre Suriane
Andre Suriane

Reputation: 13

R optim() constraint optimization does not find the first best

my problem is summarized in finding a vector X with the best solution to the problem: L is the profits, R is the restrictions, P is a constraint parameters matrix,

or max(t(L)%*%X)

or P%*%X <= R.

I find a solution for X, but not the best, which would be fb = c(.217,0,0,23,2865,0,13,427). How do I find the best solution?

code:

X<-matrix(rep(1,6),6,1)
P<-matrix(c(
    1, 1, 1, 2, 0, 0,
    0, 1, 1, 2, 1, 1,
    99.4, 37.75, 19.75, 54.40, 74.75, 53,
    2.400, 1.540, 0, 0, 0, 0,
    2.400, 1.960, 0, 0, 0, 0,
    1.800, 3.300, 5.330, 0, 0, 0,
    0, 0, 2.070, 0, 8.700, 0,
    0, 0, .436, 0, 19.100, 12.363,
    0, 3.000, .364, 0, 9.100, 26.737  ), 
9,6,1)

L <- matrix(c(83.4, 72.35, 27.3, 72.05, 217.25, 455), 6,1)
R <- matrix(c(60,60,2000,351,448,479,338,424,359),9,1)

farm<- function(par, P,R,  L){
   trues<- P%*%par<=R
 if (min(trues)==1 && min(par)>=0) {
    return(-t(L)%*%par)
     }
 else{
     return(0)
 }
}

mtds = c("Nelder-Mead", "BFGS", "CG", "L-BFGS-B", "SANN","Brent")
out <- optim(par = X,  # initial guess
             fn  = farm,
             P   = P,
             R   = R ,
             L   = L,
             method  = mtds[5]) 


# my result
t(L)%*%out$par
#A matrix: 1 × 1 of type dbl
#7419.596

# the first best
fb<- matrix(c(.217,0,0,23.2865,0,13.427),6,1)
t(L)%*%fb
#A matrix: 1 × 1 of type dbl
#7805.175

Upvotes: 1

Views: 208

Answers (2)

Enrico Schumann
Enrico Schumann

Reputation: 1493

That looks very much like a model that could be solved by a linear programme.

library("Rglpk")
Rglpk_solve_LP(obj = L,
               mat = P,
               dir = rep("<=", 9),
               rhs = R,
               max = TRUE)

Upvotes: 1

ThomasIsCoding
ThomasIsCoding

Reputation: 101024

I think you can try fmincon from package pracma

library(pracma)
objfun <- function(x) -t(L)%*%x
res <- fmincon(x0 = X,fn = objfun,A = P,b = R,lb = rep(0,length(X)))

and you will see that

> res$par
[1]  4.201711e-16 -1.239088e-15  1.863081e-17  2.310286e+01
[5]  5.566620e-01  1.323762e+01

> -res$value
         [,1]
[1,] 7808.615

Upvotes: 1

Related Questions