本文不作为教学向文章。
开篇膜拜 Pecco:算法学习笔记(84): 后缀数组 - 知乎 (zhihu.com)
有些时候,其实 \(O(n \log^2 n)\) 的排序也挺好。又短又简单。
其中 \(rk[i]\) 表示从第 \(i\) 个字符开始的后缀的排名,\(sa[i]\) 表示排名为 \(i\) 的后缀开始的位置。
#include <iostream>
#include <cstring>
#include <algorithm>
#include <vector>
const int N = 1e6 + 7;
char s[N];
int n;
int SA[N << 1], temp[2][N << 1];
int main() {
scanf("%s", s);
n = strlen(s);
for (int i(0); i < n; ++i) SA[i] = i, temp[0][i] = s[i];
for (int w(1), k(0); w < n; w <<= 1, k ^= 1) {
int *rank = temp[k], *backup = temp[k ^ 1];
std::sort(SA, SA + n, [&](const int &x, const int &y) {
return rank[x] ^ rank[y] ? rank[x] < rank[y] : rank[x + w] < rank[y + w];
});
for (int p(1), i(0); i < n; ++i) {
if (rank[SA[i]] == rank[SA[i + 1]] && rank[SA[i] + w] == rank[SA[i + 1] + w])
backup[SA[i]] = p;
else backup[SA[i]] = p++;
}
}
for (int i(0); i < n; ++i) printf("%d ", SA[i] + 1);
return 0;
}
那么接下来考虑利用基数排序优化一个 \(\log n\)。
其本质是按照上一次的 \(rk[i]\) 为第一关键字, \(rk[i + w]\) 为第二关键字排序。
而排序之后,其 \(rk\) 本身就是有序的,所以基数排序按照第二关键字排序的部分可以简化。
直接把后 \(w\) 个提到后面就是了。
然后考虑对于第一关键字桶排序,嗯,看代码。
#include <iostream>
#include <algorithm>
#include <cstring>
using namespace std;
const int N = 4e6 + 3;
char s[N];
int sa[N], tmp[2][N];
int cnt[N]; // radix_sort的计数桶
int * getSA(int n) {
int m = 128; // 值域
int *x = tmp[0], *y = tmp[1];
// 第一次排序
for (int i = 1; i <= n; ++i)
++cnt[x[i] = s[i]];
// 计数前缀和
for (int i = 1; i <= m; ++i)
cnt[i] += cnt[i - 1];
for (int i = n; i; --i)
sa[cnt[x[i]]--] = i;
// 开始之后的排序,对于 (rank[x], rank[x + k]) 进行排序。
for (int p, k = 1; k < n; k <<= 1) {
p = 0;
// 由于对于 [n - k + 1, n],其 rank[x + k] 一定为0,故会被放在最前面
for (int i = n - k + 1; i <= n; ++i) y[++p] = i;
// 很明显,rk已经是有序的,前k名一定是已经被放入的(rank[x+k] = 0)
// 所以,我们只需要将后面的 n-k 个按顺序加入即可
// (对于此次的第二关键词排序即使对上一次的第一关键词排序,已经是有序的,所以直接加进去)
for (int i = 1; i <= n; ++i) {
if (sa[i] > k) y[++p] = sa[i] - k;
}
// 现在开始对于第一关键字排序
// 清空计数桶,实际上可能不需要?
for (int i = 0; i <= m; ++i) cnt[i] = 0;
// x[i] 实际上就是上一次的rank,也就是第一关键字
// 实际上这里可以写作 ++cnt[x[y[i]]]; 效果是一样的
for (int i = 1; i <= n; ++i) ++cnt[x[i]];
for (int i = 1; i <= m; ++i) cnt[i] += cnt[i - 1];
// i外套了一层y[i], 实际上就是按照第二关键词排好的顺序放回原数组中。
for (int i = n; i; --i)
sa[cnt[x[y[i]]]--] = y[i], y[i] = 0;
// 为了避免memset,类似于滚动数组的思想
swap(x, y);
// 此时y也就是之前的rank,那么我们现在要取得此时的rank
p = 0;
for (int i = 1; i <= n; ++i) {
x[sa[i]] = (y[sa[i]] == y[sa[i - 1]] && y[sa[i] + k] == y[sa[i - 1] + k]) ? p : ++p;
}
if (p >= n) break; // 完成排序,每一个都不一样了。
m = p; // 改变值域,最终为n
}
return x; // 最终的rank
}
int main() {
cin.tie(0)->sync_with_stdio(false);
std::cin >> (s + 1);
int n = strlen(s + 1);
int * rk = getSA(n);
for (int i = 1; i <= n; ++i) {
std::cout << sa[i] << ' ';
} std::cout << '\n';
return 0;
for (int i = 1; i <= n; ++i) {
std::cout << rk[i] << ' ';
} std::cout << '\n';
return 0;
}
观察排序后的字符串:
来自 Pecco 大佬。
可以看到,相邻的后缀之间可能有一些共同前缀。
那么我们可以利用后缀数组找出 Longest Common Prefix,也就是所谓的 LCP。
假设 \(H[i] = lcp(sa[i], sa[i - 1])\)。
以及 \(h[i] = H[rk[i]]\)。
有定理 \(h[i] \ge h[i - 1] - 1\)。
其文本理解为:对于后缀 \(i\),与字典序在其前面一个的后缀的的最长公共前缀的长度,大于上一个后缀与其字典序前一个后缀的LCP减一。
还是很抽象……证明来自 OI wiki
于是我们可以 \(O(n)\) 的求出 \(H\) 数组了。
void getHt(int n) {
for (int i = 1, k = 0; i <= n; ++i) {
if (k > 0) --k;
while (s[i + k] == s[sa[rk[i] - 1] + k])
++k;
H[rk[i]] = k;
}
}
求一个字符串的本质不同子串个数就可以利用 \(H\)。
其式子为:
测试:不同子串个数 - 洛谷
求两个串 \(i, j\) 的最长公共前缀。
利用 \(H\) 转化为 RMQ
问题
\(\to lcp(i, j) = \min_{k = rk[i] + 1}^{rk[j]} H[i]\)
感性理解为字典序上相差越大,相同越小。如此而已。
其实感觉这个没啥用,完全可以被哈希二分水过,而且时空常数还小很多……
题目:最长公共前缀 - 洛谷
哈希第 \(k\) 大。
现在我还只会不同位置的相同子串算作一个的情况。
按照 \(sa[i]\) 中的东西向后,每一个子串对字典序的贡献为 \(n - sa[i] + 1 - H[i]\)。
所以依次遍历就是了。
for (int i = 1; i <= n; ++i) {
int difc = n - sa[i] + 1 - H[i];
if (k > difc) {
k -= difc; continue;
}
for (int j = sa[i]; j <= sa[i] + k + H[i] - 1; ++j)
putchar(A[j - 1]);
k = 0;
}
if (k > 0) puts("-1");