查找最长的相邻重复非重叠子串
Find longest adjacent repeating non-overlapping substring
（此问题与音乐无关，但我以音乐为例
一个用例。）
在音乐中，结构乐句的一种常见方式是作为音符序列
中间部分重复一次或多次。因此，短语
由介绍、循环部分和结尾组成。这是一个
示例：
[ E E E F G A F F G A F F G A F C D ]
我们可以"see"前奏是[E E E]重复部分是[F G A
F ] 结尾是 [ C D ]。所以拆分列表的方法是
[ [ E E E ] 3 [ F G A F ] [ C D ] ]
其中第一项是介绍，第二项是
重复部分是重复的，第三部分是结尾。
我需要一个算法来执行这样的拆分。
但有一点需要注意，可能有多种方法可以
拆分列表。例如，上面的列表可以拆分为：
[ [ E E E F G A ] 2 [ F F G A ] [ F C D ] ]
但这是一个更糟糕的拆分，因为前奏和结尾都更长。所以
该算法的标准是找到最大化
循环部分的长度和最小化的组合长度
介绍和结尾。这意味着
的正确拆分
[ A C C C C C C C C C A ]
是
[ [ A ] 9 [ C ] [ A ] ]
因为前奏和结尾的总长度是 2 并且长度
循环部分的是 9.
此外，虽然前奏和结尾可以为空，但只有 "true" 重复
允许。所以下面的分割是不允许的：
[ [ ] 1 [ E E E F G A F F G A F F G A F C D ] [ ] ]
将其视为找到最佳 "compression"
顺序。请注意，某些序列可能没有任何重复：
[ A B C D ]
对于这些退化的情况，任何合理的结果都是允许的。
这是我的算法实现：
def find_longest_repeating_non_overlapping_subseq(seq):
candidates = []
for i in range(len(seq)):
candidate_max = len(seq[i + 1:]) // 2
for j in range(1, candidate_max + 1):
candidate, remaining = seq[i:i + j], seq[i + j:]
n_reps = 1
len_candidate = len(candidate)
while remaining[:len_candidate] == candidate:
n_reps += 1
remaining = remaining[len_candidate:]
if n_reps > 1:
candidates.append((seq[:i], n_reps,
candidate, remaining))
if not candidates:
return (type(seq)(), 1, seq, type(seq)())
def score_candidate(candidate):
intro, reps, loop, outro = candidate
return reps - len(intro) - len(outro)
return sorted(candidates, key = score_candidate)[-1]
我不确定它是否正确，但它通过了我的简单测试
描述。它的问题是它的速度很慢。我看过
在后缀树上，但它们似乎不适合我的用例，因为
我之后的子串应该是非重叠和相邻的。
看起来您正在尝试做的几乎是 LZ77 压缩算法。您可以根据我链接到的维基百科文章中的参考实现检查您的代码。
这是我对你所说内容的实现。它与您的非常相似，但它会跳过已检查为重复先前子字符串的子字符串。
from collections import namedtuple
SubSequence = namedtuple('SubSequence', ['start', 'length', 'reps'])
def longest_repeating_subseq(original: str):
winner = SubSequence(start=0, length=0, reps=0)
checked = set()
subsequences = ( # Evaluates lazily during iteration
SubSequence(start=start, length=length, reps=1)
for start in range(len(original))
for length in range(1, len(original) - start)
if (start, length) not in checked)
for s in subsequences:
subseq = original[s.start : s.start + s.length]
for reps, next_start in enumerate(
range(s.start + s.length, len(original), s.length),
start=1):
if subseq != original[next_start : next_start + s.length]:
break
else:
checked.add((next_start, s.length))
s = s._replace(reps=reps)
if s.reps > 1 and (
(s.length * s.reps > winner.length * winner.reps)
or ( # When total lengths are equal, prefer the shorter substring
s.length * s.reps == winner.length * winner.reps
and s.reps > winner.reps)):
winner = s
# Check for default case with no repetitions
if winner.reps == 0:
winner = SubSequence(start=0, length=len(original), reps=1)
return (
original[ : winner.start],
winner.reps,
original[winner.start : winner.start + winner.length],
original[winner.start + winner.length * winner.reps : ])
def test(seq, *, expect):
print(f'Testing longest_repeating_subseq for {seq}')
result = longest_repeating_subseq(seq)
print(f'Expected {expect}, got {result}')
print(f'Test {"passed" if result == expect else "failed"}')
print()
if __name__ == '__main__':
test('EEEFGAFFGAFFGAFCD', expect=('EEE', 3, 'FGAF', 'CD'))
test('ACCCCCCCCCA', expect=('A', 9, 'C', 'A'))
test('ABCD', expect=('', 1, 'ABCD', ''))
为我传递你的所有三个例子。这似乎是一种可能有很多奇怪的边缘情况的事情，但鉴于它是一种优化的蛮力，它可能更多的是更新规范而不是修复代码本身中的错误。
这里有一种方法很明显quadratic-time，但是常数因子相对较低，因为它不构建除长度为1的子串对象之外的任何子串对象。结果是一个二元组，
bestlen, list_of_results
其中bestlen是重复相邻块的最长子串的长度，每个结果都是一个三元组，
start_index, width, numreps
表示被重复的子串是
the_string[start_index : start_index + width]
并且有 numreps 个相邻的。永远是那个
bestlen == width * numreps
问题描述有歧义。例如，考虑这个输出：
>>> crunch2("aaaaaabababa")
(6, [(0, 1, 6), (0, 2, 3), (5, 2, 3), (6, 2, 3), (0, 3, 2)])
所以它找到了 5 种方法来查看“最长”的长度为 6：
- 开头的“a”重复了 6 次。
- 开头的“aa”重复了3次。
- 最左边的“ab”重复了 3 次。
- 最左边的“ba”重复了 3 次。
- 开头的“aaa”重复了2次。
它没有 return 前奏或结尾，因为从它的作用中推断这些是微不足道的 return:
- 简介是
the_string[: start_index]。
- 结尾是
the_string[start_index + bestlen :]。
如果没有重复的相邻块，则returns
(0, [])
其他示例（来自您的post）：
>>> crunch2("EEEFGAFFGAFFGAFCD")
(12, [(3, 4, 3)])
>>> crunch2("ACCCCCCCCCA")
(9, [(1, 1, 9), (1, 3, 3)])
>>> crunch2("ABCD")
(0, [])
其工作原理的关键：假设您有相邻的重复块，每个块的宽度为 W。然后考虑将原始字符串与左移 W:
的字符串进行比较时会发生什么
... block1 block2 ... blockN-1 blockN ...
... block2 block3 ... blockN ... ...
然后你会在相同的位置得到 (N-1)*W 个连续的相等字符。但是 也 在另一个方向上起作用：如果你向左移动 W 并找到 (N-1)*W 个连续的相等字符，那么你可以推断：
block1 == block2
block2 == block3
...
blockN-1 == blockN
所以所有 N 个块必须是块 1 的重复。
所以代码反复将原始字符串（的副本）左移一个字符，然后从左到右遍历两个字符以确定最长的相等字符。这只需要一次比较一对字符。为了使“左移”有效（恒定时间），字符串的副本存储在 collections.deque.
中
编辑：update() 在许多情况下做了太多徒劳的工作；替换它。
def crunch2(s):
from collections import deque
# There are zcount equal characters starting
# at index starti.
def update(starti, zcount):
nonlocal bestlen
while zcount >= width:
numreps = 1 + zcount // width
count = width * numreps
if count >= bestlen:
if count > bestlen:
results.clear()
results.append((starti, width, numreps))
bestlen = count
else:
break
zcount -= 1
starti += 1
bestlen, results = 0, []
t = deque(s)
for width in range(1, len(s) // 2 + 1):
t.popleft()
zcount = 0
for i, (a, b) in enumerate(zip(s, t)):
if a == b:
if not zcount: # new run starts here
starti = i
zcount += 1
# else a != b, so equal run (if any) ended
elif zcount:
update(starti, zcount)
zcount = 0
if zcount:
update(starti, zcount)
return bestlen, results
使用正则表达式
[由于大小限制删除了这个]
使用后缀数组
这是我迄今为止发现的最快的速度，但仍然会被激发为 quadratic-time 行为。
请注意，是否找到重叠字符串并不重要。正如上面对 crunch2() 程序的解释（此处以次要方式进行了详细说明）：
- 给定字符串
s，长度为 n = len(s)。
- 给定整数
i 和 j 以及 0 <= i < j < n。
那么如果w = j-i，并且c是s[i:]和s[j:]之间共有的前导字符数，则块s[i:j]（长度w) 重复，从 s[i] 开始，共 1 + c // w 次。
下面的程序就是这样直接找到所有重复的相邻块，并记住那些最大总长度的块。 Returns 与 crunch2() 相同的结果，但有时顺序不同。
后缀数组简化了搜索，但很难消除它。后缀数组直接查找 <i, j> 最大 c 对，但仅限制搜索以最大化 w * (1 + c // w)。最坏的情况是 letter * number 形式的字符串，例如 "a" * 10000.
我没有给出下面 sa 模块的代码。它是 long-winded 并且后缀数组的任何实现都将计算相同的东西。 suffix_array() 的输出：
sa是后缀数组，range(n)的唯一排列使得对于range(1, n)中的所有i，s[sa[i-1]:] < s[sa[i]:] .
rank 这里没有使用
对于range(1, n)中的i，lcp[i]给出从sa[i-1]和[=79开始的后缀之间最长公共前缀的长度=].
为什么会赢？部分原因是它永远不必搜索以相同字母开头的后缀（后缀数组通过构造使它们相邻），并且检查重复块以及它是否是新的最佳块，无论块有多大或重复了多少次。如上所述，这只是 c 和 w.
上的简单算术
免责声明：后缀 arrays/trees 对我来说就像连分数：我可以在必要时使用它们，并且可以对结果感到惊奇，但它们让我头疼。敏感，敏感，敏感。
def crunch4(s):
from sa import suffix_array
sa, rank, lcp = suffix_array(s)
bestlen, results = 0, []
n = len(s)
for sai in range(n-1):
i = sa[sai]
c = n
for saj in range(sai + 1, n):
c = min(c, lcp[saj])
if not c:
break
j = sa[saj]
w = abs(i - j)
if c < w:
continue
numreps = 1 + c // w
assert numreps > 1
total = w * numreps
if total >= bestlen:
if total > bestlen:
results.clear()
bestlen = total
results.append((min(i, j), w, numreps))
return bestlen, results
一些时间
我将一个普通的英文单词文件读成了一个字符串，xs。每行一个字：
>>> len(xs)
209755
>>> xs.count('\n')
25481
所以大约 210K 字节中有大约 25K 字。这些是quadratic-time算法，所以我没想到它会很快，但是crunch2()下班后仍然是运行-并且仍然运行当我让它过夜时。
这导致我意识到它的 update() 函数可以做大量无用的工作，使算法整体上更像 cubic-time。所以我修复了它。那么：
>>> crunch2(xs)
(44, [(63750, 22, 2)])
>>> xs[63750 : 63750+50]
'\nelectroencephalograph\nelectroencephalography\nelec'
这花了大约 38 分钟，这在我预期的范围内。
正则表达式版本 crunch3() 用了不到十分之一秒！
>>> crunch3(xs)
(8, [(19308, 4, 2), (47240, 4, 2)])
>>> xs[19308 : 19308+10]
'beriberi\nB'
>>> xs[47240 : 47240+10]
'couscous\nc'
如前所述，正则表达式版本可能找不到最佳答案，但这里还有其他因素在起作用：默认情况下，“.”不匹配换行符，因此代码实际上进行了多次 tiny 搜索。文件中的 ~25K 换行符中的每一个都有效地结束了本地搜索范围。使用 re.DOTALL 标志编译正则表达式（因此换行符不会被特殊处理）：
>>> crunch3(xs) # with DOTALL
(44, [(63750, 22, 2)])
14 分钟多一点。
最后，
>>> crunch4(xs)
(44, [(63750, 22, 2)])
不到 9 分钟。构建后缀数组的时间是其中微不足道的一部分（不到一秒）。这实际上非常令人印象深刻，因为 not-always-right 暴力正则表达式版本速度较慢，尽管 运行 几乎完全“以 C 速度”。
但这是相对的。从绝对意义上来说，所有这些仍然是猪慢:-(
注意：下一节中的版本将其缩短到 5 秒以下（！）。
快得多
这个采用了完全不同的方法。对于上面的大字典示例，它在不到 5 秒的时间内得到了正确答案。
我为此感到非常自豪 ;-) 这是出乎意料的，我以前从未见过这种方法。它不进行任何字符串搜索，仅对索引集进行整数运算。
对于 letter * largish_integer 形式的输入，它仍然非常慢。实际上，只要存在至少两个（不一定是相邻的，甚至 non-overlapping!）子字符串（所考虑的当前长度）的副本，它就会一直增加 1。因此，例如，在
'x' * 1000000
它将尝试从 1 到 999999 的所有子字符串大小。
但是，看起来可以通过反复将当前大小加倍（而不是仅仅加 1）来大大改善，同时保存 类，进行混合形式的二进制搜索来定位存在重复的最大子串大小。
我将把它留作 reader 的毫无疑问的乏味练习。我在这里的工作完成了 ;-)
def crunch5(text):
from collections import namedtuple, defaultdict
# For all integers i and j in IxSet x.s,
# text[i : i + x.w] == text[j : j + x.w].
# That is, it's the set of all indices at which a specific
# substring of length x.w is found.
# In general, we only care about repeated substrings here,
# so weed out those that would otherwise have len(x.s) == 1.
IxSet = namedtuple("IxSet", "s w")
bestlen, results = 0, []
# Compute sets of indices for repeated (not necessarily
# adjacent!) substrings of length xs[0].w + ys[0].w, by looking
# at the cross product of the index sets in xs and ys.
def combine(xs, ys):
xw, yw = xs[0].w, ys[0].w
neww = xw + yw
result = []
for y in ys:
shifted = set(i - xw for i in y.s if i >= xw)
for x in xs:
ok = shifted & x.s
if len(ok) > 1:
result.append(IxSet(ok, neww))
return result
# Check an index set for _adjacent_ repeated substrings.
def check(s):
nonlocal bestlen
x, w = s.s.copy(), s.w
while x:
current = start = x.pop()
count = 1
while current + w in x:
count += 1
current += w
x.remove(current)
while start - w in x:
count += 1
start -= w
x.remove(start)
if count > 1:
total = count * w
if total >= bestlen:
if total > bestlen:
results.clear()
bestlen = total
results.append((start, w, count))
ch2ixs = defaultdict(set)
for i, ch in enumerate(text):
ch2ixs[ch].add(i)
size1 = [IxSet(s, 1)
for s in ch2ixs.values()
if len(s) > 1]
del ch2ixs
for x in size1:
check(x)
current_class = size1
# Repeatedly increase size by 1 until current_class becomes
# empty. At that point, there are no repeated substrings at all
# (adjacent or not) of the then-current size (or larger).
while current_class:
current_class = combine(current_class, size1)
for x in current_class:
check(x)
return bestlen, results
而且更快
crunch6() 将较大的字典示例放在我的盒子上不到 2 秒。它结合了 crunch4()（后缀和 lcp 数组）和 crunch5()（在一组索引中查找具有给定步幅的所有算术级数）的想法。
和crunch5()一样，这也循环了等于重复的最长子串的长度（重叠与否）多一倍的次数。因为如果没有长度为 n 的重复，则对于任何大于 n 的大小都有 none。这使得在不考虑重叠的情况下查找重复更容易，因为这是一个可利用的限制。当将“胜利”限制到相邻的重复时，它就会崩溃。例如，在“abcabc”中没有偶数长度为 1 的相邻重复，但有一个长度为 3。这似乎使任何形式的直接二分搜索都是徒劳的（是否存在大小为 n 没有说明任何其他大小的相邻重复的存在）。
'x' * n 形式的输入仍然很糟糕。有从 1 到 n-1.
的所有长度的重复
观察：我给出的所有程序都会生成所有可能的方法来分解最大长度的重复相邻块。例如，对于一个9个“x”的字符串，它表示可以将“x”重复9次或将“xxx”重复3次得到。所以，令人惊讶的是，它们也都可以用作因式分解算法 ;-)
def crunch6(text):
from sa import suffix_array
sa, rank, lcp = suffix_array(text)
bestlen, results = 0, []
n = len(text)
# Generate maximal sets of indices s such that for all i and j
# in s the suffixes starting at s[i] and s[j] start with a
# common prefix of at least len minc.
def genixs(minc, sa=sa, lcp=lcp, n=n):
i = 1
while i < n:
c = lcp[i]
if c < minc:
i += 1
continue
ixs = {sa[i-1], sa[i]}
i += 1
while i < n:
c = min(c, lcp[i])
if c < minc:
yield ixs
i += 1
break
else:
ixs.add(sa[i])
i += 1
else: # ran off the end of lcp
yield ixs
# Check an index set for _adjacent_ repeated substrings
# w apart. CAUTION: this empties s.
def check(s, w):
nonlocal bestlen
while s:
current = start = s.pop()
count = 1
while current + w in s:
count += 1
current += w
s.remove(current)
while start - w in s:
count += 1
start -= w
s.remove(start)
if count > 1:
total = count * w
if total >= bestlen:
if total > bestlen:
results.clear()
bestlen = total
results.append((start, w, count))
c = 0
found = True
while found:
c += 1
found = False
for s in genixs(c):
found = True
check(s, c)
return bestlen, results
总是很快，而且发布，但有时会出错
在生物信息学中，事实证明这是在“串联重复”、“串联阵列”和“简单序列重复”(SSR) 的名称下进行研究的。您可以搜索这些术语以找到相当多的学术论文，其中一些声称 worst-case linear-time 算法。
但这些似乎分为两个阵营：
- Linear-time 要描述的那种算法，实际上是错误的:-(
- 算法如此复杂，甚至试图将它们变成功能代码都需要投入精力:-(
在第一阵营中，有几篇论文归结为上面的crunch4()，但是没有它的内循环。我将在后面加上代码，crunch4a()。这是一个例子：
"SA-SSR: a suffix array-based algorithm for exhaustive and efficient SSR discovery in large genetic sequences."
Pickett et alia
crunch4a()总是很快，但有时会出错。事实上，它为此处出现的每个示例找到至少一个最大重复拉伸，在几分之一秒内解决了较大的字典示例，并且对 'x' * 1000000 形式的字符串没有问题。大部分时间花在构建上后缀和 lcp 数组。但它可能会失败：
>>> x = "bcdabcdbcd"
>>> crunch4(x) # finds repeated bcd at end
(6, [(4, 3, 2)])
>>> crunch4a(x) # finds nothing
(0, [])
问题是无法保证 相关 后缀在后缀数组中相邻。以“b”开头的后缀排序如下：
bcd
bcdabcdbcd
bcdbcd
要通过这种方法找到尾随的重复块，需要比较第一个和第三个。这就是为什么 crunch4() 有一个内部循环，尝试 all 对以一个共同的字母开头。相关对可以由后缀数组中任意数量的其他后缀分隔。但这也使算法成为二次时间。
# only look at adjacent entries - fast, but sometimes wrong
def crunch4a(s):
from sa import suffix_array
sa, rank, lcp = suffix_array(s)
bestlen, results = 0, []
n = len(s)
for sai in range(1, n):
i, j = sa[sai - 1], sa[sai]
c = lcp[sai]
w = abs(i - j)
if c >= w:
numreps = 1 + c // w
total = w * numreps
if total >= bestlen:
if total > bestlen:
results.clear()
bestlen = total
results.append((min(i, j), w, numreps))
return bestlen, results
O(n log n)
虽然我没有编码，但这篇论文对我来说很合适：
"Simple and Flexible Detection of Contiguous Repeats Using a Suffix Tree"
Jens Stoye, Dan Gusfield
不过，要达到 sub-quadratic 算法需要做出一些妥协。例如，"x" * n 有 n-1 个形式为 "x"*2 的子字符串，n-2 形式为 "x"*3，...，所以有 O(n**2)那些单独的。所以任何找到所有这些的算法也必然是最好的二次时间。
阅读论文了解详情；-)您正在寻找的一个概念是“原始”：我相信您只需要 S*n 形式的重复，其中 S 本身不能表示为重复较短的字符串。所以，例如，"x" * 10 是原始的，但 "xx" * 5 不是。
迈向 O(n log n)
的一步
crunch9()是我在评论中提到的“蛮力”算法的实现，来自：
"The enhanced suffix array and its applications to genome analysis"
Ibrahim et alia
http://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.93.2217&rep=rep1&type=pdf
那里的实现草图只发现“分支串联”重复，我在这里添加了代码来推断任意数量的重复，并且也包括 non-branching 重复。虽然它仍然是 O(n**2) 最坏的情况，但对于您在评论中指向的 seq 字符串，它比此处的任何其他内容都快得多。照原样，它重现（除了订单）与此处大多数其他程序相同的详尽说明。
这篇论文继续努力将最坏的情况减少到 O(n log n)，但这大大减慢了它的速度。因此，它会更加努力地战斗。我承认我失去了兴趣 ;-)
# Generate lcp intervals from the lcp array.
def genlcpi(lcp):
lcp.append(0)
stack = [(0, 0)]
for i in range(1, len(lcp)):
c = lcp[i]
lb = i - 1
while c < stack[-1][0]:
i_c, lb = stack.pop()
interval = i_c, lb, i - 1
yield interval
if c > stack[-1][0]:
stack.append((c, lb))
lcp.pop()
def crunch9(text):
from sa import suffix_array
sa, rank, lcp = suffix_array(text)
bestlen, results = 0, []
n = len(text)
# generate branching tandem repeats
def gen_btr(text=text, n=n, sa=sa):
for c, lb, rb in genlcpi(lcp):
i = sa[lb]
basic = text[i : i + c]
# Binary searches to find subrange beginning with
# basic+basic. A more gonzo implementation would do this
# character by character, never materialzing the common
# prefix in `basic`.
rb += 1
hi = rb
while lb < hi: # like bisect.bisect_left
mid = (lb + hi) // 2
i = sa[mid] + c
if text[i : i + c] < basic:
lb = mid + 1
else:
hi = mid
lo = lb
while lo < rb: # like bisect.bisect_right
mid = (lo + rb) // 2
i = sa[mid] + c
if basic < text[i : i + c]:
rb = mid
else:
lo = mid + 1
lead = basic[0]
for sai in range(lb, rb):
i = sa[sai]
j = i + 2*c
assert j <= n
if j < n and text[j] == lead:
continue # it's non-branching
yield (i, c, 2)
for start, c, _ in gen_btr():
# extend left
numreps = 2
for i in range(start - c, -1, -c):
if all(text[i+k] == text[start+k] for k in range(c)):
start = i
numreps += 1
else:
break
totallen = c * numreps
if totallen < bestlen:
continue
if totallen > bestlen:
bestlen = totallen
results.clear()
results.append((start, c, numreps))
# add non-branches
while start:
if text[start - 1] == text[start + c - 1]:
start -= 1
results.append((start, c, numreps))
else:
break
return bestlen, results
赚取奖励积分;-)
对于某些技术意义 ;-) crunch11() 是 worst-case O(n log n)。除了后缀和 lcp 数组，这还需要 rank 数组，sa 的逆：
assert all(rank[sa[i]] == sa[rank[i]] == i for i in range(len(sa)))
正如代码注释所指出的，它还依赖于 Python 3 的速度（range() 行为）。这很肤浅，但重写起来会很乏味。
描述此内容的论文有几个错误，所以如果此代码与您阅读的内容不完全匹配，请不要大惊小怪。相反，完全按照他们所说的执行，它将失败。
也就是说，代码变得复杂得令人不安，我不能保证没有错误。它适用于我尝试过的一切。
'x' * 1000000 形式的输入仍然不够快，但显然不再是 quadratic-time。例如，将同一个字母重复一百万次的字符串在将近 30 秒内完成。这里的大多数其他程序永远不会结束 ;-)
编辑：将 genlcpi() 更改为使用 semi-open Python 范围；对 crunch11() 进行了大部分外观更改；添加了“提前退出”，在最坏的情况下（例如 'x' * 1000000）可以节省大约三分之一的时间。
# Generate lcp intervals from the lcp array.
def genlcpi(lcp):
lcp.append(0)
stack = [(0, 0)]
for i in range(1, len(lcp)):
c = lcp[i]
lb = i - 1
while c < stack[-1][0]:
i_c, lb = stack.pop()
yield (i_c, lb, i)
if c > stack[-1][0]:
stack.append((c, lb))
lcp.pop()
def crunch11(text):
from sa import suffix_array
sa, rank, lcp = suffix_array(text)
bestlen, results = 0, []
n = len(text)
# Generate branching tandem repeats.
# (i, c, 2) is branching tandem iff
# i+c in interval with prefix text[i : i+c], and
# i+c not in subinterval with prefix text[i : i+c + 1]
# Caution: this pragmatically relies on that, in Python 3,
# `range()` returns a tiny object with O(1) membership testing.
# In Python 2 it returns a list - ahould still work, but very
# much slower.
def gen_btr(text=text, n=n, sa=sa, rank=rank):
from itertools import chain
for c, lb, rb in genlcpi(lcp):
origlb, origrb = lb, rb
origrange = range(lb, rb)
i = sa[lb]
lead = text[i]
# Binary searches to find subrange beginning with
# text[i : i+c+1]. Note we take slices of length 1
# rather than just index to avoid special-casing for
# i >= n.
# A more elaborate traversal of the lcp array could also
# give us a list of child intervals, and then we'd just
# need to pick the right one. But that would be even
# more hairy code, and unclear to me it would actually
# help the worst cases (yes, the interval can be large,
# but so can a list of child intervals).
hi = rb
while lb < hi: # like bisect.bisect_left
mid = (lb + hi) // 2
i = sa[mid] + c
if text[i : i+1] < lead:
lb = mid + 1
else:
hi = mid
lo = lb
while lo < rb: # like bisect.bisect_right
mid = (lo + rb) // 2
i = sa[mid] + c
if lead < text[i : i+1]:
rb = mid
else:
lo = mid + 1
subrange = range(lb, rb)
if 2 * len(subrange) <= len(origrange):
# Subrange is at most half the size.
# Iterate over it to find candidates i, starting
# with wa. If i+c is also in origrange, but not
# in subrange, good: then i is of the form wwx.
for sai in subrange:
i = sa[sai]
ic = i + c
if ic < n:
r = rank[ic]
if r in origrange and r not in subrange:
yield (i, c, 2, subrange)
else:
# Iterate over the parts outside subrange instead.
# Candidates i are then the trailing wx in the
# hoped-for wwx. We win if i-c is in subrange too
# (or, for that matter, if it's in origrange).
for sai in chain(range(origlb, lb),
range(rb, origrb)):
ic = sa[sai] - c
if ic >= 0 and rank[ic] in subrange:
yield (ic, c, 2, subrange)
for start, c, numreps, irange in gen_btr():
# extend left
crange = range(start - c, -1, -c)
if (numreps + len(crange)) * c < bestlen:
continue
for i in crange:
if rank[i] in irange:
start = i
numreps += 1
else:
break
# check for best
totallen = c * numreps
if totallen < bestlen:
continue
if totallen > bestlen:
bestlen = totallen
results.clear()
results.append((start, c, numreps))
# add non-branches
while start and text[start - 1] == text[start + c - 1]:
start -= 1
results.append((start, c, numreps))
return bestlen, results
（此问题与音乐无关，但我以音乐为例 一个用例。）
在音乐中，结构乐句的一种常见方式是作为音符序列 中间部分重复一次或多次。因此，短语 由介绍、循环部分和结尾组成。这是一个 示例：
[ E E E F G A F F G A F F G A F C D ]
我们可以"see"前奏是[E E E]重复部分是[F G A F ] 结尾是 [ C D ]。所以拆分列表的方法是
[ [ E E E ] 3 [ F G A F ] [ C D ] ]
其中第一项是介绍，第二项是 重复部分是重复的，第三部分是结尾。
我需要一个算法来执行这样的拆分。
但有一点需要注意，可能有多种方法可以 拆分列表。例如，上面的列表可以拆分为：
[ [ E E E F G A ] 2 [ F F G A ] [ F C D ] ]
但这是一个更糟糕的拆分，因为前奏和结尾都更长。所以 该算法的标准是找到最大化 循环部分的长度和最小化的组合长度 介绍和结尾。这意味着
的正确拆分[ A C C C C C C C C C A ]
是
[ [ A ] 9 [ C ] [ A ] ]
因为前奏和结尾的总长度是 2 并且长度 循环部分的是 9.
此外，虽然前奏和结尾可以为空，但只有 "true" 重复 允许。所以下面的分割是不允许的：
[ [ ] 1 [ E E E F G A F F G A F F G A F C D ] [ ] ]
将其视为找到最佳 "compression" 顺序。请注意，某些序列可能没有任何重复：
[ A B C D ]
对于这些退化的情况，任何合理的结果都是允许的。
这是我的算法实现：
def find_longest_repeating_non_overlapping_subseq(seq):
candidates = []
for i in range(len(seq)):
candidate_max = len(seq[i + 1:]) // 2
for j in range(1, candidate_max + 1):
candidate, remaining = seq[i:i + j], seq[i + j:]
n_reps = 1
len_candidate = len(candidate)
while remaining[:len_candidate] == candidate:
n_reps += 1
remaining = remaining[len_candidate:]
if n_reps > 1:
candidates.append((seq[:i], n_reps,
candidate, remaining))
if not candidates:
return (type(seq)(), 1, seq, type(seq)())
def score_candidate(candidate):
intro, reps, loop, outro = candidate
return reps - len(intro) - len(outro)
return sorted(candidates, key = score_candidate)[-1]
我不确定它是否正确，但它通过了我的简单测试 描述。它的问题是它的速度很慢。我看过 在后缀树上，但它们似乎不适合我的用例，因为 我之后的子串应该是非重叠和相邻的。
看起来您正在尝试做的几乎是 LZ77 压缩算法。您可以根据我链接到的维基百科文章中的参考实现检查您的代码。
这是我对你所说内容的实现。它与您的非常相似，但它会跳过已检查为重复先前子字符串的子字符串。
from collections import namedtuple
SubSequence = namedtuple('SubSequence', ['start', 'length', 'reps'])
def longest_repeating_subseq(original: str):
winner = SubSequence(start=0, length=0, reps=0)
checked = set()
subsequences = ( # Evaluates lazily during iteration
SubSequence(start=start, length=length, reps=1)
for start in range(len(original))
for length in range(1, len(original) - start)
if (start, length) not in checked)
for s in subsequences:
subseq = original[s.start : s.start + s.length]
for reps, next_start in enumerate(
range(s.start + s.length, len(original), s.length),
start=1):
if subseq != original[next_start : next_start + s.length]:
break
else:
checked.add((next_start, s.length))
s = s._replace(reps=reps)
if s.reps > 1 and (
(s.length * s.reps > winner.length * winner.reps)
or ( # When total lengths are equal, prefer the shorter substring
s.length * s.reps == winner.length * winner.reps
and s.reps > winner.reps)):
winner = s
# Check for default case with no repetitions
if winner.reps == 0:
winner = SubSequence(start=0, length=len(original), reps=1)
return (
original[ : winner.start],
winner.reps,
original[winner.start : winner.start + winner.length],
original[winner.start + winner.length * winner.reps : ])
def test(seq, *, expect):
print(f'Testing longest_repeating_subseq for {seq}')
result = longest_repeating_subseq(seq)
print(f'Expected {expect}, got {result}')
print(f'Test {"passed" if result == expect else "failed"}')
print()
if __name__ == '__main__':
test('EEEFGAFFGAFFGAFCD', expect=('EEE', 3, 'FGAF', 'CD'))
test('ACCCCCCCCCA', expect=('A', 9, 'C', 'A'))
test('ABCD', expect=('', 1, 'ABCD', ''))
为我传递你的所有三个例子。这似乎是一种可能有很多奇怪的边缘情况的事情，但鉴于它是一种优化的蛮力，它可能更多的是更新规范而不是修复代码本身中的错误。
这里有一种方法很明显quadratic-time，但是常数因子相对较低，因为它不构建除长度为1的子串对象之外的任何子串对象。结果是一个二元组，
bestlen, list_of_results
其中bestlen是重复相邻块的最长子串的长度，每个结果都是一个三元组，
start_index, width, numreps
表示被重复的子串是
the_string[start_index : start_index + width]
并且有 numreps 个相邻的。永远是那个
bestlen == width * numreps
问题描述有歧义。例如，考虑这个输出：
>>> crunch2("aaaaaabababa")
(6, [(0, 1, 6), (0, 2, 3), (5, 2, 3), (6, 2, 3), (0, 3, 2)])
所以它找到了 5 种方法来查看“最长”的长度为 6：
- 开头的“a”重复了 6 次。
- 开头的“aa”重复了3次。
- 最左边的“ab”重复了 3 次。
- 最左边的“ba”重复了 3 次。
- 开头的“aaa”重复了2次。
它没有 return 前奏或结尾，因为从它的作用中推断这些是微不足道的 return:
- 简介是
the_string[: start_index]。 - 结尾是
the_string[start_index + bestlen :]。
如果没有重复的相邻块，则returns
(0, [])
其他示例（来自您的post）：
>>> crunch2("EEEFGAFFGAFFGAFCD")
(12, [(3, 4, 3)])
>>> crunch2("ACCCCCCCCCA")
(9, [(1, 1, 9), (1, 3, 3)])
>>> crunch2("ABCD")
(0, [])
其工作原理的关键：假设您有相邻的重复块，每个块的宽度为 W。然后考虑将原始字符串与左移 W:
... block1 block2 ... blockN-1 blockN ...
... block2 block3 ... blockN ... ...
然后你会在相同的位置得到 (N-1)*W 个连续的相等字符。但是 也 在另一个方向上起作用：如果你向左移动 W 并找到 (N-1)*W 个连续的相等字符，那么你可以推断：
block1 == block2
block2 == block3
...
blockN-1 == blockN
所以所有 N 个块必须是块 1 的重复。
所以代码反复将原始字符串（的副本）左移一个字符，然后从左到右遍历两个字符以确定最长的相等字符。这只需要一次比较一对字符。为了使“左移”有效（恒定时间），字符串的副本存储在 collections.deque.
编辑：update() 在许多情况下做了太多徒劳的工作；替换它。
def crunch2(s):
from collections import deque
# There are zcount equal characters starting
# at index starti.
def update(starti, zcount):
nonlocal bestlen
while zcount >= width:
numreps = 1 + zcount // width
count = width * numreps
if count >= bestlen:
if count > bestlen:
results.clear()
results.append((starti, width, numreps))
bestlen = count
else:
break
zcount -= 1
starti += 1
bestlen, results = 0, []
t = deque(s)
for width in range(1, len(s) // 2 + 1):
t.popleft()
zcount = 0
for i, (a, b) in enumerate(zip(s, t)):
if a == b:
if not zcount: # new run starts here
starti = i
zcount += 1
# else a != b, so equal run (if any) ended
elif zcount:
update(starti, zcount)
zcount = 0
if zcount:
update(starti, zcount)
return bestlen, results
使用正则表达式
[由于大小限制删除了这个]
使用后缀数组
这是我迄今为止发现的最快的速度，但仍然会被激发为 quadratic-time 行为。
请注意，是否找到重叠字符串并不重要。正如上面对 crunch2() 程序的解释（此处以次要方式进行了详细说明）：
- 给定字符串
s，长度为n = len(s)。 - 给定整数
i和j以及0 <= i < j < n。
那么如果w = j-i，并且c是s[i:]和s[j:]之间共有的前导字符数，则块s[i:j]（长度w) 重复，从 s[i] 开始，共 1 + c // w 次。
下面的程序就是这样直接找到所有重复的相邻块，并记住那些最大总长度的块。 Returns 与 crunch2() 相同的结果，但有时顺序不同。
后缀数组简化了搜索，但很难消除它。后缀数组直接查找 <i, j> 最大 c 对，但仅限制搜索以最大化 w * (1 + c // w)。最坏的情况是 letter * number 形式的字符串，例如 "a" * 10000.
我没有给出下面 sa 模块的代码。它是 long-winded 并且后缀数组的任何实现都将计算相同的东西。 suffix_array() 的输出：
sa是后缀数组，range(n)的唯一排列使得对于range(1, n)中的所有i，s[sa[i-1]:] < s[sa[i]:].rank这里没有使用对于
range(1, n)中的i，lcp[i]给出从sa[i-1]和[=79开始的后缀之间最长公共前缀的长度=].
为什么会赢？部分原因是它永远不必搜索以相同字母开头的后缀（后缀数组通过构造使它们相邻），并且检查重复块以及它是否是新的最佳块，无论块有多大或重复了多少次。如上所述，这只是 c 和 w.
免责声明：后缀 arrays/trees 对我来说就像连分数：我可以在必要时使用它们，并且可以对结果感到惊奇，但它们让我头疼。敏感，敏感，敏感。
def crunch4(s):
from sa import suffix_array
sa, rank, lcp = suffix_array(s)
bestlen, results = 0, []
n = len(s)
for sai in range(n-1):
i = sa[sai]
c = n
for saj in range(sai + 1, n):
c = min(c, lcp[saj])
if not c:
break
j = sa[saj]
w = abs(i - j)
if c < w:
continue
numreps = 1 + c // w
assert numreps > 1
total = w * numreps
if total >= bestlen:
if total > bestlen:
results.clear()
bestlen = total
results.append((min(i, j), w, numreps))
return bestlen, results
一些时间
我将一个普通的英文单词文件读成了一个字符串，xs。每行一个字：
>>> len(xs)
209755
>>> xs.count('\n')
25481
所以大约 210K 字节中有大约 25K 字。这些是quadratic-time算法，所以我没想到它会很快，但是crunch2()下班后仍然是运行-并且仍然运行当我让它过夜时。
这导致我意识到它的 update() 函数可以做大量无用的工作，使算法整体上更像 cubic-time。所以我修复了它。那么：
>>> crunch2(xs)
(44, [(63750, 22, 2)])
>>> xs[63750 : 63750+50]
'\nelectroencephalograph\nelectroencephalography\nelec'
这花了大约 38 分钟，这在我预期的范围内。
正则表达式版本 crunch3() 用了不到十分之一秒！
>>> crunch3(xs)
(8, [(19308, 4, 2), (47240, 4, 2)])
>>> xs[19308 : 19308+10]
'beriberi\nB'
>>> xs[47240 : 47240+10]
'couscous\nc'
如前所述，正则表达式版本可能找不到最佳答案，但这里还有其他因素在起作用：默认情况下，“.”不匹配换行符，因此代码实际上进行了多次 tiny 搜索。文件中的 ~25K 换行符中的每一个都有效地结束了本地搜索范围。使用 re.DOTALL 标志编译正则表达式（因此换行符不会被特殊处理）：
>>> crunch3(xs) # with DOTALL
(44, [(63750, 22, 2)])
14 分钟多一点。
最后，
>>> crunch4(xs)
(44, [(63750, 22, 2)])
不到 9 分钟。构建后缀数组的时间是其中微不足道的一部分（不到一秒）。这实际上非常令人印象深刻，因为 not-always-right 暴力正则表达式版本速度较慢，尽管 运行 几乎完全“以 C 速度”。
但这是相对的。从绝对意义上来说，所有这些仍然是猪慢:-(
注意：下一节中的版本将其缩短到 5 秒以下（！）。
快得多
这个采用了完全不同的方法。对于上面的大字典示例，它在不到 5 秒的时间内得到了正确答案。
我为此感到非常自豪 ;-) 这是出乎意料的，我以前从未见过这种方法。它不进行任何字符串搜索，仅对索引集进行整数运算。
对于 letter * largish_integer 形式的输入，它仍然非常慢。实际上，只要存在至少两个（不一定是相邻的，甚至 non-overlapping!）子字符串（所考虑的当前长度）的副本，它就会一直增加 1。因此，例如，在
'x' * 1000000
它将尝试从 1 到 999999 的所有子字符串大小。
但是，看起来可以通过反复将当前大小加倍（而不是仅仅加 1）来大大改善，同时保存 类，进行混合形式的二进制搜索来定位存在重复的最大子串大小。
我将把它留作 reader 的毫无疑问的乏味练习。我在这里的工作完成了 ;-)
def crunch5(text):
from collections import namedtuple, defaultdict
# For all integers i and j in IxSet x.s,
# text[i : i + x.w] == text[j : j + x.w].
# That is, it's the set of all indices at which a specific
# substring of length x.w is found.
# In general, we only care about repeated substrings here,
# so weed out those that would otherwise have len(x.s) == 1.
IxSet = namedtuple("IxSet", "s w")
bestlen, results = 0, []
# Compute sets of indices for repeated (not necessarily
# adjacent!) substrings of length xs[0].w + ys[0].w, by looking
# at the cross product of the index sets in xs and ys.
def combine(xs, ys):
xw, yw = xs[0].w, ys[0].w
neww = xw + yw
result = []
for y in ys:
shifted = set(i - xw for i in y.s if i >= xw)
for x in xs:
ok = shifted & x.s
if len(ok) > 1:
result.append(IxSet(ok, neww))
return result
# Check an index set for _adjacent_ repeated substrings.
def check(s):
nonlocal bestlen
x, w = s.s.copy(), s.w
while x:
current = start = x.pop()
count = 1
while current + w in x:
count += 1
current += w
x.remove(current)
while start - w in x:
count += 1
start -= w
x.remove(start)
if count > 1:
total = count * w
if total >= bestlen:
if total > bestlen:
results.clear()
bestlen = total
results.append((start, w, count))
ch2ixs = defaultdict(set)
for i, ch in enumerate(text):
ch2ixs[ch].add(i)
size1 = [IxSet(s, 1)
for s in ch2ixs.values()
if len(s) > 1]
del ch2ixs
for x in size1:
check(x)
current_class = size1
# Repeatedly increase size by 1 until current_class becomes
# empty. At that point, there are no repeated substrings at all
# (adjacent or not) of the then-current size (or larger).
while current_class:
current_class = combine(current_class, size1)
for x in current_class:
check(x)
return bestlen, results
而且更快
crunch6() 将较大的字典示例放在我的盒子上不到 2 秒。它结合了 crunch4()（后缀和 lcp 数组）和 crunch5()（在一组索引中查找具有给定步幅的所有算术级数）的想法。
和crunch5()一样，这也循环了等于重复的最长子串的长度（重叠与否）多一倍的次数。因为如果没有长度为 n 的重复，则对于任何大于 n 的大小都有 none。这使得在不考虑重叠的情况下查找重复更容易，因为这是一个可利用的限制。当将“胜利”限制到相邻的重复时，它就会崩溃。例如，在“abcabc”中没有偶数长度为 1 的相邻重复，但有一个长度为 3。这似乎使任何形式的直接二分搜索都是徒劳的（是否存在大小为 n 没有说明任何其他大小的相邻重复的存在）。
'x' * n 形式的输入仍然很糟糕。有从 1 到 n-1.
观察：我给出的所有程序都会生成所有可能的方法来分解最大长度的重复相邻块。例如，对于一个9个“x”的字符串，它表示可以将“x”重复9次或将“xxx”重复3次得到。所以，令人惊讶的是，它们也都可以用作因式分解算法 ;-)
def crunch6(text):
from sa import suffix_array
sa, rank, lcp = suffix_array(text)
bestlen, results = 0, []
n = len(text)
# Generate maximal sets of indices s such that for all i and j
# in s the suffixes starting at s[i] and s[j] start with a
# common prefix of at least len minc.
def genixs(minc, sa=sa, lcp=lcp, n=n):
i = 1
while i < n:
c = lcp[i]
if c < minc:
i += 1
continue
ixs = {sa[i-1], sa[i]}
i += 1
while i < n:
c = min(c, lcp[i])
if c < minc:
yield ixs
i += 1
break
else:
ixs.add(sa[i])
i += 1
else: # ran off the end of lcp
yield ixs
# Check an index set for _adjacent_ repeated substrings
# w apart. CAUTION: this empties s.
def check(s, w):
nonlocal bestlen
while s:
current = start = s.pop()
count = 1
while current + w in s:
count += 1
current += w
s.remove(current)
while start - w in s:
count += 1
start -= w
s.remove(start)
if count > 1:
total = count * w
if total >= bestlen:
if total > bestlen:
results.clear()
bestlen = total
results.append((start, w, count))
c = 0
found = True
while found:
c += 1
found = False
for s in genixs(c):
found = True
check(s, c)
return bestlen, results
总是很快，而且发布，但有时会出错
在生物信息学中，事实证明这是在“串联重复”、“串联阵列”和“简单序列重复”(SSR) 的名称下进行研究的。您可以搜索这些术语以找到相当多的学术论文，其中一些声称 worst-case linear-time 算法。
但这些似乎分为两个阵营：
- Linear-time 要描述的那种算法，实际上是错误的:-(
- 算法如此复杂，甚至试图将它们变成功能代码都需要投入精力:-(
在第一阵营中，有几篇论文归结为上面的crunch4()，但是没有它的内循环。我将在后面加上代码，crunch4a()。这是一个例子：
"SA-SSR: a suffix array-based algorithm for exhaustive and efficient SSR discovery in large genetic sequences."
Pickett et alia
crunch4a()总是很快，但有时会出错。事实上，它为此处出现的每个示例找到至少一个最大重复拉伸，在几分之一秒内解决了较大的字典示例，并且对 'x' * 1000000 形式的字符串没有问题。大部分时间花在构建上后缀和 lcp 数组。但它可能会失败：
>>> x = "bcdabcdbcd"
>>> crunch4(x) # finds repeated bcd at end
(6, [(4, 3, 2)])
>>> crunch4a(x) # finds nothing
(0, [])
问题是无法保证 相关 后缀在后缀数组中相邻。以“b”开头的后缀排序如下：
bcd
bcdabcdbcd
bcdbcd
要通过这种方法找到尾随的重复块，需要比较第一个和第三个。这就是为什么 crunch4() 有一个内部循环，尝试 all 对以一个共同的字母开头。相关对可以由后缀数组中任意数量的其他后缀分隔。但这也使算法成为二次时间。
# only look at adjacent entries - fast, but sometimes wrong
def crunch4a(s):
from sa import suffix_array
sa, rank, lcp = suffix_array(s)
bestlen, results = 0, []
n = len(s)
for sai in range(1, n):
i, j = sa[sai - 1], sa[sai]
c = lcp[sai]
w = abs(i - j)
if c >= w:
numreps = 1 + c // w
total = w * numreps
if total >= bestlen:
if total > bestlen:
results.clear()
bestlen = total
results.append((min(i, j), w, numreps))
return bestlen, results
O(n log n)
虽然我没有编码，但这篇论文对我来说很合适：
"Simple and Flexible Detection of Contiguous Repeats Using a Suffix Tree"
Jens Stoye, Dan Gusfield
不过，要达到 sub-quadratic 算法需要做出一些妥协。例如，"x" * n 有 n-1 个形式为 "x"*2 的子字符串，n-2 形式为 "x"*3，...，所以有 O(n**2)那些单独的。所以任何找到所有这些的算法也必然是最好的二次时间。
阅读论文了解详情；-)您正在寻找的一个概念是“原始”：我相信您只需要 S*n 形式的重复，其中 S 本身不能表示为重复较短的字符串。所以，例如，"x" * 10 是原始的，但 "xx" * 5 不是。
迈向 O(n log n)
的一步
crunch9()是我在评论中提到的“蛮力”算法的实现，来自：
"The enhanced suffix array and its applications to genome analysis"
Ibrahim et alia
http://citeseerx.ist.psu.edu/viewdoc/download?doi=10.1.1.93.2217&rep=rep1&type=pdf
那里的实现草图只发现“分支串联”重复，我在这里添加了代码来推断任意数量的重复，并且也包括 non-branching 重复。虽然它仍然是 O(n**2) 最坏的情况，但对于您在评论中指向的 seq 字符串，它比此处的任何其他内容都快得多。照原样，它重现（除了订单）与此处大多数其他程序相同的详尽说明。
这篇论文继续努力将最坏的情况减少到 O(n log n)，但这大大减慢了它的速度。因此，它会更加努力地战斗。我承认我失去了兴趣 ;-)
# Generate lcp intervals from the lcp array.
def genlcpi(lcp):
lcp.append(0)
stack = [(0, 0)]
for i in range(1, len(lcp)):
c = lcp[i]
lb = i - 1
while c < stack[-1][0]:
i_c, lb = stack.pop()
interval = i_c, lb, i - 1
yield interval
if c > stack[-1][0]:
stack.append((c, lb))
lcp.pop()
def crunch9(text):
from sa import suffix_array
sa, rank, lcp = suffix_array(text)
bestlen, results = 0, []
n = len(text)
# generate branching tandem repeats
def gen_btr(text=text, n=n, sa=sa):
for c, lb, rb in genlcpi(lcp):
i = sa[lb]
basic = text[i : i + c]
# Binary searches to find subrange beginning with
# basic+basic. A more gonzo implementation would do this
# character by character, never materialzing the common
# prefix in `basic`.
rb += 1
hi = rb
while lb < hi: # like bisect.bisect_left
mid = (lb + hi) // 2
i = sa[mid] + c
if text[i : i + c] < basic:
lb = mid + 1
else:
hi = mid
lo = lb
while lo < rb: # like bisect.bisect_right
mid = (lo + rb) // 2
i = sa[mid] + c
if basic < text[i : i + c]:
rb = mid
else:
lo = mid + 1
lead = basic[0]
for sai in range(lb, rb):
i = sa[sai]
j = i + 2*c
assert j <= n
if j < n and text[j] == lead:
continue # it's non-branching
yield (i, c, 2)
for start, c, _ in gen_btr():
# extend left
numreps = 2
for i in range(start - c, -1, -c):
if all(text[i+k] == text[start+k] for k in range(c)):
start = i
numreps += 1
else:
break
totallen = c * numreps
if totallen < bestlen:
continue
if totallen > bestlen:
bestlen = totallen
results.clear()
results.append((start, c, numreps))
# add non-branches
while start:
if text[start - 1] == text[start + c - 1]:
start -= 1
results.append((start, c, numreps))
else:
break
return bestlen, results
赚取奖励积分;-)
对于某些技术意义 ;-) crunch11() 是 worst-case O(n log n)。除了后缀和 lcp 数组，这还需要 rank 数组，sa 的逆：
assert all(rank[sa[i]] == sa[rank[i]] == i for i in range(len(sa)))
正如代码注释所指出的，它还依赖于 Python 3 的速度（range() 行为）。这很肤浅，但重写起来会很乏味。
描述此内容的论文有几个错误，所以如果此代码与您阅读的内容不完全匹配，请不要大惊小怪。相反，完全按照他们所说的执行，它将失败。
也就是说，代码变得复杂得令人不安，我不能保证没有错误。它适用于我尝试过的一切。
'x' * 1000000 形式的输入仍然不够快，但显然不再是 quadratic-time。例如，将同一个字母重复一百万次的字符串在将近 30 秒内完成。这里的大多数其他程序永远不会结束 ;-)
编辑：将 genlcpi() 更改为使用 semi-open Python 范围；对 crunch11() 进行了大部分外观更改；添加了“提前退出”，在最坏的情况下（例如 'x' * 1000000）可以节省大约三分之一的时间。
# Generate lcp intervals from the lcp array.
def genlcpi(lcp):
lcp.append(0)
stack = [(0, 0)]
for i in range(1, len(lcp)):
c = lcp[i]
lb = i - 1
while c < stack[-1][0]:
i_c, lb = stack.pop()
yield (i_c, lb, i)
if c > stack[-1][0]:
stack.append((c, lb))
lcp.pop()
def crunch11(text):
from sa import suffix_array
sa, rank, lcp = suffix_array(text)
bestlen, results = 0, []
n = len(text)
# Generate branching tandem repeats.
# (i, c, 2) is branching tandem iff
# i+c in interval with prefix text[i : i+c], and
# i+c not in subinterval with prefix text[i : i+c + 1]
# Caution: this pragmatically relies on that, in Python 3,
# `range()` returns a tiny object with O(1) membership testing.
# In Python 2 it returns a list - ahould still work, but very
# much slower.
def gen_btr(text=text, n=n, sa=sa, rank=rank):
from itertools import chain
for c, lb, rb in genlcpi(lcp):
origlb, origrb = lb, rb
origrange = range(lb, rb)
i = sa[lb]
lead = text[i]
# Binary searches to find subrange beginning with
# text[i : i+c+1]. Note we take slices of length 1
# rather than just index to avoid special-casing for
# i >= n.
# A more elaborate traversal of the lcp array could also
# give us a list of child intervals, and then we'd just
# need to pick the right one. But that would be even
# more hairy code, and unclear to me it would actually
# help the worst cases (yes, the interval can be large,
# but so can a list of child intervals).
hi = rb
while lb < hi: # like bisect.bisect_left
mid = (lb + hi) // 2
i = sa[mid] + c
if text[i : i+1] < lead:
lb = mid + 1
else:
hi = mid
lo = lb
while lo < rb: # like bisect.bisect_right
mid = (lo + rb) // 2
i = sa[mid] + c
if lead < text[i : i+1]:
rb = mid
else:
lo = mid + 1
subrange = range(lb, rb)
if 2 * len(subrange) <= len(origrange):
# Subrange is at most half the size.
# Iterate over it to find candidates i, starting
# with wa. If i+c is also in origrange, but not
# in subrange, good: then i is of the form wwx.
for sai in subrange:
i = sa[sai]
ic = i + c
if ic < n:
r = rank[ic]
if r in origrange and r not in subrange:
yield (i, c, 2, subrange)
else:
# Iterate over the parts outside subrange instead.
# Candidates i are then the trailing wx in the
# hoped-for wwx. We win if i-c is in subrange too
# (or, for that matter, if it's in origrange).
for sai in chain(range(origlb, lb),
range(rb, origrb)):
ic = sa[sai] - c
if ic >= 0 and rank[ic] in subrange:
yield (ic, c, 2, subrange)
for start, c, numreps, irange in gen_btr():
# extend left
crange = range(start - c, -1, -c)
if (numreps + len(crange)) * c < bestlen:
continue
for i in crange:
if rank[i] in irange:
start = i
numreps += 1
else:
break
# check for best
totallen = c * numreps
if totallen < bestlen:
continue
if totallen > bestlen:
bestlen = totallen
results.clear()
results.append((start, c, numreps))
# add non-branches
while start and text[start - 1] == text[start + c - 1]:
start -= 1
results.append((start, c, numreps))
return bestlen, results