后缀三姐妹

08/10 09:07
阅读数 35

绝对不咕


写在前面

会考虑整个与标题相关的二次创作。
什么时候有能力再说


前置小碎骨

计数排序

可以参考:OI-wiki 计数排序

计数排序是一种与桶排序类似的排序方法。
将长度为 \(n\) 的数列 \(a\) 排序后放入 \(b\) 的代码如下,
其中 \(w\) 为值域,即 \(\max\{a_i\}\)

int a[kMaxn], b[kMaxn];
for (int i = 1; i <= n; ++ i) ++ cnt[a[i]];
for (int i = 1; i <= w; ++ i) cnt[i] += cnt[i - 1];
for (int i = n; i >= 1; -- i) b[cnt[a[i]] --] = a[i];

其中,在对 \(cnt\) 求前缀和后,
\(cnt_i\) 为小于 \(i\) 的数的数量,即为 \(i\) 的排名。
因此在下一步中,可以根据排名赋值。

复杂度为 \(O(n+w)\),在值域与 \(n\) 同阶时复杂度比较优秀。

基数排序

这玩意比较水,参考 OI-wiki 基数排序

个人认为基数排序只是一种思想,并不算一种排序方法。
它仅仅是将 k 个排序关键字分开依次考虑,实际每次排序还是靠计数排序实现。

请务必充分理解!


一些约定:

\(S[i:j]\):由字符串 \(S\)\(S_i\sim S_j\) 构成的子串。
\(S_1<S_2\):字符串 \(S_1\) 的字典序 \(<S_2\)
后缀:从某个位置 \(i\) 开始,到串末尾结束的子串,后缀 \(i\) 等价于子串 \(S[i:n]\)


后缀数组

网上部分题解都直接开讲优化后面目全非的代码。
*这太野蛮了*
这里参考了 OI-wiki 上的讲解。

定义

字符串 \(S\) 的后缀数组 \(A\),被定义为一个数组,内容是其所有后缀 按字典序排序后的起始下标。
有: \(S[A_{i-1}:n]<S[A_i:n]\) 成立。

举例:这里有一个可爱的字符串:\(S=\text{yuyuko}\)
\(\text{k<o<u<y}\),它的后缀数组 \(A = [5,6,4,2,3,1]\)
具体地,有:

排名 1 2 3 4 5 6
下标 \(5\) \(6\) \(4\) \(2\) \(3\) \(1\)
后缀 \(\text{ko}\) \(\text{o}\) \(\text{uko}\) \(\text{uyuko}\) \(\text{yuko}\) \(\text{yuyuko}\)

显然不同后缀的排名必然不同(因为长度不等)


倍增法构造

先将所有 \(S[i:j]\) 进行排序。
每次通过 \(S[i:i+2^{k-1}-1]\) 的大小关系,求出 \(S[i:i+2^k-1]\) 的大小关系。

对于 \(S[i:i+2^k-1]\)\(S[j:j+2^k-1]\),分别将它们裂开,分成两成长度为 \(i+2^{k-1}\) 的串。

\(A_i = S[i:i+2^{k-1}-1]\)\(B_i = S[i+2^{k-1}:i+2^k-1]\)
考虑字典序排序的过程,则 \(S[i:i+2^k-1] <S[j:j+2^k-1]\) 的条件为:

\[[A_i<A_j] \operatorname{or}\ [A_i=A_j \operatorname{and} B_i<B_j] \]

考虑每一次倍增时,都使用 sort 按双关键字 \(A_i\)\(B_i\) 进行排序。
时间复杂度显然为 \(O(n\log^2 n)\)


优化

sort 太慢啦!
发现后缀数组值域即为 \(n\),又是多关键字排序。
考虑基数排序。

上面已经给出一个用于比较的式子:

\[[A_i<A_j] \operatorname{or}\ [A_i=A_j \operatorname{and} B_i<B_j] \]

\(A_i,B_i\) 大小关系已知,直接基数排序实现即可。
先将 \(B_i\) 作为第二关键字排序。
再将 \(A_i\) 作为第一关键字排序。

单次计数排序复杂度 \(O(n + w)\)\(w\) 为值域,最大与 \(n\) 同阶)。
排序变为 \(O(n)\) 级别,时间复杂度 \(O(n\log n)\)


代码

P3809 【模板】后缀排序

这是一份没有优化过的代码,是对上述过程的直接实现。
只能获得 73 分。
可发现代码的空间复杂度为 \(O(n)\)
代码实现较为复杂,下面会进行详细讲解。


//知识点:SA
/*
By:Luckyblock
I love Marisa;
But Marisa has died;
*/
#include <cstdio>
#include <ctype.h>
#include <cstring>
#include <algorithm>
#define ll long long
const int kMaxn = 1e6 + 10;
//=============================================================
char S[kMaxn];
//sa[i]: 倍增过程中子串[i:i+2^k-1]的排名,
//rk[i] 排名为i的子串 [i:i+2^k-1],
//它们互为反函数。
//rk 和 oldrk 要开2倍空间,下面会提到原因。
int n, m, sa[kMaxn], rk[kMaxn << 1], oldrk[kMaxn << 1];
int id[kMaxn], cnt[kMaxn]; //用于计数排序的两个tmp数组
//=============================================================
inline int read() {
  int f = 1, w = 0; char ch = getchar();
  for (; !isdigit(ch); ch = getchar()) if (ch == '-') f = -1;
  for (; isdigit(ch); ch = getchar()) w = (w << 3) + (w << 1) + (ch ^ '0');
  return f * w;
}
//=============================================================
int main() {
  scanf("%s", S + 1);
  n = strlen(S + 1);
  m = std :: max(n, 300); //值域大小
  
  //初始化 rk 和 sa
  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 = n; i >= 1; -- i) sa[cnt[rk[i]] --] = i;

  //倍增过程。 
  //w = 2^{k-1},是已经推出的子串长度。
  //注意此处的 sa 数组存的并不是后缀的排名,
  //存的是指定长度子串的排名。
  for (int w = 1; w < n; w <<= 1) {
    //按照后半截 rk[i+w] 作为第二关键字排序。
    memset(cnt, 0, sizeof (cnt));
    for (int i = 1; i <= n; ++ i) id[i] = sa[i];
    for (int i = 1; i <= n; ++ i) ++ cnt[rk[id[i] + w]]; //这里有越界风险,因此开了2倍空间,否则会被卡
    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];

    //按照前半截 rk[i] 作为第一关键字排序。
    memset(cnt, 0, sizeof (cnt));
    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];

    //更新 rk 数组。
    //这里可以滚动数组一下,但是可读性会比较差(
    //卡常可以写一下。
    for (int i = 1; i <= n; ++ i) oldrk[i] = rk[i];
    for (int p = 0, i = 1; i <= n; ++ i) {
      if (oldrk[sa[i]] == oldrk[sa[i - 1]] &&  //判断两个子串是否相等。
          oldrk[sa[i] + w] == oldrk[sa[i - 1] + w]) { //这里有越界风险,因此开了2倍空间,否则会被卡
        rk[sa[i]] = p;
      } else {
        rk[sa[i]] = ++ p;
      }
    }
  }
  for (int i = 1; i <= n; ++ i) printf("%d ", sa[i]);
  return 0;
}

这里定义了两个数组:
\(sa_i\):倍增中 排名为 \(i\) 的长度为 \(2^{k-1}\) 的子串。
\(rk_i\):倍增过程中子串 \(S[i:i+2^k-1]\) 的排名,
显然它们互为反函数,\(sa_{rk_i}=i, rk_{sa_i} = i\)



首先初始化 \(rk\)\(sa\)

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 = n; i >= 1; -- i) sa[cnt[rk[i]] --] = i;

初始化 \(rk_i = S_i\),即 \(S_i\)\(\text{ASCII}\) 值。
虽然这样不满足值域在 \([1,n]\) 内,但体现了大小关系,可用于更新。 \(rk\) 的值之后还会更新。

子串长度为 \(1\),直接根据 \(rk_i\) 计数排序 \(sa\) 即可。


之后进入倍增。

每次倍增先后按照 后半截,前半截的 \(rk\) 作为关键字排序来更新 \(sa\)
\(id\) 是一个 tmp 数组,存排序前的的 \(sa\)

memset(cnt, 0, sizeof (cnt));
for (int i = 1; i <= n; ++ i) id[i] = sa[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];
memset(cnt, 0, sizeof (cnt));
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];

排后半截时 会枚举到 \(id[i]+w > n\) 怎么办?
考虑实际意义,出现此情况,表示该子串后半截为空。
空串字典序最小,考虑直接把 \(rk\) 开成两倍空间,则 \(rk[i]=0(i>n)\) 恒成立。
防止了越界,也处理了空串的字典序。



更新 \(rk\) 数组。

for (int i = 1; i <= n; ++ i) oldrk[i] = rk[i];
for (int p = 0, i = 1; 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;
  }
}

\(sa\)\(rk\) 的反函数。
这里相当于根据有序的 \(sa\),离散化并去重 \(rk\)

考虑两个子串 \(rk\) 相等的条件。
显然,当其前后两半均相等时,两子串相同,其 \(rk\) 才相同,则有上述的判断。

这里也会出现空串的情况,注意 2 倍空间。


再优化

被卡常了,排两次计数排序太慢啦!
观察对后半截排序时的特殊性质:

考虑更新前的 \(sa_i\) 的含义:排名为 \(i\) 的长度为 \(2^{k-1}\) 的子串。

在本次排序中,\(sa_i\) 是长度为 \(2^k\) 的子串 \(sa_{i}-2^{k-1}\) 的后半截。 \(sa_i\) 的排名将作为排序的关键字。

\(sa_i\) 的排名为 \(i\),则排序后 \(sa_{i}-2^{k-1}\) 的排名必为 \(i\)
考虑直接赋值,那么第一次计数排序就可以写成这样:

int p = 0;
for (int i = n; i > n - w; -- i) id[++ p] = i; //后半截为空的串
for (int i = 1; i <= n; ++ i) { //根据后半截,直接推整个串的排名
  if (sa[i] > w) id[++ p] = sa[i] - w;
}

注意后半截为空串的情况,这样的串排名相同且最小。


以及一些奇怪的常数优化:

减小值域。
发现值域大小 \(m\) 与计数排序复杂度有关。
其最小值应为 \(rk\) 的最大值,在更新 \(rk\) 时将其更新即可。

减少数组嵌套的使用,从而减少不连续内存访问。
在第二次计数排序时,将 \(rk_{id_i}\) 存下来。

用cmp函数判断两个子串是否相同。
同样是减少不连续内存访问,详见代码。


最终代码

//知识点:SA
/*
By:Luckyblock
I love Marisa;
But Marisa has died;
*/
#include <cstdio>
#include <ctype.h>
#include <cstring>
#include <algorithm>
#define ll long long
const int kMaxn = 1e6 + 10;
//=============================================================
char S[kMaxn];
int n, m, sa[kMaxn], rk[kMaxn << 1], oldrk[kMaxn << 1];
int id[kMaxn], cnt[kMaxn], rkid[kMaxn];
//=============================================================
inline int read() {
  int f = 1, w = 0; char ch = getchar();
  for (; !isdigit(ch); ch = getchar()) if (ch == '-') f = -1;
  for (; isdigit(ch); ch = getchar()) w = (w << 3) + (w << 1) + (ch ^ '0');
  return f * w;
}
bool cmp(int x, int y, int w) { //判断两个子串是否相等。
  return oldrk[x] == oldrk[y] && 
         oldrk[x + w] == oldrk[y + w]; 
}
//=============================================================
int main() {
  scanf("%s", S + 1);
  n = strlen(S + 1);
  m = std :: max(n, 300); //值域大小
  
  //初始化 sa数组
  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 = n; i >= 1; -- i) sa[cnt[rk[i]] --] = i;

  //倍增过程。 
  //此处 w = 2^{k-1},是已经推出的子串长度。
  //注意此处的 sa 数组存的并不是后缀的排名,
  //存的是指定长度子串的排名。
  for (int p, w = 1; w < n; w <<= 1) {
    //按照后半截 rk[i+w] 作为第二关键字排序。
    p = 0;
    for (int i = n; i > n - w; -- i) id[++ p] = i; //后半截为空的串
    for (int i = 1; i <= n; ++ i) { //根据后半截,直接推整个串的排名
      if (sa[i] > w) id[++ p] = sa[i] - w;
    }

    //按照前半截 rk[i] 作为第一关键字排序。
    memset(cnt, 0, sizeof (cnt));
    for (int i = 1; i <= n; ++ i) ++ cnt[(rkid[i] = rk[id[i]])];
    for (int i = 1; i <= m; ++ i) cnt[i] += cnt[i - 1];
    for (int i = n; i >= 1; -- i) sa[cnt[rkid[i]] --] = id[i];

    //更新 rk 数组。
    //这里可以滚动数组一下,但是可读性会比较差(
    //卡常可以写一下。
    std ::swap(rk, oldrk);
    m = 0; //直接更新值域 m
    for (int i = 1; i <= n; ++ i) {
      rk[sa[i]] = (m += (cmp(sa[i], sa[i - 1], w) ^ 1));
    }
  }
  for (int i = 1; i <= n; ++ i) printf("%d ", sa[i]);
  return 0;
}

后缀树


后缀自动机

涛队钦定我写一篇博客。
那我就咕了


写在最后

参考资料:

OI-wiki SA
后缀数组详解 - 自为风月马前卒
「后缀排序SA」学习笔记 - Rainy7

写句子。
I want to be a rap...
犹豫了一下。
raper。
身 败 名 裂



展开阅读全文
打赏
0
0 收藏
分享
加载中
更多评论
打赏
0 评论
0 收藏
0
分享
在线直播报名
返回顶部
顶部