From f0db1a3cc1cb834bf5345a56e9bbc7e2f71a3199 Mon Sep 17 00:00:00 2001 From: Joao-Dionisio Date: Tue, 8 Sep 2026 15:54:22 +0100 Subject: [PATCH 1/3] Add heur_lpt.py example: warm-start a scheduling MIP from a Heur plugin --- CHANGELOG.md | 1 + examples/finished/heur_lpt.py | 192 ++++++++++++++++++++++++++++++++++ 2 files changed, 193 insertions(+) create mode 100644 examples/finished/heur_lpt.py diff --git a/CHANGELOG.md b/CHANGELOG.md index ae37c3270..a3c07dfaf 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -2,6 +2,7 @@ ## Unreleased ### Added +- Added example `heur_lpt.py`: warm-starting a scheduling MIP with a `Heur` plugin that builds its solution in the original space - Added methods: `getNNodesLeft()`, `getNRuns()`, `getNReoptRuns()`, `addNNodes()` with tests - Added `addConsCumulative()` for SCIP cumulative constraints (#1222) - `Expr` and `GenExpr` support `__pos__` magic method like `+Expr` or `+GenExpr` diff --git a/examples/finished/heur_lpt.py b/examples/finished/heur_lpt.py new file mode 100644 index 000000000..1409d6687 --- /dev/null +++ b/examples/finished/heur_lpt.py @@ -0,0 +1,192 @@ +""" +Example showing a custom primal heuristic using PySCIPOpt's Heur plugin. + +The heuristic warm-starts a scheduling MIP: before SCIP processes the root +node, it builds a feasible schedule with a fast list-scheduling rule and +hands it to SCIP as an incumbent solution. + +The problem is parallel machine scheduling with release dates and makespan +objective, Pm|r_j|C_max: n jobs with processing times p_j and release dates +r_j have to be assigned to m identical machines and sequenced so that the +last job finishes as early as possible. The MIP is a disjunctive (big-M) +formulation. The heuristic is the LPT (longest processing time first) list +scheduling rule: whenever a machine becomes free, it starts the longest job +that has already been released. +""" + +from pyscipopt import Model, Heur, SCIP_RESULT, SCIP_HEURTIMING, quicksum + + +def lpt_schedule(p, r, m): + """ + LPT list scheduling for Pm|r_j|C_max. + + Whenever a machine becomes free, start the longest released job on it. + If no job is released yet, wait for the next release. + + Returns a dict job -> (machine, start time) and the makespan. + """ + unscheduled = set(range(len(p))) + free_at = [0] * m + schedule = {} + + while unscheduled: + machine = min(range(m), key=lambda k: free_at[k]) + now = free_at[machine] + + released = [j for j in unscheduled if r[j] <= now] + if not released: + now = min(r[j] for j in unscheduled) + released = [j for j in unscheduled if r[j] <= now] + + job = max(released, key=lambda j: p[j]) + schedule[job] = (machine, now) + free_at[machine] = now + p[job] + unscheduled.remove(job) + + return schedule, max(free_at) + + +class LPTHeur(Heur): + """ + Primal heuristic that offers the LPT schedule to SCIP as a solution. + + It needs the problem data and the model's variables to translate the + schedule into variable values. + """ + + def __init__(self, p, r, m, start, assign, before, makespan): + super().__init__() + self.p = p + self.r = r + self.m = m + self.start = start + self.assign = assign + self.before = before + self.makespan = makespan + self.done = False + + def heurexec(self, heurtiming, nodeinfeasible): + # The schedule does not depend on the search state, so one run is + # enough. SCIP processes the root node again after a restart, which + # would call the heuristic a second time otherwise. + if self.done: + return {"result": SCIP_RESULT.DIDNOTRUN} + self.done = True + + schedule, cmax = lpt_schedule(self.p, self.r, self.m) + + # The solution is built in the original space: after presolving, SCIP + # may have fixed or aggregated variables (here symmetry handling fixes + # the machine of some jobs), and setting a conflicting value on such a + # variable in a transformed solution is an error. Passing the heuristic + # to createOrigSol tells SCIP who found the solution, so it shows up + # under this heuristic's display character in the log. + sol = self.model.createOrigSol(self) + + sol[self.makespan] = cmax + for j, (machine, st) in schedule.items(): + sol[self.start[j]] = st + for k in range(self.m): + sol[self.assign[j, k]] = 1 if k == machine else 0 + + # Jobs on the same machine never overlap, so the one starting first + # also finishes before the other starts. + for i, (machine_i, st_i) in schedule.items(): + for j, (machine_j, st_j) in schedule.items(): + if i != j: + same_machine = machine_i == machine_j + sol[self.before[i, j]] = 1 if same_machine and st_i < st_j else 0 + + # trySol checks feasibility and stores the solution if it is accepted. + accepted = self.model.trySol(sol) + print(f"LPT heuristic: makespan {cmax} {'accepted' if accepted else 'rejected'}") + + if accepted: + return {"result": SCIP_RESULT.FOUNDSOL} + return {"result": SCIP_RESULT.DIDNOTFIND} + + +def build_model(p, r, m): + """ + Disjunctive big-M model for Pm|r_j|C_max. + + Variables: + start[j] - start time of job j + assign[j, k] - 1 if job j runs on machine k + before[i, j] - 1 if job i finishes before job j starts + makespan - completion time of the last job + + Returns the model and its variables. + """ + n = len(p) + jobs = range(n) + machines = range(m) + + # Running all jobs on one machine after the last release is always feasible, + # so no job needs to start after this horizon. + horizon = max(r) + sum(p) + + model = Model("Pm|r_j|C_max") + + start = {j: model.addVar(vtype="C", lb=r[j], ub=horizon - p[j], name=f"start_{j}") for j in jobs} + assign = {(j, k): model.addVar(vtype="B", name=f"assign_{j}_{k}") for j in jobs for k in machines} + before = {(i, j): model.addVar(vtype="B", name=f"before_{i}_{j}") for i in jobs for j in jobs if i != j} + makespan = model.addVar(vtype="C", ub=horizon, name="makespan") + + for j in jobs: + model.addCons(quicksum(assign[j, k] for k in machines) == 1) + model.addCons(start[j] + p[j] <= makespan) + + for i in jobs: + for j in jobs: + if i != j: + # if i is sequenced before j, j cannot start until i is done + model.addCons(start[i] + p[i] <= start[j] + horizon * (1 - before[i, j])) + + for i in jobs: + for j in jobs: + if i < j: + for k in machines: + # two jobs on the same machine have to be sequenced + model.addCons(assign[i, k] + assign[j, k] - 1 <= before[i, j] + before[j, i]) + + model.setObjective(makespan, "minimize") + + return model, start, assign, before, makespan + + +def print_schedule(p, schedule, m): + for k in range(m): + on_k = sorted((st, j) for j, (machine, st) in schedule.items() if machine == k) + jobs = " ".join(f"job {j} [{st}, {st + p[j]})" for st, j in on_k) + print(f" machine {k}: {jobs}") + + +if __name__ == "__main__": + # On this instance the greedy LPT schedule is not optimal, so the log shows + # SCIP starting from the heuristic's incumbent and improving on it. + p = [2, 2, 2, 9, 8, 7, 4, 5, 3, 6] + r = [0, 0, 0, 1, 1, 1, 3, 5, 2, 4] + m = 3 + + schedule, cmax = lpt_schedule(p, r, m) + print(f"LPT schedule, makespan {cmax}:") + print_schedule(p, schedule, m) + + model, start, assign, before, makespan = build_model(p, r, m) + + heur = LPTHeur(p, r, m, start, assign, before, makespan) + # freq=0: only call the heuristic at depth 0, i.e. at the root node. + model.includeHeur(heur, "lpt", "LPT list scheduling warm start", "L", + freq=0, timingmask=SCIP_HEURTIMING.BEFORENODE) + + model.optimize() + + print(f"\nSCIP status: {model.getStatus()}") + print(f"optimal schedule, makespan {model.getObjVal():g}:") + best = {} + for j in range(len(p)): + machine = next(k for k in range(m) if model.getVal(assign[j, k]) > 0.5) + best[j] = (machine, round(model.getVal(start[j]))) + print_schedule(p, best, m) From fcf00f8d4724a8e63330af98bfec670f27dd917a Mon Sep 17 00:00:00 2001 From: Joao-Dionisio Date: Tue, 8 Sep 2026 16:41:05 +0100 Subject: [PATCH 2/3] Add scheduling_logical.py example: scheduling with and/or/indicator, disjunction and cardinality constraints --- CHANGELOG.md | 1 + examples/finished/scheduling_logical.py | 184 ++++++++++++++++++++++++ 2 files changed, 185 insertions(+) create mode 100644 examples/finished/scheduling_logical.py diff --git a/CHANGELOG.md b/CHANGELOG.md index a3c07dfaf..c56c027c0 100644 --- a/CHANGELOG.md +++ b/CHANGELOG.md @@ -3,6 +3,7 @@ ## Unreleased ### Added - Added example `heur_lpt.py`: warm-starting a scheduling MIP with a `Heur` plugin that builds its solution in the original space +- Added example `scheduling_logical.py`: the same scheduling problem modelled with and/or/indicator, disjunction and cardinality constraints - Added methods: `getNNodesLeft()`, `getNRuns()`, `getNReoptRuns()`, `addNNodes()` with tests - Added `addConsCumulative()` for SCIP cumulative constraints (#1222) - `Expr` and `GenExpr` support `__pos__` magic method like `+Expr` or `+GenExpr` diff --git a/examples/finished/scheduling_logical.py b/examples/finished/scheduling_logical.py new file mode 100644 index 000000000..7ac1e06b1 --- /dev/null +++ b/examples/finished/scheduling_logical.py @@ -0,0 +1,184 @@ +""" +Example of modelling with SCIP's logical constraint handlers. + +The problem is parallel machine scheduling with release dates and makespan +objective, Pm|r_j|C_max: n jobs with processing times p_j and release dates +r_j have to be assigned to m identical machines and sequenced so that the +last job finishes as early as possible. Jobs on the same machine must not +overlap, and that condition is where the formulations differ. Each one +states it with a different constraint type instead of big-M constraints: + + and_or_indicator_model and, or and indicator constraints + disjunction_model disjunction constraints + time_indexed_model cardinality constraints + +All three are solved on the same instance and reach the same makespan. +They differ a lot in how much branching SCIP needs, because the constraint +types differ in what they contribute to the LP relaxation. The disjunction +formulation shows the extreme case: a disjunction is enforced by branching +alone, so the model needs an extra linear constraint to solve in reasonable +time. +""" + +from pyscipopt import Model, quicksum + + +def base_model(p, r, m, name): + """ + Variables and constraints shared by all formulations. + + Every job is assigned to exactly one machine, starts after its release + date and finishes before the makespan. What is missing is the condition + that jobs on the same machine do not overlap. + """ + jobs = range(len(p)) + machines = range(m) + + # Running all jobs on one machine after the last release is always feasible, + # so no job needs to start after this horizon. + horizon = max(r) + sum(p) + + model = Model(name) + + start = {j: model.addVar(vtype="C", lb=r[j], ub=horizon - p[j], name=f"start_{j}") for j in jobs} + assign = {(j, k): model.addVar(vtype="B", name=f"assign_{j}_{k}") for j in jobs for k in machines} + makespan = model.addVar(vtype="C", ub=horizon, name="makespan") + + for j in jobs: + # Exactly one machine per job. An xor constraint would be wrong here: + # it fixes the parity of the number of true variables, so with three + # machines it would also allow assigning a job to all three. + model.addCons(quicksum(assign[j, k] for k in machines) == 1) + model.addCons(start[j] + p[j] <= makespan) + + model.setObjective(makespan, "minimize") + + return model, start, assign, makespan + + +def and_or_indicator_model(p, r, m): + """ + before[i, j] = 1 means job i finishes before job j starts, enforced by an + indicator constraint. For every pair of jobs and every machine, an and + constraint detects that both jobs run on that machine, an or constraint + detects that the pair is sequenced one way or the other, and an indicator + constraint requires the second whenever the first holds. + """ + model, start, assign, makespan = base_model(p, r, m, "and-or-indicator") + jobs = range(len(p)) + machines = range(m) + + before = {(i, j): model.addVar(vtype="B", name=f"before_{i}_{j}") for i in jobs for j in jobs if i != j} + + for i in jobs: + for j in jobs: + if i != j: + model.addConsIndicator(start[i] + p[i] <= start[j], binvar=before[i, j], name=f"seq_{i}_{j}") + + for i in jobs: + for j in jobs: + if i < j: + # ordered = before[i, j] or before[j, i] + ordered = model.addVar(vtype="B", name=f"ordered_{i}_{j}") + model.addConsOr([before[i, j], before[j, i]], ordered) + for k in machines: + # same = assign[i, k] and assign[j, k] + same = model.addVar(vtype="B", name=f"same_{i}_{j}_{k}") + model.addConsAnd([assign[i, k], assign[j, k]], same) + # same machine => the pair is ordered + model.addConsIndicator(ordered >= 1, binvar=same, name=f"ordered_if_same_{i}_{j}_{k}") + + return model, start, assign, makespan + + +def disjunction_model(p, r, m): + """ + For every pair of jobs and every machine, a disjunction constraint states + that one of the following holds: job i is not on the machine, job j is not + on the machine, i finishes before j starts, or j finishes before i starts. + + No auxiliary binary variables are needed for the sequencing, but SCIP + enforces a disjunction by branching on its members only, so the LP + relaxation knows nothing about the non-overlap condition. On its own the + model makes SCIP enumerate schedules. The load constraint at the end is a + weaker linear consequence of non-overlap that gives the LP a useful bound; + remove it to see the difference. + """ + model, start, assign, makespan = base_model(p, r, m, "disjunction") + jobs = range(len(p)) + machines = range(m) + + for i in jobs: + for j in jobs: + if i < j: + for k in machines: + model.addConsDisjunction( + [assign[i, k] <= 0, assign[j, k] <= 0, + start[i] + p[i] <= start[j], start[j] + p[j] <= start[i]], + name=f"no_overlap_{i}_{j}_{k}", + ) + + for k in machines: + # jobs on one machine run one after the other, so their total + # processing time is a lower bound on the makespan + model.addCons(quicksum(p[j] * assign[j, k] for j in jobs) <= makespan, name=f"load_{k}") + + return model, start, assign, makespan + + +def time_indexed_model(p, r, m): + """ + x[j, k, t] = 1 means job j starts on machine k at time t. Each job starts + exactly once, and a cardinality constraint per machine and time step + allows at most one of the jobs that would be running then to be started. + Start times and machine assignments are recovered linearly from x. + """ + model, start, assign, makespan = base_model(p, r, m, "time-indexed") + jobs = range(len(p)) + machines = range(m) + horizon = max(r) + sum(p) + + # possible start times of each job + slots = {j: range(r[j], horizon - p[j] + 1) for j in jobs} + + x = {} + for j in jobs: + for k in machines: + for t in slots[j]: + x[j, k, t] = model.addVar(vtype="B", name=f"x_{j}_{k}_{t}") + + for j in jobs: + model.addCons(quicksum(x[j, k, t] for k in machines for t in slots[j]) == 1) + model.addCons(start[j] == quicksum(t * x[j, k, t] for k in machines for t in slots[j])) + for k in machines: + model.addCons(assign[j, k] == quicksum(x[j, k, t] for t in slots[j])) + + for k in machines: + for t in range(horizon): + # jobs that occupy machine k at time t if started at s + running = [x[j, k, s] for j in jobs for s in range(t - p[j] + 1, t + 1) if (j, k, s) in x] + if len(running) > 1: + model.addConsCardinality(running, 1, name=f"one_job_{k}_{t}") + + return model, start, assign, makespan + + +def print_schedule(p, m, start, assign, model): + for k in range(m): + on_k = sorted((round(model.getVal(start[j])), j) for j in range(len(p)) if model.getVal(assign[j, k]) > 0.5) + jobs = " ".join(f"job {j} [{st}, {st + p[j]})" for st, j in on_k) + print(f" machine {k}: {jobs}") + + +if __name__ == "__main__": + p = [2, 2, 2, 9, 8, 7, 4, 5, 3, 6] + r = [0, 0, 0, 1, 1, 1, 3, 5, 2, 4] + m = 3 + + for build in (and_or_indicator_model, disjunction_model, time_indexed_model): + model, start, assign, makespan = build(p, r, m) + model.hideOutput() + model.optimize() + print(f"{model.getProbName()}: {model.getStatus()}, makespan {model.getObjVal():g}, " + f"{model.getNNodes()} nodes, {model.getSolvingTime():.2f}s") + print_schedule(p, m, start, assign, model) From 1642a61c6d6ff5656d0ded1608c45641fe8006e8 Mon Sep 17 00:00:00 2001 From: Joao-Dionisio Date: Tue, 8 Sep 2026 16:48:16 +0100 Subject: [PATCH 3/3] Shorten comments in the scheduling examples --- examples/finished/heur_lpt.py | 79 +++++++------------------ examples/finished/scheduling_logical.py | 78 +++++++----------------- 2 files changed, 44 insertions(+), 113 deletions(-) diff --git a/examples/finished/heur_lpt.py b/examples/finished/heur_lpt.py index 1409d6687..30e66c01b 100644 --- a/examples/finished/heur_lpt.py +++ b/examples/finished/heur_lpt.py @@ -1,17 +1,10 @@ """ -Example showing a custom primal heuristic using PySCIPOpt's Heur plugin. - -The heuristic warm-starts a scheduling MIP: before SCIP processes the root -node, it builds a feasible schedule with a fast list-scheduling rule and -hands it to SCIP as an incumbent solution. - -The problem is parallel machine scheduling with release dates and makespan -objective, Pm|r_j|C_max: n jobs with processing times p_j and release dates -r_j have to be assigned to m identical machines and sequenced so that the -last job finishes as early as possible. The MIP is a disjunctive (big-M) -formulation. The heuristic is the LPT (longest processing time first) list -scheduling rule: whenever a machine becomes free, it starts the longest job -that has already been released. +Parallel machine scheduling with release dates (Pm|r_j|C_max), solved with a +big-M formulation that is warm-started by a custom heuristic. + +The heuristic is a Heur plugin that runs once before the root node. It builds +an LPT (longest processing time first) list schedule and hands it to SCIP as +an incumbent. """ from pyscipopt import Model, Heur, SCIP_RESULT, SCIP_HEURTIMING, quicksum @@ -19,12 +12,8 @@ def lpt_schedule(p, r, m): """ - LPT list scheduling for Pm|r_j|C_max. - - Whenever a machine becomes free, start the longest released job on it. - If no job is released yet, wait for the next release. - - Returns a dict job -> (machine, start time) and the makespan. + LPT list scheduling: whenever a machine becomes free, start the longest + released job on it. Returns {job: (machine, start)} and the makespan. """ unscheduled = set(range(len(p))) free_at = [0] * m @@ -48,12 +37,7 @@ def lpt_schedule(p, r, m): class LPTHeur(Heur): - """ - Primal heuristic that offers the LPT schedule to SCIP as a solution. - - It needs the problem data and the model's variables to translate the - schedule into variable values. - """ + """Offers the LPT schedule to SCIP as a solution.""" def __init__(self, p, r, m, start, assign, before, makespan): super().__init__() @@ -67,21 +51,16 @@ def __init__(self, p, r, m, start, assign, before, makespan): self.done = False def heurexec(self, heurtiming, nodeinfeasible): - # The schedule does not depend on the search state, so one run is - # enough. SCIP processes the root node again after a restart, which - # would call the heuristic a second time otherwise. + # run once; SCIP processes the root again after a restart if self.done: return {"result": SCIP_RESULT.DIDNOTRUN} self.done = True schedule, cmax = lpt_schedule(self.p, self.r, self.m) - # The solution is built in the original space: after presolving, SCIP - # may have fixed or aggregated variables (here symmetry handling fixes - # the machine of some jobs), and setting a conflicting value on such a - # variable in a transformed solution is an error. Passing the heuristic - # to createOrigSol tells SCIP who found the solution, so it shows up - # under this heuristic's display character in the log. + # Build the solution in the original space. Presolving may have fixed + # some variables (symmetry handling does so here), and setting a + # conflicting value on a fixed variable of a transformed solution fails. sol = self.model.createOrigSol(self) sol[self.makespan] = cmax @@ -90,15 +69,12 @@ def heurexec(self, heurtiming, nodeinfeasible): for k in range(self.m): sol[self.assign[j, k]] = 1 if k == machine else 0 - # Jobs on the same machine never overlap, so the one starting first - # also finishes before the other starts. for i, (machine_i, st_i) in schedule.items(): for j, (machine_j, st_j) in schedule.items(): if i != j: same_machine = machine_i == machine_j sol[self.before[i, j]] = 1 if same_machine and st_i < st_j else 0 - # trySol checks feasibility and stores the solution if it is accepted. accepted = self.model.trySol(sol) print(f"LPT heuristic: makespan {cmax} {'accepted' if accepted else 'rejected'}") @@ -109,23 +85,16 @@ def heurexec(self, heurtiming, nodeinfeasible): def build_model(p, r, m): """ - Disjunctive big-M model for Pm|r_j|C_max. + Disjunctive big-M model. - Variables: - start[j] - start time of job j - assign[j, k] - 1 if job j runs on machine k - before[i, j] - 1 if job i finishes before job j starts - makespan - completion time of the last job - - Returns the model and its variables. + start[j] - start time of job j + assign[j, k] - 1 if job j runs on machine k + before[i, j] - 1 if job i finishes before job j starts + makespan - completion time of the last job """ - n = len(p) - jobs = range(n) + jobs = range(len(p)) machines = range(m) - - # Running all jobs on one machine after the last release is always feasible, - # so no job needs to start after this horizon. - horizon = max(r) + sum(p) + horizon = max(r) + sum(p) # no job needs to start later than this model = Model("Pm|r_j|C_max") @@ -141,14 +110,13 @@ def build_model(p, r, m): for i in jobs: for j in jobs: if i != j: - # if i is sequenced before j, j cannot start until i is done model.addCons(start[i] + p[i] <= start[j] + horizon * (1 - before[i, j])) + # jobs on the same machine have to be sequenced for i in jobs: for j in jobs: if i < j: for k in machines: - # two jobs on the same machine have to be sequenced model.addCons(assign[i, k] + assign[j, k] - 1 <= before[i, j] + before[j, i]) model.setObjective(makespan, "minimize") @@ -164,8 +132,6 @@ def print_schedule(p, schedule, m): if __name__ == "__main__": - # On this instance the greedy LPT schedule is not optimal, so the log shows - # SCIP starting from the heuristic's incumbent and improving on it. p = [2, 2, 2, 9, 8, 7, 4, 5, 3, 6] r = [0, 0, 0, 1, 1, 1, 3, 5, 2, 4] m = 3 @@ -177,9 +143,8 @@ def print_schedule(p, schedule, m): model, start, assign, before, makespan = build_model(p, r, m) heur = LPTHeur(p, r, m, start, assign, before, makespan) - # freq=0: only call the heuristic at depth 0, i.e. at the root node. model.includeHeur(heur, "lpt", "LPT list scheduling warm start", "L", - freq=0, timingmask=SCIP_HEURTIMING.BEFORENODE) + freq=0, timingmask=SCIP_HEURTIMING.BEFORENODE) # freq=0: root node only model.optimize() diff --git a/examples/finished/scheduling_logical.py b/examples/finished/scheduling_logical.py index 7ac1e06b1..29e3a08c4 100644 --- a/examples/finished/scheduling_logical.py +++ b/examples/finished/scheduling_logical.py @@ -1,42 +1,23 @@ """ -Example of modelling with SCIP's logical constraint handlers. - -The problem is parallel machine scheduling with release dates and makespan -objective, Pm|r_j|C_max: n jobs with processing times p_j and release dates -r_j have to be assigned to m identical machines and sequenced so that the -last job finishes as early as possible. Jobs on the same machine must not -overlap, and that condition is where the formulations differ. Each one -states it with a different constraint type instead of big-M constraints: +Parallel machine scheduling with release dates (Pm|r_j|C_max), modelled three +ways with SCIP's logical constraints instead of big-M constraints: and_or_indicator_model and, or and indicator constraints disjunction_model disjunction constraints time_indexed_model cardinality constraints -All three are solved on the same instance and reach the same makespan. -They differ a lot in how much branching SCIP needs, because the constraint -types differ in what they contribute to the LP relaxation. The disjunction -formulation shows the extreme case: a disjunction is enforced by branching -alone, so the model needs an extra linear constraint to solve in reasonable -time. +The formulations differ in what they contribute to the LP relaxation, which +shows in the number of nodes SCIP needs. """ from pyscipopt import Model, quicksum def base_model(p, r, m, name): - """ - Variables and constraints shared by all formulations. - - Every job is assigned to exactly one machine, starts after its release - date and finishes before the makespan. What is missing is the condition - that jobs on the same machine do not overlap. - """ + """Everything except the condition that jobs on one machine do not overlap.""" jobs = range(len(p)) machines = range(m) - - # Running all jobs on one machine after the last release is always feasible, - # so no job needs to start after this horizon. - horizon = max(r) + sum(p) + horizon = max(r) + sum(p) # no job needs to start later than this model = Model(name) @@ -45,9 +26,8 @@ def base_model(p, r, m, name): makespan = model.addVar(vtype="C", ub=horizon, name="makespan") for j in jobs: - # Exactly one machine per job. An xor constraint would be wrong here: - # it fixes the parity of the number of true variables, so with three - # machines it would also allow assigning a job to all three. + # not an xor constraint: xor fixes the parity, so with three machines + # a job could be assigned to all of them model.addCons(quicksum(assign[j, k] for k in machines) == 1) model.addCons(start[j] + p[j] <= makespan) @@ -58,11 +38,9 @@ def base_model(p, r, m, name): def and_or_indicator_model(p, r, m): """ - before[i, j] = 1 means job i finishes before job j starts, enforced by an - indicator constraint. For every pair of jobs and every machine, an and - constraint detects that both jobs run on that machine, an or constraint - detects that the pair is sequenced one way or the other, and an indicator - constraint requires the second whenever the first holds. + before[i, j] = 1 forces j to start after i is done (indicator). For each + pair and machine: same = assign[i, k] and assign[j, k], ordered = + before[i, j] or before[j, i], and same implies ordered (indicator). """ model, start, assign, makespan = base_model(p, r, m, "and-or-indicator") jobs = range(len(p)) @@ -78,14 +56,11 @@ def and_or_indicator_model(p, r, m): for i in jobs: for j in jobs: if i < j: - # ordered = before[i, j] or before[j, i] ordered = model.addVar(vtype="B", name=f"ordered_{i}_{j}") model.addConsOr([before[i, j], before[j, i]], ordered) for k in machines: - # same = assign[i, k] and assign[j, k] same = model.addVar(vtype="B", name=f"same_{i}_{j}_{k}") model.addConsAnd([assign[i, k], assign[j, k]], same) - # same machine => the pair is ordered model.addConsIndicator(ordered >= 1, binvar=same, name=f"ordered_if_same_{i}_{j}_{k}") return model, start, assign, makespan @@ -93,16 +68,12 @@ def and_or_indicator_model(p, r, m): def disjunction_model(p, r, m): """ - For every pair of jobs and every machine, a disjunction constraint states - that one of the following holds: job i is not on the machine, job j is not - on the machine, i finishes before j starts, or j finishes before i starts. - - No auxiliary binary variables are needed for the sequencing, but SCIP - enforces a disjunction by branching on its members only, so the LP - relaxation knows nothing about the non-overlap condition. On its own the - model makes SCIP enumerate schedules. The load constraint at the end is a - weaker linear consequence of non-overlap that gives the LP a useful bound; - remove it to see the difference. + For each pair and machine: i is not on the machine, or j is not, or i + finishes before j starts, or j finishes before i starts. + + A disjunction is enforced by branching only and adds nothing to the LP + relaxation, so without the load constraint at the end SCIP enumerates + schedules. Try removing it. """ model, start, assign, makespan = base_model(p, r, m, "disjunction") jobs = range(len(p)) @@ -119,8 +90,7 @@ def disjunction_model(p, r, m): ) for k in machines: - # jobs on one machine run one after the other, so their total - # processing time is a lower bound on the makespan + # total processing time on a machine is a lower bound on the makespan model.addCons(quicksum(p[j] * assign[j, k] for j in jobs) <= makespan, name=f"load_{k}") return model, start, assign, makespan @@ -128,18 +98,14 @@ def disjunction_model(p, r, m): def time_indexed_model(p, r, m): """ - x[j, k, t] = 1 means job j starts on machine k at time t. Each job starts - exactly once, and a cardinality constraint per machine and time step - allows at most one of the jobs that would be running then to be started. - Start times and machine assignments are recovered linearly from x. + x[j, k, t] = 1 if job j starts on machine k at time t. A cardinality + constraint per machine and time step allows at most one running job. """ model, start, assign, makespan = base_model(p, r, m, "time-indexed") jobs = range(len(p)) machines = range(m) horizon = max(r) + sum(p) - - # possible start times of each job - slots = {j: range(r[j], horizon - p[j] + 1) for j in jobs} + slots = {j: range(r[j], horizon - p[j] + 1) for j in jobs} # possible start times x = {} for j in jobs: @@ -155,7 +121,7 @@ def time_indexed_model(p, r, m): for k in machines: for t in range(horizon): - # jobs that occupy machine k at time t if started at s + # jobs that would be running on machine k at time t running = [x[j, k, s] for j in jobs for s in range(t - p[j] + 1, t + 1) if (j, k, s) in x] if len(running) > 1: model.addConsCardinality(running, 1, name=f"one_job_{k}_{t}")