Masterproblem setup + solve subproblem with dual values

src/Hospitalplanning.jl

@@ -14,7 +14,9 @@ include("Patient.jl")
 include("GenerateSampleData.jl")
 
 include("MIP/datastructures.jl")
+include("MIP/SupportFunctions.jl")
 include("MIP/subproblem.jl")
+include("MIP/masterproblem.jl")
 
 
 

src/MIP/SupportFunctions.jl

@@ -0,0 +1,23 @@
+getDays(resource,subcal) = (x->x.date[1]).(filter(day ->day.date in subcal,resource.calendar.workdays))
+
+function getTimeSlot(resource,dayID)
+    workdays = filter(x->x.date[1]==dayID,resource.calendar.workdays)
+    length(workdays) > 1 ? error("Multiple days with same dayID") :
+    length(workdays) == 1 ? (x->x.intID).(filter(x-> x.status == free ,workdays[1].timeslots)) : Int64[]
+end
+getTSendtime(resource,dayID,slotID) = Dates.value(filter(x->x.date[1]==dayID,resource.calendar.workdays)[1].timeslots[slotID].endTime)/60000000000
+function getTSendtime(resource,dayID,slotID,refDate)
+    workday = filter(x->x.date[1]==dayID,resource.calendar.workdays)[1]
+    datetime = workday.timeslots[slotID].endTime + workday.date[2]
+    Dates.value(datetime-refDate)/60000000000
+end
+function getTSstarttime(resource,dayID,slotID,refDate)
+    workday = filter(x->x.date[1]==dayID,resource.calendar.workdays)[1]
+    test = length(workday.timeslots)
+    slotID > test ? error("slotID out of range") :
+    datetime = workday.timeslots[slotID].startTime + workday.date[2]
+    Dates.value(datetime-refDate)/60000000000
+end
+getTSstarttime(resource,dayID,slotID) = Dates.value(filter(x->x.date[1]==dayID,resource.calendar.workdays)[1].timeslots[slotID].startTime)/60000000000
+getResource(resources,visits) = findqualifiedResourceIDs(resources,visits)
+getTdelta(v1,v2) = v1 == 1 ? 1440 - 8*60 : 0  #TODO create this funtion correctly

src/MIP/datastructures.jl

@@ -1,21 +1,38 @@
+mutable struct Masterproblem
+    model::JuMP.Model
+    consref_offtime
+    consref_onepatient
+    convexitycons
+    lambda
+    closingtime
+end
+
+mutable struct Subproblem
+    model::JuMP.Model
+    xvars; yvars; tvars; kvars
+    V; D; D_v; J; J_d; I; Ts; Te
+    Subproblem() = new()
+end
 
 
-mutable struct Subproblems
-    models::Array{JuMP.Model,1}
-    xvars
-    yvars
-    tvars
-    kvars
 
-    Subproblems(nsub) = new(Vector{JuMP.Model}(undef,nsub),[],[],[],[])
-end
+function addsubproblem(subproblems,nsub,sub::JuMP.Model,xvars,yvars,tvars,kvars, V, D ,D_v, J, J_d, I, Ts, Te)
+        cursub = Subproblem()
+        cursub.model = sub
+        cursub.xvars = xvars
+        cursub.yvars = yvars
+        cursub.tvars = tvars
+        cursub.kvars = kvars
+        cursub.V = V
+        cursub.D = D
+        cursub.D_v = D_v
+        cursub.J = J
+        cursub.J_d = J_d
+        cursub.I = I
+        cursub.Ts = Ts
+        cursub.Te = Te
 
-function addsubproblem(subproblems::Subproblems,sub::JuMP.Model,xvars,yvars,tvars,kvars,nsub)
-        subproblems.models[nsub] = sub
-        push!(subproblems.xvars, xvars)
-        push!(subproblems.yvars, yvars)
-        push!(subproblems.tvars, tvars)
-        push!(subproblems.kvars, kvars)
+        subproblems[nsub] = cursub
 end
 
 function getsubproblem(subproblems,patient)

src/MIP/masterproblem.jl

@@ -0,0 +1,79 @@
+
+
+function setupmaster(patients,resources,subMastercalendar)
+
+
+    master = Model(with_optimizer(Gurobi.Optimizer,OutputFlag=1))
+
+    K = length(patients)
+
+    P = (x->x.intID).(patients)
+    D = (x->x.intID).(resources)
+    J = keys(subMastercalendar)
+    J_d = Dict(d => getDays(resources[d],subMastercalendar) for d in D)
+    I = Dict(d => Dict(j => getTimeSlot(resources[d],j) for j in J) for d in D )
+    @variable(master,lambda[1:K] >= 0)
+
+    @variable(master,closingtime[d in D, j in J_d[d]] >= 0)
+
+    @objective(master, Min, sum(closingtime[d,j] for d in D, j in J_d[d]) + sum(100000*lambda[m] for p in P, m in 1:K))
+
+    @constraint(master,consref_offtime[d in D, j in J_d[d] ,i in I[d][j]], sum(lambda[m]*10000 for m in 1:K) <= closingtime[d,j])
+
+    @constraint(master,convexitycons[p in P],
+       sum(lambda[m] for m in 1:K if m == p) == 1 )
+
+    println("test1")
+    @constraint(master,consref_onepatient[ d in D, j in J_d[d], i in I[d][j]],
+        sum(lambda[m]*0 for p in P, m in 1:K) <=1 )
+        println("test2")
+
+    return Masterproblem(master,consref_offtime,consref_onepatient,convexitycons,lambda,closingtime)
+    # return consref_offtime, consref_onepatient, convexitycons,lambda,closingtime
+end
+
+
+function addcolumntomaster(masterproblem::Masterproblem,yvalues,xvalues,tvalues,patient,param,iteration,plandata,srproblem)
+    P, I, J, D, R, DR, w1,w2,activities,M,timeofday = get(param)
+    touchedconstraints = ConstraintRef[]
+    constraint_coefficients = Float64[]
+
+    mptemp = zeros(length(P))
+    mptemp[patient] = 1
+    Atemp = zeros(length(P),length(DR),length(I),length(J))
+    ttemp = zeros(length(DR),length(I),length(J))
+
+    for i in I, j in J, d in D
+        # if tvalues[j,d] == xtimes[i,j,d]
+        if xvalues[i,j,d] != 0
+            push!(touchedconstraints,masterproblem.consref_offtime[i,j,d])
+            push!(constraint_coefficients,tvalues[j,d])
+            ttemp[d,i,j] = tvalues[j,d]
+        end
+    end
+    for i in I, j in J, d in DR
+        if xvalues[i,j,d] != 0
+            push!(touchedconstraints,masterproblem.consref_onepatient[i,j,d])
+            push!(constraint_coefficients,1)
+            Atemp[patient,d,i,j] = 1
+        end
+    end
+
+    push!(touchedconstraints,masterproblem.convexitycons[patient])
+    push!(constraint_coefficients,1)
+
+    push!(masterproblem.lambda,@variable(
+        masterproblem.model,
+        lower_bound = 0,
+        base_name = "lambda_new[$(patient),$(iteration)]" ,
+        # objective = sum(yvalues)*w2[patient],
+        #inconstraints= touchedconstraints,
+        #coefficients = constraint_coefficients
+    ))
+JuMP.set_objective_coefficient(masterproblem.model,masterproblem.lambda[end],sum(yvalues)*w2[patient])
+JuMP.set_coefficient.(touchedconstraints,masterproblem.lambda[end],constraint_coefficients)
+
+
+    addplan(plandata,Atemp,ttemp,sum(yvalues)*w2[patient],mptemp)
+
+end

src/MIP/subproblem.jl

@@ -1,28 +1,5 @@
-getDays(resource,subcal) = (x->x.date[1]).(filter(day ->day.date in subcal,resource.calendar.workdays))
 
-function getTimeSlot(resource,dayID)
-    workdays = filter(x->x.date[1]==dayID,resource.calendar.workdays)
-    length(workdays) > 1 ? error("Multiple days with same dayID") :
-    length(workdays) == 1 ? (x->x.intID).(workdays[1].timeslots) : Int64[]
-end
-getTSendtime(resource,dayID,slotID) = Dates.value(filter(x->x.date[1]==dayID,resource.calendar.workdays)[1].timeslots[slotID].endTime)/60000000000
-function getTSendtime(resource,dayID,slotID,refDate)
-    workday = filter(x->x.date[1]==dayID,resource.calendar.workdays)[1]
-    datetime = workday.timeslots[slotID].endTime + workday.date[2]
-    Dates.value(datetime-refDate)/60000000000
-end
-function getTSstarttime(resource,dayID,slotID,refDate)
-    workday = filter(x->x.date[1]==dayID,resource.calendar.workdays)[1]
-    test = length(workday.timeslots)
-    slotID > test ? error("slotID out of range") :
-    datetime = workday.timeslots[slotID].startTime + workday.date[2]
-    Dates.value(datetime-refDate)/60000000000
-end
-getTSstarttime(resource,dayID,slotID) = Dates.value(filter(x->x.date[1]==dayID,resource.calendar.workdays)[1].timeslots[slotID].startTime)/60000000000
-getResource(resources,visits) = findqualifiedResourceIDs(resources,visits)
-getTdelta(v1,v2) = v1 == 1 ? 1440 - 8*60 : 0  #TODO create this funtion correctly
-
-function setup_sub!(subproblems::Subproblems,patient::Patient,resources::Array{Resource},subMastercalendar::Dict{Int,Date})
+function setup_sub!(subproblems,patient::Patient,resources::Array{Resource},subMastercalendar::Dict{Int,Date})
     M = 10000
     sub = Model(with_optimizer(Gurobi.Optimizer,OutputFlag=1))
     refDate = minimum(subMastercalendar)[2]+Time(0)
@@ -62,24 +39,24 @@ function setup_sub!(subproblems::Subproblems,patient::Patient,resources::Array{R
             -(1-kvars[v2,v1])*M) #TODO Te, Ts, Tdelta
 
 
-    addsubproblem(subproblems,sub,xvars,yvars,tvars,kvars,patient.intID)
+    addsubproblem(subproblems,patient.intID,sub,xvars,yvars,tvars,kvars,V,D,D_v,J,J_d,I,Ts,Te)
 end
 
 
 
 
-# function solveSub(subproblems,phi,pi,kappa,patient,param::Parameters)
-#     P, I, J, D, R, DR, w1,w2,activities,M,timeofday = get(param)
-#     sub,xvars_sub, yvars_sub, tvars_sub = getsubproblem(subproblems,patient)
-#
-#     @objective(sub, Min, sum(yvars_sub[j]*w2[patient] for j in J) - sum(tvars_sub[j,d]*phi[i,j,d] for i in I, j in J, d in D ) - sum(xvars_sub[i,j,d] * pi[i,j,d] for i in I, j in J, d in DR)-kappa[patient] )
-#
-#
-#     optimize!(sub)
-#     status = termination_status(sub)
-#     if status != MOI.TerminationStatusCode(1)
-#         throw("Error: Non optimal sub-problem")
-#     end
-#
-#     return objective_value(sub), value.(xvars_sub), value.(yvars_sub),value.(tvars_sub)
-# end
+function solveSub(sub,phi,pi,kappa,patient)
+
+
+
+    @objective(sub.model, Min, sum(sub.yvars[j] for j in sub.J) - sum(sub.tvars[d,j]*phi[d,j,i] for d in sub.D,j in sub.J_d[d],i in sub.I[d][j] ) - sum(sub.xvars[d,j,i] * pi[d,j,i] for d in sub.D,j in sub.J_d[d],i in sub.I[d][j])-kappa[patient] )
+
+
+    optimize!(sub.model)
+    status = termination_status(sub.model)
+    if status != MOI.TerminationStatusCode(1)
+        throw("Error: Non optimal sub-problem")
+    end
+
+    return objective_value(sub.model), value.(sub.xvars), value.(sub.yvars),value.(sub.tvars)
+end

test/current.jl

@@ -1,6 +1,7 @@
 
 using Dates
 using Debugger
+using JuMP
 using Hospitalplanning
 HP = Hospitalplanning
 
@@ -14,7 +15,7 @@ mastercalendar = HP.MasterCalendar(startdate,enddate)
 mastercalendar[1]
 columns = Dict(:Consultation => "Consultation" , :Telefon =>"Telephone" , :TTE => "TTE", :AEKG => "AEKG",  :MR => "MR", :Holter=>"Holter")
 patients = HP.readPatientTable(path,sheet,columns,mastercalendar)
-
+patients
 
 path_resourceOverview = "C:/Users/hebb/.julia/dev/Hospitalplanning/test/Sample data/GUCHamb_Timeslot_test.xlsx"
 sheet_amb = "GUCH AMB"
@@ -41,9 +42,15 @@ for i in 1:length(patients)
     termination_status(subproblems.models[i])
 end
 
-using JuMP
-JuMP.optimize!(subproblems.models[1])
-termination_status(subproblems.models[1])
+mp = HP.setupmaster(patients,GUCHAmb_resources,submastercalendar)
+
+optimize!(mp.model)
+
+phi = dual.(mp.consref_offtime)
+pi = dual.(mp.consref_onepatient)
+kappa = dual.(mp.convexitycons)
+
+
 #--------------------------------------------------------------------------
 
 HP.getDays(GUCHAmb_resources[2],submastercalendar)