Secure your code as it's written. Use Snyk Code to scan source code in minutes - no build needed - and fix issues immediately.
best_z_ub = None
best_z_lb = None
for tr in range(tries):
inst = gen_inst(nj)
J = inst["J"]
dur = inst["d"]
c = inst["c"]
r = inst["r"]
S = inst["S"]
EST = inst["EST"]
mip = create_mip("gurobi", J, dur, S, c, r, EST)
mip.verbose = 0
mip.relax()
mip.optimize()
assert mip.status == OptimizationStatus.OPTIMAL
z_relax = mip.objective_value
mip = create_mip("gurobi", J, dur, S, c, r, EST)
mip.verbose = 0
mip.max_nodes = 512
mip.optimize()
z_lb = mip.objective_bound
z_ub = mip.objective_value
assert mip.status in [
OptimizationStatus.OPTIMAL,
OptimizationStatus.FEASIBLE,
]
sum_dur = sum(dur[i] for i in J)
score = (
(z_ub - z_relax) / sum_dur
+ 2 * (z_ub - z_relax)
+ 5 * (z_ub - z_lb)
def create_mip(solver, J, dur, S, c, r, EST, relax=False, sense=MINIMIZE):
"""Creates a mip model to solve the RCPSP"""
NR = len(c)
mip = Model(solver_name=solver)
sd = sum(dur[j] for j in J)
vt = CONTINUOUS if relax else BINARY
x = [
{
t: mip.add_var("x(%d,%d)" % (j, t), var_type=vt)
for t in range(EST[j], sd + 1)
}
for j in J
]
TJ = [set(x[j].keys()) for j in J]
T = set()
for j in J:
T = T.union(TJ[j])
if sense == MINIMIZE:
mip.objective = minimize(xsum(t * x[J[-1]][t] for t in TJ[-1]))
def create_mip(solver, J, dur, S, c, r, EST, relax=False, sense=MINIMIZE):
"""Creates a mip model to solve the RCPSP"""
NR = len(c)
mip = Model(solver_name=solver)
sd = sum(dur[j] for j in J)
vt = CONTINUOUS if relax else BINARY
x = [
{
t: mip.add_var("x(%d,%d)" % (j, t), var_type=vt)
for t in range(EST[j], sd + 1)
}
for j in J
]
TJ = [set(x[j].keys()) for j in J]
T = set()
for j in J:
T = T.union(TJ[j])
EST = inst["EST"]
mip = create_mip("gurobi", J, dur, S, c, r, EST)
mip.verbose = 0
mip.relax()
mip.optimize()
assert mip.status == OptimizationStatus.OPTIMAL
z_relax = mip.objective_value
mip = create_mip("gurobi", J, dur, S, c, r, EST)
mip.verbose = 0
mip.max_nodes = 512
mip.optimize()
z_lb = mip.objective_bound
z_ub = mip.objective_value
assert mip.status in [
OptimizationStatus.OPTIMAL,
OptimizationStatus.FEASIBLE,
]
sum_dur = sum(dur[i] for i in J)
score = (
(z_ub - z_relax) / sum_dur
+ 2 * (z_ub - z_relax)
+ 5 * (z_ub - z_lb)
)
assert score >= -1e-10
if score > best_score:
best_score = score
best_inst = inst
best_z_relax = z_relax
best_z_lb = z_lb
best_z_ub = z_ub
inst = best_inst
with open("./data/rcpsp-%d-%d.json" % (nj, ii + 1), "w") as outfile:
sd = sum(dur[j] for j in J)
vt = CONTINUOUS if relax else BINARY
x = [
{
t: mip.add_var("x(%d,%d)" % (j, t), var_type=vt)
for t in range(EST[j], sd + 1)
}
for j in J
]
TJ = [set(x[j].keys()) for j in J]
T = set()
for j in J:
T = T.union(TJ[j])
if sense == MINIMIZE:
mip.objective = minimize(xsum(t * x[J[-1]][t] for t in TJ[-1]))
else:
mip.objective = maximize(xsum(t * x[J[-1]][t] for t in TJ[-1]))
# one time per job
for j in J:
mip += xsum(x[j][t] for t in TJ[j]) == 1, "selTime(%d)" % j
# precedences
for (u, v) in S:
mip += (
xsum(t * x[v][t] for t in TJ[v])
>= xsum(t * x[u][t] for t in TJ[u]) + dur[u],
"prec(%d,%d)" % (u, v),
)
# resource usage
def create_mip(solver, J, dur, S, c, r, EST, relax=False, sense=MINIMIZE):
"""Creates a mip model to solve the RCPSP"""
NR = len(c)
mip = Model(solver_name=solver)
sd = sum(dur[j] for j in J)
vt = CONTINUOUS if relax else BINARY
x = [
{
t: mip.add_var("x(%d,%d)" % (j, t), var_type=vt)
for t in range(EST[j], sd + 1)
}
for j in J
]
TJ = [set(x[j].keys()) for j in J]
T = set()
for j in J:
T = T.union(TJ[j])
if sense == MINIMIZE:
mip.objective = minimize(xsum(t * x[J[-1]][t] for t in TJ[-1]))
else:
mip.objective = maximize(xsum(t * x[J[-1]][t] for t in TJ[-1]))
if v.name.startswith('x(') and v.x >= 0.99]
U = [int(v.name.split('(')[1].split(',')[0]) for v, f in r]
V = [int(v.name.split(')')[0].split(',')[1]) for v, f in r]
N, cp = set(U+V), CutPool()
# output nodes for each node
out = defaultdict(lambda: list())
for i in range(len(U)):
out[U[i]].append(V[i])
for n in N:
S = set(subtour(N, out, n))
if S:
arcsInS = [(v, f) for i, (v, f) in enumerate(r)
if U[i] in S and V[i] in S]
if sum(f for v, f in arcsInS) >= (len(S)-1)+1e-4:
cut = xsum(1.0*v for v, fm in arcsInS) <= len(S)-1
cp.add(cut)
for cut in cp.cuts:
model += cut
model.start = [(x[seq[i]][seq[i+1]], 1) for i in range(len(seq)-1)]
model.max_seconds = timeLimit
model.threads = threads
model.optimize()
end = time()
print(model.status)
print(model.solver_name)
objv = 1e20
gap = 1e20
n_nodes = 0
if model.status in [OptimizationStatus.OPTIMAL, OptimizationStatus.FEASIBLE]:
out.write('route with total distance %g found: 0'
% (model.objective_value))
nc = 0
while True:
n_nodes += 1
nc = [i for i in V if x[nc][i].x >= 0.99][0]
out.write(' -> %s' % nc)
if nc == 0:
break
out.write('\n')
objv = model.objective_value
gap = model.gap
if n_nodes != n:
err.write('incomplete route (%d from %d) generated.\n' % (n_nodes, n))
exit(1)
def get_int_param(self, name: str) -> int:
res = ffi.new("int *")
env = GRBgetenv(self._model)
error = GRBgetintparam(env, name.encode("utf-8"), res)
if error != 0:
raise ParameterNotAvailable(
"Error getting gurobi integer parameter {}".format(name)
)
return res[0]
def get_objective(self) -> LinExpr:
obj = cbclib.Cbc_getObjCoefficients(self._model)
if obj == ffi.NULL:
raise ParameterNotAvailable("Error getting objective function coefficients")
return (
xsum(
obj[j] * self.model.vars[j]
for j in range(self.num_cols())
if abs(obj[j]) >= 1e-15
)
+ self._objconst
)