A while ago I started trying to make a fair performance comparison of Pyomo and JuMP. I created a GitHub repository to share the code I use and I am happy if anyone wants to contribute. As of now, I am trying to find a more efficient implementation for the Pyomo models than the intuitive one i came up so far:
IJKLM
def pyomo(I, IJK, JKL, KLM, solve):
model = pyo.ConcreteModel()
model.I = pyo.Set(initialize=I)
x_list = [
(i, j, k, l, m) for (i, j, k) in IJK for l in JKL[j, k] for m in KLM[k, l]
]
constraint_dict_i = {
ii: ((i, j, k, l, m) for (i, j, k, l, m) in x_list if i == ii) for ii in I
}
model.x_list = pyo.Set(initialize=x_list)
model.c_dict_i = pyo.Set(model.I, initialize=constraint_dict_i)
model.z = pyo.Param(default=1)
model.x = pyo.Var(model.x_list, domain=pyo.NonNegativeReals)
model.OBJ = pyo.Objective(expr=model.z)
model.ei = pyo.Constraint(model.I, rule=ei_rule)
if solve:
opt = pyo.SolverFactory("gurobi")
opt.solve(model)
def ei_rule(model, i):
return sum(model.x[i, j, k, l, m] for i, j, k, l, m in model.c_dict_i[i]) >= 0
Supply Chain
def intuitive_pyomo(I, L, M, IJ, JK, IK, KL, LM, D, solve):
model = pyo.ConcreteModel()
model.I = pyo.Set(initialize=I)
model.L = pyo.Set(initialize=L)
model.M = pyo.Set(initialize=M)
model.IJ = pyo.Set(initialize=IJ)
model.JK = pyo.Set(initialize=JK)
model.IK = pyo.Set(initialize=IK)
model.KL = pyo.Set(initialize=KL)
model.LM = pyo.Set(initialize=LM)
model.f = pyo.Param(default=1)
model.d = pyo.Param(model.I, model.M, initialize=D)
model.x = pyo.Var(
[
(i, j, k)
for (i, j) in model.IJ
for (jj, k) in model.JK
if jj == j
for (ii, kk) in model.IK
if (ii == i) and (kk == k)
],
domain=pyo.NonNegativeReals,
)
model.y = pyo.Var(
[(i, k, l) for i in model.I for (k, l) in model.KL], domain=pyo.NonNegativeReals
)
model.z = pyo.Var(
[(i, l, m) for i in model.I for (l, m) in model.LM], domain=pyo.NonNegativeReals
)
model.OBJ = pyo.Objective(expr=model.f)
model.production = pyo.Constraint(model.IK, rule=intuitive_production_rule)
model.transport = pyo.Constraint(model.I, model.L, rule=intuitive_transport_rule)
model.demand = pyo.Constraint(model.I, model.M, rule=intuitive_demand_rule)
# model.write("int.lp")
if solve:
opt = pyo.SolverFactory("gurobi")
opt.solve(model)
def intuitive_production_rule(model, i, k):
lhs = [
model.x[i, j, k]
for (ii, j) in model.IJ
if ii == i
for (jj, kk) in model.JK
if (jj == j) and (kk == k)
]
rhs = [model.y[i, k, l] for (kk, l) in model.KL if kk == k]
if lhs or rhs:
return sum(lhs) >= sum(rhs)
else:
return pyo.Constraint.Skip
def intuitive_transport_rule(model, i, l):
lhs = [model.y[i, k, l] for (k, ll) in model.KL if ll == l]
rhs = [model.z[i, l, m] for (lll, m) in model.LM if lll == l]
if lhs or rhs:
return sum(lhs) >= sum(rhs)
else:
return pyo.Constraint.Skip
def intuitive_demand_rule(model, i, m):
return sum(model.z[i, l, m] for (l, mm) in model.LM if mm == m) >= model.d[i, m]
I measure performance in model generation time for two exemplary models I refer to as IJKLM and Supply Chain.
Models:
These are the results for increasing instance sizes:
Can anyone help to improve Pyomos performance?










pyomoside, but one of the driving issues here is that your data is somewhat nonsensical. At least on the Supply Chain Model. (I didn't look at the other). If you inspect the sizes of some of the sets you are generating, when the size ofIis 4600, the domain ofxis only ~230. Meaning less than 5% of your items are even capable of being manufactured/delivered. I would expect the domain ofxto me many times larger than|I|, somewhere closer to|I| x |J| x |K|.python main_supply_chain.py. I had tomkdir supply_chain/databefore running though.