#!/usr/bin/env python3# Copyright 2025, Gurobi Optimization, LLC# Assign workers to shifts; each worker may or may not be available on a# particular day. We use lexicographic optimization to solve the model:# first, we minimize the linear sum of the slacks. Then, we constrain# the sum of the slacks, and we minimize a quadratic objective that# tries to balance the workload among the workers.importgurobipyasgpfromgurobipyimportGRBimportsys# Number of workers required for each shiftshifts,shiftRequirements=gp.multidict({"Mon1":3,"Tue2":2,"Wed3":4,"Thu4":4,"Fri5":5,"Sat6":6,"Sun7":5,"Mon8":2,"Tue9":2,"Wed10":3,"Thu11":4,"Fri12":6,"Sat13":7,"Sun14":5,})# Amount each worker is paid to work one shiftworkers,pay=gp.multidict({"Amy":10,"Bob":12,"Cathy":10,"Dan":8,"Ed":8,"Fred":9,"Gu":11,})# Worker availabilityavailability=gp.tuplelist([("Amy","Tue2"),("Amy","Wed3"),("Amy","Fri5"),("Amy","Sun7"),("Amy","Tue9"),("Amy","Wed10"),("Amy","Thu11"),("Amy","Fri12"),("Amy","Sat13"),("Amy","Sun14"),("Bob","Mon1"),("Bob","Tue2"),("Bob","Fri5"),("Bob","Sat6"),("Bob","Mon8"),("Bob","Thu11"),("Bob","Sat13"),("Cathy","Wed3"),("Cathy","Thu4"),("Cathy","Fri5"),("Cathy","Sun7"),("Cathy","Mon8"),("Cathy","Tue9"),("Cathy","Wed10"),("Cathy","Thu11"),("Cathy","Fri12"),("Cathy","Sat13"),("Cathy","Sun14"),("Dan","Tue2"),("Dan","Wed3"),("Dan","Fri5"),("Dan","Sat6"),("Dan","Mon8"),("Dan","Tue9"),("Dan","Wed10"),("Dan","Thu11"),("Dan","Fri12"),("Dan","Sat13"),("Dan","Sun14"),("Ed","Mon1"),("Ed","Tue2"),("Ed","Wed3"),("Ed","Thu4"),("Ed","Fri5"),("Ed","Sun7"),("Ed","Mon8"),("Ed","Tue9"),("Ed","Thu11"),("Ed","Sat13"),("Ed","Sun14"),("Fred","Mon1"),("Fred","Tue2"),("Fred","Wed3"),("Fred","Sat6"),("Fred","Mon8"),("Fred","Tue9"),("Fred","Fri12"),("Fred","Sat13"),("Fred","Sun14"),("Gu","Mon1"),("Gu","Tue2"),("Gu","Wed3"),("Gu","Fri5"),("Gu","Sat6"),("Gu","Sun7"),("Gu","Mon8"),("Gu","Tue9"),("Gu","Wed10"),("Gu","Thu11"),("Gu","Fri12"),("Gu","Sat13"),("Gu","Sun14"),])# Modelm=gp.Model("assignment")# Assignment variables: x[w,s] == 1 if worker w is assigned to shift s.# This is no longer a pure assignment model, so we must use binary variables.x=m.addVars(availability,vtype=GRB.BINARY,name="x")# Slack variables for each shift constraint so that the shifts can# be satisfiedslacks=m.addVars(shifts,name="Slack")# Variable to represent the total slacktotSlack=m.addVar(name="totSlack")# Variables to count the total shifts worked by each workertotShifts=m.addVars(workers,name="TotShifts")# Constraint: assign exactly shiftRequirements[s] workers to each shift s,# plus the slackreqCts=m.addConstrs((slacks[s]+x.sum("*",s)==shiftRequirements[s]forsinshifts),"_")# Constraint: set totSlack equal to the total slackm.addConstr(totSlack==slacks.sum(),"totSlack")# Constraint: compute the total number of shifts for each workerm.addConstrs((totShifts[w]==x.sum(w)forwinworkers),"totShifts")# Objective: minimize the total slack# Note that this replaces the previous 'pay' objective coefficientsm.setObjective(totSlack)# OptimizedefsolveAndPrint():m.optimize()status=m.statusifstatusin(GRB.INF_OR_UNBD,GRB.INFEASIBLE,GRB.UNBOUNDED):print("The model cannot be solved because it is infeasible or\               unbounded")sys.exit(1)ifstatus!=GRB.OPTIMAL:print(f"Optimization was stopped with status{status}")sys.exit(0)# Print total slack and the number of shifts worked for each workerprint("")print(f"Total slack required:{totSlack.X:g}")forwinworkers:print(f"{w} worked{totShifts[w].X:g} shifts")print("")solveAndPrint()# Constrain the slack by setting its upper and lower boundstotSlack.UB=totSlack.XtotSlack.LB=totSlack.X# Variable to count the average number of shifts workedavgShifts=m.addVar(name="avgShifts")# Variables to count the difference from average for each worker;# note that these variables can take negative values.diffShifts=m.addVars(workers,lb=-GRB.INFINITY,name="Diff")# Constraint: compute the average number of shifts workedm.addConstr(len(workers)*avgShifts==totShifts.sum(),"avgShifts")# Constraint: compute the difference from the average number of shiftsm.addConstrs((diffShifts[w]==totShifts[w]-avgShiftsforwinworkers),"Diff")# Objective: minimize the sum of the square of the difference from the# average number of shifts workedm.setObjective(gp.quicksum(diffShifts[w]*diffShifts[w]forwinworkers))# OptimizesolveAndPrint()