如何在Python中实现一个针对TSP的动态规划算法？

Question

我希望使用Python中的动态规划算法来解决TSP问题。该问题如下：

• 输入：表示为点列表的城市。例如，[(1,2), (0.3, 4.5), (9, 3)...]。城市之间的距离定义为欧几里得距离。
• 输出：此实例的旅行商最小成本，向下舍入到最近的整数。

伪代码如下：

Let A = 2-D array, indexed by subsets of {1, 2, ,3, ..., n} that contains 1 and destinations j belongs to {1, 2, 3,...n}
1. Base case:
2.          if S = {0}, then A[S, 1] = 0;
3.          else, A[S, 1] = Infinity.
4.for m = 2, 3, ..., n:   // m = subproblem size
5.    for each subset of {1, 2,...,n} of size m that contains 1:
6.        for each j belongs to S and j != 1:
7.            A[S, j] = the least value of A[S-{j},k]+the distance of k and j for every k belongs to S that doesn't equal to j
8.Return the least value of A[{1,2..n},j]+the distance between j and 1 for every j = 2, 3,...n.

我的困惑是：

如何使用子集索引列表，即如何有效地实现伪代码中的第5行。

- Pavel

---

重新阅读您的问题。itertools.combinations应该可以完成工作。 - Konstantin

@Alik 但是我该如何使用集合来索引列表？ - Pavel

@Alik 是的，我必须在Python中使用一个集合S来索引一个列表A。但我不知道该怎么做。现在我有一个想法。使用n位整数表示n个城市。如果第二位是1，则第二个城市在子集中。然后使用这个n位整数来索引列表。但我困惑如何从另一个n位整数表示的集合中获取表示为n位整数的子集。例如，1001表示第一个和第四个城市在此子集中，以便我可以使用9来索引列表。但是如何从0000中获取1001、1100、1010、0110、0101、0011作为所有包含2个城市的子集？ - Pavel

另外，我建议您使用字典而不是列表。 - Konstantin

另一个想法：您可以在字典中使用元组作为键。行 d[(1, 2, 3)] = somevalue 是完全有效的。 - Konstantin

Answer 1

您可以将集合编码为整数：整数的第i个比特位表示第i个城市的状态（即我们是否将其包含在子集中）。

例如，35₁₀ = 100011₂ 将表示城市 {1, 2, 6}。这里我从最右边的比特位开始计数，它代表城市1。

为了使用这种子集表示法索引列表，您应创建长度为 2ⁿ 的二维数组：

# Assuming n is given.
A = [[0 for i in xrange(n)] for j in xrange(2 ** n)]

这是因为使用n位整数可以表示{1,2,...,n}的每个子集（记住，每个位对应于恰好一个城市）。

这种表示方法给您带来了许多不错的可能性：

# Check whether some city (1-indexed) is inside subset.
if (1 << (i - 1)) & x:
    print 'city %d is inside subset!' % i

# In particular, checking for city #1 is super-easy:
if x & 1:
    print 'city 1 is inside subset!'

# Iterate over subsets with increasing cardinality:
subsets = range(1, 2 ** n)
for subset in sorted(subsets, key=lambda x: bin(x).count('1')):
    print subset, 
# For n=4 prints "1 2 4 8 3 5 6 9 10 12 7 11 13 14 15"

# Obtain a subset y, which is the same as x, 
# except city #j (1-indexed) is removed:
y = x ^ (1 << (j - 1))  # Note that city #j must be inside x.

这是我如何实现你的伪代码（警告：未进行测试）：

# INFINITY and n are defined somewhere above.
A = [[INFINITY for i in xrange(n)] for j in xrange(2 ** n)]
# Base case (I guess it should read "if S = {1}, then A[S, 1] = 0",
because otherwise S = {0} is not a valid index to A, according to line #1)
A[1][1] = 0
# Iterate over all subsets:
subsets = range(1, 2 ** n)
for subset in sorted(subsets, key=lambda x: bin(x).count('1')):
    if not subset & 1:
        # City #1 is not presented.
        continue
    for j in xrange(2, n + 1):
        if not (1 << (j - 1)) & subset:
            # City #j is not presented.
            continue
        for k in xrange(1, n + 1):
            if k == j or not (1 << (k - 1)) & subset:
                continue
            A[subset][j] = min(A[subset][j], A[subset ^ (1 << (j - 1))][k] + get_dist(j, k))

除了具备实现您的伪代码所需的所有功能之外，这种方法比使用元组\字典更快。

- Mike Koltsov