0
$\begingroup$

My scheduling problem has a specific pattern of 2 days of work followed by 2 days off. Additionally, there is a week off after every 4 weeks of work. The scheduling matrix (P) has five different schedules, where 0 represents not working, and 1 represents working. The constraints of my problem include pattern assignment constraints and worker efficiency constraints.

P = [[0, 0, 0, 0, 0, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0],
     [1, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1],
     [0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1],
     [0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 0, 0, 1, 1, 0],
     [1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0]]

My objective function aims to maximize the days spent together by married couples while minimizing the days they both work, to make it easier to take care of the children. I am trying to find a compromise between the needs of the employees and the company by assuring at least a certain percentage of the workforce (e.g., 40%) is present every day throughout the year while providing generous rest days for the employees. The goal of my project is to increase the workforce proportion to its highest while satisfying the objective functions for 100 employees. However, my current algorithm works well only for up to 30 people. I have noticed that for some reason, certain values of N (such as 14 and 19) result in the highest availability (K can increase up to 4), but I am not sure why. Currently I believe it might be related to how I distribute the worker’s efficiency.

NOTATION:

*𝑐_π‘šπ‘›: π‘€π‘Žπ‘‘π‘Ÿπ‘–π‘₯ π‘€β„Žπ‘’π‘Ÿπ‘’ π‘π‘Žπ‘‘π‘‘π‘’π‘Ÿπ‘› π‘š π‘Žπ‘›π‘‘ 𝑛 π‘π‘œπ‘‘β„Ž π‘Žπ‘Ÿπ‘’ π‘œπ‘“π‘“ π‘€π‘œπ‘Ÿπ‘˜

𝑀_π‘šπ‘›: π‘€π‘Žπ‘‘π‘Ÿπ‘–π‘₯ π‘€β„Žπ‘’π‘Ÿπ‘’ π‘π‘Žπ‘‘π‘‘π‘’π‘Ÿπ‘› π‘š π‘Žπ‘›π‘‘ 𝑛 π‘π‘œπ‘‘β„Ž π‘Žπ‘Ÿπ‘’ π‘Žπ‘‘ π‘€π‘œπ‘Ÿπ‘˜

π‘π‘Ž_𝑖𝑗: π΅π‘–π‘›π‘Žπ‘Ÿπ‘¦ π‘π‘Žπ‘Ÿπ‘Žπ‘šπ‘’π‘‘π‘’π‘Ÿ π‘€β„Žπ‘–π‘β„Ž π‘ β„Žπ‘œπ‘€π‘  π‘€β„Žπ‘’π‘‘β„Žπ‘’π‘Ÿ π‘π‘’π‘Ÿπ‘ π‘œπ‘› 𝑖 π‘Žπ‘›π‘‘ 𝑗 π‘Žπ‘Ÿπ‘’ π‘π‘œπ‘’π‘π‘™π‘’π‘  (Assumed half of the workforce was married)

π‘Ž_𝑖: 𝐸𝑓𝑓𝑖𝑐𝑖𝑒𝑛𝑐𝑦 π‘œπ‘“ π‘π‘’π‘Ÿπ‘ π‘œπ‘› 𝑖(πΆπ‘’π‘Ÿπ‘Ÿπ‘’π‘›π‘‘π‘™π‘¦ π‘šπ‘œπ‘‘π‘’π‘™π‘™π‘’π‘‘ π‘Žπ‘  π‘›π‘œπ‘Ÿπ‘šπ‘Žπ‘™ π‘‘π‘–π‘ π‘‘π‘Ÿπ‘–π‘π‘’π‘‘π‘–π‘œπ‘›)

π‘₯_π‘–π‘š: π΅π‘–π‘›π‘Žπ‘Ÿπ‘¦ π·π‘’π‘π‘–π‘ π‘–π‘œπ‘› π‘‰π‘Žπ‘Ÿπ‘–π‘Žπ‘π‘™π‘’ π‘€β„Žπ‘’π‘‘β„Žπ‘’π‘Ÿ π‘π‘’π‘Ÿπ‘ π‘œπ‘› 𝑖 𝑔𝑒𝑑𝑠 π‘Žπ‘ π‘ π‘–π‘”π‘›π‘’π‘‘ π‘π‘Žπ‘‘π‘‘π‘’π‘Ÿπ‘› π‘š

𝑦_𝑖𝑑: π‘Šβ„Žπ‘’π‘‘β„Žπ‘’π‘Ÿ π‘π‘’π‘Ÿπ‘ π‘œπ‘› 𝑖 π‘€π‘œπ‘Ÿπ‘˜π‘  π‘œπ‘› π‘‘π‘Žπ‘¦ 𝑑

𝑧_π‘šπ‘‘: Binary parameter whether pattern m works on day d

N: number of employees *

MODEL

Max Ξ£_𝑖 Ξ£_𝑗 Ξ£_π‘š Ξ£_𝑛 (𝛼(π‘π‘Ž_𝑖𝑗 𝑐_π‘šπ‘› π‘₯_π‘–π‘š π‘₯_𝑗𝑛)-(1βˆ’π›Ό)(π‘π‘Ž_𝑖𝑗 𝑀_π‘šπ‘› π‘₯_π‘–π‘š π‘₯_𝑗𝑛 ))

s.t.

  1. sum(π‘₯_π‘–π‘š)==1≫Person can only have one schedule

  2. 𝑦_𝑖𝑑 == π‘ π‘’π‘š(𝑧_π‘šπ‘‘ π‘₯_π‘–π‘š )≫ πΉπ‘œπ‘Ÿ π‘π‘’π‘Ÿπ‘ π‘œπ‘› 𝑖 π‘‘π‘œ π‘€π‘œπ‘Ÿπ‘˜ π‘œπ‘› π‘‘π‘Žπ‘¦ 𝑑

  3. sum(𝑦_𝑖𝑑 π‘Ž_𝑖) β‰₯π‘˜Γ—π‘ π‘’π‘šπ‘ π‘’π‘š π‘Šβ„Žπ‘’π‘Ÿπ‘’ π‘˜ 𝑖𝑠 π‘‘β„Žπ‘’ π‘π‘Ÿπ‘œπ‘π‘œπ‘Ÿπ‘‘π‘–π‘œπ‘› π‘Žπ‘›π‘‘ π‘ π‘’π‘šπ‘ π‘’π‘š 𝑖𝑠 π‘‘β„Žπ‘’ π‘ π‘’π‘š π‘œπ‘“ 𝑒𝑓𝑓𝑖𝑐𝑖𝑛𝑐𝑖𝑒𝑠 π‘œπ‘“ π‘’π‘šπ‘π‘™π‘œπ‘¦π‘’π‘’π‘  k initially set as 0.3 and improved

  4. 𝑒_π‘–π‘—π‘šπ‘›β‰€π‘₯_π‘–π‘š Linearise 𝑒_π‘–π‘—π‘šπ‘›β‰€π‘₯_𝑗𝑛 𝑒_π‘–π‘—π‘šπ‘›β‰₯π‘₯_π‘–π‘š+ π‘₯_π‘—π‘›βˆ’1

  5. 𝛼 is currently set as 0.7

My questions are:

1. Is there a better way to model the distribution of worker efficiency? Currently, I am using a normal distribution

a = [0] * N
smart = round(N * 0.021)
ok = round(smart + N * 0.136)
notbad = round(ok + N * 0.682)
bad = round(notbad + N * 0.136)
dumb = round(bad + N * 0.021)
for j in range(N):
    if j <= smart:
        a[j] = 5
    elif smart < j <= ok:
        a[j] = 4
    elif ok < j <= notbad:
        a[j] = 3
    elif notbad < j <= bad:
        a[j] = 2
    elif bad < j <= dumb:
        a[j] = 1

2. Why do certain values of N, such as 14 and 19, result in better outcomes in terms of availability (K)? >> the maximum for 14 and 19 was 0.4 while it was only 0.34 for 11 employees

Here is a full version of the code.

from ortools.linear_solver import pywraplp
import numpy as np
import time
import random
import signal
import multiprocessing as mp





class Me:

    def __init__(self, kkk, N):

        self.hihi(kkk, N)

    def hihi(self, kkk, N):

        startTime = time.time()

        print("When N is ", N)
        # D = 35 # five weeks

        # while()
        D = 35
        
        I = [i for i in range(N)]  # set of workers : 0~39, 40λͺ…
        # k = 0.3
        k = kkk
        print("k is ", k)
        # Patterns
        P = []

        #5 DIFFERENT PATTERNS EMPLOYEES CAN BE ASSIGNED TO
        P = [[0, 0, 0, 0, 0, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0],
             [1, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1],
             [0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1],
             [0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0, 1, 1, 0, 0, 1, 1, 0],
             [1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, 0]]

        z = P.copy()

        # pattern-pattern (c, w)
        c = [[0] * len(P)] * len(P)  # When both patterns are off work
        w = [[0] * len(P)] * len(P)  # When both patterns are on work

        for m in range(len(P)):  # mth pattern
            c[m] = [sum([(1 - a) * (1 - b) for a, b in zip(P[m], P[n])]) for n in range(len(P))]
            w[m] = [sum([a * b for a, b in zip(P[m], P[n])]) for n in range(len(P))]

        #OMIT C MATRIX LESS THAN 9 AND W MATRIX BIGGER THAN 5;
        #AS EVEN WITH A CONVENTIONAL 5 DAY WEEK, THERE WOULD BE 10 DAYS WHEN THEY ARE BOTH OFF WORK
        for m in range(len(P)):
            for n in range(len(c[m])):
                if c[m][n] < 9 or w[m][n] > 5:

                    c[m][n] = 0
                    w[m][n] = 0

        print("C matrix", c)
        print("W matrix", w)

        #EFFICIENCY OF WORKERS OF NORMAL DIST
        a = [0] * N
        smart = round(N * 0.021)
        ok = round(smart + N * 0.136)
        notbad = round(ok + N * 0.682)
        bad = round(notbad + N * 0.136)
        dumb = round(bad + N * 0.021)
        for j in range(N):
            if j <= smart:
                a[j] = 5
            elif smart < j <= ok:
                a[j] = 4
            elif ok < j <= notbad:
                a[j] = 3
            elif notbad < j <= bad:
                a[j] = 2
            elif bad < j <= dumb:
                a[j] = 1
        print("original a", a)
        for kk in range(N):
            b = random.randint(0, N - 1)
            tmp = a[kk]
            a[kk] = a[b]
            a[b] = tmp
        print("a", a)
        sumsum = 0
        for aa in a:
            sumsum = sumsum + aa
        print(sumsum)

        # ASSUMED HALF OF THE WORKFORCE IS MARRIED
        pa = [[0] * len(I)] * len(I)
        for i in [2 * j for j in range(int(len(I) / 4))]:
            pa[i] = [0] * len(I)
            pa[i + 1] = [0] * len(I)
            pa[i][i + 1] = 1
            pa[i + 1][i] = 1
            # p[2*i+1][2*i] = 1
        print("pa", pa)

        # Solver
        # Create the mip solver with the CBC backend.
        solver = pywraplp.Solver.CreateSolver('CBC')

        x = []
        y = []
        u = []

        for i in I:
            t1 = []
            t2 = []
            for m in range(len(P)):
                # for m in range(20):
                t1.append(solver.IntVar(0, 1, ''))
            for d in range(len(P[0])):
                # for d in range(35):
                t2.append(solver.IntVar(0, 1, ''))
            x.append(t1)
            y.append(t2)

        for i in I:
            t1 = []
            for j in I:
                t2 = []
                for m in range(len(P)):
                    # for m in range(20):
                    t3 = []
                    for n in range(len(P)):
                        # for n in range(20):
                        t3.append(solver.IntVar(0, 1, ''))
                    t2.append(t3)
                t1.append(t2)
            u.append(t1)

        #print('Number of variables =', model.NumVari)

        # CONSTRAINTS
        # pattern assignment constraint
        for i in I:
            solver.Add(sum(x[i]) == 1)

        # worker efficiency
        for d in range(len(P[0])):

            solver.Add(sum(y[i][d] * a[i] for i in I) >= k * sumsum)


        # defining y
        for i in I:
            for d in range(len(P[0])):

                solver.Add(y[i][d] == sum(z[m][d] * x[i][m] for m in range(len(P))))
                #THIS SHOULD BE A SOFT CONSTRAINT INSTEAD OF A HARD CONSTRAINT BUT NEED MULTIPLE FEASIBLE SOLUTIONS

        # defining u : linearize!
        for i in I:
            for j in I:
                for m in range(len(P)):
                    # for m in range(20):
                    for n in range(len(P)):
                        # for n in range(20):
                        solver.Add(u[i][j][m][n] <= x[i][m])
                        solver.Add(u[i][j][m][n] <= x[j][n])
                        solver.Add(u[i][j][m][n] >= x[i][m] + x[j][n] - 1)

        # OBJECTIVE FUNCTION
        objective_terms = []
        alpha = 0.7
        for i in I:
            for j in I:
                for m in range(len(P)):
                    # for m in range(20):
                    for n in range(len(P)):
                        # for n in range(20):
                        # objective_terms.append(0)
                        objective_terms.append(alpha * pa[i][j] * c[m][n] * u[i][j][m][n])
                        objective_terms.append(- (1 - alpha) * pa[i][j] * w[m][n] * u[i][j][m][n])

        solver.Maximize(solver.Sum(objective_terms))  # MAXIMISING THE OBJECTIVE FUNCTION

        # Solve
        solver.SetTimeLimit(100000)#SETTING A TIME LIMIT OF 100 SECONDS
        status = solver.Solve()

        print("hae")

        # Print solution.
        if status == pywraplp.Solver.OPTIMAL or status == pywraplp.Solver.FEASIBLE:

            print('Optimal objective function = ', solver.Objective().Value(), '\n')
            print('Optimal assignment\n')
            for i in I:
                print('Schedule to worker ', i, ': ', int(sum(m * x[i][m].solution_value() for m in range(len(P)))))
                for m in range(len(P)):
                    if x[i][m].solution_value() == 1:
                        print(P[m])
                        runTime = time.time() - startTime
                        print("--- %s seconds ---" % (runTime))
            print("finding better solution")
            # if N<11 or k<0.33:
            # if N!=11 and N!=12:
            if N != 0:
                kiki = k + 0.01
                kfork = round(kiki, 3)
                self.hihi(kfork, N)


        else:
            print("Sorry we couldn't find any feasible or optimal solution.")




if __name__ == '__main__':

    for i in range(10):
        for N in range(5, 45):
            a = Me(0.3, N)
Improve this question
$\endgroup$
3
  • $\begingroup$ It might help if you explained what $K$ and $N$ are and add tags for the language and solver you are using. $\endgroup$
    – prubin ♦
    Commented Jul 1, 2023 at 22:00
  • $\begingroup$ Please clarify your specific problem or provide additional details to highlight exactly what you need. As it's currently written, it's hard to tell exactly what you're asking. $\endgroup$
    – Community Bot
    Commented Jul 1, 2023 at 22:59
  • $\begingroup$ @prubin Thanks for the suggestion! I have changed it according to your suggestion $\endgroup$
    – abdula smith
    Commented Jul 3, 2023 at 4:41

1 Answer 1

Reset to default
1
$\begingroup$

I don't see anything wrong following normal dist for distributing efficiency level. As for the outcomes, more employees with higher efficiency means more choices with more slack on constraints. As for half of them being married it should mean married pairs, otherwise all would be married- fewer choices.

Just confirming if you are using a binary variable $u_{ij}^d$ to count married couple assignment with constraints

$ p_dx_{i,p}+p'_{d}x_{j,p'} \le u_{ij}^d + 1 \ \ \forall p,p' \in P$
$ u_{ij}^d \le p_dx_{i,p}$
$ u_{ij}^d \le p_dx_{j,p} \ \ \forall p \ \ \forall d \ \ \forall (i,j) \in $ Couples

And then minimizing $\sum_d \sum_{(ij)}u_{ij}^d$

Improve this answer
$\endgroup$
1
  • $\begingroup$ I have edited my post . Models and notations are in the above post. Thanks for the help! $\endgroup$
    – abdula smith
    Commented Jul 3, 2023 at 4:45

Not the answer you're looking for? Browse other questions tagged or ask your own question.