This file shows how to solve the minimum cost flow problem for the "Reliable Construction Co." prototype example from Hillier and Lieberman, Section 9.8. To write an explicit and self-contained program, we use a concrete Pyomo model. Note, however, that an abstract model with external .dat
file would be a natural choice here.
from pyomo.environ import *
from pyomo.opt import *
opt = solvers.SolverFactory("glpk")
The nodes of the project network correspond to tasks. The regular completion time of task $i$ is $T_i$. However, up to $R_i$ units of time can be saved by "crashing" at a cost of $C_i$ per unit time.
T = {'Start':0,
'A':2,
'B':4,
'C':10,
'D':6,
'E':4,
'F':5,
'G':7,
'H':9,
'I':7,
'J':8,
'K':4,
'L':5,
'M':2,
'N':6,
'Finish':0}
R = {'Start':0,
'A':1,
'B':2,
'C':3,
'D':2,
'E':1,
'F':2,
'G':3,
'H':3,
'I':2,
'J':2,
'K':1,
'L':2,
'M':1,
'N':3,
'Finish':0}
C = {'Start':0,
'A':100,
'B':50,
'C':80,
'D':40,
'E':160,
'F':40,
'G':40,
'H':60,
'I':30,
'J':30,
'K':40,
'L':50,
'M':100,
'N':60,
'Finish':0}
We now define the network arc. For this problem, it is most natural to specify the list of dependencies for each task. A list of all arcs is then easily constructed by one line of Python code.
P = {'Start':[],
'A':['Start'],
'B':['A'],
'C':['B'],
'D':['C'],
'E':['C'],
'F':['E'],
'G':['D'],
'H':['G','E'],
'I':['C'],
'J':['F','I'],
'K':['J'],
'L':['J'],
'M':['H'],
'N':['K','L'],
'Finish':['M','N']}
N = list(T.keys())
A = [(j,i) for i in N for j in P[i]]
The decision variables are $t_i$, the starting time of task $i$ and the crash time $x_i$, i.e., the time saved on task $i$ by crashing.
model = ConcreteModel()
model.t = Var(N, within=NonNegativeReals)
model.x = Var(N, within=NonNegativeReals)
Whenever task $j$ depends on task $i$, task $j$ must not start before $i$ is completed, i.e. $$ t_j \geq t_i+ T_i - x_i \,. $$ The project starts at time $0$ and must finish before time $40$. Further, for each task $i$ the crash time must satisfy $x_i \leq R_i$.
def time_rule(model, i, j):
return model.t[j] >= model.t[i] + T[i] - model.x[i]
model.time = Constraint(A, rule=time_rule)
model.start = Constraint(expr = model.t['Start'] == 0)
model.finish = Constraint(expr = model.t['Finish'] <= 40)
def crash_rule(model, n):
return model.x[n] <= R[n]
model.crash = Constraint(N, rule=crash_rule)
Our objective is minimizing the cost of crashing.
model.cost = Objective(expr = sum(model.x[n]*C[n] for n in N), sense=minimize)
We now solve and inspect the solution.
results = opt.solve(model)
model.x.get_values()
model.t.get_values()
model.cost.expr()
In the setting of the example, this means that we can reduce the project time to 40 weeks if we crash task F and task J by two weeks each at a total additional cost of $140000.