连接所有岛屿的最低成本是多少？

有一个大小为 N x M的网格。一些细胞是由“0”表示的 岛屿，另一些是 水。每个水池上都有一个数字，表示在这个水池上建桥的成本。你必须找到所有岛屿可以连接的最低成本。如果一个单元共享一个边或一个顶点，那么它就连接到另一个单元。

什么算法可以用来解决这个问题？如果 N，M 的值非常小，比如说 NxM < = 100，那么什么可以用作蛮力方法？

示例 : 在给定的图像中，绿色单元格表示岛屿，蓝色单元格表示水，浅蓝色单元格表示应该在其上建立桥梁的单元格。因此，对于下面的图像，答案将是 17。

http://i.imgur.com/ClcboBy.png

起初，我想把所有的岛屿都标记为节点，然后用一座最短的桥把每一对岛屿连接起来。这样问题就可以简化为最小生成树，但在这种方法中，我忽略了边缘重叠的情况。在下面的图片中，任意两个岛屿之间的最短距离是 7(用黄色标记) ，所以通过使用最小生成树，答案是 14，但是答案应该是 11(用浅蓝色标记)。

image2

12982 次浏览
小开

一种用伪代码表示的暴力方法:

start with a horrible "best" answer
given an nxm map,
try all 2^(n*m) combinations of bridge/no-bridge for each cell
if the result is connected, and better than previous best, store it


return best

在 C + + 中，这可以写成

// map = linearized map; map[x*n + y] is the equivalent of map2d[y][x]
// nm = n*m
// bridged = true if bridge there, false if not. Also linearized
// nBridged = depth of recursion (= current bridge being considered)
// cost = total cost of bridges in 'bridged'
// best, bestCost = best answer so far. Initialized to "horrible"
void findBestBridges(char map[], int nm,
bool bridged[], int nBridged, int cost, bool best[], int &bestCost) {
if (nBridged == nm) {
if (connected(map, nm, bridged) && cost < bestCost) {
memcpy(best, bridged, nBridged);
bestCost = best;
}
return;
}
if (map[nBridged] != 0) {
// try with a bridge there
bridged[nBridged] = true;
cost += map[nBridged];


// see how it turns out
findBestBridges(map, nm, bridged, nBridged+1, cost, best, bestCost);


// remove bridge for further recursion
bridged[nBridged] = false;
cost -= map[nBridged];
}
// and try without a bridge there
findBestBridges(map, nm, bridged, nBridged+1, cost, best, bestCost);
}

在进行第一次调用之后(我假设您正在将您的2d 映射转换为1d 数组以便于进行复制) ，bestCost将包含最佳答案的成本，而 best将包含产生最佳答案的桥接模式。这是，无论多么缓慢。

优化措施:

  • 通过使用“桥限制”，并运行增加最大桥数量的算法，您可以在不探索整个树的情况下找到最小的答案。找到一个单桥的答案，如果它存在，将是 O (纳米)而不是 O (2 ^ 纳米)-一个巨大的改进。
  • 一旦超过 bestCost，就可以避免搜索(通过停止递归; 这也称为“修剪”) ，因为事后继续查找是没有意义的。如果情况不能好转，就别再挖了。
  • 如果在查看“坏”候选者之前先查看“好”候选者，那么上面的修剪工作会更好(实际上，单元格都是按照从左到右、从上到下的顺序查看的)。一个很好的启发式方法是将接近几个未连接组件的单元视为优先级高于未连接的单元。但是，一旦添加了启发式搜索，搜索就开始类似于 A * (并且还需要某种优先级队列)。
  • 避免重复的桥梁和不通往任何地方的桥梁。任何不断开岛屿网络的桥梁，如果删除是多余的。

像 A * 这样的通用搜索算法允许更快的搜索，尽管寻找更好的启发式算法并不是一个简单的任务。对于一种更加针对具体问题的方法，使用 斯坦纳树上的现有结果，正如@Gassa 所建议的那样，是可行的方法。然而，请注意，根据这个 加里和约翰逊的论文，在正交网格上建立 Steiner 树的问题是 NP 完全的。

如果“足够好”就足够了，那么遗传算法可能很快就能找到可以接受的解决方案，只要您添加一些关键的启发式算法来选择首选的桥梁布局。

小开

这个问题是 Steiner 树的一个变体，称为 节点加权 Steiner 树，专门化了一类图。给出一个节点加权无向图，其中一些节点是终端，节点加权 Steiner 树可以找到包含所有终端的最便宜的节点集，从而得到一个连通子图。遗憾的是，我在一些粗略的搜索中似乎找不到任何解决方案。

为了表示为一个整数规划，为每个非终端节点设置一个0-1变量，然后对于从起始图中移除连接两个终端的所有非终端节点子集，要求子集中的变量之和至少为1。这会产生太多的约束，因此您必须延迟地强制它们，使用一个有效的节点连接算法(基本上是最大流)来检测最大违反的约束。很抱歉没有提供详细信息，但是即使你已经熟悉了整数规划，实现起来也会很痛苦。

小开

考虑到这个问题发生在网格中，而且你有明确的参数，我会通过创建一个最小生成树系统地消除问题空间来处理这个问题。在这样做的时候，如果你用普里姆算法来解决这个问题，对我来说是有意义的。

不幸的是，您现在遇到了抽象网格以创建一组节点和边的问题... 因此本文的真正问题是 如何将 n × m 网格转换为{ V }和{ E } ？

这个转换过程，一目了然，可能是 NP 难的，因为组合的绝对数量可能(假设所有的水路费用是相同的)。要处理路径重叠的实例，应该考虑使用 虚拟岛屿。

当这样做时，运行 Prim 的算法，你应该达到最佳的解决方案。

我不相信动态编程可以有效地运行在这里，因为没有可观察的最优性原则。如果我们发现两个岛屿之间的最小费用，这并不一定意味着我们可以找到这两个岛屿和第三个岛屿之间的最小费用，或另一个岛屿的子集将(根据我的定义，找到 MST 通过 Prim)连接。

如果您希望代码(伪代码或其他代码)将您的网格转换为一组{ V }和{ E } ，请给我一个私有消息，我将研究拼接在一起的实现。

小开
最佳答案

为了解决这个问题，我将使用一个整数规划框架并定义三组决策变量:

  • X _ ij : 一个二进制指示变量，用于判断我们是否在水域建造了一座桥梁(i，j)。
  • Y _ ijbcn : 表示水位(i，j)是否是连接 b 岛和 c 岛的第 n 个位置的二进制指示器。
  • L _ bc : 表示岛屿 b 和 c 是否直接相连的二进制指示变量(也就是说，你只能在从 b 到 c 的桥上行走)。

对于桥梁建设成本 (西班牙语)，最小化的目标值是 sum_ij c_ij * x_ij。我们需要在模型中加入以下约束:

  • 我们需要确保 Y _ ijbcn变量是有效的。我们总是只能达到一个水广场，如果我们建立一个桥梁那里，所以 y_ijbcn <= x_ij为每个水位置(i，j)。此外，如果(i，j)不是边界岛 b，则 y_ijbc1必须等于0。最后，对于 n > 1，y_ijbcn只能在步骤 n-1中使用邻近的水位置时使用。将 N(i, j)定义为与(i，j)相邻的水方块，这相当于 y_ijbcn <= sum_{(l, m) in N(i, j)} y_lmbc(n-1)。
  • 我们需要确保 我不知道变量只有在 b 和 c 是连接的情况下才被设置。如果我们将 I(c)定义为与 c 岛接壤的位置，这可以用 l_bc <= sum_{(i, j) in I(c), n} y_ijbcn来实现。
  • 我们需要确保所有岛屿直接或间接地联系在一起。这可以通过以下方式实现: 对于岛屿的每个非空真子集 S，要求 S 中至少有一个岛与 S 补语中的至少一个岛相连，我们称之为 S’。在约束中，我们可以通过为大小 < = K/2(其中 K 是岛屿的数目)的每个非空集 S 添加约束来实现这一点，sum_{b in S} sum_{c in S'} l_bc >= 1。

对于有 k 岛、 W 水方和指定的最大路径长度 N 的问题，这是一个带有 O(K^2WN)变量和 O(K^2WN + 2^K)约束的混合整数规划模型。显然，当问题的大小变得越来越大时，这将变得难以处理，但是对于您所关心的大小，这可能是可以解决的。为了了解可伸缩性，我将使用 Pulp 包在 python 中实现它。让我们先从小一点的7 x 9地图开始，在问题的底部有3个岛屿:

import itertools
import pulp
water = {(0, 2): 2.0, (0, 3): 1.0, (0, 4): 1.0, (0, 5): 1.0, (0, 6): 2.0,
(1, 0): 2.0, (1, 1): 9.0, (1, 2): 1.0, (1, 3): 9.0, (1, 4): 9.0,
(1, 5): 9.0, (1, 6): 1.0, (1, 7): 9.0, (1, 8): 2.0,
(2, 0): 1.0, (2, 1): 9.0, (2, 2): 9.0, (2, 3): 1.0, (2, 4): 9.0,
(2, 5): 1.0, (2, 6): 9.0, (2, 7): 9.0, (2, 8): 1.0,
(3, 0): 9.0, (3, 1): 1.0, (3, 2): 9.0, (3, 3): 9.0, (3, 4): 5.0,
(3, 5): 9.0, (3, 6): 9.0, (3, 7): 1.0, (3, 8): 9.0,
(4, 0): 9.0, (4, 1): 9.0, (4, 2): 1.0, (4, 3): 9.0, (4, 4): 1.0,
(4, 5): 9.0, (4, 6): 1.0, (4, 7): 9.0, (4, 8): 9.0,
(5, 0): 9.0, (5, 1): 9.0, (5, 2): 9.0, (5, 3): 2.0, (5, 4): 1.0,
(5, 5): 2.0, (5, 6): 9.0, (5, 7): 9.0, (5, 8): 9.0,
(6, 0): 9.0, (6, 1): 9.0, (6, 2): 9.0, (6, 6): 9.0, (6, 7): 9.0,
(6, 8): 9.0}
islands = {0: [(0, 0), (0, 1)], 1: [(0, 7), (0, 8)], 2: [(6, 3), (6, 4), (6, 5)]}
N = 6


# Island borders
iborders = {}
for k in islands:
iborders[k] = {}
for i, j in islands[k]:
for dx in [-1, 0, 1]:
for dy in [-1, 0, 1]:
if (i+dx, j+dy) in water:
iborders[k][(i+dx, j+dy)] = True


# Create models with specified variables
x = pulp.LpVariable.dicts("x", water.keys(), lowBound=0, upBound=1, cat=pulp.LpInteger)
pairs = [(b, c) for b in islands for c in islands if b < c]
yvals = []
for i, j in water:
for b, c in pairs:
for n in range(N):
yvals.append((i, j, b, c, n))


y = pulp.LpVariable.dicts("y", yvals, lowBound=0, upBound=1)
l = pulp.LpVariable.dicts("l", pairs, lowBound=0, upBound=1)
mod = pulp.LpProblem("Islands", pulp.LpMinimize)


# Objective
mod += sum([water[k] * x[k] for k in water])


# Valid y
for k in yvals:
i, j, b, c, n = k
mod += y[k] <= x[(i, j)]
if n == 0 and not (i, j) in iborders[b]:
mod += y[k] == 0
elif n > 0:
mod += y[k] <= sum([y[(i+dx, j+dy, b, c, n-1)] for dx in [-1, 0, 1] for dy in [-1, 0, 1] if (i+dx, j+dy) in water])


# Valid l
for b, c in pairs:
mod += l[(b, c)] <= sum([y[(i, j, B, C, n)] for i, j, B, C, n in yvals if (i, j) in iborders[c] and B==b and C==c])


# All islands connected (directly or indirectly)
ikeys = islands.keys()
for size in range(1, len(ikeys)/2+1):
for S in itertools.combinations(ikeys, size):
thisSubset = {m: True for m in S}
Sprime = [m for m in ikeys if not m in thisSubset]
mod += sum([l[(min(b, c), max(b, c))] for b in S for c in Sprime]) >= 1


# Solve and output
mod.solve()
for row in range(min([m[0] for m in water]), max([m[0] for m in water])+1):
for col in range(min([m[1] for m in water]), max([m[1] for m in water])+1):
if (row, col) in water:
if x[(row, col)].value() > 0.999:
print "B",
else:
print "-",
else:
print "I",
print ""

这需要1.4秒运行使用默认的解决方案从纸浆包(CBC 解决方案)和输出正确的解决方案:

I I - - - - - I I
- - B - - - B - -
- - - B - B - - -
- - - - B - - - -
- - - - B - - - -
- - - - B - - - -
- - - I I I - - -

接下来，考虑问题顶部的完整问题，这是一个13 x 14的网格，包含7个岛屿:

water = {(i, j): 1.0 for i in range(13) for j in range(14)}
islands = {0: [(0, 0), (0, 1), (1, 0), (1, 1), (2, 0), (2, 1)],
1: [(9, 0), (9, 1), (10, 0), (10, 1), (10, 2), (11, 0), (11, 1),
(11, 2), (12, 0)],
2: [(0, 7), (0, 8), (1, 7), (1, 8), (2, 7)],
3: [(7, 7), (8, 6), (8, 7), (8, 8), (9, 7)],
4: [(0, 11), (0, 12), (0, 13), (1, 12)],
5: [(4, 10), (4, 11), (5, 10), (5, 11)],
6: [(11, 8), (11, 9), (11, 13), (12, 8), (12, 9), (12, 10), (12, 11),
(12, 12), (12, 13)]}
for k in islands:
for i, j in islands[k]:
del water[(i, j)]


for i, j in [(10, 7), (10, 8), (10, 9), (10, 10), (10, 11), (10, 12),
(11, 7), (12, 7)]:
water[(i, j)] = 20.0


N = 7

MIP 求解者通常相对较快地获得好的解，然后花费大量的时间试图证明解的最优性。使用与上面相同的求解代码，程序不会在30分钟内完成。但是，您可以为求解程序提供一个超时以得到一个近似解:

mod.solve(pulp.solvers.PULP_CBC_CMD(maxSeconds=120))

这产生了一个目标价值17的解决方案:

I I - - - - - I I - - I I I
I I - - - - - I I - - - I -
I I - - - - - I - B - B - -
- - B - - - B - - - B - - -
- - - B - B - - - - I I - -
- - - - B - - - - - I I - -
- - - - - B - - - - - B - -
- - - - - B - I - - - - B -
- - - - B - I I I - - B - -
I I - B - - - I - - - - B -
I I I - - - - - - - - - - B
I I I - - - - - I I - - - I
I - - - - - - - I I I I I I

为了提高您获得的解决方案的质量，您可以使用商业 MIP 解决方案(如果您在学术机构，这是免费的，否则可能不免费)。例如，下面是 Gurobi 6.0.4的性能，同样有2分钟的时间限制(尽管我们从解决方案日志中读到，解决方案的人在7秒内找到了当前的最佳解决方案) :

mod.solve(pulp.solvers.GUROBI(timeLimit=120))

这实际上找到了目标值16的解决方案，比 OP 手工找到的解决方案要好！

I I - - - - - I I - - I I I
I I - - - - - I I - - - I -
I I - - - - - I - B - B - -
- - B - - - - - - - B - - -
- - - B - - - - - - I I - -
- - - - B - - - - - I I - -
- - - - - B - - B B - - - -
- - - - - B - I - - B - - -
- - - - B - I I I - - B - -
I I - B - - - I - - - - B -
I I I - - - - - - - - - - B
I I I - - - - - I I - - - I
I - - - - - - - I I I I I I

提交答案