[这与Minimum set cover有关]
我想用计算机解决以下 n 的小难题。考虑所有长度为 n 的 2^n 个二进制向量。对于每一个,你恰好删除了 n/3 个位,留下一个二进制向量长度 2n/3(假设 n 是 3 的整数倍)。目标是选择要删除的位,以尽量减少最后保留的长度为 2n/3 的不同二进制向量的数量。
例如,对于 n = 3,最佳答案是 2 个不同的向量 11 和 00。对于 n = 6,它是 4,对于 n = 9,它是 6,对于 n = 12,它是 10。
我之前曾尝试将此问题作为以下类型的最小集覆盖问题来解决。所有列表仅包含 1 和 0。
我说列表 A
覆盖列表 B
如果你可以通过插入从 A
生成 B
恰好是 x
个符号。
考虑所有长度为 n
的 1 和 0 的 2^n 列表,并设置 x = n/3
。我想计算一组最小长度的列表 2n/3
来覆盖它们。 David Eisenstat 提供了将这个最小集合覆盖问题转换为混合整数规划问题的代码,该问题可以输入 CPLEX(或开源的 http://scip.zib.de/)。
from collections import defaultdict
from itertools import product, combinations
def all_fill(source, num):
output_len = (len(source) + num)
for where in combinations(range(output_len), len(source)):
poss = ([[0, 1]] * output_len)
for (w, s) in zip(where, source):
poss[w] = [s]
for tup in product(*poss):
(yield tup)
def variable_name(seq):
return ('x' + ''.join((str(s) for s in seq)))
n = 12
shortn = ((2 * n) // 3)
x = (n // 3)
all_seqs = list(product([0, 1], repeat=shortn))
hit_sets = defaultdict(set)
for seq in all_seqs:
for fill in all_fill(seq, x):
hit_sets[fill].add(seq)
print('Minimize')
print(' + '.join((variable_name(seq) for seq in all_seqs)))
print('Subject To')
for (fill, seqs) in hit_sets.items():
print(' + '.join((variable_name(seq) for seq in seqs)), '>=', 1)
print('Binary')
for seq in all_seqs:
print(variable_name(seq))
print('End')
问题是,如果您设置 n=15,那么它输出的实例对于我能找到的任何求解器来说都太大了。有没有更有效的方法来解决这个问题,以便我可以解决 n=15 甚至 n = 18?
最佳答案
这并不能解决您的问题(好吧,不够快),但是您没有得到很多想法,其他人可能会在这里找到有用的东西。
这是一个简短的纯 Python 3 程序,使用带有一些贪婪排序启发式的回溯搜索。它非常快速地解决了 N = 3、6 和 9 个实例。它也很快找到了 N=12 大小为 10 的封面,但显然需要更长的时间来耗尽搜索空间(我没时间了,它仍在运行)。对于 N=15,初始化时间已经很慢了。
比特串在这里由普通的 N 位整数表示,因此占用的存储空间很小。这是为了简化用更快的语言重新编码的过程。它确实大量使用整数集,但没有使用其他“高级”数据结构。
希望这对某人有帮助!但很明显,随着 N 的增加,可能性的组合爆炸式增长确保了如果不深入研究问题的数学原理,任何事情都不会“足够快”。
def dump(cover):
for s in sorted(cover):
print(" {:0{width}b}".format(s, width=I))
def new_best(cover):
global best_cover, best_size
assert len(cover) < best_size
best_size = len(cover)
best_cover = cover.copy()
print("N =", N, "new best cover, size", best_size)
dump(best_cover)
def initialize(N, X, I):
from itertools import combinations
# Map a "wide" (length N) bitstring to the set of all
# "narrow" (length I) bitstrings that generate it.
w2n = [set() for _ in range(2**N)]
# Map a narrow bitstring to all the wide bitstrings
# it generates.
n2w = [set() for _ in range(2**I)]
for wide, wset in enumerate(w2n):
for t in combinations(range(N), X):
narrow = wide
for i in reversed(t): # largest i to smallest
hi, lo = divmod(narrow, 1 << i)
narrow = ((hi >> 1) << i) | lo
wset.add(narrow)
n2w[narrow].add(wide)
return w2n, n2w
def solve(needed, cover):
if len(cover) >= best_size:
return
if not needed:
new_best(cover)
return
# Find something needed with minimal generating set.
_, winner = min((len(w2n[g]), g) for g in needed)
# And order its generators by how much reduction they make
# to `needed`.
for g in sorted(w2n[winner],
key=lambda g: len(needed & n2w[g]),
reverse=True):
cover.add(g)
solve(needed - n2w[g], cover)
cover.remove(g)
N = 9 # CHANGE THIS TO WHAT YOU WANT
assert N % 3 == 0
X = N // 3 # number of bits to exclude
I = N - X # number of bits to include
print("initializing")
w2n, n2w = initialize(N, X, I)
best_cover = None
best_size = 2**I + 1 # "infinity"
print("solving")
solve(set(range(2**N)), set())
N=9 的示例输出:
initializing
solving
N = 9 new best cover, size 6
000000
000111
001100
110011
111000
111111
跟进
对于 N=12 这最终完成,确认最小覆盖集包含 10 个元素(它在开始时很快就找到了)。我没有计时,但至少用了 5 个小时。
这是为什么呢?因为它接近脑死亡 ;-) 完全 天真的搜索将尝试 256 个 8 位短字符串的所有子集。有 2**256 个这样的子集,大约 1.2e77 - 它不会在宇宙的预期生命周期内完成 ;-)
这里的排序技巧首先检测到“全 0”和“全 1”短字符串必须在任何覆盖集中,因此选择它们。这让我们“只”看到剩下的 254 个短字符串。然后贪婪的“选择覆盖最多的元素”策略很快找到一个总共有 11 个元素的覆盖集,不久之后找到一个有 10 个元素的覆盖集。这恰好是最优的,但需要很长时间才能穷尽所有其他可能性。
在这一点上,任何对达到 10 个元素的覆盖集的尝试都将中止(那时它不可能小于 10 个元素!)。如果这也完全天真地完成,则需要尝试添加(到“全 0”和“全 1”字符串)剩余 254 个的所有 8 元素子集,并且 254-choose-8 大约是 3.8e14。比 1.2e77 小得多 - 但仍然太大而不实用。了解代码如何设法做得比这好得多,这是一个有趣的练习。提示:与本题数据有很大关系。
工业级求解器无比复杂。我对这个简单的小程序在较小的问题实例上的表现如此出色感到惊喜!它很幸运。
但是对于 N=15,这种简单的方法是没有希望的。它很快找到了一个包含 18 个元素的封面,但至少在几个小时内没有任何明显的进展。在内部,它仍在使用包含数百(甚至数千)个元素的 needed
集,这使得 solve()
的主体非常昂贵。它还有 2**10 - 2 = 1022 个短字符串需要考虑,1022-choose-16 大约是 6e34。我不认为即使将这段代码速度提高一百万倍,它也不会有明显的帮助。
虽然尝试很有趣 :-)
和一个小的重写
此版本在 N=12 的完整运行中至少快 6 倍,只需提前一级切断无用的搜索。还可以加快初始化速度,并通过将 2**N w2n
集合更改为列表(不使用集合操作)来减少内存使用。 N=15 仍然没有希望,尽管 :-(
def dump(cover):
for s in sorted(cover):
print(" {:0{width}b}".format(s, width=I))
def new_best(cover):
global best_cover, best_size
assert len(cover) < best_size
best_size = len(cover)
best_cover = cover.copy()
print("N =", N, "new best cover, size", best_size)
dump(best_cover)
def initialize(N, X, I):
from itertools import combinations
# Map a "wide" (length N) bitstring to the set of all
# "narrow" (length I) bitstrings that generate it.
w2n = [set() for _ in range(2**N)]
# Map a narrow bitstring to all the wide bitstrings
# it generates.
n2w = [set() for _ in range(2**I)]
# mask[i] is a string of i 1-bits
mask = [2**i - 1 for i in range(N)]
for t in combinations(range(N), X):
t = t[::-1] # largest i to smallest
for wide, wset in enumerate(w2n):
narrow = wide
for i in t: # delete bit 2**i
narrow = ((narrow >> (i+1)) << i) | (narrow & mask[i])
wset.add(narrow)
n2w[narrow].add(wide)
# release some space
for i, s in enumerate(w2n):
w2n[i] = list(s)
return w2n, n2w
def solve(needed, cover):
if not needed:
if len(cover) < best_size:
new_best(cover)
return
if len(cover) >= best_size - 1:
# can't possibly be extended to a cover < best_size
return
# Find something needed with minimal generating set.
_, winner = min((len(w2n[g]), g) for g in needed)
# And order its generators by how much reduction they make
# to `needed`.
for g in sorted(w2n[winner],
key=lambda g: len(needed & n2w[g]),
reverse=True):
cover.add(g)
solve(needed - n2w[g], cover)
cover.remove(g)
N = 9 # CHANGE THIS TO WHAT YOU WANT
assert N % 3 == 0
X = N // 3 # number of bits to exclude
I = N - X # number of bits to include
print("initializing")
w2n, n2w = initialize(N, X, I)
best_cover = None
best_size = 2**I + 1 # "infinity"
print("solving")
solve(set(range(2**N)), set())
print("best for N =", N, "has size", best_size)
dump(best_cover)
关于python - 如何加速代码来解决位删除难题,我们在Stack Overflow上找到一个类似的问题: https://stackoverflow.com/questions/19314465/