diff --git a/data/cases.csv b/data/cases.csv new file mode 100644 index 0000000..e893c11 --- /dev/null +++ b/data/cases.csv @@ -0,0 +1,31 @@ +CaseID,ConsultantID,Procedure,Speciality,Median Duration,TargetDeadline +1,11,Cataract Surgery,Ophthalmology,45,07/07/2020 +2,11,Vitrectomy ,Ophthalmology,70,17/07/2020 +3,11,Cataract Surgery,Ophthalmology,45,05/06/2020 +4,11,Cataract Surgery,Ophthalmology,45,28/06/2020 +5,11,Cataract Surgery,Ophthalmology,45,22/07/2020 +6,11,Cataract Surgery,Ophthalmology,45,03/07/2020 +7,11,Cataract Surgery,Ophthalmology,45,03/07/2020 +8,11,Eyelid Lesion Excision,Ophthalmology,30,11/06/2020 +9,11,Vitrectomy,Ophthalmology,70,04/07/2020 +10,11,Cataract Surgery,Ophthalmology,45,14/07/2020 +11,11,Cataract Surgery,Ophthalmology,45,06/07/2020 +12,11,Cataract Surgery,Ophthalmology,45,16/07/2020 +13,11,Trabeculectomy,Ophthalmology,100,01/08/2020 +14,11,Vitrectomy,Ophthalmology,70,15/06/2020 +15,11,Vitrectomy,Ophthalmology,70,20/06/2020 +16,11,Cataract Surgery,Ophthalmology,45,18/06/2020 +17,11,Cataract Surgery,Ophthalmology,45,20/07/2020 +18,11,Trabeculectomy,Ophthalmology,100,04/07/2020 +19,11,Eyelid Skin Graft,Ophthalmology,120,15/06/2020 +20,11,Cataract Surgery,Ophthalmology,45,13/07/2020 +21,11,Cataract Surgery,Ophthalmology,45,04/07/2020 +22,11,Cataract Surgery,Ophthalmology,45,28/06/2020 +23,11,Cataract Surgery,Ophthalmology,45,22/07/2020 +24,11,Cataract Surgery,Ophthalmology,45,06/07/2020 +25,11,Vitrectomy,Ophthalmology,70,08/06/2020 +26,11,Cataract Surgery,Ophthalmology,45,14/07/2020 +27,11,Cataract Surgery,Ophthalmology,45,06/07/2020 +28,11,Cataract Surgery,Ophthalmology,45,16/07/2020 +29,11,Trabeculectomy,Ophthalmology,100,11/07/2020 +30,11,Vitrectomy,Ophthalmology,70,27/06/2020 diff --git a/data/sessions.csv b/data/sessions.csv new file mode 100644 index 0000000..d84a04e --- /dev/null +++ b/data/sessions.csv @@ -0,0 +1,5 @@ +SessionID,Date,Start,End,Duration,ConsultantID,Specialty +1001,03/06/2020,08:30:00,18:00:00,570,11,Ophthalmology +1002,10/06/2020,08:30:00,18:00:00,570,11,Ophthalmology +1003,17/06/2020,08:30:00,18:00:00,570,11,Ophthalmology +1004,25/06/2020,08:30:00,13:00:00,270,11,Ophthalmology diff --git a/modelling/scheduler.py b/modelling/scheduler.py deleted file mode 100644 index 55f68b6..0000000 --- a/modelling/scheduler.py +++ /dev/null @@ -1,100 +0,0 @@ -import pandas as pd -import numpy as np -import pyomo.environ as pe -import pyomo.gdp as pyogdp -from pyomo.core.base.set_types import Any -import os - - -# check with Taha if code is too similar to Alstom? -class TheatreScheduler: - - def __init__(self, case_file_path, session_file_path): - """ - Read case and session data into Pandas DataFrames - Args: - case_file_path (str): path to case data in CSV format - session_file_path (str): path to theatre session data in CSV format - """ - try: - self.df_cases = pd.read_csv(case_file_path) - except FileNotFoundError: - print("Case data not found.") - - try: - self.df_sessions = pd.read_csv(session_file_path) - except FileNotFoundError: - print("Session data not found") - - def _generate_case_durations(self): - """ - Generate mapping of cases IDs to median case time for the procedure - Returns: - (dict): dictionary with CaseID as key and median case time (mins) for procedure as value - """ - return pd.Series(self.df_cases["Median Duration"].values, index=self.df_cases["CaseID"]).to_dict() - - def _generate_session_durations(self): - """ - Generate mapping of all theatre sessions IDs to session duration in minutes - Returns: - (dict): dictionary with SessionID as key and session duration as value - """ - return pd.Series(self.df_sessions["Duration"].values, index=self.df_sessions["SessionID"]).to_dict() - - def create_model(self): - model = pe.ConcreteModel() - - # Model Data - model.CASES = pe.Set(initialize=self.df_cases["CaseID"].tolist()) - model.SESSIONS = pe.Set(initialize=self.df_sessions["SessionID"].tolist()) - model.TASKS = pe.Set(initialize=model.CASES * model.SESSIONS, dimen=2) - model.CASE_DURATIONS = pe.Param(model.CASES, initialize=self._generate_case_durations()) - model.SESSION_DURATIONS = pe.param(model.SESSIONS, initalize=self._generate_session_durations()) - model.M = pe.Param(initialize=1e7) # big M - max_util = 0.85 - - # Decision Variables - model.SESSION_ASSIGNED = pe.Var() - model.CASE_START_TIME = pe.Var() - model.UTILISATION = pe.Var() - - # Objective - def objective_function(model): - pass - model.OBJECTIVE = pe.Objective() - - # Constraints - # Case start time must be after start time of assigned theatre session - def case_start_time(): - pass - model.CASE_START = pe.Constraint() - - # Case end time must be before end time of assigned theatre session - def case_end_time(): - pass - model.CASE_END_TIME = pe.Constraint() - - # Cases can be assigned to a maximum of one session - def session_assignment(): - pass - model.SESSION_ASSIGNMENT = pe.Constraint() - - def disjunctions(): - pass - model.DISJUNCTIONS = pyogdp.Disjunction() - - def theatre_util(): - pass - model.THEATRE_UTIL = pe.Constraint() - - return model - - def solve(self): - pass - - -if __name__ == "__main__": - case_path = os.path.join(os.path.dirname(os.getcwd()), "data", "case_data.csv") - session_path = os.path.join(os.path.dirname(os.getcwd()), "data", "session_data.csv") - scheduler = TheatreScheduler(case_file_path=case_path, session_file_path=session_path) diff --git a/modelling/scheduler_median.py b/modelling/scheduler_median.py new file mode 100644 index 0000000..a1d7c24 --- /dev/null +++ b/modelling/scheduler_median.py @@ -0,0 +1,241 @@ +import pandas as pd +import numpy as np +import pyomo.environ as pe +import pyomo.gdp as pyogdp +import os +import matplotlib.pyplot as plt +import matplotlib.cm as cm +from itertools import product + + +# check with Taha if code is too similar to Alstom? +class TheatreScheduler: + + def __init__(self, case_file_path, session_file_path): + """ + Read case and session data into Pandas DataFrames + Args: + case_file_path (str): path to case data in CSV format + session_file_path (str): path to theatre session data in CSV format + """ + try: + self.df_cases = pd.read_csv(case_file_path) + except FileNotFoundError: + print("Case data not found.") + + try: + self.df_sessions = pd.read_csv(session_file_path) + except FileNotFoundError: + print("Session data not found") + + self.model = self.create_model() + + def _generate_case_durations(self): + """ + Generate mapping of cases IDs to median case time for the procedure + Returns: + (dict): dictionary with CaseID as key and median case time (mins) for procedure as value + """ + return pd.Series(self.df_cases["Median Duration"].values, index=self.df_cases["CaseID"]).to_dict() + + def _generate_session_durations(self): + """ + Generate mapping of all theatre sessions IDs to session duration in minutes + Returns: + (dict): dictionary with SessionID as key and session duration as value + """ + return pd.Series(self.df_sessions["Duration"].values, index=self.df_sessions["SessionID"]).to_dict() + + def _generate_session_start_times(self): + """ + Generate mapping from SessionID to session start time + Returns: + (dict): dictionary with SessionID as key and start time in minutes since midnight as value + """ + # Convert session start time from HH:MM:SS format into seconds elapsed since midnight + self.df_sessions.loc[:, "Start"] = pd.to_timedelta(self.df_sessions["Start"]) + self.df_sessions.loc[:, "Start"] = self.df_sessions["Start"].dt.total_seconds() / 60 + return pd.Series(self.df_sessions["Start"].values, index=self.df_sessions["SessionID"]).to_dict() + + def _generate_disjunctions(self): + """ + #TODO + Returns: + disjunctions (list): list of tuples containing disjunctions + """ + cases = self.df_cases["CaseID"].to_list() + sessions = self.df_sessions["SessionID"].to_list() + disjunctions = [] + for (case1, case2, session) in product(cases, cases, sessions): + if (case1 != case2) and (case2, case1, session) not in disjunctions: + disjunctions.append((case1, case2, session)) + + return disjunctions + + def create_model(self): + model = pe.ConcreteModel() + + # Model Data + model.CASES = pe.Set(initialize=self.df_cases["CaseID"].tolist()) + model.SESSIONS = pe.Set(initialize=self.df_sessions["SessionID"].tolist()) + model.TASKS = pe.Set(initialize=model.CASES * model.SESSIONS, dimen=2) + model.CASE_DURATION = pe.Param(model.CASES, initialize=self._generate_case_durations()) + model.SESSION_DURATION = pe.Param(model.SESSIONS, initialize=self._generate_session_durations()) + model.SESSION_START_TIME = pe.Param(model.SESSIONS, initialize=self._generate_session_start_times()) + model.DISJUNCTIONS = pe.Set(initialize=self._generate_disjunctions(), dimen=3) + + ub = 1440 # seconds in a day + model.M = pe.Param(initialize=1e3*ub) # big M + max_util = 0.85 + num_sessions = self.df_sessions.shape[0] + + # Decision Variables + model.SESSION_ASSIGNED = pe.Var(model.TASKS, domain=pe.Binary) + model.CASE_START_TIME = pe.Var(model.TASKS, bounds=(0, ub), within=pe.PositiveReals) + model.UTILISATION = pe.Var(model.SESSIONS, bounds=(0, 1), within=pe.PositiveReals) + model.MEDIAN_UTIL = pe.Var(bounds=(0, ub), within=pe.PositiveReals) + model.DUMMY_BINARY = pe.Var(model.SESSIONS, domain=pe.Binary) + model.CANCEL_SESSION = pe.Var(model.SESSIONS, domain=pe.Binary, within=pe.PositiveReals) + + # Objective + def objective_function(model): + #return pe.summation(model.UTILISATION) + return model.MEDIAN_UTIL + model.OBJECTIVE = pe.Objective(rule=objective_function, sense=pe.maximize) + + # Constraints + + # TODO add constraint to complete before deadline if it is assigned + # TODO add constraint to make tasks follow each other without gaps? + + # Case start time must be after start time of assigned theatre session + def case_start_time(model, case, session): + return model.CASE_START_TIME[case, session] >= model.SESSION_START_TIME[session] - \ + ((1 - model.SESSION_ASSIGNED[(case, session)])*model.M) + model.CASE_START = pe.Constraint(model.TASKS, rule=case_start_time) + + # Case end time must be before end time of assigned theatre session + def case_end_time(model, case, session): + return model.CASE_START_TIME[case, session] + model.CASE_DURATION[case] <= model.SESSION_START_TIME[session] + \ + model.SESSION_DURATION[session]*max_util + ((1 - model.SESSION_ASSIGNED[(case, session)]) * model.M) + model.CASE_END_TIME = pe.Constraint(model.TASKS, rule=case_end_time) + + # Cases can be assigned to a maximum of one session + def session_assignment(model, case): + return sum([model.SESSION_ASSIGNED[(case, session)] for session in model.SESSIONS]) <= 1 + model.SESSION_ASSIGNMENT = pe.Constraint(model.CASES, rule=session_assignment) + + def no_case_overlap(model, case1, case2, session): + return [model.CASE_START_TIME[case1, session] + model.CASE_DURATION[case1] <= model.CASE_START_TIME[case2, session] + \ + ((2 - model.SESSION_ASSIGNED[case1, session] - model.SESSION_ASSIGNED[case2, session])*model.M), + model.CASE_START_TIME[case2, session] + model.CASE_DURATION[case2] <= model.CASE_START_TIME[case1, session] + \ + ((2 - model.SESSION_ASSIGNED[case1, session] - model.SESSION_ASSIGNED[case2, session])*model.M)] + + model.DISJUNCTIONS_RULE = pyogdp.Disjunction(model.DISJUNCTIONS, rule=no_case_overlap) + + def theatre_util(model, session): + return model.UTILISATION[session] == (1 / model.SESSION_DURATION[session]) * \ + sum([model.SESSION_ASSIGNED[case, session]*model.CASE_DURATION[case] for case in model.CASES]) + model.THEATRE_UTIL = pe.Constraint(model.SESSIONS, rule=theatre_util) + + def cancel_sessions(model, session): # TODO + return model.CANCEL_SESSION[session] <= 1 - model.M*sum([model.SESSION_ASSIGNED[case, session] for case in model.CASES]) + model.SET_CANCEL_SESSIONS = pe.Constraint(model.SESSIONS, rule=cancel_sessions) + + def force_cancel_sessions(model, session): + return sum([model.SESSION_ASSIGNED[case, session] for case in model.CASES]) <= 0 + model.M*(1-model.CANCEL_SESSION[session]) + #model.FORCE_CANCEL_SESSIONS = pe.Constraint(model.SESSIONS, rule=force_cancel_sessions) + + def set_dummy_variable(model): + return sum([model.DUMMY_BINARY[session] for session in model.SESSIONS]) == np.floor(num_sessions/2) + model.FLOOR = pe.Constraint(rule=set_dummy_variable) + + def set_median_util(model, session): + return model.MEDIAN_UTIL <= model.UTILISATION[session] + model.DUMMY_BINARY[session]*model.M + model.SET_MEDIAN_UTIL = pe.Constraint(model.SESSIONS, rule=set_median_util) + + pe.TransformationFactory("gdp.bigm").apply_to(model) + + return model + + def solve(self, solver_name, options=None, solver_path=None): + + if solver_path is not None: + solver = pe.SolverFactory(solver_name, executable=solver_path) + else: + solver = pe.SolverFactory(solver_name) + + # TODO remove - too similar to alstom + if options is not None: + for key, value in options.items(): + solver.options[key] = value + + solver_results = solver.solve(self.model, tee=True) + + results = [{"Case": case, + "Session": session, + "Start": self.model.CASE_START_TIME[case, session](), + "Assignment": self.model.SESSION_ASSIGNED[case, session]()} + for (case, session) in self.model.TASKS] + + self.df_times = pd.DataFrame(results) + + all_cases = self.model.CASES.value_list + cases_assigned = [] + cases_missed = [] + for (case, session) in self.model.SESSION_ASSIGNED: + if self.model.SESSION_ASSIGNED[case, session] == 1: + cases_assigned.append(case) + + cases_missed = list(set(all_cases).difference(cases_assigned)) + print("Number of cases assigned = {} out of {}:".format(len(cases_assigned), len(all_cases))) + print("Cases assigned: ", cases_assigned) + print("Number of cases missed = {} out of {}:".format(len(cases_missed), len(all_cases))) + print("Cases missed: ", cases_missed) + self.model.UTILISATION.pprint() + print("Total Utilisation = {}".format(sum(self.model.UTILISATION.get_values().values()))) + print("Number of constraints = {}".format(solver_results["Problem"].__getitem__(0)["Number of constraints"])) + #self.model.SESSION_ASSIGNED.pprint() + #print(self.df_times.to_string()) + self.draw_gantt() + + def draw_gantt(self): + + df = self.df_times[self.df_times["Assignment"] == 1] + cases = sorted(list(df['Case'].unique())) + sessions = sorted(list(df['Session'].unique())) + + bar_style = {'alpha': 1.0, 'lw': 25, 'solid_capstyle': 'butt'} + text_style = {'color': 'white', 'weight': 'bold', 'ha': 'center', 'va': 'center'} + colors = cm.Dark2.colors + + df.sort_values(by=['Case', 'Session']) + df.set_index(['Case', 'Session'], inplace=True) + + fig, ax = plt.subplots(1, 1) + for c_ix, c in enumerate(cases, 1): + for s_ix, s in enumerate(sessions, 1): + if (c, s) in df.index: + xs = df.loc[(c, s), 'Start'] + xf = df.loc[(c, s), 'Start'] + \ + self.df_cases[self.df_cases["CaseID"] == c]["Median Duration"] + ax.plot([xs, xf], [s] * 2, c=colors[c_ix % 7], **bar_style) + ax.text((xs + xf) / 2, s, c, **text_style) + + ax.set_title('Session Schedule') + ax.set_xlabel('Time') + ax.set_ylabel('Sessions') + ax.grid(True) + + fig.tight_layout() + plt.show() + + +if __name__ == "__main__": + case_path = os.path.join(os.path.dirname(os.getcwd()), "data", "case_data_long.csv") + session_path = os.path.join(os.path.dirname(os.getcwd()), "data", "session_data.csv") + cbc_path = "C:\\Users\\LONLW15\\Documents\\Linear Programming\\Solvers\\cbc.exe" + + options = {"seconds": 30} + scheduler = TheatreScheduler(case_file_path=case_path, session_file_path=session_path) + scheduler.solve(solver_name="cbc", solver_path=cbc_path, options=options)