I have updated this question to specify the problem more clearly.
I've got a dataset of the form:
#sample dataset
dt <- read.table(text = 'Plant Type Ownership Capacity VC TargetPLF
1 T Base Pub 1000 1.50 0.80
2 S Base Pub 500 1.62 0.80
3 Y Base Pub 500 2.43 0.75
4 K Base Pub 500 1.93 0.80
5 J Base Pub 500 3.40 0.80', header =T)
demand <- c(1500,2000,2200,3000,2800)
My objective function looks like this:
eval_f <- function(x) #objective function
{
return((1.50*1000*x[1] + 1.62*500*x[2] + 2.43*500*x[3] + 1.93*500*x[4] + 3.40*500*x[5]) / (1000*x[1] + 500*x[2] + 500*x[3] + 500*x[4] + 500*x[5]))
}
Other arguments:
x0 <- rep(0.45,5) #initialisation values
local_opts <- list( "algorithm" = "NLOPT_LD_MMA", "xtol_rel" = 1.0e-15 ) #options for optimisation
opts = list( "algorithm"= "NLOPT_GN_ISRES","xtol_rel"= 1.0e-15,"maxeval"= 160000,"local_opts" = local_opts,"print_level" = 0 ) #optimisation algorithm
# constraint function
eval_g0_eq <- function(x) {
return(1000*x[1] + 500*x[2] + 500*x[3] + 500*x[4] + 500*x[5] - 1500)
}
lb = c(rep(0.40,5)) #lower bounds of controls
ub = dt$Target #upper bounds of controls
nloptr(x0,eval_f,lb=lb,ub=ub,eval_g_eq = eval_g0_eq,opts=opts) #optimisation
How can I conduct this optimisation 5 times, with a different constraint function each time?
So in the second iteration, the constraint function should take the second element from the vector 'demand' and become this:
# constraint function for second iteration
eval_g0_eq <- function(x) {
return(1000*x[1] + 500*x[2] + 500*x[3] + 500*x[4] + 500*x[5] - 2000)
}
nloptr(x0,eval_f,lb=lb,ub=ub,eval_g_eq = eval_g0_eq,opts=opts) #optimisation
Further, is there a way to construct eval_f and eval_g0_eq using the values in dt, and not manually as done now?
I am not aware of this nlopt
function or library where it comes from, but I think you do not need to create objective function
or constraint function
repeatedly. You may use the R's vectorised calculation potential to do this.
dt <- read.table(text = 'Plant Type Ownership Capacity VC TargetPLF
1 T Base Pub 1000 1.50 0.80
2 S Base Pub 500 1.62 0.80
3 Y Base Pub 500 2.43 0.75
4 K Base Pub 500 1.93 0.80
5 J Base Pub 500 3.40 0.80', header =T)
demand <- c(1500,2000,2200,3000,2800)
library(nloptr)
dt
#> Plant Type Ownership Capacity VC TargetPLF
#> 1 T Base Pub 1000 1.50 0.80
#> 2 S Base Pub 500 1.62 0.80
#> 3 Y Base Pub 500 2.43 0.75
#> 4 K Base Pub 500 1.93 0.80
#> 5 J Base Pub 500 3.40 0.80
The below code will do iterations equal to number of rows in dt
lapply(seq(nrow(dt)), function(.x) nloptr(x0 = rep(0.45,5),
eval_f = function(x) sum(dt$VC * dt$Capacity * x)/sum(dt$Capacity * x),
lb=c(rep(0.40,5)),
ub = dt$Target,
eval_g_eq = function(x) sum(dt$Capacity * x) - demand[.x] ,
opts = list( "algorithm"= "NLOPT_GN_ISRES",
"xtol_rel"= 1.0e-15,
"maxeval"= 160000,
"local_opts" = list("algorithm" = "NLOPT_LD_MMA",
"xtol_rel" = 1.0e-15 ),
"print_level" = 0 )
)
)
Its output will be
#> [[1]]
#>
#> Call:
#> nloptr(x0 = rep(0.45, 5), eval_f = function(x) sum(dt$VC * dt$Capacity *
#> x)/sum(dt$Capacity * x), lb = c(rep(0.4, 5)), ub = dt$Target,
#> eval_g_eq = function(x) sum(dt$Capacity * x) - demand[.x],
#> opts = list(algorithm = "NLOPT_GN_ISRES", xtol_rel = 1e-15,
#> maxeval = 160000, local_opts = list(algorithm = "NLOPT_LD_MMA",
#> xtol_rel = 1e-15), print_level = 0))
#>
#>
#> Minimization using NLopt version 2.4.2
#>
#> NLopt solver status: 5 ( NLOPT_MAXEVAL_REACHED: Optimization stopped because
#> maxeval (above) was reached. )
#>
#> Number of Iterations....: 160000
#> Termination conditions: xtol_rel: 1e-15 maxeval: 160000
#> Number of inequality constraints: 0
#> Number of equality constraints: 1
#> Current value of objective function: 1.9506666666645
#> Current value of controls: 0.7 0.4 0.4 0.4 0.4
#>
#>
#>
#> [[2]]
#>
#> Call:
#> nloptr(x0 = rep(0.45, 5), eval_f = function(x) sum(dt$VC * dt$Capacity *
#> x)/sum(dt$Capacity * x), lb = c(rep(0.4, 5)), ub = dt$Target,
#> eval_g_eq = function(x) sum(dt$Capacity * x) - demand[.x],
#> opts = list(algorithm = "NLOPT_GN_ISRES", xtol_rel = 1e-15,
#> maxeval = 160000, local_opts = list(algorithm = "NLOPT_LD_MMA",
#> xtol_rel = 1e-15), print_level = 0))
#>
#>
#> Minimization using NLopt version 2.4.2
#>
#> NLopt solver status: 5 ( NLOPT_MAXEVAL_REACHED: Optimization stopped because
#> maxeval (above) was reached. )
#>
#> Number of Iterations....: 160000
#> Termination conditions: xtol_rel: 1e-15 maxeval: 160000
#> Number of inequality constraints: 0
#> Number of equality constraints: 1
#> Current value of objective function: 1.893
#> Current value of controls: 0.8 0.8 0.4 0.8 0.4
#>
#>
#>
#> [[3]]
#>
#> Call:
#> nloptr(x0 = rep(0.45, 5), eval_f = function(x) sum(dt$VC * dt$Capacity *
#> x)/sum(dt$Capacity * x), lb = c(rep(0.4, 5)), ub = dt$Target,
#> eval_g_eq = function(x) sum(dt$Capacity * x) - demand[.x],
#> opts = list(algorithm = "NLOPT_GN_ISRES", xtol_rel = 1e-15,
#> maxeval = 160000, local_opts = list(algorithm = "NLOPT_LD_MMA",
#> xtol_rel = 1e-15), print_level = 0))
#>
#>
#> Minimization using NLopt version 2.4.2
#>
#> NLopt solver status: 5 ( NLOPT_MAXEVAL_REACHED: Optimization stopped because
#> maxeval (above) was reached. )
#>
#> Number of Iterations....: 160000
#> Termination conditions: xtol_rel: 1e-15 maxeval: 160000
#> Number of inequality constraints: 0
#> Number of equality constraints: 1
#> Current value of objective function: 1.92181376479112
#> Current value of controls: 0.7883007 0.7815353 0.583928 0.7845948 0.4043761
#>
#>
#>
#> [[4]]
#>
#> Call:
#> nloptr(x0 = rep(0.45, 5), eval_f = function(x) sum(dt$VC * dt$Capacity *
#> x)/sum(dt$Capacity * x), lb = c(rep(0.4, 5)), ub = dt$Target,
#> eval_g_eq = function(x) sum(dt$Capacity * x) - demand[.x],
#> opts = list(algorithm = "NLOPT_GN_ISRES", xtol_rel = 1e-15,
#> maxeval = 160000, local_opts = list(algorithm = "NLOPT_LD_MMA",
#> xtol_rel = 1e-15), print_level = 0))
#>
#>
#> Minimization using NLopt version 2.4.2
#>
#> NLopt solver status: 5 ( NLOPT_MAXEVAL_REACHED: Optimization stopped because
#> maxeval (above) was reached. )
#>
#> Number of Iterations....: 160000
#> Termination conditions: xtol_rel: 1e-15 maxeval: 160000
#> Number of inequality constraints: 0
#> Number of equality constraints: 1
#> Current value of objective function: 1.93672658072513
#> Current value of controls: 0.7870349 0.7534322 0.6660887 0.770566 0.4118776
#>
#>
#>
#> [[5]]
#>
#> Call:
#> nloptr(x0 = rep(0.45, 5), eval_f = function(x) sum(dt$VC * dt$Capacity *
#> x)/sum(dt$Capacity * x), lb = c(rep(0.4, 5)), ub = dt$Target,
#> eval_g_eq = function(x) sum(dt$Capacity * x) - demand[.x],
#> opts = list(algorithm = "NLOPT_GN_ISRES", xtol_rel = 1e-15,
#> maxeval = 160000, local_opts = list(algorithm = "NLOPT_LD_MMA",
#> xtol_rel = 1e-15), print_level = 0))
#>
#>
#> Minimization using NLopt version 2.4.2
#>
#> NLopt solver status: 5 ( NLOPT_MAXEVAL_REACHED: Optimization stopped because
#> maxeval (above) was reached. )
#>
#> Number of Iterations....: 160000
#> Termination conditions: xtol_rel: 1e-15 maxeval: 160000
#> Number of inequality constraints: 0
#> Number of equality constraints: 1
#> Current value of objective function: 1.93152444911471
#> Current value of controls: 0.795248 0.7687759 0.6712354 0.7692782 0.4034175