Problem RCPSP
From CSPLib: ``The Resource-Constrained Project Scheduling Problem is a classical problem in operations research. A number of activities are to be scheduled. Each activity has a duration and cannot be interrupted. There are a set of precedence relations between pairs of activities which state that the second activity must start after the first has finished. There are a set of renewable resources. Each resource has a maximum capacity and at any given time slot no more than this amount can be in use. Each activity has a demand (possibly zero) on each resource. The problem is usually stated as an optimization problem where the makespan (i.e., the completion time of the last activity) is minimized.’’ See CSPLib (Problem 61) for more information.
To build a COP (Constraint Optimization Problem) model, we need first to import the library PyCSP$^3$:
from pycsp3 import *
Then, we need some data. Here the data correspond to those in file ‘Rcpsp_j30-01-01’ that can be found in PSBLIB. We have 158 as horizon, four resources with specified quantities, and 30+2 jobs with specified durations, successors and resource requirements.
horizon = 158
capacities = [12, 13, 4, 12]
# below, job durations, successors and quantities (requirements)
durations = [0, 8, 4, 6, 3, 8, 5, 9, 2, 7, 9, 2, 6, 3, 9, 10, 6, 5, 3, 7, 2, 7, 2, 3, 3, 7, 8, 3, 7, 2, 2, 0]
successors = [[1, 2, 3], [5, 10, 14], [6, 7, 12], [4, 8, 9], [19], [29], [26], [11, 18, 26], [13],
[15, 24], [19, 25], [13], [16, 17], [16], [24], [20, 21], [21], [19, 21], [23, 28], [22, 24], [27],
[22], [23], [29], [29], [30], [27], [30], [31], [31], [31], []]
quantities = [[0, 0, 0, 0], [4, 0, 0, 0], [10, 0, 0, 0], [0, 0, 0, 3], [3, 0, 0, 0], [0, 0, 0, 8],
[4, 0, 0, 0], [0, 1, 0, 0], [6, 0, 0, 0], [0, 0, 0, 1], [0, 5, 0, 0], [0, 7, 0, 0], [4, 0, 0, 0],
[0, 8, 0, 0], [3, 0, 0, 0], [0, 0, 0, 5], [0, 0, 0, 8], [0, 0, 0, 7], [0, 1, 0, 0], [0, 10, 0, 0],
[0, 0, 0, 6], [2, 0, 0, 0], [3, 0, 0, 0], [0, 9, 0, 0], [4, 0, 0, 0], [0, 0, 4, 0], [0, 0, 0, 7],
[0, 8, 0, 0], [0, 7, 0, 0], [0, 7, 0, 0], [0, 0, 2, 0], [0, 0, 0, 0]]
nJobs = len(durations)
We start our COP model by introducing an array $x$ of variables, one per job. Note that the first job must start at 0.
# x[i] is the starting time of the ith job
x = VarArray(size=nJobs, dom=lambda i: {0} if i == 0 else range(horizon))
We can display the structure of the array, as well as the domain of the 2 first variables (note that all variables have the same domain, except for the first one).
print("Array x: ", x)
print("Domain of the first variable: ", x[0].dom)
print("Domain of the other variables: ", x[1].dom)
Array x: [x[0], x[1], x[2], x[3], x[4], x[5], x[6], x[7], x[8], x[9], x[10], x[11], x[12], x[13], x[14], x[15], x[16], x[17], x[18], x[19], x[20], x[21], x[22], x[23], x[24], x[25], x[26], x[27], x[28], x[29], x[30], x[31]]
Domain of the first variable: 0
Domain of the other variables: 0..157
Concerning the constraints, we first post a group of constraints Intension to ensure that precedence relations are respected.
satisfy(
# precedence constraints
x[i] + durations[i] <= x[j] for i in range(nJobs) for j in successors[i]
);
We can display the 8 first posted constraints, to control that everything is ok. Although a little bit technical, we can see that the expressions are correct (note that ‘le’ stands for ‘less than or equal to’).
print(posted()[:8])
intension(function:le(x[0],x[1]))
intension(function:le(x[0],x[2]))
intension(function:le(x[0],x[3]))
intension(function:le(add(x[1],8),x[5]))
intension(function:le(add(x[1],8),x[10]))
intension(function:le(add(x[1],8),x[14]))
intension(function:le(add(x[2],4),x[6]))
intension(function:le(add(x[2],4),x[7]))
To minimize the makespan, it suffices to post the following objective:
minimize(
x[-1]
);
We can run the solver, with this optimization task. Note that we need to check that the status returned by the solver is OPTIMUM.
if solve() is OPTIMUM:
print("Starting times of jobs: ", values(x))
print("Makespan: ", bound())
Starting times of jobs: [0, 0, 0, 0, 6, 8, 4, 4, 6, 6, 8, 13, 4, 15, 8, 13, 18, 10, 13, 17, 23, 24, 31, 33, 24, 17, 13, 25, 16, 36, 28, 38]
Makespan: 38
We obtain a very small value for the makespan, but the problem is that we have not taken resource capacities into account.
For each ressource, we need to post a constraint Cumulative that imposes a strict limit of consumption:
satisfy(
# resource constraints
Cumulative(
Task(origin=x[i], length=durations[i], height=quantities[i][k]) for i in range(nJobs) if quantities[i][k] > 0
) <= capacity for k, capacity in enumerate(capacities)
);
By calling posted(-1) we can get the last posted constraints, i.e., those that have been posted during the last call to satisfy(). One can then check that they are correctly defined.
print(posted(-1))
cumulative(origins:[x[1], x[2], x[4], x[6], x[8], x[12], x[14], x[21], x[22], x[24]], lengths:[8, 4, 3, 5, 2, 6, 9, 7, 2, 3], heights:[4, 10, 3, 4, 6, 4, 3, 2, 3, 4], condition:(le,12))
cumulative(origins:[x[7], x[10], x[11], x[13], x[18], x[19], x[23], x[27], x[28], x[29]], lengths:[9, 9, 2, 3, 3, 7, 3, 3, 7, 2], heights:[1, 5, 7, 8, 1, 10, 9, 8, 7, 7], condition:(le,13))
cumulative(origins:[x[25], x[30]], lengths:[7, 2], heights:[4, 2], condition:(le,4))
cumulative(origins:[x[3], x[5], x[9], x[15], x[16], x[17], x[20], x[26]], lengths:[6, 8, 7, 10, 6, 5, 2, 8], heights:[3, 8, 1, 5, 8, 7, 6, 7], condition:(le,12))
We can run again the solver.
if solve() is OPTIMUM:
print("Starting times of jobs: ", values(x))
print("Makespan: ", bound())
Starting times of jobs: [0, 4, 0, 0, 12, 31, 4, 4, 10, 6, 12, 13, 4, 16, 15, 13, 23, 10, 13, 26, 29, 29, 36, 38, 33, 21, 15, 34, 19, 41, 37, 43]
Makespan: 43
This time, we obtain a valid solution respecting both precedence and capacity constraints.
Finally, we give below the model in one piece. Here, the data is expected to be given by the user (in a command line).
from pycsp3 import *
horizon, capacities, jobs = data
durations, successors, quantities = zip(*jobs)
nJobs = len(jobs)
# x[i] is the starting time of the ith job
x = VarArray(size=nJobs, dom=lambda i: {0} if i == 0 else range(horizon))
satisfy(
# precedence constraints
[x[i] + durations[i] <= x[j] for i in range(nJobs) for j in successors[i]],
# resource constraints
[
Cumulative(
Task(origin=x[i], length=durations[i], height=quantities[i][k]) for i in range(nJobs) if quantities[i][k] > 0
) <= capacity for k, capacity in enumerate(capacities)
]
)
minimize(
x[-1]
)