后缀数组
后缀数组
概念
后缀，对于一个长为 \(n\) 的字符串 \(s\)，\(s[i,n]\) 即为 \(s\) 的一个后缀 \((i \in [1,n])\)。
将 \(s\) 的所有后缀进行排序。考虑两个数组 \(sa_i\) 和 \(rk_i\)，\(sa_i\) 表示排名第 \(i\) 的后缀的编号为 \(sa_i\)，\(rk_i\) 表示编号为 \(i\) 的后缀的排名。
举个例子：对于 \(s=\)aabaaab，对它的后缀进行排序：
- \(s[4,7]\) 即
aaab； - \(s[5,7]\) 即
aab； - \(s[1,7]\) 即
aabaaab； - \(s[6,7]\) 即
ab； - \(s[2,7]\) 即
abaaab； - \(s[7,7]\) 即
b； - \(s[3,7]\) 即
baaab。
则有：
| \(i\) | \(1\) | \(2\) | \(3\) | \(4\) | \(5\) | \(6\) | \(7\) |
|---|---|---|---|---|---|---|---|
| \(sa_i\) | \(4\) | \(5\) | \(1\) | \(6\) | \(2\) | \(7\) | \(3\) |
| \(rk_i\) | \(3\) | \(5\) | \(7\) | \(1\) | \(2\) | \(4\) | \(6\) |
于是有性质：\(sa[rk[i]] = rk[sa[i]]\)。
求解
倍增法求解
常用倍增法求解。
首先考虑这样一个问题：
考虑对于两个长度均为 \(2n\) 的字符串 \(s,t\) ，若 \(s>t\) 当且仅当：
- \(s[1,n]>t[1,n]\)；
- \(s[1,n]=t[1,n]\) 且 \(s[n + 1,2n]>t[n+1,2n]\)。
于是，对于许多长度为 \(2n\) 的字符串排序，可以以每个字符串的前 \(n\) 位大小为第一关键字、后 \(n\) 位大小为第二关键字进行排序，即可以将每一个 \(s\) 的前 \(n\) 位和后 \(n\) 位作为一个二元组进行排序。
考虑运用类似的方式对于后缀进行排序。考虑倍增模型，从局部到全局逐层排序。
每一轮排序枚举 \(k \in[0,\lfloor\log_2 n\rfloor]\)，对于每个 \(k\)，将 \(rk_{sa_i}\) 作为第一关键字，\(rk_{sa_i+2^k}\) 作为第二关键字，对所有后缀进行这一轮的排序，其中 \(rk\) 和 \(sa\) 是上一轮排序得出的数组。规定若 \(sa_i+2^k>n\)，则 \(rk_{sa_i+2^k}=0\)。于是在进行 \(k\) 轮倍增后即可求出 \(sa\) 数组，时间复杂度 \(O(n \log ^2 n)\)。
时间复杂度瓶颈在于 \(O(n \log n)\) 的排序。由于是二元组排序，考虑基数排序。
基数排序优化
基数排序
基数排序是一种排序。考虑对于 \(a_i\) 进行排序，考虑从低到高枚举每一位（个位、十位、百位……）即枚举 \(k(k \in\{1,2,\dots\})\)，对于 \(a_i\) 由低至高的第 \(k\) 位进行桶排，即以本轮排序前 \(a_i\) 的大小为第二关键字，以第 \(k\) 位大小为第一关键字。如此进行 \(m=\log_{10}(\max a_i)\) 轮排序后，即可得到结果。另 \(V\) 为一次桶排的值域（若 \(a_i\) 为整数则 \(V\) 为 10），则时间复杂度为 \(O(nm+\sum V)\)。
于是，对于每一轮排序，先对第二关键字即 \(rk_{sa_i+2^k}\) 进行桶排，再对第二关键字即 \(rk_{sa_i}\) 进行桶排，排序的时间复杂度从 \(O(n \log n)\) 降为 \(O(n)\)。总时间复杂度为 \(O(n \log n)\)。
给出长为 \(n\) 的仅由英文大小写字母组成的字符串 \(s\)，输出 \(s\) 的 \(sa\) 数组。\(1 \le n \le 10^6\)。
点击查看代码
const int N = 1e6 + 10;
int n;
char s[N];
int sa[N], lsa[N], rk[N << 1], lrk[N << 1], cnt[N];
// lsa 和 lrk 是拷贝的 sa 和 rk
// rk 开 2 倍是为了让超出 n 的 i+2^k 的 rk 值为 0
signed main() {
cin >> (s + 1);
n = strlen(s + 1);
int m = 128; // 初始值域要大于 'z' 的 ASCII 码
for (int i = 1; i <= n; i++) rk[i] = s[i];
for (int i = 1; i <= n; i++) cnt[rk[i]]++;
for (int i = 1; i <= m; i++) cnt[i] += cnt[i - 1];
for (int i = n; i >= 1; i--) sa[cnt[rk[i]]--] = i; // 对第一个字符进行桶排
for (int k = 1; k < n; k <<= 1, m = n) {
// 进行两次桶排实现对于二元组的基数排序
// 对第二关键字进行桶排
memset(cnt, 0, sizeof(cnt));
for (int i = 1; i <= n; i++) lsa[i] = sa[i];
for (int i = 1; i <= n; i++) cnt[rk[lsa[i] + k]]++;
for (int i = 1; i <= m; i++) cnt[i] += cnt[i - 1];
for (int i = n; i >= 1; i--) sa[cnt[rk[lsa[i] + k]]--] = lsa[i];
// 对第一关键字进行桶排
memset(cnt, 0, sizeof(cnt));
for (int i = 1; i <= n; i++) lsa[i] = sa[i];
for (int i = 1; i <= n; i++) cnt[rk[lsa[i]]]++;
for (int i = 1; i <= m; i++) cnt[i] += cnt[i - 1];
for (int i = n; i >= 1; i--) sa[cnt[rk[lsa[i]]]--] = lsa[i];
// 用求得的 sa 数组求解 rk 数组
for (int i = 1; i <= n; i++) lrk[i] = rk[i];
for (int p = 0, i = 1; i <= n; i++) {
if (lrk[sa[i]] == lrk[sa[i - 1]] && lrk[sa[i] + k] == lrk[sa[i - 1] + k]) {
rk[sa[i]] = p; // 若排序后两者完全相同，则排名也相同
} else rk[sa[i]] = ++p; // 否则排名不同
}
}
for (int i = 1; i <= n; i++)
cout << sa[i] << " ";
return 0;
}
常数优化
- 第二关键字无需桶排
- 对于 \(sa+2^k\) 大于 \(n\) 的 \(rk\) 必然为 \(0\) 就直接放到最前面。
- 否则按 \(rk\) 从低到高枚举，如果当前的 \(sa_i>k\)，说明它能与前面的字符串拼成长度为 \(2k\) 的字符串。由于是按照排名遍历，直接插入 \(sa_i-k\) 就是正确排序。
- 值域限制，每轮结束后将值域改为最大的 \(rk\) 即可，不必改成 \(n\)。
- 若 \(rk\) 的最大值等于 \(n\) 说明此时排名已经各不相同，直接退出即可。
点击查看代码
const int N = 1e6 + 10;
int n;
char s[N];
int sa[N], lsa[N], rk[N << 1], lrk[N << 1], cnt[N];
signed main() {
cin >> (s + 1);
n = strlen(s + 1);
int m = 128;
for (int i = 1; i <= n; i++) rk[i] = s[i];
for (int i = 1; i <= n; i++) cnt[rk[i]]++;
for (int i = 1; i <= m; i++) cnt[i] += cnt[i - 1];
for (int i = n; i >= 1; i--) sa[cnt[rk[i]]--] = i;
for (int k = 1, p = 0; k < n && p < n; k <<= 1, m = p) {
int cur = 0;
for (int i = n - k + 1; i <= n; i++) lsa[++cur] = i;
for (int i = 1; i <= n; i++)
if (sa[i] > k) lsa[++cur] = sa[i] - k;
memset(cnt, 0, sizeof(cnt));
for (int i = 1; i <= n; i++) cnt[rk[lsa[i]]]++;
for (int i = 1; i <= m; i++) cnt[i] += cnt[i - 1];
for (int i = n; i >= 1; i--) sa[cnt[rk[lsa[i]]]--] = lsa[i];
for (int i = 1; i <= n; i++) lrk[i] = rk[i];
p = 0;
for (int i = 1; i <= n; i++) {
if (lrk[sa[i]] == lrk[sa[i - 1]] && lrk[sa[i] + k] == lrk[sa[i - 1] + k]) {
rk[sa[i]] = p;
} else rk[sa[i]] = ++p;
}
}
for (int i = 1; i <= n; i++)
cout << sa[i] << " ";
return 0;
}
应用
我们程两个字符串的最长公共前缀（LCP）为满足 \(s[1,x]=t[1,x]\) 的最大的 \(x\)。
\(height\) 数组
定义 \(height_i=\text{LCP}(sa_i,sa_{i-1})\)，即排名第 \(i\) 的后缀和它前一位的后缀的最长公共前缀。为方便书写记作 \(ht_i\)，规定 \(ht_1=0\)。
考虑求解。
考虑证明 \(ht_{rk_i} \geq ht_{rk_{i-1}}-1\)。
当 \(ht_{rk_{i-1}} \le 1\) 时，右边小于等于 \(0\)，上式显然成立。
否则，有 \(ht_{rk_i-1}=\text{LCP}(sa_{rk_{i-1}},sa_{rk_{i-1}-1})>1\)，即后缀 \(i-1\) 与后缀 \(sa_{rk_{i-1}-1}\) 有长度为 \(ht_{rk_i-1}\) 的 LCP。考虑用 \(aA\) 表示这个 LCP，其中 \(a\) 为一个字符，\(A\) 为一个非空字符串，\(A\) 的长度为 \(ht_{rk_i-1}-1\)。那么后缀 \(i-1\) 可以表达为 \(aAB\)，后缀 \(sa_{rk_{i-1}-1}\) 可以表示为 \(aAC\)，其中 \(B>C\)，\(C\) 可能为空，\(B\) 一定非空。
于是，后缀 \(i\) 可以表示为 \(AB\)，且存在后缀 \(sa_{rk_{i-1}-1}+1\) 为 \(AC\)。因为后缀 \(sa_{rk_i-1}\) 在大小排名上只比后缀 \(sa_{rk_i}\) 小一位，而 \(AB>AC\)，则 \(AC \le\) 后缀 \(sa_{rk_i-1}< AB\)。此时显然后缀 \(i\) 即 \(AB\) 与后缀 \(sa_{rk_i-1}\) 有共同前缀 \(A\)。
可得 \(ht_{rk_i}=\text{LCP}(sa_{rk_i},sa_{rk_i-1})=\text{LCP}(i,sa_{rk_i-1}) \geq |A|=ht_{rk_{i-1}}-1\)。得证。
于是可以依此求出 \(height\) 数组。
点击查看代码
for (int k = 0, i = 1; i <= n; i++) {
if (rk[i] == 0) continue;
if (k) k--; // k 就是 ht[rk[i-1]]，于是有 ht[rk[i]] >= k-1
while (s[i + k] == s[sa[rk[i] - 1] + k]) k++; // 暴力扩展
ht[rk[i]] = k;
}
求解 LCP
考虑证明：\(\text{LCP}(sa_i,sa_j)=\min(\text{LCP}(sa_i,sa_k),\text{LCP}(sa_k,sa_j))(i \le k \le j)\)。
考虑反证法。
假设 \(p=\min(\text{LCP}(sa_i,sa_k),\text{LCP}(sa_k,sa_j))\)，则 \(\text{LCP}(i,k)\geq p,\text{LCP}(sa_k,sa_j) \geq p\)。设后缀 \(sa_i\) 为 \(A\)，后缀 \(sa_j\) 为 \(B\)，后缀 \(sa_k\) 为 \(C\)。于是 \(A,C\) 前 \(p\) 个字符相同，\(B,C\) 前 \(p\) 个字符相同。所以 \(A,B\) 前 \(p\) 个字符相同，即 \(\text{LCP}(sa_i,sa_j) \geq p\)。
设 \(\text{LCP}(sa_i,sa_j)=q\) 且 \(q>p\)。则 \(q \geq p+1\)，且 \(A_{p+1}=B_{p+1}\)。
由于 \(p=\min(\text{LCP}(sa_i,sa_k),\text{LCP}(sa_k,sa_j))\)，所以 \(A_{p+1} \neq C_{p+1}\) 或 \(B_{p+1} \neq C_{p+1}\)。于是 \(A_{p+1} \neq B_{p+1}\)，矛盾。则假设不成立，则 \(q \le p\)，即 \(\text{LCP}(sa_i,sa_j) \le p\)。
所以 \(\text{LCP}(sa_i,sa_j)=p=\min(\text{LCP}(sa_i,sa_k),\text{LCP}(sa_k,sa_j))\)，得证。
依此，易得：
至此 LCP 问题转化为求 \(ht\) 数组的区间最小值。可利用 ST 表求解。
点击查看代码
int sa[N], lsa[N], rk[N << 1], lrk[N << 1], cnt[N], ht[N];
il void SA(int n) {
int m = 128;
for (int i = 1; i <= n; i++) rk[i] = s[i];
for (int i = 1; i <= n; i++) cnt[rk[i]]++;
for (int i = 1; i <= m; i++) cnt[i] += cnt[i - 1];
for (int i = n; i >= 1; i--) sa[cnt[rk[i]]--] = i;
for (int k = 1, p = 0; k < n && p < n; k <<= 1, m = p) {
int cur = 0;
for (int i = n - k + 1; i <= n; i++) lsa[++cur] = i;
for (int i = 1; i <= n; i++)
if (sa[i] > k) lsa[++cur] = sa[i] - k;
memset(cnt, 0, sizeof(cnt));
for (int i = 1; i <= n; i++) cnt[rk[lsa[i]]]++;
for (int i = 1; i <= m; i++) cnt[i] += cnt[i - 1];
for (int i = n; i >= 1; i--) sa[cnt[rk[lsa[i]]]--] = lsa[i];
for (int i = 1; i <= n; i++) lrk[i] = rk[i];
p = 0;
for (int i = 1; i <= n; i++) {
if (lrk[sa[i]] == lrk[sa[i - 1]] && lrk[sa[i] + k] == lrk[sa[i - 1] + k]) {
rk[sa[i]] = p;
} else rk[sa[i]] = ++p;
}
}
for (int k = 0, i = 1; i <= n; i++) {
if (rk[i] == 0) continue;
if (k) k--;
while (s[i + k] == s[sa[rk[i] - 1] + k]) k++;
ht[rk[i]] = k;
}
}
int mn[N][M], Log2[N];
il void init(int n) {
for (int i = 2; i <= n; i++) Log2[i] = Log2[i / 2] + 1;
for (int i = 1; i <= n; i++) mn[i][0] = ht[i];
for (int j = 1; j <= Log2[n]; j++)
for (int i = 1; i + (1 << j) - 1 <= n; i++)
mn[i][j] = min(mn[i][j - 1], mn[i + (1 << (j - 1))][j - 1]);
}
il int getmn(int l, int r) {
if (l > r) return INF;
int s = Log2[r - l + 1];
return min(mn[l][s], mn[r - (1 << s) + 1][s]);
}
il int get\text{LCP}(int x, int y) {
if (x > y) swap(x, y);
return getmn(x + 1, y);
}
浙公网安备 33010602011771号