Extension on examples

[abb-omidi]

1. An extension on the heuristic callback

from pyscipopt import Model, quicksum, SCIP_PARAMSETTING, Heur, SCIP_RESULT, SCIP_HEURTIMING
import matplotlib.pyplot as plt
import matplotlib.patches as patches
import seaborn as sns

# =============================================================================
# Model data
# =============================================================================

def data_():
    p = [4, 9, 7, 5, 5, 4, 
         5, 9, 3, 7, 10, 7, 
         10, 9, 4, 8, 4, 5, 
         8, 6, 5, 5, 4, 4, 7, 
         9, 4, 8, 9, 5, 3, 7, 
        #  4, 9, 5, 
        #  5, 7, 8, 8, 6, 6, 3, 
        #  5, 3, 5, 4, 8, 7, 9, 5
         ]

    r = [61, 59, 0, 21, 6, 0, 
         57, 41, 0, 67, 0, 9, 
         40, 58, 0, 0, 0, 0, 
         0, 25, 0, 0, 5, 0, 22, 
         60, 0, 0, 61, 4, 0, 0, 
        #  24, 39, 0, 
        #  0, 24, 0, 
        #  0, 0, 0, 33, 28, 64, 
        #  31, 0, 2, 59, 0, 51
         ]

    assert(len(p) == len(r))
    m = 4
    nb_m = range(m)
    tasks = range(len(p))
    BigM = sum(p[j] for j in tasks) + max(r[j] for j in tasks)
    mapping = [(i, j) for i in tasks for j in tasks if i < j]
    mapping_2 = [(i, j) for i in tasks for j in tasks if i != j]
    mapping_3 = [(i, j, k) for i in tasks for j in tasks for k in tasks if i < j and j < k and i < k]
    T_horizen = range(sum(p[j] for j in tasks) + max(r[j] for j in tasks))
    return p, r, m, nb_m, tasks, BigM, mapping, mapping_2, mapping_3, T_horizen

p, r, m, nb_m, tasks, BigM, mapping, mapping_2, mapping_3, T_horizen = data_()

# =============================================================================
# The variables declaration 
# =============================================================================

def models_variables():
    model = Model()
    
    finish = {j: model.addVar(vtype="C", lb=0.0, name=f"finish[{j}]") for j in tasks}
    start  = {j: model.addVar(vtype="C", lb=0.0, name=f"start[{j}]")  for j in tasks}
    w_1 = {(i,j): model.addVar(vtype="B", name=f"w_1[{i}_{j}]") for i in tasks for j in tasks} 
    w_2 = {(i,j): model.addVar(vtype="B", name=f"w_2[{i}_{j}]") for i in tasks for j in tasks}  
    z = {(j,m): model.addVar(vtype="B", name=f"x[{j}_{m}]") for j in tasks for m in nb_m} 
    sigma = {(i,j): model.addVar(vtype="B", name=f"sigma[{i}_{j}]") for i in tasks for j in tasks}
    d = {(i,j): model.addVar(vtype="B", name=f"d[{i}_{j}]") for i in tasks for j in tasks}
    z_ = {(j,m,t): model.addVar(vtype="B", name=f"z[{j}_{m}_{t}]") for j in tasks for m in nb_m for t in T_horizen}  
    makespan = model.addVar(vtype="C", lb=0.0, name="makespan")

    return model, finish, start, w_1, w_2, z, sigma, d, z_, makespan
    
# =============================================================================
# Disjunctive model
# =============================================================================

def disjunctive_model(model, p, r, nb_m, tasks, BigM, mapping, finish, start, z, w_1, w_2, makespan):
    model.setObjective(makespan, "minimize")

    for j in tasks:
        model.addCons(quicksum(z[j,m] for m in nb_m) == 1)
        model.addCons(start[j] >= max(0, r[j]))
        model.addCons(start[j] + p[j] <= makespan)
        model.addCons(start[j] + p[j] == finish[j])

    for (i,j) in mapping:
        model.addCons(start[i] + p[i] <= start[j] + BigM * (1-w_1[i,j]))
        model.addCons(start[j] + p[j] <= start[i] + BigM * (1-w_2[i,j]))
        for m in nb_m:
            model.addCons((1-z[i,m]) + (1-z[j,m]) + w_1[i,j] + w_2[i,j] >= 1)

    return model
    
# =============================================================================
# LRPT heuristic for Pm|rj|Cmax
# =============================================================================

def lrpt_parallel_machine(p, r, m):
    """
    LRPT heuristic for Pm-rj-Cmax
    """

    n = len(p)
    unscheduled = list(range(n))
    machine_available = [0] * m

    schedule = []
    while unscheduled:

        machine = 0
        for k in range(1, m):
            if machine_available[k] < machine_available[machine]:
                machine = k

        current_time = machine_available[machine]
        available_jobs = [
            j for j in unscheduled
            if r[j] <= current_time
        ]

        if not available_jobs:
            next_job = unscheduled[0]
            for j in unscheduled:
                if r[j] < r[next_job]:
                    next_job = j

            current_time = r[next_job]
            available_jobs = [
                j for j in unscheduled
                if r[j] <= current_time
            ]

        selected_job = available_jobs[0]
        for j in available_jobs[1:]:
            if p[j] > p[selected_job]:
                selected_job = j

        start = max(current_time, r[selected_job])
        finish = start + p[selected_job]
        schedule.append(
            (selected_job, machine, start, finish)
        )

        machine_available[machine] = finish
        unscheduled.remove(selected_job)

    makespan = max(machine_available)

    return schedule, makespan
    
# =============================================================================
# SCIP Heuristic callback
# =============================================================================

class LRPTHeur(Heur):

    def __init__(self, p, r, m):
        super().__init__()
        self.p = p
        self.r = r
        self.m = m
        self.executed = False

    def heurexec(self, heurtiming, nodeinfeasible):

        if nodeinfeasible or self.executed:
            return {"result": SCIP_RESULT.DIDNOTRUN}

        self.executed = True

        schedule, cmax = lrpt_parallel_machine(self.p, self.r, self.m)
        n = len(self.p)

        start_vals = {}
        finish_vals = {}
        z_vals = {
            (j, k): 0.0
            for j in range(n)
            for k in range(self.m)
        }

        job_machine = {}
        for job, machine, st, ft in schedule:
            start_vals[job] = float(st)
            finish_vals[job] = float(ft)
            z_vals[(job, machine)] = 1.0
            job_machine[job] = machine

        w1_vals = {}
        w2_vals = {}
        for i in range(n):
            for j in range(i+1, n):
                same_machine = (
                    job_machine[i] == job_machine[j]
                )

                if not same_machine:
                    w1_vals[(i,j)] = 0.0
                    w2_vals[(i,j)] = 0.0
                elif finish_vals[i] <= start_vals[j]:
                    w1_vals[(i,j)] = 1.0
                    w2_vals[(i,j)] = 0.0
                else:
                    w1_vals[(i,j)] = 0.0
                    w2_vals[(i,j)] = 1.0

        warm_start = self.model.createSol()

        self.model.setSolVal(
            warm_start,
            makespan,
            float(cmax)
        )

        for j in range(n):
            self.model.setSolVal(
                warm_start,
                start[j],
                start_vals[j]
            )
            self.model.setSolVal(
                warm_start,
                finish[j],
                finish_vals[j]
            )

        for j in range(n):
            for k in range(self.m):
                self.model.setSolVal(
                    warm_start,
                    z[j, k],
                    z_vals[(j, k)]
                )

        for (i, j), val in w1_vals.items():
            self.model.setSolVal(
                warm_start,
                w_1[i, j],
                val
            )
            self.model.setSolVal(
                warm_start,
                w_2[i, j],
                w2_vals[(i, j)]
            )

        accepted = self.model.trySol(warm_start)
        if accepted:
            print("FOUND Solution")
            return {"result": SCIP_RESULT.FOUNDSOL}
        else:
            print("DID NOT FIND Solution")
            return {"result": SCIP_RESULT.DIDNOTRUN}
            
# =============================================================================
# Reconstructing the model
# =============================================================================
if __name__ == "__main__":

    mdl, finish, start, w_1, w_2, z, sigma, d, z_, makespan = models_variables()
    mdl = disjunctive_model(mdl, p, r, nb_m, tasks, BigM, mapping, finish, start, z, w_1, w_2, makespan)

# =============================================================================
# Opt. configuration
# =============================================================================

    mdl.setPresolve(SCIP_PARAMSETTING.AGGRESSIVE)
    mdl.setHeuristics(SCIP_PARAMSETTING.AGGRESSIVE)
    mdl.setSeparating(SCIP_PARAMSETTING.AGGRESSIVE)
    mdl.setIntParam("propagating/rootredcost/freq", -1)
    mdl.setRealParam("limits/time", 30.0)
    mdl.setRealParam("limits/gap", 0.05)

# =============================================================================
# Build and solve the models
# =============================================================================

    mdl.redirectOutput()
    mdl.printVersion()

    LRPT_Hur = LRPTHeur(p, r, m)
    mdl.includeHeur(LRPT_Hur, "LRPT_Heur", 
                    "callback_heuristics", "S",
                    timingmask = SCIP_HEURTIMING.BEFORENODE
                    )

    mdl.optimize()

# =============================================================================
# Stor the output
# =============================================================================

    start_vals = {j: mdl.getVal(start[j]) for j in tasks}
    finish_vals = {j: mdl.getVal(finish[j]) for j in tasks}
    z_flat = [mdl.getVal(z[j,m]) for j in tasks for m in nb_m]
    w1_vals = {(i,j): mdl.getVal(w_1[i,j]) for (i, j) in mapping}
    w2_vals = {(i,j): mdl.getVal(w_2[i,j]) for (i, j) in mapping}

# =============================================================================
# Build and solve the models
# =============================================================================

    results = []
    model_names = [
        "Disjunctive (Linear)"
    ]
    results.append((start_vals, finish_vals, z_flat, mdl.getVal(makespan)))

# =============================================================================
# Plotting function
# =============================================================================

def plot_multi_gantt(tasks, nb_m, all_results, model_names):
    n_models = len(all_results)
    if n_models <= 3:
        ncols = n_models
        nrows = 1
    else:
        ncols = 3
        nrows = (n_models + 2) // 3

    sns.set_theme(style="whitegrid")
    n_tasks = len(tasks)
    colors = sns.color_palette("husl", n_tasks)

    fig, axes = plt.subplots(nrows, ncols, figsize=(18, 10),
                             squeeze=False, sharex=False, sharey=False)

    legend_patches = [patches.Patch(facecolor=colors[i], edgecolor='black',
                                   label=f'Task {tasks[i]}') for i in range(n_tasks)]

    for idx, (start_dict, finish_dict, z_flat, ms) in enumerate(all_results):
        ax = axes.flat[idx]
        z_dict = {}
        flat_idx = 0
        for j in tasks:
            for m in nb_m:
                z_dict[(j,m)] = z_flat[flat_idx]
                flat_idx += 1
        task_machine = {}
        for j in tasks:
            for m in nb_m:
                if z_dict[(j,m)] > 0.5:
                    task_machine[j] = m
                    break

        for i, j in enumerate(tasks):
            m = task_machine.get(j)
            if m is None:
                continue
            start = start_dict[j]
            dur = finish_dict[j] - start
            ax.barh(m, dur, left=start, height=0.6,
                    color=colors[i], edgecolor='black', linewidth=0.8, alpha=0.85)
            label = f'Task {j}' if isinstance(j, int) else str(j)
            if dur > 0.5:
                ax.text(start + dur/2, m, label, ha='center', va='center',
                        fontsize=7, fontweight='bold', color='white', rotation=90)
            else:
                ax.text(start + dur + 0.1, m, label, ha='left', va='center',
                        fontsize=7, fontweight='bold', color='black', rotation=90)

        ax.set_yticks(list(nb_m))
        ax.set_yticklabels([f'M{m}' for m in nb_m])
        ax.set_xlabel('Time')
        ax.set_title(f"{model_names[idx]} (makespan={ms:.1f})", fontweight='bold')
        ax.set_ylim(-0.5, len(nb_m)-0.5)
        ax.invert_yaxis()
        ax.grid(axis='x', linestyle=':', alpha=0.4)

    for i in range(n_models, nrows*ncols):
        axes.flat[i].axis('off')

    fig.legend(handles=legend_patches, loc='center right',
               bbox_to_anchor=(1.01, 0.5), fontsize=8, title='Tasks')
    fig.suptitle("Comparison of Scheduling Models", fontsize=16, fontweight='bold')
    sns.despine()
    plt.tight_layout()
    plt.subplots_adjust(right=0.85)
    plt.show()

plot_multi_gantt(tasks, nb_m, results, model_names)

2. An extension on the pure logical form from the modeling aspects

from pyscipopt import Model, quicksum, SCIP_PARAMSETTING, Heur, SCIP_RESULT, SCIP_HEURTIMING
import matplotlib.pyplot as plt
import matplotlib.patches as patches
import seaborn as sns

# =============================================================================
# Model data
# =============================================================================

def data_():
    p = [4, 9, 7, 5, 5, 4,
         5, 9, 3, 7, 10, 7,
         10, 9, 4, 8, 4, 5
         ]

    r = [61, 59, 0, 21, 6, 0,
         57, 41, 0, 67, 0, 9,
         40, 58, 0, 0, 0, 0
         ]

    assert(len(p) == len(r))
    m = 4
    nb_m = range(m)
    tasks = range(len(p))
    BigM = sum(p[j] for j in tasks) + max(r[j] for j in tasks)
    mapping = [(i, j) for i in tasks for j in tasks if i < j]
    mapping_2 = [(i, j) for i in tasks for j in tasks if i != j]
    mapping_3 = [(i, j, k) for i in tasks for j in tasks for k in tasks if i < j and j < k and i < k]
    T_horizen = range(sum(p[j] for j in tasks) + max(r[j] for j in tasks))
    return p, r, m, nb_m, tasks, BigM, mapping, mapping_2, mapping_3, T_horizen

p, r, m, nb_m, tasks, BigM, mapping, mapping_2, mapping_3, T_horizen = data_()

# =============================================================================
# The variables declaration
# =============================================================================

def models_variables():
    model = Model()

    finish = {j: model.addVar(vtype="C", lb=0.0, name=f"finish[{j}]") for j in tasks}
    start  = {j: model.addVar(vtype="C", lb=0.0, name=f"start[{j}]")  for j in tasks}
    w_1 = {(i,j): model.addVar(vtype="B", name=f"w_1[{i}_{j}]") for i in tasks for j in tasks}
    w_2 = {(i,j): model.addVar(vtype="B", name=f"w_2[{i}_{j}]") for i in tasks for j in tasks}
    z = {(j,m): model.addVar(vtype="B", name=f"x[{j}_{m}]") for j in tasks for m in nb_m}
    sigma = {(i,j): model.addVar(vtype="B", name=f"sigma[{i}_{j}]") for i in tasks for j in tasks}
    d = {(i,j): model.addVar(vtype="B", name=f"d[{i}_{j}]") for i in tasks for j in tasks}
    z_ = {(j,m,t): model.addVar(vtype="B", name=f"z[{j}_{m}_{t}]") for j in tasks for m in nb_m for t in T_horizen}
    makespan = model.addVar(vtype="C", lb=0.0, name="makespan")

    return model, finish, start, w_1, w_2, z, sigma, d, z_, makespan

# =============================================================================
# Disjunctive model logical form
# =============================================================================

def disjunctive_model_logical(model, p, r, nb_m, tasks, BigM, mapping, finish, start, z, makespan):
    model.setObjective(makespan, "minimize")

    rhs = True
    for j in tasks:
        vars_j = [z[j, m] for m in nb_m]
        model.addConsXor(vars_j, rhs)

    for (i, j) in mapping:
        for m in nb_m:
            and_ = model.addVar(vtype="B", name=f"and_{i}_{j}_{m}")
            model.addConsAnd([z[i, m], z[j, m]], and_)
            or_ = model.addVar(vtype="B", name=f"or_{i}_{j}_{m}")
            model.addConsOr([w_1[i, j], w_2[i, j]], or_)
            model.addConsIndicator(or_ >= 1, and_, name=f"impl_{i}_{j}_{m}")

    for (i,j) in mapping:
        model.addConsIndicator(finish[i] - start[j] <= 0, binvar=w_1[i,j], name="indicator_1")
        model.addConsIndicator(finish[j] - start[i] <= 0, binvar=w_2[i,j], name="indicator_2")

    for j in tasks:
        model.addCons(start[j] >= max(0, r[j]))
        model.addCons(start[j] + p[j] <= makespan)
        model.addCons(start[j] + p[j] == finish[j])

    return model

# =============================================================================
# Ordering linear model logical form
# =============================================================================

def ordering_linear_model_logical(model, p, r, nb_m, tasks, BigM, mapping, mapping_2, mapping_3, finish, start, z, sigma, d, makespan):
    model.setObjective(makespan, "minimize")

    rhs = True
    for j in tasks:
        vars_j = [z[j, m] for m in nb_m]
        model.addConsXor(vars_j, rhs)

    for (i,j) in mapping_2:
        model.addConsIndicator(start[j] - finish[i] >= 0, sigma[i, j], name="indicator_1")

    for (i,j) in mapping:
        model.addConsXor([sigma[i,j], sigma[j,i], d[i,j]], True)
        
    for (i,j,k) in mapping_3:
        model.addConsDisjunction([
            (sigma[i,j] == 0),
            (sigma[j,k] == 0),
            (sigma[k,i] == 0)
        ], name=f"disj_{i}_{j}_{k}")

    for (i,j) in mapping:
        for m in nb_m:
            model.addConsDisjunction([
                (z[i,m] == 0),
                (z[j,m] == 0),
                (d[i,j] == 0)
            ], name=f"disj_{i}_{j}_{m}")

    for j in tasks:
        model.addCons(start[j] >= max(0, r[j]))
        model.addCons(start[j] + p[j] <= makespan)
        model.addCons(start[j] + p[j] == finish[j])
    
    for (i,j) in mapping:
        for m in nb_m:
            model.addCons( z[i,m] + z[j,m] + d[i,j] <= 2 )

    return model

# =============================================================================
# Time_indexed model
# =============================================================================

def time_indexed_model(model, p, r, nb_m, tasks, T_horizen, finish, start, z_, makespan):
    model.setObjective(makespan, "minimize")

    for m in nb_m:
        for t in T_horizen:
            vars_mt = [
                z_[j, m, tt]
                for j in tasks
                for tt in T_horizen
                if tt >= t - p[j] + 1 and tt <= t
            ]
            model.addConsCardinality(vars_mt, 1)

    for j in tasks:
        vars_j = [z_[j,m,t] for m in nb_m for t in T_horizen]
        model.addConsXor(vars_j, True)

    for j in tasks:
        for m in nb_m:
            for t in T_horizen:
                model.addConsIndicator(start[j] <= t, binvar=z_[j,m,t])
                model.addConsIndicator(start[j] >= t, binvar=z_[j,m,t])


        for j in tasks:
            model.addCons(start[j] >= max(0, r[j]))
            model.addCons(start[j] + p[j] <= makespan)
            model.addCons(start[j] + p[j] == finish[j])

    return model

# =============================================================================
# Build and solve the model
# =============================================================================

if __name__ == "__main__":

    results = []
    model_names = [
        "Disjunctive (Logical)",
        "Ordering (Logical)",
        "Time_Indexed"
    ]
# =============================================================================
# Disjunctive Logical
# =============================================================================

    mdl, finish, start, w_1, w_2, z, sigma, d, z_, makespan = models_variables()
    mdl = disjunctive_model_logical(mdl, p, r, nb_m, tasks, BigM, mapping, 
    finish, start, z, makespan)
    mdl.setPresolve(SCIP_PARAMSETTING.AGGRESSIVE)
    mdl.setHeuristics(SCIP_PARAMSETTING.AGGRESSIVE)
    mdl.setSeparating(SCIP_PARAMSETTING.AGGRESSIVE)
    mdl.setIntParam("propagating/rootredcost/freq", -1)
    mdl.setRealParam("limits/time", 60.0)
    mdl.setRealParam("limits/gap", 0.001)
    mdl.redirectOutput()
    mdl.printVersion()
    mdl.optimize()
    start_vals = {j: mdl.getVal(start[j]) for j in tasks}
    finish_vals = {j: mdl.getVal(finish[j]) for j in tasks}
    z_flat = [mdl.getVal(z[j,m]) for j in tasks for m in nb_m]
    results.append((start_vals, finish_vals, z_flat, mdl.getVal(makespan)))

# =============================================================================
# Linear Ordering Logical
# =============================================================================

    mdl, finish, start, w_1, w_2, z, sigma, d, z_, makespan = models_variables()
    mdl = ordering_linear_model_logical(mdl, p, r, nb_m, tasks, BigM, mapping, 
    mapping_2, mapping_3, finish, start, z, sigma, d, makespan)
    mdl.setPresolve(SCIP_PARAMSETTING.AGGRESSIVE)
    mdl.setHeuristics(SCIP_PARAMSETTING.AGGRESSIVE)
    mdl.setSeparating(SCIP_PARAMSETTING.AGGRESSIVE)
    mdl.setIntParam("propagating/rootredcost/freq", -1)
    mdl.setRealParam("limits/time", 60.0)
    mdl.setRealParam("limits/gap", 0.001)
    mdl.redirectOutput()
    mdl.printVersion()
    mdl.optimize()
    start_vals = {j: mdl.getVal(start[j]) for j in tasks}
    finish_vals = {j: mdl.getVal(finish[j]) for j in tasks}
    z_flat = [mdl.getVal(z[j,m]) for j in tasks for m in nb_m]
    results.append((start_vals, finish_vals, z_flat, mdl.getVal(makespan)))

# =============================================================================
# Time_indexed model
# =============================================================================

    mdl, finish, start, w_1, w_2, z, sigma, d, z_, makespan = models_variables()
    mdl = time_indexed_model(mdl, p, r, nb_m, tasks, T_horizen, finish, start, z_, makespan)
    mdl.setPresolve(SCIP_PARAMSETTING.AGGRESSIVE)
    mdl.setHeuristics(SCIP_PARAMSETTING.AGGRESSIVE)
    mdl.setSeparating(SCIP_PARAMSETTING.FAST)
    mdl.setIntParam("propagating/rootredcost/freq", -1)
    mdl.setRealParam("limits/time", 300.0)
    mdl.setRealParam("limits/gap", 0.001)
    mdl.redirectOutput()
    mdl.printVersion()
    mdl.optimize()
    start_vals = {j: mdl.getVal(start[j]) for j in tasks}
    finish_vals = {j: mdl.getVal(finish[j]) for j in tasks}
    z_flat = []
    for j in tasks:
        assigned_m = None
        for m in nb_m:
            for t in T_horizen:
                if mdl.getVal(z_[j,m,t]) > 0.5:
                    assigned_m = m
                    break
            if assigned_m is not None:
                break
        for m in nb_m:
            z_flat.append(1.0 if m == assigned_m else 0.0)
    results.append((start_vals, finish_vals, z_flat, mdl.getVal(makespan)))

# =============================================================================
# Plotting function
# =============================================================================

def plot_multi_gantt(tasks, nb_m, all_results, model_names):
    n_models = len(all_results)
    if n_models <= 3:
        ncols = n_models
        nrows = 1
    else:
        ncols = 3
        nrows = (n_models + 2) // 3

    sns.set_theme(style="whitegrid")
    n_tasks = len(tasks)
    colors = sns.color_palette("husl", n_tasks)

    fig, axes = plt.subplots(nrows, ncols, figsize=(18, 10),
                             squeeze=False, sharex=False, sharey=False)

    legend_patches = [patches.Patch(facecolor=colors[i], edgecolor='black',
                                   label=f'Task {tasks[i]}') for i in range(n_tasks)]

    for idx, (start_dict, finish_dict, z_flat, ms) in enumerate(all_results):
        ax = axes.flat[idx]
        z_dict = {}
        flat_idx = 0
        for j in tasks:
            for m in nb_m:
                z_dict[(j,m)] = z_flat[flat_idx]
                flat_idx += 1
        task_machine = {}
        for j in tasks:
            for m in nb_m:
                if z_dict[(j,m)] > 0.5:
                    task_machine[j] = m
                    break

        for i, j in enumerate(tasks):
            m = task_machine.get(j)
            if m is None:
                continue
            start = start_dict[j]
            dur = finish_dict[j] - start
            ax.barh(m, dur, left=start, height=0.6,
                    color=colors[i], edgecolor='black', linewidth=0.8, alpha=0.85)
            label = f'Task {j}' if isinstance(j, int) else str(j)
            if dur > 0.5:
                ax.text(start + dur/2, m, label, ha='center', va='center',
                        fontsize=7, fontweight='bold', color='white', rotation=90)
            else:
                ax.text(start + dur + 0.1, m, label, ha='left', va='center',
                        fontsize=7, fontweight='bold', color='black', rotation=90)

        ax.set_yticks(list(nb_m))
        ax.set_yticklabels([f'M{m}' for m in nb_m])
        ax.set_xlabel('Time')
        ax.set_title(f"{model_names[idx]} (makespan={ms:.1f})", fontweight='bold')
        ax.set_ylim(-0.5, len(nb_m)-0.5)
        ax.invert_yaxis()
        ax.grid(axis='x', linestyle=':', alpha=0.4)

    for i in range(n_models, nrows*ncols):
        axes.flat[i].axis('off')

    fig.legend(handles=legend_patches, loc='center right',
               bbox_to_anchor=(1.01, 0.5), fontsize=8, title='Tasks')
    fig.suptitle("Comparison of Scheduling Models", fontsize=16, fontweight='bold')
    sns.despine()
    plt.tight_layout()
    plt.subplots_adjust(right=0.85)
    plt.show()

plot_multi_gantt(tasks, nb_m, results, model_names)

[abb-omidi]

@Joao-Dionisio,

Thanks for your comments on the previous PR. Actually, I did what you suggested, however, the error was remaining. I really do not know whats happened in this one that passed all the checks!!!

Regards