4.5 KiB
4.5 KiB
<html lang="en">
<head>
</head>
</html>
In [1]:
from typing import Sequence
def find_subtours(x: Sequence[Sequence[int]]) -> Sequence[Sequence[int]]:
"""
Extracts subtours from the given 2D-array.
Args:
x: A two-dimensional array corresponding to the x variable in the TSP formulation, where
x[i][j] is 1 if arc (i, j) is used.
Returns:
A list of subtours, where each subtour is a list.
"""
N = len(x)
marked = [False] * N
subtours: list[list[int]] = []
while not all(marked):
# Index of the first non-marked city:
istart = min(range(N), key=lambda j: marked[j])
# Create the subtour:
subtour = [istart]
while True:
i = istart
for i, b in enumerate(x[subtour[-1]]):
if b:
break
if i == istart:
break
subtour.append(i)
for i in subtour:
marked[i] = True
subtours.append(subtour)
return subtours
In [44]:
from docplex.mp.model import Model
from docplex.mp.callbacks.cb_mixin import ConstraintCallbackMixin
from cplex.callbacks import LazyConstraintCallback
from docplex.mp.solution import SolveSolution
import numpy as np
from numpy.typing import NDArray
import tsp.data
class SubtourConstraintCallback(ConstraintCallbackMixin, LazyConstraintCallback):
x: NDArray
def __init__(self, env):
LazyConstraintCallback.__init__(self, env)
ConstraintCallbackMixin.__init__(self)
self.nb_callback = 0
self.nb_subtours = 0
self.single_subtour_per_call = False
def __call__(self):
n = len(self.x)
x_list = self.x.flatten().tolist()
sol: SolveSolution = self.make_solution_from_vars(x_list)
xvals = np.array(sol.get_values(x_list)).reshape(self.x.shape).astype(int)
subtours = find_subtours(xvals)
if len(subtours) == 1:
return
self.nb_callback += 1
for subtour in subtours:
xs = [
self.x[i, j].index
for i in subtour
for j in range(n)
if j not in subtour
]
self.add([xs, [1] * len(xs)], "G", 1)
self.nb_subtours += 1
if self.single_subtour_per_call:
break
with Model(name="TSP") as m:
dist = tsp.data.grid42
n = len(dist)
x = np.array([m.binary_var_list(n, name=f"x_{i}") for i in range(n)])
for i in range(n):
m.add_constraint(x[i, :].sum() == 1)
m.add_constraint(x[:, i].sum() == 1)
m.add_constraint(x[i, i] == 0)
m.set_objective("min", (x * dist).sum())
cb: SubtourConstraintCallback = m.register_callback(SubtourConstraintCallback)
cb.x = x
s = m.solve()
print(cb.nb_subtours)
print(s.get_objective_value())
699.0