nick-example-ompr.R

library(dplyr)
library(ROI)
library(ROI.plugin.glpk)
library(ompr)
library(ompr.roi)

rb_mat <- function(r,c,prob) matrix(rbinom(r*c,1,prob),r,c)

my_A <- rb_mat(r = 1000, c = 200) # 1000 heart attacks, 200 potential locations

# number of OHCA incidents
J <- nrow(my_A)
I <- ncol(my_A)
# number of AEDs
N <- 5

MIPModel() %>%
    add_variable(x[j], j = 1:J, type = "binary") %>%
    add_variable(y[i], i = 1:I, type = "binary") %>%
    # add_variable(a[i,j], i = 1:I, j = 1:J, type = "binary") %>%
    add_variable(my_A[i,j], i = 1:I, j = 1:J, type = "binary") %>%
    set_objective(sum_exp(x[j], j = 1:J), "max") %>%
    add_constraint(sum_exp(y[i], i = 1:I), "==", N) %>%
    add_constraint(sum_exp(a[i,j] * y[i], i = 1:I, j = 1:J), ">=", x[j])

# > Error in on_element(push, inplace_update_ast, get_ast_value, element) : 
# > The expression contains a variable, that is not part of the model.

Hey @njtierney, what are the values of J and I?

ompr is currently not optimised for speed and thus comes with a quite high cost of abstraction. I have some ideas on how to make ompr as fast as building a matrix, but at the moment, it is quite slow.

For maxcovr your could use ROI for example and build the coefficient matrix directly.

Hi @dirkschumacher,

I'm using some "real data" and trying to use ompr to solve my simple optimization problem, I get a message warning of "Error in eval(x$expr, x$env, emptyenv()) : object 'i' not found". Can you pls help me what is wrong. Initially it did run but had many warnings. Now for some reason, I'm getting this message.

Hotspots<- Ranking[1:1453, ] # Ranking is a dataset with estimated risk values
CMF<- data.frame(id= 1:3, cmf.name= c("CMF.Passive_FLB","CMF.Passive_FLBG","CMF.FLB_FLBG"), cmf.value= c(0.39, 0.297, 0.647))
cost<- data.frame(id= 1:3, c.name= c("c.Passive_FLB","c.Passive_FLBG","c.FLB_FLBG"), c.value= c(30000, 180000, 150000))

library(ompr)
library(magrittr)

#p means crash risk after the treatment application
p<- function(i, j) {
crash.before.treatment <- Hotspots[i, ] #from /hotspots dataset
CMFj <- CMF[j, ]
output<- crash.before.treatment$WeightedEstimate*CMFj$cmf.value
return(output)
}
#e.g.
p(1, 1)

#For excluding those site-countermeasure pairs not feasible for upgrade
#It is based on the logic that:
#If crossing type is passive, treatment options are: FLB, FLBG
#If crossing type is FLB, treatment option is FLBG
#If crossing type is FLBG, no option for upgrade

y<- function(i, j){
sitei <- Hotspots[i, ] #get from hotspots dataset
CrossingType<- sitei$ProtectionType

  if (j == 1 && CrossingType == "SRCS") {
      output<- 1
      } else {
      output<- 0
      }
  
 if (j == 2 && CrossingType == "SRCS") {
    output<- 1
  } else {
    output<- 0
  }
    
  if (j == 3 && CrossingType == "FLB") {
    output<- 1
  } else {
    output<- 0
  }
  
  return(output)

}

y(1,2) # testing

cost_var<- function(j){
cj <- cost[j, ] #get from cost dataset
costj<- cj$c.value
return (costj)
}
cost_var(3)# test

n<- 1453
m<- 3
b<- 1000000

model <- MIPModel() %>%

  add_variable(x[i, j], i = 1:n, j= 1:m, type = "binary") %>%
  
  set_objective(sum_expr(p(i,j) * x[i,j], i = 1:n, j= 1:m), "max") %>%
  
  add_constraint(sum_expr(cost_var(j) * x[i,j], i = 1:n, j= 1:m)<= b) %>%
  add_constraint(sum_expr(x[i,j], j = 1:m)<= 1, i= 1:n) %>% # at most one countermeasure used
  add_constraint(x[i,j] <= y(i,j), i= 1:n, j = 1:m) # only selective countermeasures can be implemented

model

#Thanks in advance for your feedback.

Hi @dirkschumacher,
Im trying to solve a problem on optimization .We have an list of 68 raw materials in inventory stock in a vector called Inv stock .And we have the required amount of items in an array qty_req. We have a matrix of raw materials which has list of raw materials as column and rows contains the quantity of raw material req to complete a unit of item per order. That means complete raw materials required to complete each order can be determined by multiplying the qty_req vector with the raw material matrix .

Problem: We want to find qtysold as a output which should give us the appropriate qty to be manufactured instead of qty req so that we complete atleast 80% of each order and the max revenue can be obtained. We have revenue/unit for each order .

Eg of Raw material matrix is as below:

OrderNo.	A	B	C
1	0	1	1
2	2	3	0
3	1	2	3

Eg. of req_qty vector

OrderNo.	Qty req
1	2
2	1
3	3

Eg. of Inv_stock

A	B	C
4	0	4

Eg of Revenue per unit

OrderNo.	Revenue/unit
1	2400
2	4500
3	450

Total material required for manufacturing the req quantity will become req_qty* Raw_material.

Please find below the code which is not working. Please provide your suggestions .
model <- MIPModel() %>%
add_variable(Qtysold[i],i=1:n,type = "integer" ,lb=0 ) %>%
add_variable(Total_mat[i,j],i=1:n,j=1:68, type = "continuous", lb = 0) %>%
set_objective(sum_expr(Qtysold[i]*Revenueperunit[i],i=1:n), "max")%>%
add_constraint(Qtysold[i] <= qty_req[i], i = 1:n)%>%
add_constraint(Qtysold[i]*Raw_mat==Total_mat[i,j],i=1:n,j=1:68)%>%
add_constraint(sum_expr(Total_mat[i,j],i=1:n)<=Inv_Stock[1,j],j=1:68)
result <- solve_model(model, with_ROI(solver = "glpk"))

On adding the last two lines im getting errors like unexpected operator '('.

Thanks in advance.

Hi @dirkschumacher, @njtierney,
I am solving an optimization problem using the followoing code.
library(Rglpk)
library(magrittr)
library(dplyr)
library(ROI)
library(ROI.plugin.glpk)
library(ompr)
library(ompr.roi)
data<-read.csv('input.csv',header=TRUE,sep=",")
Inv_Stock<-read.csv('INV_DATA.csv',header = TRUE,sep = ",")
Inv_Stock<- sapply(Inv_Stock,as.matrix)
qty_req<- data[,5]
qty_req<- sapply(qty_req,as.double)

Revenueperunit<-data[,4]
Revenueperunit<-sapply(Revenueperunit,as.double)

Raw_mat<-data[,-c(1,2,3,4,5)]
Raw_mat <- sapply(Raw_mat,as.matrix)

model <- MIPModel() %>%
add_variable(Qtysold[i],i=1:n,type = "integer" ,lb=0 ) %>%
add_variable(Total_mat[i,j],i=1:n,j=1:68, type = "continuous", lb = 0) %>%
set_objective(sum_expr(Qtysold[i]*Revenueperunit[i],i=1:n), "max")%>%
add_constraint(Qtysold[i] <= qty_req[i], i = 1:n)%>%
add_constraint(rawmat[j,i]*Qtysold[i]==Total_mat[j,i])%>%
add_constraint(sum_expr(Total_mat[i,j],i=1:n)<=Inv_Stock[1,j],j=1:68)
result <- solve_model(model, with_ROI(solver = "glpk"))

Context of the problem:
A vector holding inventory of total raw materials available by colums containing the raw material name and a single row with total qty available per column.
A matrix raw material with column names of raw material and each row indicating per order number, the quantity if raw material req for unit item to be manufactured.
A vector qty req indicating per order number the number of quantity asked by a customer to manufacture.
A vector revenue/unit indicating the revenue from each unit of material per order.
Problem statement:
We want to maximize total revenue
we want to determine qtysold such that atleast 80% of the qty req per order is delivered.
when we multiply qtysold with req raw material matrix , total req raw material columnwise should not exceed the total available in inventory stock.

Error: while adding last two constraints, its not giving answer, it says'(' unexpected operator found.
Error in acc[[key]] : no such index at level 1

Hi, I am now facing issue in only last line, can anyone please tell if this is making it a quadratic constraint.

model <- MIPModel() %>%
add_variable(Qtysold[i],i=1:n,type = "integer" ,lb=0 ) %>%
add_variable(Total_mat[i,j],i=1:n,j=1:68, type = "continuous", lb = 0) %>%
set_objective(sum_expr(Qtysold[i]*Revenueperunit[i],i=1:n), "max")%>%
add_constraint(Qtysold[i] <= qty_req[i], i = 1:n)%>%
add_constraint(Qtysold[i]*Raw_mat[i,j]==Total_mat[i,j],i=1:n,j=1:68)%>%
add_constraint(colSums(Total_mat[i,j],i=1:n)<=Inv_Stock[j],j=1:68)

Error:rror in check_for_unknown_vars_impl(model, the_ast) :
The expression contains a variable that is not part of the model.