「笔记」后缀数组
一些约定：
- \(\mid \sum \mid\) ：字符集大小
- \(S[i:j]\)：字符串 \(S\) 中 \(S_i \sim S_j\) 构成的子串
- \(S_1 < S_2\) ：字符串 \(S_1\) 的字典序 \(< S_2\) 。
- 后缀 \(i\) ：从 \(i\) 位置开始到字符串末尾的字串，也就是 \(S[i:n]\) 。
核心有两个数组 \(sa\) 和 \(rk\) 。
\(sa_i\) 表示将所有后缀 排序后第 \(i\) 小的后缀的编号。
\(rk_i\) 表示后缀 \(i\) 的排名
显然 \(sa\) 与 \(rk\) 互为逆数组。
不同后缀的排名必然不同（因为长度不等），所以 \(rk\) 中不会有重复的值出现。
如何求后缀数组？
- \(\mathcal O(n^2 \log n)\) ：
直接用 string + sort 排序。
- \(\mathcal O(n \log^2 n)\) ：
倍增法构造。我们以模板为例：P3809 【模板】后缀排序。
我们可以先比较出 \(2^k\) 长度的字串的大小，然后通过和相邻子串合并得到 \(2^{k+1}\) 长度的字串的大小。
对于每个长度为 \(p = 2^{k+1}\) 的串 \(S[i,i+p), S[j,j+p)\) 都能分裂成两个长度为 \(q = 2^k\) 的串 \(S[i,i+q), S[i+q,i+p)\) ，那我们在比较的时候可以先比较 \(S[i,i+q), S[j,j+q]\) 的字典序大小，在比较 \(S[i+q,i+p), S[j+q,j+p)\) 的字典序大小。
注意，在理解上，这些串的后面一部分可能是空的（并不是所有的后缀都满足 \(len \ge p\)。
#include<bits/stdc++.h>
#define LL long long
#define orz cout<<"lkp AK IOI!"<<endl
using namespace std;
const int MAXN = 1e6+5;
const int INF = 1e9+7;
const int mod = 1e9+7;
int n, w;
int rk[MAXN << 1], sa[MAXN], oldrk[MAXN << 1];
char s[MAXN];
// 对于两部分，先比较前面那一部分，再比较后面那一部分。
bool cmp(int x, int y) { return rk[x] == rk[y] ? rk[x + w] < rk[y + w] : rk[x] < rk[y]; }
int main()
{
cin >> s + 1;
n = strlen(s + 1);
// 因为每次倍增时都会根据 rk 来更新 sa, 因此仅需保证 sa 数组是一个 1~n 的排列即可。
for(int i = 1; i <= n; ++i) sa[i] = i, rk[i] = s[i];
// w 表示当前已经排完序的长度
for(w = 1; w < n; w <<= 1) {
// 按照我们规定的 cmp 排序，先比较前一部分，在比较后一部分
sort(sa + 1, sa + n + 1, cmp);
// 先复制一遍以便于后面的更新
for(int j = 1; j <= n; ++j) oldrk[j] = rk[j];
for(int j = 1, p = 0; j <= n; ++j) {
// 判断一下相邻两个**排名**的子串是否相等。
//虽然这里相等赋了相同的值，但最后一定会不一样的（因为长度不同）
// 注意这里有越界风险，所以要开两倍空间
if(oldrk[sa[j]] == oldrk[sa[j - 1]] && oldrk[sa[j] + w] == oldrk[sa[j - 1] + w]) rk[sa[j]] = p;
else rk[sa[j]] = ++p;
}
}
for(int i = 1; i <= n; ++i) printf("%d ", sa[i]); puts("");
return 0;
}
这份代码在洛谷上的测试结果是 $5.48s $。
- \(\mathcal O(n \log n)\) ：
我再来回顾一下 \(sa\) 和 \(rk\) 的定义。
\(sa_i\) 表示将所有后缀 排序后第 \(i\) 小的后缀的编号。
\(rk_i\) 表示后缀 \(i\) 的排名（别忘了，别搞反）
这个优化主要借助了 计数排序 和 基数排序 的思想。
因为计数排序的复杂度是 \(O(n+w)\)，这里我们可以认为 \(w\) 与 \(n\) 同阶。
基数排序的思想是先排其他关键字最后排第一关键字，而这个倍增排序的过程恰好可以看作双关键字排序。
所以实现的时候就是把排序换成了这两个排序即可。注意初始化。
这里依旧有一个小问题，就是排后半截时回枚举到 \(id_i + w > n\) 怎么办？我们还是把它看作空字符串，并且这也符合实际情况。（空字符串字典序最小）
/*
Work by: Suzt_ilymtics
Problem: 不知名屑题
Knowledge: 垃圾算法
Time: O(能过)
*/
#include<bits/stdc++.h>
#define LL long long
#define orz cout<<"lkp AK IOI!"<<endl
using namespace std;
const int MAXN = 1e6+5;
const int INF = 1e9+7;
const int mod = 1e9+7;
int n, m;
char s[MAXN];
int sa[MAXN << 1], rk[MAXN << 1], oldrk[MAXN << 1];
int cnt[MAXN], id[MAXN];
int main()
{
cin >> s + 1;
n = strlen(s + 1);
m = max(n, 300);
//初始化
for(int i = 1; i <= n; ++i) cnt[rk[i] = s[i]]++;
for(int i = 1; i <= m; ++i) cnt[i] += cnt[i - 1];
for(int i = 1; i <= n; ++i) sa[cnt[rk[i]]--] = i;
for(int w = 1; w < n; w <<= 1) {
// 先按照第二关键字排序
for(int i = 1; i <= m; ++i) cnt[i] = false;
for(int i = 1; i <= n; ++i) id[i] = i;
for(int i = 1; i <= n; ++i) cnt[rk[id[i] + w]] ++;
for(int i = 1; i <= m; ++i) cnt[i] += cnt[i - 1];
for(int i = n; i >= 1; --i) sa[cnt[rk[id[i] + w]]--] = id[i];
// 再按照第一关键字排序
for(int i = 1; i <= m; ++i) cnt[i] = false;
for(int i = 1; i <= n; ++i) id[i] = sa[i];
for(int i = 1; i <= n; ++i) cnt[rk[id[i]]] ++;
for(int i = 1; i <= m; ++i) cnt[i] += cnt[i - 1];
for(int i = n; i >= 1; --i) sa[cnt[rk[id[i]]]--] = id[i];
for(int i = 1; i <= n; ++i) oldrk[i] = rk[i];
for(int i = 1, p = 0; i <= n; ++i) {
if(oldrk[sa[i]] == oldrk[sa[i - 1]] && oldrk[sa[i] + w] == oldrk[sa[i - 1] + w]) rk[sa[i]] = p;
else rk[sa[i]] = ++ p;
}
}
for(int i = 1; i <= n; ++i) printf("%d ", sa[i]); puts("");
return 0;
}
然后这份代码就被优化到了 \(3.24s\)。
- 再优化
这个主要是在常数上的优化。
观察对后半截排序时的特殊性质。
考虑更新前的 \(sa_i\) 的含义：排名为 \(i\) 的长度为 \(2^{k-1}\) 的子串 \(S[sa_i : sa_i + 2^{k-1})\) 。
在本次排序中，它是长度为 \(2^k\) 的字串 \(S[sa_i - 2^{k-1}:sa_i + 2^{k-1})\) 的后半截，\(sa_i\) 的排名将作为排序的关键字。
\(S[sa_i : sa_i + 2^{k-1})\) 的排名为 \(i\) ，则第一次计排后 \(S[sa_i - 2^{k-1}:sa_i + 2^{k+1})\) 的排名必为 \(i\) 。考虑直接赋值，那么原来的第一次计排就可以写成这样：
for(int i = n; i > n - w; --i) id[++sc] = i; // 后半截为空的串
for(int i = 1; i <= n; ++i) if(sa[i] > w) id[++sc] = sa[i] - w; //根据后半截取退整个串的排名
现在被我们优化到了 \(2.83s\) 。
注意后半截为空串的情况，这样的串排名相同且最小。
-
其他常数优化
-
减小值域，值域大小 \(m\) 与复杂度有关，其最小值应为 \(rk\) 的最大值，更新 \(rk\) 时更新 \(m\) 即可。
-
减少数组嵌套的使用，减少不连续的内存访问。在第二次计数排序时，将 \(rk_{id_i}\) 存下来。
-
用 cmp 函数判断两个字串是否相同。同样是减少不连续内存的访问。
-
-
线性构建
大多数题目中，上面的倍增已经够用了。但是某些特殊题目中可以使用 DC3/SA-IS 算法实现i小奶牛相关构建后缀数组。
具体做法可以参考：诱导排序与 SA-IS 算法 与 [DC3：2009]后缀数组——处理字符串的有力工具 by. 罗穗骞。
LCP 问题
- 一些定义
\(\text{lcp(S,T)}\) 定义为字符串 \(S\) 和 \(T\) 的最长公共前缀（Longest common prefix），即最大的 \(p \le \min \{ \left\vert S \right\vert, \left\vert T \right\vert\}\)，满足 \(S_i = T_i (1 \le i \le l)\)。
下文以 \(\text{lcp(i,j)}\) 表示后缀 \(i,j\) 的最大公共前缀，并延续上文的一些概念：\(sa_i\) 排名为 \(i\) 的后缀，\(rk_i\) 后缀为 \(i\) 的排名。
这里，定义 \(\text{height}_i\) 表示排名为 \(i-1\) 和 \(i\) 的两后缀的最长公共前缀。即，
定义 \(h_i\) 表示后缀 \(i\) 和排名在 \(i\) 之前一位的后缀的最长公共后缀。即，
- 引理：LCP Lemma
此引理是证明其他引理的基础。
证明：设 \(p = \min \{ \text{lcp}(sa_i, sa_k), \text{lcp}(sa_j, sa_k) \}\) ，则有：
则 \(sa_i[1:p] = sa_j[1:p] = sa_k[1:p]\)，可得 \(\text{lcp}(sa_i, sa_k) \ge p\) 。
- 定理：LCP Theorem
由 LCP Lemma 可知显然成立。
根据这个优美的式子，求解任意两个后缀的 \(\text{lcp}\) 变为求解 \(\text{height}\) 的区间最值问题。
可通过 st 表 实现 \(\mathcal O(n \log n)\) 预处理，\(\mathcal O(1)\) 查询。
现在的问题是如何快速求 \(\text{hetght}\)。
- 推论：LCP Corollary
表示排名不相邻的两个后缀的 \(\text{lcp}\) 不超过它们之间任何相邻元素的 \(\text{lcp}\) 。
由引理 LCP lemma 显然可得。
- 引理
一个用来快速计算 \(\text{height}\) 的引理。
- 快速求 \(\text{height}\) ：
由上面的定义可知 \(h_i = \text{height}_{rk_i}\) ，秩序快速求出 \(h\)，便可 \(\mathcal O(n)\) 的获得 \(\text{height}\) 。
有上述引理可以考虑正序枚举 \(i\) 进行递推。
当 \(rk_i = 1\) 时，\(sa_{rk_i - 1}\) 不存在，特判 \(h_i = 0\)。
当 \(i=1\) 时，暴力比较出 \(\text{lcp}(i, sa_{rk_i - 1})\) ，比较次数 \(< n\) 。
若上述情况均不满足，由引理知，\(h_i = \text{lcp}(i, sa_{rk_i - 1}) \ge h_{i-1} - 1\)，两后缀前 \(h_{i-1} - 1\) 位相同，那么我们可以从第 \(h_{i-1}\) 位开始比较两后缀计算出 \(h_i\)，比较次数 \(= h_i - h_{i-1} + 2\)。
代码如下，其中 \(h_i = k\)：
for(int i = 1, k = 0; i <= n; ++i) {
if(rk[i] == 1) k = 0;
else {
if(k) --k;
int j = sa[rk[i] - 1];
while(i + k <= n && j + k <= n && S[i + k] == S[j + k]) ++k;
}
height[rk[i]] = k;
}
\(k \le n\)，且最多减 \(n\) 次，则最多会在比较中加 \(2n\) 次。因为总复杂度为 \(\mathcal O(n)\) 级别。
这里有另外一个写法，本质上一个意思：
for(int i = 1, k = 0; i <= n; ++i) {
if(k) --k;
int j = sa[rk[i] - 1];
while(i + k <= n && j + k <= n && S[i + k] == S[j + k]) ++k;
height[rk[i]] = k;
}
杂题选做
SP1811 LCS - Longest Common Substring
把 \(A\) 和 \(B\) 拼起来，中间加一个 \(\$\) ，然后跑 SA 即可，预处理 \(\text{height}\) 数组来，答案就是
其中 \(sa_i\) 和 \(sa_{i-1}\) 必须处在两个不同的字符串内。
P3181 [HAOI2016]找相同字符
如果按照上一题的思路把这两个字符串拼在一起。那么答案就是
根据上面的定理，$ \text{lcp}(sa_i, sa_j) = \displaystyle\min_{k=i+1}^{j} { \text{height}_{k} }$ 。
所以预处理一个 ST 表就可以做到 \(\mathcal O(n^2)\)。
但是这远远不够，
我们可以根据上面的引理，用单调栈来优化，实现细节看代码。时间复杂度瓶在于颈求 \(sa\)。
#include<bits/stdc++.h>
#define int long long
using namespace std;
const int MAXN = 4e5 + 10;
const int INF = 1e9 + 7;
int n, m, w, ans = 0, N;
char s[MAXN << 1];
int sa[MAXN], rk[MAXN << 1], oldrk[MAXN << 1], id[MAXN];
int cnt[MAXN];
int height[MAXN];
int sum[MAXN], val[MAXN], stc[MAXN], sc = 0;
bool cmp(int x, int y) { return rk[x] == rk[y] ? rk[x + w] < rk[y + w] : rk[x] < rk[y]; }
signed main() {
cin >> s + 1;
m = strlen(s + 1);
s[m + 1] = 'z' + 1;
cin >> s + m + 2;
n = strlen(s + 1);
N = max(300ll, n);
for(int i = 1; i <= n; ++i) cnt[rk[i] = s[i]] ++;
for(int i = 1; i <= N; ++i) cnt[i] += cnt[i - 1];
for(int i = 1; i <= n; ++i) sa[cnt[rk[i]]--] = i;
for(w = 1; w < n; w <<= 1) {
int top = 0;
for(int i = n; i > n - w; --i) id[++top] = i;
for(int i = 1; i <= n; ++i) if(sa[i] > w) id[++top] = sa[i] - w;
for(int i = 1; i <= N; ++i) cnt[i] = false;
for(int i = 1; i <= n; ++i) cnt[rk[id[i]]] ++;
for(int i = 1; i <= N; ++i) cnt[i] += cnt[i - 1];
for(int i = n; i >= 1; --i) sa[cnt[rk[id[i]]]--] = id[i];
for(int i = 1; i <= n; ++i) oldrk[i] = rk[i];
for(int i = 1, p = 0; i <= n; ++i) {
if(oldrk[sa[i]] == oldrk[sa[i - 1]] && oldrk[sa[i] + w] == oldrk[sa[i - 1] + w]) rk[sa[i]] = p;
else rk[sa[i]] = ++p;
}
}
for(int i = 1, k = 0; i <= n; ++i) {
if(k) --k;
while(i + k <= n && sa[rk[i] - 1] + k <= n && s[i + k] == s[sa[rk[i] - 1] + k]) ++k;
height[rk[i]] = k;
}
++ m;
stc[0] = 1, val[0] = 0;
for(int i = 1; i <= n; ++i) sum[i] = sum[i - 1] + (sa[i] < m);
for(int i = 1; i <= n; ++i) {
while(sc && height[i] < height[stc[sc]]) --sc;
stc[++sc] = i;
val[sc] = val[sc - 1] + (sum[i - 1] - sum[stc[sc - 1] - 1]) * height[i];
if(sa[i] > m) ans += val[sc];
}
stc[0] = 1, val[0] = 0;
for(int i = 1; i <= n; ++i) sum[i] = sum[i - 1] + (sa[i] > m);
for(int i = 1; i <= n; ++i) {
while(sc && height[i] < height[stc[sc]]) --sc;
stc[++sc] = i;
val[sc] = val[sc - 1] + (sum[i - 1] - sum[stc[sc - 1] - 1]) * height[i];
if(sa[i] < m) ans += val[sc];
}
printf("%lld\n", ans);
return 0;
}
