Wicelo
Wicelo

Reputation: 2426

constrained optimization of a complicated function

I need to optimize the following to find the maximum value for r1:

ad = 0.95*M_D + 0.28*G_D + 0.43*S_D + 2.25*Q_D
as = 0.017*M_A + 0.0064*G_A + 0.0076*S_A + 0.034*Q_A
ccb = 0.0093*M_CC+ 0.0028*G_CC + 0.0042*S_CC + 0.0186*Q_CC
ccd = 0.0223*M_CD + 0.0056*G_CD + 0.0078*S_CD + 0.0446*Q_CD
apb = 1.28*M_P + 2.56*Q_P 
r1=(1+ccb*(1+ccd))*ad*as*100/(130-apb)

subject to the following constraints:

0 <= M_CD <= M_CC <= M_A <= M_D <= M_P <= 9
0 <= G_CD <= G_CC <= G_A <= G_D <= 9
0 <= S_CD <= S_CC <= S_A <= S_D <= 9
0 <= Q_CD <= Q_CC <= Q_A <= Q_D <= Q_P <= 3

The approach I've tried before doesn't work very well and I'm hoping to find a better solution.

Upvotes: 1

Views: 260

Answers (2)

Wicelo
Wicelo

Reputation: 2426

I finally solved my problem not with vectorization but with C. My program containing 14nested loops executes 100 to 1000 times faster with C than with R ! Which is sad because I didn't learn anything new from that and it prooves that R can be pretty useless on some problems but what can we do.

Upvotes: 0

user1609452
user1609452

Reputation: 4444

Once the problem is stated correctly you maybe able to firstly map the parameters to lower and upper bounds of [0,1]. You can then implement the inequalities inside your function and optimise using an algorithm which accepts basic lower and upper bound constraints. nlminb could be used but the vignette suggests the algorithm used may not be the best.

UPDATE:

With OP revised function

dumFun <- function(p){
    p[1] -> M_CD; p[2] -> M_CC; p[3] -> M_A; p[4] -> M_D; p[5] -> M_P;
    M_P <- 9*M_P; M_D <- M_P*M_D; M_A <- M_D*M_A; M_CC <- M_A*M_CC; M_CD <- M_CC*M_CD; 
    p[6] -> G_CD; p[7] -> G_CC; p[8] -> G_A; p[9] -> G_D;
    G_D <- 9*G_D; G_A <- G_D*G_A; G_CC <- G_A*G_CC; G_CD <- G_CC*G_CD; 
    p[10] -> S_CD; p[11] -> S_CC; p[12] -> S_A; p[13] -> S_D;
    S_D <- 9*S_D; S_A <- S_D*S_A; S_CC <- S_A*S_CC; S_CD <- S_CC*S_CD; 
    p[14] -> Q_CD; p[15] -> Q_CC; p[16] -> Q_A; p[17] -> Q_D; p[18] -> Q_P;
    Q_P <- 3*Q_P; Q_D <- Q_P*Q_D; Q_A <- Q_D*Q_A; Q_CC <- Q_A*Q_CC; Q_CD <- Q_CC*Q_CD; 

    ad = 0.95*M_D + 0.28*G_D + 0.43*S_D + 2.25*Q_D
    as = 0.017*M_A + 0.0064*G_A + 0.0076*S_A + 0.034*Q_A
    ccb = 0.0093*M_CC+ 0.0028*G_CC + 0.0042*S_CC + 0.0186*Q_CC
    ccd = 0.0223*M_CD + 0.0056*G_CD + 0.0078*S_CD + 0.0446*Q_CD
    apb = 1.28*M_P + 2.56*Q_P 
    r1=(1+ccb*(1+ccd))*ad*as*100/(130-apb)
    -r1
}
require(minqa)
p <- rep(.1, 18)
bobyqa(p, dumFun, lower = rep(0, length(p)), upper = rep(1, length(p)))
parameter estimates: 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1, 1 
objective: -9.65605526502482 
number of function evaluations: 97 
> 

Upvotes: 2

Related Questions