I am solving a transport optimization problem given a set of constraints. The following are the three key data sets that I have
#demand file demand - has demand(DEMAND) across 4821(DPP) sale points(D)
head(demand)
D PP DEMAND DPP
1 ADILABAD (V) - T:11001 OPC:PACK 131.00 ADILABAD (V) - T:11001:OPC:PACK
2 ADILABAD (V) - T:13003 OPC:PACK 235.00 ADILABAD (V) - T:13003:OPC:PACK
3 ADILABAD (V) - T:2006 PPC:PACK 30.00 ADILABAD (V) - T:2006:PPC:PACK
4 ADILABAD (V) - T:4001 OPC:PACK 30.00 ADILABAD (V) - T:4001:OPC:PACK
5 ADILABAD (V) - T:7006 OPC:NPACK 34.84 ADILABAD (V) - T:7006:OPC:NPACK
6 AHMEDABAD:1001 OPC:PACK 442.10 AHMEDABAD:1001:OPC:PACK
#Capacity file cc - has capacity constraint (MaxP, MinP) across 1823 sources(SOURCE)
head(cc,4)
SOURCE MinP MaxP
1 CHILAMKUR:P:OPC:NPACK:0:R 900 10806
2 CHILAMKUR:P:OPC:NPACK:0:W 900 10806
3 CHILAMKUR:P:OPC:PACK:0:R 5628 67536
4 CHILAMKUR:P:OPC:PACK:0:W 5628 67536
#LandingCost file LCMat - This is a matrix with the landing cost to deliver the product across the demand location (DPP) from a given source(SOURCE). This is an 1823 x 4821 matrix. Since the landing costs to all locations do not exist from a given source, I have replace that with a huge cost (10^6) to such DPPs.
I am using the OMPR package in R to optimize shipping material to meet the demand. This is potentially a very simple transport problem but it is taking a lot of time. I am using a 16GB ram machine
The following is the code. Could anyone guide me on what I should do better?
a = Sys.time()
grid = expand.grid(i = 1:nrow(LCMat),j = 1:ncol(LCMat))
grid_solve = grid[which(LCMat < 10^6),]
grid_notsolve = grid[which(LCMat >= 10^6),]
model <- MILPModel() %>%
add_variable(x[grid$i, grid$j],lb = 0, type = "continuous") %>%
add_constraint(x[grid_notsolve$i, grid_notsolve$j] == 0) %>%
add_constraint(sum_over(x[i,j], i = 1:nrow(LCMat)) <= demand$DEMAND[j], j = 1:ncol(LCMat)) %>%
add_constraint(sum_over(x[i,j], j = 1:ncol(LCMat)) <= cc$MaxP[i], i = 1:nrow(LCMat)) %>%
add_constraint(sum_over(x[i,j], j = 1:ncol(LCMat)) >= cc$MinP[i], i = 1:nrow(LCMat)) %>%
set_objective(sum_expr(LCMat[grid_solve$i,grid_solve$j]*x[grid_solve$i,grid_solve$j]),"min")
solution = model %>% solve_model(with_ROI(solver = "glpk", verbose = TRUE))
Sys.time() - a