[Python] Spreading out shift assignments in constraint solver (ortools)

I used the Google OR-Tools Employee Scheduling script (thanks by the way) to make a on-call scheduler. Everything works fine and it is doing what it is supposed to. It makes sure each person works about the same amount of "shifts" (two week periods), it lets certain shifts be requested and I added a constraint where it won't let someone work shifts back to back.

Everyone gets the same amount of shifts (as much as possible):

df['Name'].value_counts()
Out[42]:
Jeff 7
Bubba 7
Sarah 6
Scott 6
Name: Name, dtype: int64


One thing I notice is that it will use up a person as much as it can before moving on to the next. E.g it will go 1-2-1-2-1-3...3-4-3-2-3-4. As opposed to 1-2-3-4-1-2-3-4...

print(df)
Name Date Shift
0 Sarah 01-09-2022 On Call
1 Scott 01-23-2022 On Call
2 Sarah 02-06-2022 On Call
3 Scott 02-20-2022 On Call
4 Sarah 03-06-2022 On Call
5 Jeff 03-20-2022 On Call
6 Sarah 04-03-2022 On Call
7 Jeff 04-17-2022 On Call
8 Sarah 05-01-2022 On Call
9 Jeff 05-15-2022 On Call
10 Sarah 05-29-2022 On Call
11 Jeff 06-12-2022 On Call
12 Bubba 06-26-2022 On Call
13 Jeff 07-10-2022 On Call
14 Bubba 07-24-2022 On Call
15 Jeff 08-07-2022 On Call
16 Scott 08-21-2022 On Call
17 Bubba 09-04-2022 On Call
18 Jeff 09-18-2022 On Call
19 Bubba 10-02-2022 On Call
20 Scott 10-16-2022 On Call
21 Bubba 10-30-2022 On Call
22 Scott 11-13-2022 On Call
23 Bubba 11-27-2022 On Call
24 Scott 12-11-2022 On Call
25 Bubba 12-25-2022 On Call


(See how it kind of burns up one person--like Sarah at the start. I would expect Scott to be a lot like this, because of the -1 in request to not be on call during a large stretch--but a more even spread amongst everyone else would be ideal).

So I have two questions:

1. Is there a way to make it distribute the people more evenly?
2. Also, can I add another level of constraint in here by identifying certain time periods that contain holidays and then equally distribute those as well (kind of like a shift inside a shift)?

Here is my script:

# %% imports
import pandas as pd
from ortools.sat.python import cp_model

# %% Data for the model
num_employees = 4
num_shifts = 1
num_oncall_shifts = 26
all_employees = range(num_employees)
all_shifts = range(num_shifts)
all_oncall_shifts = range(num_oncall_shifts)

dict_shift_name = {0: 'On Call'}

dict_emp_name = {0: 'Bubba', 1: 'Scott', 2: 'Jeff', 3: 'Sarah'}

dict_dates = {
0: '01-09-2022',
1: '01-23-2022',
2: '02-06-2022',
3: '02-20-2022',
4: '03-06-2022',
5: '03-20-2022',
6: '04-03-2022',
7: '04-17-2022',
8: '05-01-2022',
9: '05-15-2022',
10: '05-29-2022',
11: '06-12-2022',
12: '06-26-2022',
13: '07-10-2022',
14: '07-24-2022',
15: '08-07-2022',
16: '08-21-2022',
17: '09-04-2022',
18: '09-18-2022',
19: '10-02-2022',
20: '10-16-2022',
21: '10-30-2022',
22: '11-13-2022',
23: '11-27-2022',
24: '12-11-2022',
25: '12-25-2022'
}

shift_requests = [
[
#Employee 0 Bubba
#1/09 1/23 2/06 2/20 3/06 3/20 4/03 4/17 5/01 5/15 5/29 6/12
[0], [0], [0], [0], [0], [0], [0], [0], [0], [0], [0], [0],
#6/26 7/10 7/24 8/07 8/21 9/04 9/18 10/02 10/16 10/30 11/13 11/27
[0], [0], [0], [0], [-4], [0], [0], [0], [0], [0], [0], [0],
#12/11 12/25
[0], [0]
],

[
#Employee 1 Scott
#1/09 1/23 2/06 2/20 3/06 3/20 4/03 4/17 5/01 5/15 5/29 6/12
[0],[0],[0],[0],[0],[0],[0],[0],[0],[0],[-1],[-1],
#6/26 7/10 7/24 8/07 8/21 9/04 9/18 10/02 10/16 10/30 11/13 11/27
[-1],[-1],[-1],[-1],[-1],[-1],[-1],[0],[0],[0],[0],[0],
#12/11 12/25
[0],[0]
],

[
#Employee 2 Jeff
#1/09 1/23 2/06 2/20 3/06 3/20 4/03 4/17 5/01 5/15 5/29 6/12
[0],[0],[0],[0],[0],[0],[0],[0],[0],[0],[0],[0],
#6/26 7/10 7/24 8/07 8/21 9/04 9/18 10/02 10/16 10/30 11/13 11/27
[0],[0],[0],[0],[-2],[0],[0],[0],[0],[0],[0],[0],
#12/11 12/25
[0],[0]
],

[
#Employee 3 Sarah
#1/09 1/23 2/06 2/20 3/06 3/20 4/03 4/17 5/01 5/15 5/29 6/12
[0],[0],[0],[0],[0],[0],[0],[0],[0],[0],[0],[0],
#6/26 7/10 7/24 8/07 8/21 9/04 9/18 10/02 10/16 10/30 11/13 11/27
[0],[0],[0],[0],[-3],[0],[0],[0],[0],[0],[0],[0],
#12/11 12/25
[0],[0]
],
]

# dataframe
df = pd.DataFrame(columns=['Name', 'Date', 'Shift'])

# %% Create the Model
model = cp_model.CpModel()

# %% Create the variables

# Shift variables# Creates shift variables.
# shifts[(n, d, s)]: nurse 'n' works shift 's' on day 'd'.
shifts = {}
for n in all_employees:
for d in all_oncall_shifts:
for s in all_shifts:
shifts[(n, d, s)] = model.NewBoolVar('shift_n%id%is%i' % (n, d, s))

# %% Add constraints
# Each shift is assigned to exactly one employee in .
for d in all_oncall_shifts :
for s in all_shifts:
model.AddExactlyOne(shifts[(n, d, s)] for n in all_employees)

# Each employee works at most one shift per oncall_shifts.
for n in all_employees:
for d in all_oncall_shifts:
model.AddAtMostOne(shifts[(n, d, s)] for s in all_shifts)

# Try to distribute the shifts evenly, so that each employee works
# min_shifts_per_employee shifts. If this is not possible, because the total
# number of shifts is not divisible by the number of employee, some employees will
# be assigned one more shift.
min_shifts_per_employee = (num_shifts * num_oncall_shifts) // num_employees
if num_shifts * num_oncall_shifts % num_employees == 0:
max_shifts_per_employee = min_shifts_per_employee
else:
max_shifts_per_employee = min_shifts_per_employee + 1
for n in all_employees:
num_shifts_worked = 0
for d in all_oncall_shifts:
for s in all_shifts:
num_shifts_worked += shifts[(n, d, s)]
model.Add(min_shifts_per_employee <= num_shifts_worked)
model.Add(num_shifts_worked <= max_shifts_per_employee)

# "penalize" working shift back to back
for d in all_employees:
for b in all_oncall_shifts[:-1]:
for r in all_shifts:
for r1 in all_shifts:
model.AddImplication(shifts[(d, b, r)], shifts[(d, b+1, r1)].Not())

# %% Objective
model.Maximize(
sum(shift_requests[n][d] * shifts[(n, d, s)] for n in all_employees
for d in all_oncall_shifts for s in all_shifts))

# %% Solve
# Creates the solver and solve.
solver = cp_model.CpSolver()
status = solver.Solve(model)
if status == cp_model.OPTIMAL:
print('Solution:')
for d in all_oncall_shifts:
print('On Call Starts: ', dict_dates[d])
for n in all_employees:
for s in all_shifts:
if solver.Value(shifts[(n, d, s)]) == 1:
if shift_requests[n][d] == 1:
print(dict_emp_name[n], ' is ', dict_shift_name, '(requested).')
else:
print(dict_emp_name[n], ' is ', dict_shift_name,
'(not requested).')

list_append = [dict_emp_name[n], dict_dates[d], dict_shift_name]

df.loc[len(df)] = list_append
print()
print(f'Number of shift requests met = {solver.ObjectiveValue()}',
f'(out of {num_employees * min_shifts_per_employee})')
else:
print('No optimal solution found !')

# %% Stats
print('\nStatistics')
print(' - conflicts: %i' % solver.NumConflicts())
print(' - branches : %i' % solver.NumBranches())
print(' - wall time: %f s' % solver.WallTime())


I get that it is difficult to program "fairness" into something, so any help would be greatly appreciated.

*** edit--added examples ***

Continue reading...