All,
I am working on the following code. Takes some inputs from two files (portfolio and blend) and does some sorting/ functions, then runs a monte carlo-esque simulation into a linear program. I Think everything is set up correctly, but I am getting some reactivity errors that I can't pin down. Any help would be nice.
library(shiny)
library(triangle)
library(lpSolveAPI)
# External CSV data
start_mat<-matrix(rep("A",150),ncol=15,nrow=10)
family_vector<-c("A","B","C","C","D",rep("E",5))
widgetweight_vect<-rep(5,10)
comp_weight_vector<-rep(4,10)
widgets_vector<-rep(100,10)
portfolio<-cbind.data.frame(start_mat,family_vector,widgetweight_vect,comp_weight_vector,widgets_vector)
start_mat<-matrix(rep("A",150)ncol=15,rnow=10)
family_vector<-c("A","B","C","C","D",rep("E",5))
widgetweight_vect<-rep(5,10)
comp_weight_vector<-rep(4,10)
widgets_vector<-rep(100,10)
portolfio<-c(start_mat,family_vector,widgetweight_vect,comp_weight_vectorwidgets_vector)
blend<-matrix(rep(.20,25),nrow=5, ncol=5)
colnames(blend)<-c("a","b","c","d","e")
rownames(blend)<-c("a","b","c","d","e")
ui <- shinyUI(fluidPage(
br(),
actionButton("numb", "generate "),
column(3,numericInput("plan", label = h3("plan"), value = 50)),
column(3,numericInput("recweeks", label = h3("recweeks"), value = 50)),
column(3,numericInput("targetgen", label = h3("targetgen"), min= .1, max = .2, value = .17)),
column(3,numericInput("process_min", label = h3("process_min"), value = 5000)),
column(3,numericInput("process_target_trsp", label = h3("process_target_trsp"), min= .1, max = 1, value = .95)),
column(3,numericInput("process_max_trsp", label = h3("process_max_trsp"), min= .1, max = 1, value = 1)),
column(3,numericInput("planning_cycle", label = h3("planning_cycle"), min= 1, max = 7, value = 2)),
column(3,numericInput("rec_starting_inventory", label = h3("rec_starting_inventory"), min= 1, max = 10, value = 5)),
column(3,numericInput("process_min_campaign", label = h3("process_min_campaign"), min= 1, max = 2000, value = 50)),
column(3,numericInput("widgets_target", label = h3("widgets_target"), min= 5000, max = 7000, value = 5000)),
column(3,numericInput("widgets_variation_perc", label = h3("widgets_variation_perc"), min= .01, max = .2, value = .05)),
column(3,numericInput("process_recoup_max", label = h3("process_recoup_max"), min= .20, max = .50, value = .4)),
br(),
br(),
plotOutput("Plot")
)
)
server <- shinyServer(function(input, output) {
plan<-reactive({input$plan})
recweeks<-reactive({input$recweeks})
targetgen<-reactive({input$targetgen})
process_min<-reactive({input$process_min})
process_target_trsp<-reactive({input$process_target_trsp})
process_max_trsp<-reactive({input$process_max_trsp})
planning_cycle<-reactive({input$planning_cycle})
rec_starting_inventory<-reactive({input$rec_starting_inventory})
process_min_campaign<-reactive({input$process_min_campaign})
widgets_target<-reactive({input$widgets_target})
widgets_variation_perc<-reactive({input$widgets_variation_perc})
process_recoup_max<-reactive({input$process_recoup_max})
model<- eventReactive(input$numb, {
totalsim<-as.numeric(plan()*recweeks())
blend_length<-length(blend)
nrow<-function(x) dim(x)[1]
globalinventory<-vector(mode="numeric", length=(length(blend)-1))
dims<-nrow(portfolio)
rec_stats<-matrix(c(0,0,0,0,0,0,0),ncol=7,nrow=1)
inventory_byBL<-matrix(c(rep(0,length(blend)-1)),ncol=(length(blend)-1),nrow=1)
#Set the random number seed
set.seed(100)
#Develop table of process performances
proc_location<-0.0
proc_myscale<-85.16507
proc_myshape<-13.35404
process_trs<-matrix((rweibull(totalsim, shape=proc_myshape, scale=proc_myscale)+proc_location)/100, nrow=totalsim, ncol=1)
process_trs<-replace(process_trs,process_trs>process_max_trsp(),process_max_trsp())
process_perf<-process_trs/process_target_trsp()
process_index<-1
#write.csv(process_perf,file="MEZprocess_perf.csv")
#Develop tables of proced rec generations
location<-5.5646+4
myscale<-8.07516
myshape<-3.76987
recgen_perc<-matrix((rweibull(totalsim*dims, shape=myshape, scale=myscale)+location)/100, nrow=dims, ncol=totalsim)
#Develop the tables of plan Variations
trimin<-0.1
trimax<-2.2
trimode<-.95
tri<-rtriangle(totalsim*dims, a=trimin,b=trimax,c=trimode)
trimat<-matrix(tri,nrow=dims,ncol=totalsim)
#uniform distribution for some additional variation
unif_range<-as.numeric(widgets_variation_perc())
do_var<-matrix(runif(totalsim,min=(1-unif_range),max=(1+unif_range)),nrow=totalsim,ncol=1)
#Create the rec generation vector
widgetweight<-portfolio[,17]
comp_weight<-portfolio[,18]
family<-as.factor(portfolio[,16])
widgets<-portfolio[,19]
for(iii in 1:plan){
local({
widgets_base<-widgets*trimat[,iii]
total_widgets<-sum(widgets_base)
widgets_base<-widgets_base*widgets_target/total_widgets
for (ii in 1:recweeks){ ##Iterates through different rec generations
local({
widgets_var<-widgets_base*do_var[process_index]
recgen_gr<-widgetweight*(recgen_perc[,ii])*widgets_var
parentBL_gr<-comp_weight*(1+targetgen)*widgets_base*process_perf[process_index]
newportfolio<-data.frame(family,recgen_gr,parentBL_gr)
#Now create the table of the aggregate values by family for the rec generated and the parent proc available for consumption
rectable<-aggregate(. ~ newportfolio$family, newportfolio, sum)
rectable
rectable<-rectable[,-2]
colnames(rectable)<-c("Family","recgen_gr","parentBL_gr")
rectable$parentBL_name<-substr(rectable$Family,1,5)
rectable
#Combine the rec generation information with the blend tables. This will be used to generated the constraints for the LP
#The two constraints that must be generated are the rec generation and the parent proc quantity
combinedtable<-merge(rectable,blend,all = TRUE)
combinedtable[is.na(combinedtable)] <- 0
combinedtable
parent.table<-data.frame(combinedtable$parentBL_name,combinedtable$parentBL_gr)
parent.table
colnames(parent.table)<-c("parentBL_name","parentBL_gr")
#rec generation constraint created for later use in the LP
generation.constraint<-combinedtable$recgen_gr
#Now parent proc constraint
blend_parent<-colnames(combinedtable)
blend_parent<-unlist(colnames(combinedtable[5:ncol(combinedtable)]))
blend_parent<-data.frame(blend_parent)
dim(blend_parent)
colnames(blend_parent)<-("parentBL_name")
parent_gen<-merge(blend_parent,parent.table)
parent_gen$gen<-as.factor(match(parent_gen$parentBL_name,blend_parent$parentBL_name))
parent_gen<-parent_gen[order(parent_gen$gen),]
parent.constraint<-parent_gen[,2]
rnd <- function(x) trunc(x+sign(x)*0.5+.5)
#Create the production wheel based on the minimum run requirement for Z
parent.wheel<-parent.constraint
v<-length(parent.wheel)
parent.wheel<-ifelse(parent.wheel/process_min_campaign>1,1,rnd(process_min_campaign/parent.wheel))
parent.wheel.adj<-ifelse(process_index%%parent.wheel==0,parent.wheel,0)
#Now multiple the results of the production wheel to the parent proc constraint to make the proces available to parent proced rec
parent.constraint_adj<-parent.constraint*parent.wheel.adj
##################################################################################################################################################################
#Create LP with decision variables based on the number of constraints
n<-length(generation.constraint)
#Add a factor for "overage" dummy proc if there is rec that cannot be consumed
m<-length(parent.constraint_adj)+1
model<-make.lp(0,n*m)
rhs.gen<-generation.constraint+globalinventory+ifelse(process_index==1,generation.constraint*rec_starting_inventory,0)
rhs.proc<-parent.constraint_adj
blend_con<-cbind(1/blend[,2:(length(blend))],over=1000)
rownames(blend_con)<-blend[,1]
####Row Constraints###########################
add.constraint(model,rep(1,m), "=", rhs.gen[1],indices=c(1:m))
for (i in 1:(n-1)){
start<-i*m+1
end<-i*m+m
add.constraint(model,rep(1,m), "=", rhs.gen[i+1],indices=c(start:end))
}
for (i in 1:n) {
col<-c()
for (ii in seq(i,n*m,m)){
col<-cbind(col,ii)
}
col
add.constraint(model,blend_con[,i],"<=",rhs.proc[i],indices=c(col))
}
lp.control(model,sense='max')
objective<-1:(n*m)
remove<-seq(m,(n*m),m)
objective<-objective[-remove]
length(objective)
set.objfn(model,rep(1,length(objective)),indices=(objective))
solve(model)
write.lp(model,filename="rec_linearprogram.lp")
get.objective(model)
length(get.variables(model))
table<-get.variables(model)
tablemat<-t(matrix(table,ncol=n,nrow=m))
sum(rhs.gen)
get.primal.solution(model, orig = FALSE)
#get.variables(model)
values<-get.variables(model)
values<-values[-objective]
sum(values)+get.objective(model)
#print(model)
#Actual consumption based on rules
tablemat_con<- ifelse(tablemat[,1:(m-1)]<process_min,0,tablemat[,1:(m-1)])
real_consumption<-sum(tablemat_con) ##Write this value to a table "rec consumed"
#Remainders that Cannot be consumed because of run rule
table_mat_remainder<-ifelse(tablemat[,1:(m-1)]<process_min,tablemat[,1:(m-1)],0)
table_mat_remainder<-rowSums(table_mat_remainder)
#table_mat_remainder ##Write this value to a table "Lost Opportunity"
missed_op<-sum(table_mat_remainder)
overage_next_day<-tablemat[,m]
globalinventory<-overage_next_day+ table_mat_remainder ##Write this value to a table "New Inventory"
globalinventory_sum<-sum(globalinventory)
real_con_perc<-real_consumption/(sum(parent_gen[,2])*do_var[process_index]) ##achieved rec consumptions
max_con_perc<-(real_consumption+missed_op)/(sum(parent_gen[,2])*do_var[process_index]) ##Unconstrained rec consumption
actual_gen_perc<-sum(recgen_gr)/(sum(parent_gen[,2])*do_var[process_index])
colnames(rec_stats)<-c("RealCon","GlobalInven","RealConPer","MaxPotent","ActrecGen","process Performance","tires per day")
rec_stats<<-rbind(rec_stats,c(real_consumption,globalinventory_sum,real_con_perc,max_con_perc,actual_gen_perc,process_perf[process_index],sum(widgets_var)))
inventory_byBL<<-rbind(inventory_byBL,c(globalinventory))
process_index<-process_index+1
})
}
#Plot 1
plot1<- matplot(rec_stats[,3:5]*100,col=c("springgreen4","grey50","red"),lwd=.5,lty=3,type = c("b"),pch=1,main="Evolution of rec Consumption/Generation",ylab="%", xlab="Days",ylim=c(0,30))
list(plot1=plot1)
})
}
})
output$Plot<-renderPlot({model()$plot1})
})
shinyApp(ui = ui, server = server)
Errors:
Warning: Error in .getReactiveEnvironment()$currentContext: Operation not allowed without an active reactive context. (You tried to do something that can only be done from inside a reactive expression or observer.)
Stack trace (innermost first):
48: .getReactiveEnvironment()$currentContext
47: .dependents$register
46: plan
45: server [#14]
4: <Anonymous>
3: do.call
2: print.shiny.appobj
1: print
Error in .getReactiveEnvironment()$currentContext() :
Operation not allowed without an active reactive context. (You tried to do something that can only be done from inside a reactive expression or observer.)
PS I know I could do a lot of this coding better with plyr and some apply vs loops, but it's what I have.
There were miscellaneous errors. There were a few off-by one errors, and a few reactive functions that were missing the ()
(like plan
needed to be plan()
).
A few things to point out here:
browser()
statement) and then stopping and inspecting the variables there.plan
) where it needs a numerical value, it will probably give you a rather confusing error message - like the one you got. So if you see confusing error message, look around for missing brackets like these ()
.df
or sum
. Type in sum = sum+1
for example and try and make sense of the message.I had to comment out a constraint block, and change a few off-by-one mismatches to get this to work. It probably will be easy for you to fix these though.
Anyway here is the code:
library(shiny)
library(triangle)
library(lpSolveAPI)
# External CSV data
start_mat <- matrix(rep("A",150),ncol = 15,nrow = 10)
family_vector <- c("A","B","C","C","D",rep("E",5))
widgetweight_vect <- rep(5,10)
comp_weight_vector <- rep(4,10)
widgets_vector <- rep(100,10)
portfolio <- cbind.data.frame(start_mat,family_vector,widgetweight_vect,comp_weight_vector,widgets_vector)
start_mat <- matrix(rep("A",150),ncol = 15,nrow = 10)
family_vector <- c("A","B","C","C","D",rep("E",5))
widgetweight_vect <- rep(5,10)
comp_weight_vector <- rep(4,10)
widgets_vector <- rep(100,10)
portolfio <- c(start_mat,family_vector,widgetweight_vect,comp_weight_vector,widgets_vector)
blend <- matrix(rep(.20,25),nrow = 5,ncol = 5)
colnames(blend) <- c("a","b","c","d","e")
rownames(blend) <- c("a","b","c","d","e")
ui <- shinyUI(fluidPage(
br(),
actionButton("numb","generate "),
column(3,numericInput("plan",label = h3("plan"),value = 50)),
column(3,numericInput("recweeks",label = h3("recweeks"),value = 50)),
column(3,numericInput("targetgen",label = h3("targetgen"),min = .1,max = .2,value = .17)),
column(3,numericInput("process_min",label = h3("process_min"),value = 5000)),
column(3,numericInput("process_target_trsp",label = h3("process_target_trsp"),min = .1,max = 1,value = .95)),
column(3,numericInput("process_max_trsp",label = h3("process_max_trsp"),min = .1,max = 1,value = 1)),
column(3,numericInput("planning_cycle",label = h3("planning_cycle"),min = 1,max = 7,value = 2)),
column(3,numericInput("rec_starting_inventory",label = h3("rec_starting_inventory"),min = 1,max = 10,value = 5)),
column(3,numericInput("process_min_campaign",label = h3("process_min_campaign"),min = 1,max = 2000,value = 50)),
column(3,numericInput("widgets_target",label = h3("widgets_target"),min = 5000,max = 7000,value = 5000)),
column(3,numericInput("widgets_variation_perc",label = h3("widgets_variation_perc"),min = .01,max = .2,value = .05)),
column(3,numericInput("process_recoup_max",label = h3("process_recoup_max"),min = .20,max = .50,value = .4)),
br(),
br(),
plotOutput("Plot")
)
)
server <- shinyServer(function(input,output) {
plan <- reactive({ input$plan })
recweeks <- reactive({ input$recweeks })
targetgen <- reactive({ input$targetgen })
process_min <- reactive({ input$process_min })
process_target_trsp <- reactive({ input$process_target_trsp })
process_max_trsp <- reactive({ input$process_max_trsp })
planning_cycle <- reactive({ input$planning_cycle })
rec_starting_inventory <- reactive({ input$rec_starting_inventory })
process_min_campaign <- reactive({ input$process_min_campaign })
widgets_target <- reactive({ input$widgets_target })
widgets_variation_perc <- reactive({ input$widgets_variation_perc })
process_recoup_max <- reactive({ input$process_recoup_max })
model <- eventReactive(input$numb,{
totalsim <- as.numeric(plan() * recweeks())
blend_length <- length(blend)
nrow <- function(x) dim(x)[1]
globalinventory <- vector(mode = "numeric",length = (length(blend))) # -1 ?
dims <- nrow(portfolio)
rec_stats <- matrix(c(0,0,0,0,0,0,0),ncol = 7,nrow = 1)
inventory_byBL <- matrix(c(rep(0,length(blend) - 0)),ncol = (length(blend) - 0 ),nrow = 1) # changed 1 to 0 twice
#Set the random number seed
set.seed(100)
#Develop table of process performances
proc_location <- 0.0
proc_myscale <- 85.16507
proc_myshape <- 13.35404
process_trs <- matrix((rweibull(totalsim,shape = proc_myshape,scale = proc_myscale) + proc_location) / 100,nrow = totalsim,ncol = 1)
process_trs <- replace(process_trs,process_trs > process_max_trsp(),process_max_trsp())
process_perf <- process_trs / process_target_trsp()
process_index <- 1
#write.csv(process_perf,file="MEZprocess_perf.csv")
#Develop tables of proced rec generations
location <- 5.5646 + 4
myscale <- 8.07516
myshape <- 3.76987
recgen_perc <- matrix((rweibull(totalsim * dims,shape = myshape,scale = myscale) + location) / 100,nrow = dims,ncol = totalsim)
#Develop the tables of plan Variations
trimin <- 0.1
trimax <- 2.2
trimode <- .95
tri <- rtriangle(totalsim * dims,a = trimin,b = trimax,c = trimode)
trimat <- matrix(tri,nrow = dims,ncol = totalsim)
#uniform distribution for some additional variation
unif_range <- as.numeric(widgets_variation_perc())
do_var <- matrix(runif(totalsim,min = (1 - unif_range),max = (1 + unif_range)),nrow = totalsim,ncol = 1)
#Create the rec generation vector
widgetweight <- portfolio[,17]
comp_weight <- portfolio[,18]
family <- as.factor(portfolio[,16])
widgets <- portfolio[,19]
for (iii in 1:plan()) {
print(sprintf("iii:%d",iii))
local({
widgets_base <- widgets * trimat[,iii]
total_widgets <- sum(widgets_base)
widgets_base <- widgets_base * widgets_target() / total_widgets
for (ii in 1:recweeks()) {
##Iterates through different rec generations
# print(sprintf(" ii:%d",ii))
local( {
widgets_var <- widgets_base * do_var[process_index]
recgen_gr <- widgetweight * (recgen_perc[,ii]) * widgets_var
parentBL_gr <- comp_weight * (1 + targetgen()) * widgets_base * process_perf[process_index]
newportfolio <- data.frame(family,recgen_gr,parentBL_gr)
#Now create the table of the aggregate values by family for the rec generated and the parent proc available for consumption
rectable <- aggregate(. ~ newportfolio$family,newportfolio,sum)
rectable
rectable <- rectable[,-2]
colnames(rectable) <- c("Family","recgen_gr","parentBL_gr")
rectable$parentBL_name <- substr(rectable$Family,1,5)
rectable
#Combine the rec generation information with the blend tables. This will be used to generated the constraints for the LP
#The two constraints that must be generated are the rec generation and the parent proc quantity
combinedtable <- merge(rectable,blend,all = TRUE)
combinedtable[is.na(combinedtable)] <- 0
combinedtable
parent.table <- data.frame(combinedtable$parentBL_name,combinedtable$parentBL_gr)
parent.table
colnames(parent.table) <- c("parentBL_name","parentBL_gr")
#rec generation constraint created for later use in the LP
generation.constraint <- combinedtable$recgen_gr
#Now parent proc constraint
blend_parent <- colnames(combinedtable)
blend_parent <- unlist(colnames(combinedtable[5:ncol(combinedtable)]))
blend_parent <- data.frame(blend_parent)
dim(blend_parent)
colnames(blend_parent) <- ("parentBL_name")
parent_gen <- merge(blend_parent,parent.table)
parent_gen$gen <- as.factor(match(parent_gen$parentBL_name,blend_parent$parentBL_name))
parent_gen <- parent_gen[order(parent_gen$gen),]
parent.constraint <- parent_gen[,2]
rnd <- function(x) trunc(x + sign(x) * 0.5 + .5)
#Create the production wheel based on the minimum run requirement for Z
parent.wheel <- parent.constraint
v <- length(parent.wheel)
parent.wheel <- ifelse(parent.wheel / process_min_campaign() > 1,1,rnd(process_min_campaign() / parent.wheel))
parent.wheel.adj <- ifelse(process_index %% parent.wheel == 0,parent.wheel,0)
#Now multiple the results of the production wheel to the parent proc constraint to make the proces available to parent proced rec
parent.constraint_adj <- parent.constraint * parent.wheel.adj
##################################################################################################################################################################
#Create LP with decision variables based on the number of constraints
n <- length(generation.constraint)
#Add a factor for "overage" dummy proc if there is rec that cannot be consumed
m <- length(parent.constraint_adj) + 1
model <- make.lp(0,n * m)
rhs.gen <- generation.constraint + globalinventory + ifelse(process_index == 1,generation.constraint * rec_starting_inventory(),0)
rhs.proc <- parent.constraint_adj
blend_con <- cbind(1 / blend[,2:(dim(blend)[2])],over = 1000)
rownames(blend_con) <- blend[,1]
####Row Constraints###########################
add.constraint(model,rep(1,m),"=",rhs.gen[1],indices = c(1:m))
for (i in 1:(n - 1)) {
start <- i * m + 1
end <- i * m + m
add.constraint(model,rep(1,m),"=",rhs.gen[i + 1],indices = c(start:end))
}
#for (i in 1:n) {
#col <- c()
#for (ii in seq(i,n * m,m)) {
#col <- cbind(col,ii)
#}
#col
#add.constraint(model,blend_con[,i],"<=",rhs.proc[i],indices = c(col))
#}
lp.control(model,sense = 'max')
objective <- 1:(n * m)
remove <- seq(m,(n * m),m)
objective <- objective[-remove]
length(objective)
set.objfn(model,rep(1,length(objective)),indices = (objective))
solve(model)
write.lp(model,filename = "rec_linearprogram.lp")
get.objective(model)
length(get.variables(model))
table <- get.variables(model)
tablemat <- t(matrix(table,ncol = n,nrow = m))
sum(rhs.gen)
get.primal.solution(model,orig = FALSE)
#get.variables(model)
values <- get.variables(model)
values <- values[-objective]
sum(values) + get.objective(model)
#print(model)
#Actual consumption based on rules
tablemat_con <- ifelse(tablemat[,1:(m - 1)] < process_min(),0,tablemat[,1:(m - 1)])
real_consumption <- sum(tablemat_con) ##Write this value to a table "rec consumed"
#Remainders that Cannot be consumed because of run rule
table_mat_remainder <- ifelse(tablemat[,1:(m - 1)] < process_min(),tablemat[,1:(m - 1)],0)
#table_mat_remainder <- rowSums(table_mat_remainder)
table_mat_remainder <- sum(table_mat_remainder)
#table_mat_remainder ##Write this value to a table "Lost Opportunity"
missed_op <- sum(table_mat_remainder)
overage_next_day <- tablemat[,m]
globalinventory <- overage_next_day + table_mat_remainder ##Write this value to a table "New Inventory"
globalinventory_sum <- sum(globalinventory)
real_con_perc <- real_consumption / (sum(parent_gen[,2]) * do_var[process_index]) ##achieved rec consumptions
max_con_perc <- (real_consumption + missed_op) / (sum(parent_gen[,2]) * do_var[process_index]) ##Unconstrained rec consumption
actual_gen_perc <- sum(recgen_gr) / (sum(parent_gen[,2]) * do_var[process_index])
colnames(rec_stats) <- c("RealCon","GlobalInven","RealConPer","MaxPotent","ActrecGen","process Performance","tires per day")
rec_stats <<- rbind(rec_stats,c(real_consumption,globalinventory_sum,real_con_perc,max_con_perc,actual_gen_perc,process_perf[process_index],sum(widgets_var)))
inventory_byBL <<- rbind(inventory_byBL,c(globalinventory))
process_index <- process_index + 1
})
}
#Plot 1
plot1 <- matplot(rec_stats[,c(1,2,6,7)] * 100,col = c("springgreen4","grey50","red"),lwd = .5,lty = 3,type = c("b"),pch = 1,main = "Evolution of rec Consumption/Generation",ylab = "%",xlab = "Days",ylim = c(0,30))
list(plot1 = plot1)
})
}
})
output$Plot <- renderPlot({ model()$plot1 })
})
shinyApp(ui = ui,server = server)
And this is what it outputs: