Комментарии 42
Оставлю для интересующихся SOTA.
https://en.m.wikipedia.org/wiki/Concorde_TSP_Solver
Решает инстансы до несколько десятков тысяч вершин.
Под капотом integer programming подход с branch and cut.
Если я правильно понимаю, подход применим только для симметричной задачи и только на плоскости. Я же предлагаю решать в общем виде. Хотя да результаты интересные.
Да симметричный, но есть трюк с конвертацией асимметричной в симметричную. Не только на плоскости, веса не должны удовлетворять неравенство треугольника.
неравенство треугольника
Это то которое a+b ≤ c ?
Скорее D(a, b) <= D(a, c) + D(b, c)
https://en.wikipedia.org/wiki/Travelling_salesman_problem#Metric
Почитал подробнее, это прорывной метод для огромных матриц с заранее известной метрикой пространства, я в восхищении!
Забавно, что в нём так же используется алгоритм ветвления и границ, для выбора лучшего решения.
На самом деле формулировка через integer programming позволяет решать многие задачи дискретной оптимизации в том числе и NP так, что это довольно стандартный и старый подход. Но дьявол как обычно в деталях реализации.
Если у нас не каждый город имеет соединение с другим городом напрямую, а расстояние между городами непосредственно связанными между собой одинаковое (1 прыжок), мы сможем воспользоваться этим методом?
Вы проводили оценки эффективности своих оптимизаций? Избавление от рекурсии не выглядит чем-то эффективным. Особенно по сравнению с постоянными выделениями и освобождениями памяти, которых как раз можно было бы избежать, перемещаясь по статическим массивам.
Эффект есть, и довольно значительный, иначе не было бы смысла об этом писать. Конкретные значения я не приведу, это требует отдельного обзора. Избавление от рекурсии как таковой не самоцель, это требуется для того чтоб не хранить отдельные матрицы на каждом шаге ветвления.
Очень грубый пример: матрица 32х32 требует в среднем где-то 10млн рабочих шагов, возьмём за среднее значение хранимой матрицы шага как 20х20 ~ 1кб. Получается нам нужно хранить в памяти более 9 ГБ информации. А если матрица ещё больше?
Дублирование данных вообще никак не связано с рекурсией.
Я возможно вас не совсем правильно понимаю, поясните свою мысль.
Я имел в виду рекурсию в широком смысле, как возможность не вызывать функцию дополнительный раз, когда нужно передавать что-то в виде параметров, которые будут хранится в стеке.
При стандартном алгоритме ветвей и границ мы в любом случае обязаны хранить обе ветви на каждом ветвлении. В моём варианте мы копируем только одну ветвь, а вторую обходим итеративно в цикле не выделяя память отдельно.
Тут надо начать с того, что хвостовая рекурсия алгоритмически строго эквивалентна циклу, это просто способ записи управляющих конструкций. Некоторые компиляторы прямо преобразуют хвостовую рекурсию в цикл.
А писал я выше о том, что выделение памяти вообще с рекурсией никак не связано. Вы можете весь свой алгоритм переписать на статических переменных и вообще не выделять память, в том числе и в рекурсивных функциях. И экономить вам надо выделения и освобождения памяти, а не рекурсивные вызовы, которые сами по себе – просто передача управления. А на стеке в таком случае будут передаваться просто индексы в массивах.
Метод ветвей и границ не является точным тк его скорость и качество зависят от размерности задачи. Говорить что метод ветвей и границ решает NP problem задачи, значит говорить о том, что они соответствуют P problem и все эти квантовые отжиги и тп не нужны. Генетика и пр метаэвристики тоже решают задачи NP до определленой размерности К числа городов точно, но потом находят приближённое к оптимальному решению, но гораздо быстрее чем ветви и границы.
Вы не правы, метод ветвей и границ с полным перебором — это точный метод. Другое дело что его практически всегда используют с эвристиками, значительно уменьшающими число вычислений, но тогда да алгоритм перестаёт быть точным.
Зависит от эвристик. Иногда они не ломают точность метода. Можно, например, выходить из ветви, если оценка снизу скажет, что все равно получится не лучше оптимального решения. Если же оценка снизу не строгая, то — да, есть шанс потерять оптимальное решение.
Можно, например, выходить из ветви, если оценка снизу скажет, что все равно получится не лучше оптимального решения
Так тут и используется этот подход, но это не является эвристикой.
Ну так оценку можно по разному делать. Вот тут и есть пространство для всяких эвристик.
Ну какие могут быть эвристики в точном решении? Эта реализация как раз и создавалась для поиска точного решения максимально быстро. Для эвристики существуют другие более выигрышные подходы.
Допустим города на плоскости. Тогда можно выходить, если минимальное расстояние, нужное для замыкания всех висячих концов сделает текущее решение хуже известного лучшего (высчитанное как тупо эвклидово расстояние, не смотря на то, какие ребра там есть). При чем не факт, что именно с такой стоимостью вы сможете все замкнуть, да и вообще у вас будет один большой цикл в конце, а не куча мелких. Эта эвристика считается легко, в некоторых случаях может ускорить решение, но при этом не ломает поиск точного решения.
Вот это вот вычисление оценки через эвклидовы расстояния и есть эвристика, в моем понимании. Оно не ищиет минимальное требуемое расстояние для замыкания цикла точно — это было бы полноценное решение задачи коммивояжора. Оно делает лишь какую-то грубую оценку, забивая на все подряд, лишь бы просто считалось. Но, поскольку эта оценка снизу, то метод ветвей и границ все еще находит оптимальное решение.
Допустим города на плоскости
Тут не всё так однозначно, для плоскости есть более эффективные решения. А где гарантия что задающие вершины располагаются не в трёхмерном пространстве, n- мерном, или расчёт строения белков. Мне хотелось создать универсальный алгоритм.
оценки через эвклидовы расстояния
Над этим нужно думать.
Пробовал играться с минимальными остовыми деревьями как барьерами на плоскости, для уменьшения числа расчётов. Но и там есть исключения, не стал добавлять в текст, ибо остались вопросы.
Чем более конкретная задача, тем лучше работают эвристики. Но можно придумать и что-то совсем общее. Например, брать для каждой висячей вершины мнимальное исходящее ребро (даже если оно замыкает мелкий цикл). Работать это будет не так эффективно, как пример выше.
Но неважно. Оценка эвристическая? Да. Метод все еще точный? Да.
А вообще, есть и общие алгоритмы основанные на эвристиках фактически. Намриер A*. Там надо иметь некоторые оценки оставшегося пути снизу. А уж как их считать — это уже в алгоритме не указано. Но эти оценки всегда будут эвристическими, ведь точно оставшееся расстояние не подсчитать.
Вы поймите, я тут не критикую вашу реализацию, а лишь спорю с вашим утверждением, что эвристики взаимосиключительны с алгоритмами, ищущими точное решение.
NP задачи очевидно решаются брутфорсом имеющим в худшем случаи экспоненциальную сложность. То что брутфорс оптимизирован и решает какие-то реальные инстансы гораздо быстрее, с теоретической точки зрения ничего не значит и не делает его принадлежащим P.
Но так ли важно качество кода, если он медленный?
Примерно такой же вопрос задают авторы Pony
unsigned int l = 0;
проклято.
unsigned int line_cache[n + 1];
я правильно понял, что кэши хранят первым элементом сумму строки/столбца?
Вообще, есть ощущение, что сишный алгоритм можно ещё поускорять, т.к. вроде некоторые куски имеют достаточно четкие ограничения для распараллеливания.
А суперперестановки пробовали генерировать?
Про проклятие не уловил, можно подробнее?
Кэши хранят в [i] - элементе минимум2 по [i] строке.
line_cache[0] - нулевой элемент не используется, так как нулевая строка матрицы содержит индексы. Мне думается что в этом нагруженном куске делать пересчёт индекса на – 1 на каждом элементе цикла хуже неиспользуемых 4 байт.
Алгоритм можно ещё ускорить тут вы правы.
Чем тут помогут суперперестановки?
А вот тут +1 и вправду лишнее:
unsigned int line_cache[n + 1];
Суперперестановки это уже про применение коммивояжера, а не про оптимизацию вашего алгоритма. Задача генерации последовательности суперперестановок сводилось к задаче коммивояжёра (wiki). Правда, в постановке статьи упоминается эвристический алгоритм и ассиметричный случай TSP. Потому и спросил - пробовали ли написать что-то подобное для вашего решателя.
Спасибо за интересную статью. Не думали выложить код на github?
Муравьи как роевой интеллект и эмерджентная система решают задачу коммивояжера на раз-два. Хотелось бы увидеть сравнение эффективности таких систем в подобной статье.
Ещё добавлю пять копеек к вопросу эвристик. Для симметричной версии задачи на плоскости можно несколько сократить число расчётов убрав заведомо неоптимальные связи. Например, построив для графа минимальное остовое древо, ветви которого будут своеобразными стенами, отсекающими части графа.
Это простая идея может в некоторых случаях сильно ускорить расчёт, но и она не идеальна. Вот пример того, когда эвристика не позволяет найти оптимальное решение.
![Зелёным цветом - минимальное остовое древо Красным - минимальный путь с учётом стен из ветвей остова (длина 4799) Синим пунктиром - фактический минимальный путь (длина 4796) Зелёным цветом - минимальное остовое древо Красным - минимальный путь с учётом стен из ветвей остова (длина 4799) Синим пунктиром - фактический минимальный путь (длина 4796)](https://habrastorage.org/getpro/habr/upload_files/ae3/b57/c63/ae3b57c633f1e9998da24f5ce4bbd3c0.png)
Ну да время сократили, но это не минимальное решение.
Код
# Генерация случайного массива и сравнение у метода ветвей и границ с остовным беревом и без
import random
import matplotlib.pyplot as plt
import networkx as nx
import math as mt
import numpy as np
from ctypes import *
from datetime import datetime
import json
#-----------------------------------------------------------------
def distance(x1, y1, x2, y2):
return mt.sqrt((x2 - x1) ** 2 + (y2 - y1) ** 2)
#-----------------------------------------------------------------
def Intersection(ax1, ay1, ax2, ay2, bx1, by1, bx2, by2):
v1 = (bx2 - bx1) * (ay1 - by1) - (by2 - by1) * (ax1 - bx1)
v2 = (bx2 - bx1) * (ay2 - by1) - (by2 - by1) * (ax2 - bx1)
v3 = (ax2 - ax1) * (by1 - ay1) - (ay2 - ay1) * (bx1 - ax1)
v4 = (ax2 - ax1) * (by2 - ay1) - (ay2 - ay1) * (bx2 - ax1)
if (v1 * v2 < 0) and (v3 * v4 < 0):
return True
else:
return False
#-----------------------------------------------------------------
def tsp_branch(n, py_arr, lib):
if n < 2:
return {}
flatten_arr = list(np.concatenate(py_arr).flat)
l = [-1] * (n * n - len(flatten_arr))
flatten_arr = (flatten_arr + l)[:n * n:]
int_arr = (c_int * (n * n))(*flatten_arr)
res = lib.tsp_branch(n, byref(int_arr))
if res > 0:
l = list(int_arr)[:res:]
return {'len' : l.pop(0), 'steps' : l.pop(0), 'path' : l}
else:
return {}
#-----------------------------------------------------------------
random.seed(40)
n = 35
INF = -1
#-----------------------------------------------------------------
v1 = []
points = {}
for i in range(n):
points[i] = (random.randint(1,1000), random.randint(1,1000))
input_matrix = []
for i, vi in points.items():
m1 = []
for j, vj in points.items():
if i==j:
m1.append(INF)
else:
m1.append(round(distance(vi[0], vi[1], vj[0], vj[1])))
v1.append([i,j,round(distance(vi[0], vi[1], vj[0], vj[1]))])
input_matrix.append(m1.copy())
#-----------------------------------------------------------------
lib_c = cdll.LoadLibrary(r"tsp_branch.dll")
lib_c.tsp_branch.argtypes = [c_int, c_void_p]
lib_c.tsp_branch.restype = c_int
#-----------------------------------------------------------------
start_time = datetime.now()
res1 = tsp_branch(n, input_matrix, lib_c)
date_diff = (datetime.now() - start_time).total_seconds()
print(date_diff)
print('Без остова:', res1['len'], res1['steps'], res1['path'])
if 'path' in res1:
d1 = []
for i, v in enumerate(res1['path']):
d1.append([res1['path'][i-1], res1['path'][i]])
#-----------------------------------------------------------------
plt.figure(figsize=(5, 5))
graph = nx.Graph()
graph.add_nodes_from(points)
# Добавляем дуги в граф
for i in v1:
graph.add_edge(i[0], i[1], weight=i[2])
mst = nx.minimum_spanning_tree(graph)
d2 = mst.edges()
for k in mst.edges():
for i in range(len(points)):
for j in range(len(points)):
if (i < j) and (k[0] != i) and (k[0] != j) and (k[1] != i) and (k[1] != j):
if Intersection(points[k[0]][0], points[k[0]][1], points[k[1]][0], points[k[1]][1], points[i][0], points[i][1], points[j][0], points[j][1]):
input_matrix[i][j] = INF
input_matrix[j][i] = INF
# Убираем удалённые рёбра в оригинальном графе
if graph.has_edge(i, j):
graph.remove_edge(i, j)
start_time = datetime.now()
res2 = tsp_branch(n, input_matrix, lib_c)
date_diff = (datetime.now() - start_time).total_seconds()
print(date_diff)
print('С остовом:', res2['len'], res2['steps'], res2['path'])
if 'path' in res2:
d3 = []
for i, v in enumerate(res2['path']):
d3.append([res2['path'][i-1], res2['path'][i]])
# -----------------------------------------------------------------
nx.draw(graph, points, width=1, edge_color="#C0C0C0", with_labels=True, node_size=200, font_size=10)
nx.draw(graph, points, width=7, edge_color="green", edgelist=d2, style="-", node_size=0, alpha=1)
nx.draw(graph, points, width=4, edge_color="red", edgelist=d3, style="-", node_size=0, alpha=1)
nx.draw(graph, points, width=3, edge_color="blue", edgelist=d1, style=":", node_size=0)
Далее по ходу текста я буду выдвигать утверждения, которые кому-то могут показаться спорными, просто примите за факт
хах. Типа всех не согласных прошу заткнуться.
Из всей статьи меня смутила только эта формулировка) в целом я не на столько умный что б критиковать вашу статью ?
Задача коммивояжера (TSP) точное решение — метод ветвей и границ