再学串串(七):哈希,倍增 & 诱导排序与 SA-IS 算法

· · 算法·理论

过去的章节

  1. 再学串串(一):我真的学会 KMP 了吗?评「clx201022 式的傲慢」
  2. 再学串串(二):AC 自动机不是自动 AC 机
  3. 再学串串(三):被构造自动机上计数 DP 题吓晕
  4. 再学串串(四):后缀是后缀的后缀是后缀的后缀
  5. 再学串串(五):谁会不喜欢可爱的小马(拉车)呢?
  6. 再学串串(六):回文自动机动自文回:)六(串串学再

目前已同步在博客园更新。

如果觉得写得好麻烦动手点个赞谢谢喵,本喵会很开心的喵。

第七站:后缀数组(Suffix Array)

后缀数组包含两个数组:sark

给定一个串 T=c_1 c_2 c_3\dots c_n,对其所有后缀排序,每个后缀用其起始点标记。rk_i 即为 T[i,n]T 所有后缀中的大小关系,sa_i 表示第 i 小的后缀。两者均为 1\sim n 的排列且互逆,即 sa_{rk_i}=i=rk_{sa_i}

因为按字典序排序所以也算作是后缀排序。

一般会在一些和字典序关联较强的题目中出现。

朴素求法

简单的问题复杂化,归约到串 T 的任意 n 个子串排序。

预处理哈希可以做到 \mathcal{O}(\log n) 的字符串比较,总时间复杂度 \mathcal{O}(n\log^2 n)

:::info[63 pts]

#include<bits/stdc++.h>
using namespace std;
typedef long long ll;
constexpr int maxn = 1e6, sigma = 'z';
constexpr ll mod = 1e9 + 7, base = 13331;
ll pw[maxn + 1], h[maxn + 1];
string s;
int sa[maxn + 1], rk[maxn + 1];
int main()
{
    cin.tie(nullptr)->sync_with_stdio(false);
    cin >> s;
    int n = s.size();
    s = " " + s;
    pw[0] = 1;
    for (int i = 1; i <= n; i++) pw[i] = pw[i - 1] * base % mod;
    for (int i = 1; i <= n; i++) h[i] = (h[i - 1] * base % mod + s[i]) % mod;
    auto geth = [&](int l, int r) -> ll
    {
        return (mod + h[r] - (h[l - 1] * pw[r - l + 1]) % mod) % mod;
    };
    auto cmp = [&](int x, int y) -> bool // 后缀 x < 后缀 y 吗?
    {
        if (x == y) return false;
        int len = 0, r = n - max(x, y) + 1; // LCP 长度
        while (len < r)
        {
            int mid = (len + r + 1) >> 1;
            if(geth(x, x + mid - 1) == geth(y, y + mid - 1)) len = mid;
            else r = mid - 1;
        }
        // len == r == LCP 长度
        if (max(x, y) + len == n + 1)
        {
            // 在相同的部分内全部相等
            // n - x + 1 < n - y + 1
            return x > y; // 那就看后缀长度了
        }
        return s[x + len] < s[y + len]; // 直接比第一个不同的字符大小
    };
    for (int i = 1; i <= n; i++) sa[i] = i, rk[i] = 1;
    sort(sa + 1, sa + n + 1, cmp);
    for (int i = 1; i <= n; i++) rk[sa[i]] = i;
    for (int i = 1; i <= n; i++) cout << sa[i] << " \n"[i == n];
    return 0;
}

:::

很不幸,这就是朴素做法的极限了,这是因为我们把它归约为了一个更复杂的问题。(如果字符串 T 任意 n 个子串排序有更低复杂度不要打我,这是我会的最低复杂度了)

能不能利用一些后缀的性质呢?

倍增做法

从哈希 \mathcal{O}(\log n) 比较中可以发现一个要点:字典序的比较只考虑相同段^lcp后的第一个字符。

引理 1.1

A<B$ 当且仅当 $A_{\operatorname{LCP}(A,B)+1} < B_{\operatorname{LCP}(A,B)+1}

比较两个后缀时,可以将其都分成长度相同的两段比较(超出的部分用 \epsilon 补齐,设 \epsilon 的字典序小于任何其他字符)。

那么只需要先比较前面一段子串,前一段子串相等再比较后一段子串即可。

但怎么求更小子串之间的大小关系呢?很简单,因为主串是同一个所以只要求每个起点开头的相同长度的所有子串之间的大小关系即可,而这个问题是可以递归求直到只剩单个字符的。

考虑怎样劈成更小的段。显然从中间劈开会使得要比较大小的每一段长度相等,这样会非常好处理。

经典递归转迭代,从单个字符不停倍增直到涵盖整个主串为止。

核心是上一个长度为 l 的迭代求得的大小关系在求下一个长度为 2l 迭代比较时是可复用的。

:::error[快速幂错误示范] ```cpp typedef unsigned long long ull; ull qpow(ull b, ull e, ull p) { if (!e) return 1; ull a = qpow(b, e / 2, p) * qpow(b, e / 2, p) % p; if (e & 1) (a *= b) %= p; return a; } ``` ::: :::success[快速幂正确示范] ```cpp typedef unsigned long long ull; ull qpow(ull b, ull e, ull p) { if (!e) return 1; ull a = qpow(b, e / 2, p); (a *= a) %= p; if (e & 1) (a *= b) %= p; return a; } ``` ::: 以 $\texttt{qwqwqaqwq}$ 为例。 ![](https://cdn.luogu.com.cn/upload/image_hosting/ax3qfpji.png) (排序方式为从小到大,虚线表示同上相等,省略) 左边是 $rk$,右边是 $sa$。 长度倍增到 $l\ge|T|$ 时必然可以涵盖 $T$ 的所有后缀所以倍增次数是 $\mathcal{O}(\log n)$ 的。 图示用到了一个小优化:当 $rk$ 为一个 $1\sim n$ 的排列时,终止倍增。因为此时已经可以根据现有信息判断各个后缀间的大小了。 注意到每次排序只需要用到前一次的 $rk$,而 $rk$ 因为是排列所以值域是 $\mathcal{O}(n)$ 的,可以使用复杂度与值域有关的计数排序。 为什么不说基数排序呢?事实上对于第二关键字的排序可以用一些神秘的手段简化掉,每次倍增进行一次计数排序即可。 时间复杂度 $\mathcal{O}(n\log n)$。 :::success[[倍增 + 计数排序](https://www.luogu.com.cn/record/279839506)] ```cpp #include<bits/stdc++.h> using namespace std; typedef long long ll; constexpr int maxn = 1e6, sigma = 'z'; string s; array<int, maxn + 1> sa, cnt, id; array<int, maxn + 1> rk, oldrk; int main() { cin.tie(nullptr)->sync_with_stdio(false); cin >> s; int n = s.length(); s = " " + s; for (int i = 1; i <= n; ++i) ++cnt[rk[i] = s[i]]; // 字符集比较小不用离散化 for (int i = 1; i <= sigma; ++i) cnt[i] += cnt[i - 1]; // 前缀和 for (int i = n; i >= 1; --i) sa[cnt[rk[i]]--] = i; // 初始计数排序 for (int delta = 1, p = 0; delta < n; delta <<= 1) // 倍增步长 { int cur = 0; for (int i = n - delta + 1; i <= n; i++) id[++cur] = i; // 超出的直接在前面 for (int i = 1; i <= n; ++i) if (sa[i] > delta) id[++cur] = sa[i] - delta; // 剩下的按原顺序放入 // 先提前预处理对第二关键字的排序 for (int i = 1; i <= p; ++i) cnt[i] = 0; // 清空 for (int i = 1; i <= n; ++i) ++cnt[rk[i]]; // 计数 for (int i = 1; i <= p; ++i) cnt[i] += cnt[i - 1]; // 前缀和 for (int i = n; i >= 1; --i) sa[cnt[rk[id[i]]]--] = id[i]; // 放(因为是 cnt 递减所以倒着枚举) p = 0; oldrk = rk; auto key_equal = [&](int i, int j) -> bool { if (oldrk[sa[i]] != oldrk[sa[j]]) return false; if ((sa[i] + delta > n) && (sa[j] + delta > n)) return true; if ((sa[i] + delta > n) || (sa[j] + delta > n)) return false; if (oldrk[sa[i] + delta] != oldrk[sa[j] + delta]) return false; return true; }; rk[sa[1]] = ++p; for (int i = 2; i <= n; ++i) { if (!key_equal(i, i - 1)) p++; rk[sa[i]] = p; } if (p == n) break; } for (int i = 1; i <= n; i++) cout << sa[i] << " \n"[i == n]; return 0; } ``` ::: 这个思路还可以拓展到一些图上求字典序最大的问题,具体可见《[钻石笑传之 Collect Collect 璧](https://www.luogu.com.cn/article/qemb738u)》。 ### 性质分析 有没有可能通过一些性质的分析,优化到 $\Theta(N)$ 呢? 观察前面两种做法,核心都是在字典序上做文章。字典序的本质是什么呢?找到两个串的 LCP[^lcp],然后比较后一位。 ![](https://cdn.luogu.com.cn/upload/image_hosting/dmm8c6l3.png) 如果主串 $T$ 的每一位都互不相同,那计数排序做法不用倍增一开始就可以得到结果,时间复杂度 $\Theta(N)$。 如同 KMP 匹配失败不耗时反而是匹配成功耗时,求 SA 也是后缀的前缀不同不耗时反而后缀的前缀相同耗时。每一位都相同的情况下倍增做法有最坏复杂度 $\Theta(N\log N)$。 若读者阅读过本系列的 SAM 部分,「后缀的前缀」这个子串对您来说一定不陌生。没错,后缀自动机和后缀数组的关系不止是有一个「后缀」的前缀。事实上,您可以使用后缀自动机 $\Theta(N)$ 地辅助求后缀数组。这一部分会放到 SP 里细说,敬请期待。 相邻后缀的关系是最近的,因此考虑先尝试得出相邻后缀的大小关系,再做下一步打算。 ![](https://cdn.luogu.com.cn/upload/image_hosting/gmefkz6u.png) > **引理 1.2** > > 若 $A'=c+A,B'=c+B$,则 $A'<B'$ 当且仅当 $A<B$。 发现相邻的 $(i,i+1)$ 后缀对之间的大小关系可以直接 $\Theta(1)$ 求出或递归地转移到 $(i+1,i+2)$ 之间的大小关系。 将 $T[i,n]<T[i+1,n]$ 的 $i$ 叫做 $S$ 型后缀,$T[i,n]>T[i+1,n]$ 的 $i$ 叫做 $L$ 型后缀。(显然长度不同两者不可能相等) 设 $n+1$ 为 $S$ 型后缀,$n$ 为 $L$ 型后缀。 对于 $i\in[1, n)$,有: > **引理 2.1**(后缀类型递推性质) > > $i$ 是 $S$ 型后缀当且仅当以下任意一项成立 > > 1. $T_i<T_{i+1}

否则,也即以下任意一项成立

  1. T_i>T_{i+1}

此时 iL 型后缀。

引理 2.2(后缀类型指导排序)

若后缀 A,B 满足 A_1=B_1AS 型,BL 型,则 A>B

可由引理 1.2引理 2.1 推出。

LMS 子串

LMS 的意思是 Leftmost S-type。

称满足 i-1L 型后缀,iS 型后缀的 i 为一个 * 型后缀。称 T[i,j](i< j) 是一个 LMS 子串当且仅当 ij 是两个相邻的 * 型后缀。

注意 * 型后缀同时也是 S 型后缀。

引理 2.3

作为一个例外的边界情况。 **引理 2.4** 对于 LMS 子串 $K\ne \epsilon$,总有 $|K|>2$。 由 $*$ 型的定义可知。 **引理 2.5**(原串折半) 任意串 $T$ 中 LMS 子串的数量不超过 $\lceil \frac{n}{2}\rceil$。 根据**引理 2.4** 可知。 **引理 2.6** 一个字符串 $T$ 的所有 LMS 子串长度之和不超过 $2n$。 所有 LMS 子串间除了 $*$ 型后缀开头的第一个字符外均不交。 **引理 2.7** 若令每个字符的字典序包含它的后缀类型作为第二关键字,其中 $\text{S-type} >\text{L-type}$。那么不存在一个 LMS 子串是另一个 LMS 子串的真前缀。 设 $A$ 是 $B$ 的真前缀,因为 $A$ 以 $*$ 型结尾,根据 $*$ 型的定义则 $B$ 中出现的 $A$ 最后必然包含一个 $*$ 型在中间,这与 LMS 子串的定义矛盾。 ``` acbbccbbccbab# SLSSLLSSLLLSLS __*___*____*_* ``` (`#` 在此指代 $\epsilon$) 虽然 $\texttt{bbccb}$ 是 $\texttt{bbccba}$ 的真前缀,但是前者的最后一个 $\texttt{b}$ 是 $S$ 型,后者的最后一个 $\texttt{b}$ 是 $L$ 型。

称将 T 中所有 LMS 子串单独提出,离散化后作为字符组成的新字符串为 T'

T =
mmiissiissiippii#
LLSSLLSSLLSSLLLLS
__*___*___*_____*

2 = 
iissi
SSLLS

1 =
iippii#
SSLLLLS

0 =
#
S

T' =
2210

引理 2.8(问题缩减)

比较 $T'$ 中字符相当于比较 $T$ 中 LMS 子串。根据**引理 2.7**,不存在一个 LMS 子串是另一个 LMS 子串的真前缀所以 $T'$ 中两个任意后缀比较相当于 $T$ 中任意两个 $*$ 型后缀比较。 比如 $T$ 中 $*$ 型后缀 $i$ 和 $j$ 比较,按 LMS 子串划分成若干段。决定大小关系的即为相同段后一个个不同段。选用 LMS 子串的核心是不存在任意一个 LMS 子串为另一个 LMS 子串的真前缀:这样能够保证要么就在一段内的比较要么决出结果,要么扫过的前缀全等,不会出现虽然未得出结果但是扫过的前缀不全等的问题。 ![](https://cdn.luogu.com.cn/upload/image_hosting/lekn8td9.png)

需要注意的是这里只是 型后缀的字典序关系,与其它后缀无关。

诱导排序与 SA - IS 算法

sa' 型后缀的后缀数组,rk' 为其对应的 sa'_{rk'_i}=i=rk'_{sa'_i}

自己手玩几把可以发现后缀排序时总是首字母相同的后缀相邻(显然)。

#
aaaab#
aaab#
aab#
aabaaaab#
ab#
abaaaab#
b#
baaaab#

引理 2.2 可知首字母相同的情况下总是类型相同的后缀相邻。将每种字符和两种类型看作关键字就可以计数排序了。

但你会发现首字母与类型均相同的后缀之间大小不好判,现在就要利用 sa' 的信息了。

sa' 诱导 sa 排序,就是诱导排序了。诱导排序的流程如下:

开始前记得初始化 sa 数组以及确定 S 型桶 L 型桶的起始位置,其中 S 型桶逆序排放。(计数排序基操)

  1. 按顺序扫一遍 sa'* 型后缀,依次放入 sa 中。
  2. 顺序遍历 sa,然后对于所有非空的 sa_isa_i-1L 型的依次放入 sa 中。
  3. 重新确定 S 型起始位置,准备重排 * 型后缀。逆序遍历 sa 同时将 sa_i-1S 型的依次放入 sa 中。

整体流程有点像基数排序的说。

对于 sa'T' 中每个字符都不同的情况见下文,否则递归求总是可以转化成前者。

对 LMS 子串的排序,诱导排序流程中 1. 改为将 * 型后缀的首个字符按任意顺序放入桶中即可。

定义后缀 i 的 LMS 前缀为从后缀 i 的首个字符开始一直到一个 * 型后缀的开头。

接下来证明改动后每一步的正确性。

  1. 桶的排列按照字典序保证了放入后的有序和首个字符、后缀类型不同时有序。
  2. 假设放入第 k+1L 型 LMS 前缀 i 时前面放入的 kL 型 LMS 前缀保持有序。此时只需要考虑首字母相同的其他 L 型 LMS 前缀即可。假设存在 j 的 LMS 前缀比 i 的 LMS 前缀大,因为首字母相同,所以 j+1 的 LMS 前缀比 i+1 的 LMS 前缀大。与假设矛盾。根据数学归纳法可知在 k=0 时情况成立,推广到任意 k
  3. 类似,但这时 i 应理解为从其右边第一个 * 型后缀开头开始,一步一步向左加字符得到的 LMS 子串。换句话说此步后所有的 LMS 子串就已经完成排序了。

时间复杂度是 T(n)=T(\lceil\frac{n}{2}\rceil)+\Theta(n),根据主定理容易得出 T(n)=\Theta(n)

:::success[SA-IS]

#include <bits/stdc++.h>
using namespace std;

typedef long long ll;
template<typename sigmatype = char>
vector<int> SAIS(const basic_string<sigmatype> &T, int Sigma) 
{
    int n = T.size(); // 包含末尾空字符的长度
    if (n <= 1) return {0};

    vector<bool> isStype(n + 1);
    vector<int> pos(n + 1), name(n + 1, -1), sa(n, -1);
    vector<int> bucket(Sigma + 1, 0), bucketL(Sigma + 1), bucketS(Sigma + 1);

    auto isLMSchar = [&](int x) -> bool 
    {
        return x > 0 && isStype[x] && !isStype[x - 1];
    }; 

    auto equal_substring = [&](int x, int y) -> bool 
    {
        if (x == y) return true;
        while (true) {
            if (T[x] != T[y]) return false;
            ++x, ++y;
            bool x_is_lms = isLMSchar(x);
            bool y_is_lms = isLMSchar(y);
            if (x_is_lms || y_is_lms) 
            {
                return x_is_lms && y_is_lms;
            }
        }
        return false;
    }; 

    auto induced_sort = [&]() -> void 
    {
        for (int i = 0; i < n; i++) 
        {
            if (sa[i] > 0 && !isStype[sa[i] - 1]) 
            {
                sa[bucketL[T[sa[i] - 1]]++] = sa[i] - 1;
            }
        }

        bucketS[0] = 0;
        for (int i = 1; i <= Sigma; i++) bucketS[i] = bucket[i] - 1;
        for (int i = n - 1; i >= 0; i--) 
        {
            if (sa[i] > 0 && isStype[sa[i] - 1]) 
            {
                sa[bucketS[T[sa[i] - 1]]--] = sa[i] - 1;
            }
        }
    }; 

    // 1. 确定后缀类型
    isStype[n] = 1; // 空字符是 S 型
    for (int i = n - 1; i >= 0; i--) 
    {
        if (T[i] < T[i + 1]) isStype[i] = 1; 
        else if (T[i] > T[i + 1]) isStype[i] = 0; 
        else isStype[i] = isStype[i + 1];
    }

    // 2. 寻找 LMS 字符并放入桶中初步排序
    int cnt = 0; 
    for (int i = 1; i <= n; i++) 
    {
        if (isStype[i] && !isStype[i - 1]) pos[cnt++] = i;
    }

    // 计算桶边界
    for (int i = 0; i < n; i++) bucket[T[i]]++;
    bucketL[0] = bucketS[0] = 0;
    for (int i = 1; i <= Sigma; i++) 
    {
        bucket[i] += bucket[i - 1];
        bucketL[i] = bucket[i - 1]; 
        bucketS[i] = bucket[i] - 1; 
    }

    for (int i = 0; i < cnt; i++) {
        sa[bucketS[T[pos[i]]]--] = pos[i];
    }

    // 3. 第一次诱导排序
    induced_sort();

    // 4. 离散化命名
    int lstLMS = -1, namecnt = 0;
    bool flag = false; 
    for (int i = 0; i < n; i++) {
        if (isLMSchar(sa[i])) {
            if (lstLMS != -1) {
                if (!equal_substring(sa[i], lstLMS)) ++namecnt;
                else flag = true; 
            } else {
                namecnt = 0;
            }
            name[lstLMS = sa[i]] = namecnt;
        }
    }
    name[n] = 0; 

    // 5. 生成 T1 并递归求解
    basic_string<int> T1(cnt, 0);
    int post = 0;
    for (int i = 0; i <= n; i++) {
        if (name[i] != -1) {
            T1[post++] = name[i];
        }
    }

    vector<int> sa1(cnt);
    if (!flag) {
        for (int i = 0; i < cnt; i++) sa1[T1[i]] = i;
    } else {
        sa1 = SAIS<int>(T1, namecnt + 1); 
    }

    // 6. 最终诱导排序
    fill(bucket.begin(), bucket.end(), 0);
    for (int i = 0; i < n; i++) bucket[T[i]]++;
    bucketL[0] = bucketS[0] = 0;
    for (int i = 1; i <= Sigma; i++) 
    {
        bucket[i] += bucket[i - 1];
        bucketL[i] = bucket[i - 1];
        bucketS[i] = bucket[i] - 1;
    }

    fill(sa.begin(), sa.end(), -1);
    for (int i = cnt - 1; i >= 0; i--) sa[bucketS[T[pos[sa1[i]]]]--] = pos[sa1[i]];
    induced_sort();

    return sa;
}

int main() 
{
    ios::sync_with_stdio(false);
    cin.tie(nullptr);

    string s;
    cin >> s;
    // 关键修复:手动追加一个值为 0 的空字符,确保最后一个后缀是 S 型
    s.push_back(0); 

    // 'z' 的 ASCII 为 122,加上空字符 0,值域设为 123 足够安全
    auto sa = SAIS(s, 123); 

    // 输出结果,跳过索引 0 (即空字符的位置),并将 0-indexed 转为 1-indexed
    for (size_t i = 1; i < sa.size(); i++) 
    {
        cout << sa[i] + 1;
        if (i != sa.size() - 1) cout << " ";
    }
    cout << "\n";
    return 0;
}

/*

AI 使用声明:此段代码使用 Qwen3.7 - 千问辅助调试

但文章原有文字内容和代码主体均为纯人工写就

神必 Lambda 码风也是人工的,AI 只是辅助调试用

*/

:::

例题什么的后面会再开一篇专题,毕竟到这里已经超过一万字符(算上代码)了。

参考资料

创作声明

本文遵循 CC BY-NC-SA 4.0 协议。

保证本文文字说明部分未使用 GenAI 工具辅助创作。

转载例如下:

[原文](https://www.luogu.com.cn/article/73yg70al)作者为 [clx201022](https://www.luogu.com.cn/user/552688),转载人保证会遵循 [CC BY-NC-SA 4.0](https://creativecommons.org/licenses/by-nc-sa/4.0/deed.zh-hans) 协议:以适当的方式署名,不以盈利为目的进行转载,并以相同方式共享。