felixm/discrete_optimization
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Eod. Add clustering, but did not make anything better besides that.

Browse Source
This commit is contained in:
Felix Martin
2020-01-12 23:06:01 -05:00
parent a2f9761517
commit 65fc139e65
1 changed files with 225 additions and 111 deletions
Download Patch File Download Diff File Expand all files Collapse all files

336
tsp/tsp.py
View File

@@ -1,10 +1,33 @@
import math import math
from functools import lru_cache from functools import lru_cache
from collections import namedtuple from collections import namedtuple
from geometry import intersect from random import shuffle
import time import time
Point = namedtuple("P", ['name', 'x', 'y'])
class Point(object):
def __init__(self, index, x, y):
self.index = index
self.x = x
self.y = y
self.cluster_x = None
self.cluster_y = None
# A list tuples. The first field is the distance
# and the second is another point (i.e. a neighbor).
self.neighbors = []
def add_neighbors(self, neighbors):
neighbors = [(n, distance(self, n)) for n in neighbors]
self.neighbors = sorted(neighbors, key=lambda t: t[1])
def __str__(self):
# m = "P_{}({}, {})".format(self.index, self.x, self.y)
# m = "P_{}({}, {})".format(self.index, self.cluster_x, self.cluster_y)
m = "P({})".format(self.index)
return m
def __repr__(self):
return self.__str__()
def parse_input_data(input_data): def parse_input_data(input_data):
@@ -20,44 +43,13 @@ def float_is_equal(a, b):
return False return False
def plot_graph(points):
try:
import matplotlib.pyplot as plt
except ModuleNotFoundError:
return
def plot_arrows():
for i in range(len(points)):
p1 = points[i - 1]
p2 = points[i]
plot_arrow(p1, p2)
def plot_arrow(p1, p2):
x = p1.x
y = p1.y
dx = p2.x - x
dy = p2.y - y
opt = {'head_width': 0.4, 'head_length': 0.4, 'width': 0.05,
'length_includes_head': True}
plt.arrow(x, y, dx, dy, **opt)
def plot_points():
for p in points:
plt.plot(p.x, p.y, '')
plt.text(p.x, p.y, ' ' + p.name)
plot_points()
plot_arrows()
plt.show()
def prepare_output_data(points): def prepare_output_data(points):
# Basic plausibility checks # Basic plausibility checks
assert(len(set(points)) == len(points)) assert(len(set(points)) == len(points))
assert(len(points) > 4) assert(len(points) > 4)
obj = total_distance(points) obj = total_distance(points)
output_data = '%.2f' % obj + ' ' + str(0) + '\n' output_data = '%.2f' % obj + ' ' + str(0) + '\n'
output_data += ' '.join(map(lambda p: p.name, points)) output_data += ' '.join(map(lambda p: str(p.index), points))
return output_data return output_data
@@ -125,63 +117,38 @@ def swap_edges(i, j, points, current_distance=0):
return current_distance return current_distance
def local_search_2_opt(points):
current_total = total_distance(points)
ignore_set = set()
while True:
pi, i = longest_distance(points, ignore_set)
ignore_set.add(pi)
if not pi:
break
best_new_total = current_total
best_points = None
swap = None
for j in range(len(points)):
if j in [i, i + 1, i + 2]:
continue
new_points = list(points)
swap_edges(i, j - 1, new_points)
new_total = total_distance(new_points)
if new_total < best_new_total:
swap = (points[i], points[j - 1])
best_new_total = new_total
best_points = new_points
if best_new_total < current_total:
current_total = best_new_total
points = best_points
ignore_set = set()
return points
def reorder_points_greedy(points): def reorder_points_greedy(points):
current_point = points[0] best_length = float("inf")
solution = [current_point] best_solution = None
points = points[1:]
while points: for i in range(1000):
min_length = 999999 shuffle(points)
min_point = None current_point, points = points[0], points[1:]
for next_point in points: solution = [current_point]
new_length = distance(current_point, next_point)
if new_length < min_length:
min_length = new_length
min_point = next_point
current_point = min_point
solution.append(current_point)
points.remove(current_point)
return solution while points:
next_point = None
# Select the closest point as the following one.
for neighbor, _ in current_point.neighbors:
if neighbor in points:
next_point = neighbor
points.remove(next_point)
break
# If none of the neighbors could be selected use any point.
if next_point is None:
next_point = points.pop()
def print_swap(i, j, points): solution.append(next_point)
print("Swap:", points[i].name, " <-> ", points[j].name) current_point = next_point
total_length = total_distance(solution)
points = solution
if total_length < best_length:
best_length = total_length
best_solution = solution.copy()
def get_indices(current_index, points): return best_solution
for i in range(len(points)):
yield i
def k_opt(p1_index, points, steps): def k_opt(p1_index, points, steps):
@@ -194,8 +161,8 @@ def k_opt(p1_index, points, steps):
ignore_set.add(p2) ignore_set.add(p2)
p4_index = None p4_index = None
#for p3_index in range(len(points)): # TODO(felixm): Keep track of current indices and then make this more efficient.
for p3_index in get_indices(p2_index, points): for p3_index in range(len(points)):
p3 = points[p3_index] p3 = points[p3_index]
p4 = points[p3_index - 1] p4 = points[p3_index - 1]
if p4 in ignore_set or p4 is p1: if p4 in ignore_set or p4 is p1:
@@ -246,36 +213,154 @@ def local_search_k_opt(points):
return points return points
def split_into_sections(points):
x_min, x_max, y_min, y_max = float("inf"), 0, float("inf"), 0 class Map(object):
for p in points: # Create Map. Cluster points into regions. Calculate distances only to own
if p.x < x_min: x_min = p.x # and neighbor regions. We can actually cluster in O(n) when we know how
if p.x > x_max: x_max = p.x # high and wide the clusters are. Once we have that working we go from
if p.y < y_min: y_min = p.y # there
if p.y > y_max: y_max = p.y CLUSTER_SIZE = 3 # How many points we want per cluster.
return
def __init__(self, points):
self.points = points
self.num_points = len(points)
self.calc_corners()
self.calc_cluster_dim()
self.sort_points_into_clusters()
self.add_neighbors_to_points()
def calc_cluster_dim(self):
clusters = self.num_points // self.CLUSTER_SIZE
# Calculate number of clusters to have a square
self.clusters_x = math.ceil(math.sqrt(clusters))
self.clusters_y = self.clusters_x
self.clusters_total = self.clusters_x ** 2
self.cluster_x_dim = (self.x_max - self.x_min) / self.clusters_x
self.cluster_y_dim = (self.y_max - self.y_min) / self.clusters_y
def add_neighbors_to_points(self):
""" Add all points from the surrounding clusters to each point. """
for p in self.points:
clusters_x = [p.cluster_x]
clusters_y = [p.cluster_y]
if p.cluster_x - 1 >= 0:
clusters_x.append(p.cluster_x - 1)
if p.cluster_x + 1 < self.clusters_x:
clusters_x.append(p.cluster_x + 1)
if p.cluster_y - 1 >= 0:
clusters_y.append(p.cluster_y - 1)
if p.cluster_y + 1 < self.clusters_y:
clusters_y.append(p.cluster_y + 1)
clusters = [(x, y)
for x in clusters_x
for y in clusters_y]
neighbors = []
for x, y in clusters:
for p2 in self.clusters[x][y]:
if p is not p2:
neighbors.append(p2)
p.add_neighbors(neighbors)
def sort_points_into_clusters(self):
self.clusters = [[[]
for x in range(self.clusters_y)]
for y in range(self.clusters_y)]
for p in self.points:
cluster_x = int((p.x - self.x_min) // self.cluster_x_dim)
cluster_y = int((p.y - self.y_min) // self.cluster_y_dim)
# If the point is on the outer edge of the highest cluster
# the index will be outside the correct range. We put it
# into the closes cluster.
if cluster_x == self.clusters_x:
cluster_x -= 1
if cluster_y == self.clusters_y:
cluster_y -= 1
self.clusters[cluster_x][cluster_y].append(p)
p.cluster_x = cluster_x
p.cluster_y = cluster_y
def calc_corners(self):
x_min, x_max = float("inf"), float("-inf")
y_min, y_max = float("inf"), float("-inf")
for p in self.points:
if p.x < x_min:
x_min = p.x
if p.x > x_max:
x_max = p.x
if p.y < y_min:
y_min = p.y
if p.y > y_max:
y_max = p.y
self.x_min = x_min
self.x_max = x_max
self.y_min = y_min
self.y_max = y_max
def plot(self):
try:
import matplotlib.pyplot as plt
except ModuleNotFoundError:
return
def plot_grid():
for x_i in range(self.clusters_x + 1):
x_1 = self.x_min + x_i * self.cluster_x_dim
x_2 = x_1
y_1 = self.y_min
y_2 = self.y_max
plt.plot([x_1, x_2], [y_1, y_2], 'b:')
for y_i in range(self.clusters_y + 1):
x_1 = self.x_min
x_2 = self.x_max
y_1 = self.y_min + y_i * self.cluster_y_dim
y_2 = y_1
plt.plot([x_1, x_2], [y_1, y_2], 'b:')
def plot_arrows():
for i in range(self.num_points):
p1 = self.points[i - 1]
p2 = self.points[i]
plot_arrow(p1, p2)
def plot_arrow(p1, p2):
x = p1.x
y = p1.y
dx = p2.x - x
dy = p2.y - y
opt = {'head_width': 0.4, 'head_length': 0.4, 'width': 0.05,
'length_includes_head': True}
plt.arrow(x, y, dx, dy, **opt)
def plot_points():
for i, p in enumerate(self.points):
plt.plot(p.x, p.y, '')
plt.text(p.x, p.y, ' ' + str(p))
for nb, _ in p.neighbors:
# plt.plot([p.x, nb.x], [p.y, nb.y], 'r--')
pass
plot_points()
plot_grid()
plot_arrows()
plt.show()
def solve_it(input_data): def solve_it(input_data):
points = parse_input_data(input_data) points = parse_input_data(input_data)
num_points = len(points) # Initialiaze map before algorithm because it clusters the points
# and adds the neighbors to each point.
m = Map(points)
m.points = reorder_points_greedy(points)
# FIXME(felixm): Don't do this here.
m.points = local_search_k_opt(m.points)
m.plot()
if num_points == 51: return prepare_output_data(m.points)
return """428.98 0
47 26 6 36 12 30 23 35 13 7 19 40 11 42 18 16 44 14 15 38 50 39 43 29 21 37 20 25 1 31 22 48 49 17 32 0 33 5 2 28 10 9 45 3 46 8 4 34 24 41 27"""
elif num_points == 100:
return """21930.64 0
5 21 99 11 32 20 87 88 77 37 47 7 83 39 74 66 57 71 24 3 55 96 80 14 16 4 91 13 69 28 62 64 76 34 2 50 89 61 95 73 81 56 31 58 27 75 10 86 78 67 98 65 0 12 93 15 97 33 60 1 45 36 46 30 94 82 49 23 6 85 63 48 68 41 59 42 53 9 18 52 22 8 90 38 70 17 79 26 29 51 84 72 19 25 40 43 44 35 54 92
"""
elif num_points < 2000:
points = reorder_points_greedy(points)
points = local_search_k_opt(points)
#sections = split_into_sections(points)
#points = local_search_2_opt(points)
# plot_graph(points)
return prepare_output_data(points)
if __name__ == "__main__": if __name__ == "__main__":
@@ -284,3 +369,32 @@ if __name__ == "__main__":
input_data = input_data_file.read() input_data = input_data_file.read()
print(solve_it(input_data)) print(solve_it(input_data))
def local_search_2_opt(points):
current_total = total_distance(points)
ignore_set = set()
while True:
pi, i = longest_distance(points, ignore_set)
ignore_set.add(pi)
if not pi:
break
best_new_total = current_total
best_points = None
swap = None
for j in range(len(points)):
if j in [i, i + 1, i + 2]:
continue
new_points = list(points)
swap_edges(i, j - 1, new_points)
new_total = total_distance(new_points)
if new_total < best_new_total:
swap = (points[i], points[j - 1])
best_new_total = new_total
best_points = new_points
if best_new_total < current_total:
current_total = best_new_total
points = best_points
ignore_set = set()
return points