【模板】后缀排序 | DC3算法构造后缀数组

· · 题解

2021/02/22 更新:修改了排版使其符合题解规范。

CSP 之前写篇题解加加 RP

DC3 算法构造后缀数组

题目:P3809 【模板】后缀排序

看这篇题解之前建议学习一下倍增法求后缀数组,不然可能有点难以理解。

切入正题,DC3 是一种 O(n) 时空复杂度构造后缀数组的算法,然而因为它常数大(相较于 SA-IS),而且比倍增法难打,所以几乎没人打。

算法主要流程如下:

  1. 把所有的后缀分成两部分,记后缀起始的下标为 i,第一部分为 i\bmod 3=0 的后缀(以下称为 A 类后缀),第二部分为 i\bmod 3\neq0 的后缀(以下称为 B 类后缀)。

  2. 利用基数排序,对所有 B 类后缀按照前三个字符排序。

  3. 如果前三个字符无法将所有 B 类后缀排好序,把每三个相邻的字符看作一个字符,构造新的字符串,递归对 B 类后缀排序。

  4. 利用 B 类后缀的顺序可以利用基数排序很快的求出 A 类后缀的顺序。

  5. 对 A 类后缀和 B 类后缀进行归并排序,然后 SA 数组就求出来了。

下面用一组样例演示一下 DC3 算法的执行过程。

这里我们要排序的字符串是 aaababaaca,下面这张图标注出了所有 B 类后缀的前三个字符。

对这些前三个字符进行排序,方法是三次基数排序,从后往前每次排一个字符,用上一次排序的结果做第二关键字。图中的数字表示排序后的排名。

发现前三个字符不能确定 B 类后缀的顺序(有两个 aba,排名都是 2),于是构造新字符串,递归执行 DC3 算法,如图,构造出的新串是 1230245,下一行表示新串每个后缀的排名。

从图中可以看出,新串的构造方式为:所有 i\bmod 3=1 的后缀按顺序排列在左侧,i\bmod 3=2 的后缀按顺序排列在右侧,中间用一个 0 字符(就是一个空的 A 类后缀)分隔。这样构造可以保证新串中相邻两个字符在原串中表示相邻六个字符。同时中间用 0 字符分隔,可以保证两类后缀的排序相互独立,不会互相干扰。(代码中可以看到,其实只有 n\bmod 3=1 的时候需要加这个 0 字符,因为如果不是这样,代表最后一个 A 类后缀的三个字符末尾本来就是 0 字符)。于是,我们可以从新串后缀的排序中得到所有 B 类后缀的顺序。

无论是仅排序三个字符还是递归排序,我们现在已经得到所有 B 类后缀的排名了(下图中第三行)。注意到对于所有 A 类后缀,都相当于一个字符加上一个 B 类后缀,而 B 类后缀的顺序已经求出。所以用第一个字符做第一关键字,其后的 B 类后缀排名做第二关键字,一次基数排序即可求出 A 类后缀的排名(图中第四行)。图中第五、六行写出了按顺序排列的 A 类后缀和 B 类后缀起始下标。

现在两类后缀都已经有序,只需要一次归并就能将两类后缀的顺序合并。归并就需要比较两类后缀的大小。记需要比较大小的 A 类后缀起始下标为 i,A 类后缀起始下标为 j

如果位置 i,j 上的字符不同,直接比较这两个字符即可。否则需要比较后缀 i+1j+1,有两种情况:

然后就可以求出最终答案了。

由于本题下标是从 1 开始计算,所以程序输出 10 1 2 7 5 3 8 6 4 9

代码实现细节非常多,全都写在注释里了。

无注释版代码

#include <cstdio>
#define N 1000001
int rk[3*N+100], sa[3*N], bu[N], x[N], y[N];
void sort(int *rk, int *a, int *b, int n, int m) {
    //rk是要排序的串,a是要排序的字符的下标,其顺序是第二关键字
    for (int i = 0; i <= m; i++) bu[i] = 0;
    for (int i = 0; i < n; i++) bu[rk[a[i]]]++;
    for (int i = 1; i <= m; i++) bu[i] += bu[i-1];
    for (int i = n-1; i >= 0; i--) b[--bu[rk[a[i]]]] = a[i];
}
bool cmp3(int *r, int x, int y) {
    return r[x] == r[y] && r[x+1] == r[y+1] && r[x+2] == r[y+2];
}
bool cmpt(int* r, int x, int y) { //比较两类后缀的大小
    if (r[x] != r[y]) return r[x] < r[y];
    if (x%3 == 1) return bu[x+1] < bu[y+1];
    else return !cmpt(r, y+1, x+1);
}
void DC3(int *rk, int *sa, int n, int m) {
    bool h = (n%3 == 1); if (h) rk[n++] = 0;
    //如果n%3==1,在原串末尾增加一个空的A类后缀,以防递归排序越界
    int *rn = rk+n+2, *san = sa+n, ln = 0, p;
    //ln: B类后缀的数量
    for (int i = 0; i < n; i++)
        if (i % 3) x[ln++] = i;
    rk[n] = rk[n+1] = 0; //在原串后增加两个空字符,方便处理
    sort(rk+2, x, y, ln, m);
    sort(rk+1, y, x, ln, m);
    sort(rk, x, y, ln, m);
    int ta = 0, td = (n+1)/3;
    //ta: A类后缀的数量
    //td: i%3==1的后缀的数量
    #define F(x) (x/3)+(x%3==1?0:td)
    //利用三个字符构造新串
    //F(x): 原串的后缀x在新串中的位置
    rn[F(y[0])] = p = 1;
    for (int i = 1; i < ln; i++) {
        if (!cmp3(rk, y[i], y[i-1])) p++;
        rn[F(y[i])] = p;
    }
    if (p < ln) DC3(rn, san, ln, p); //递归对B类后缀进行排序
    else for (int i = 0; i < ln; i++)
        if (rn[i]) san[rn[i]-1] = i;
    for (int i = 0; i < ln; i++)
        if (san[i] < td) y[ta++] = san[i]*3;
    //对于所有i%3==1的后缀,其排名是之前一个A类后缀排序用的第二关键字
    sort(rk, y, x, ta, m);
    #define G(x) (x>=td?(x-td)*3+2:x*3+1)
    //G(x): 新串的后缀x在原串中的位置
    for (int i = 0; i < ln; i++)
        bu[y[i] = G(san[i])] = i; //这里的bu[x]就是B类后缀x的排名
    bu[n] = -1;
    int i = 0, j = h; p = 0; //如果h等于1,则第一个A类后缀是为了防止越界加入的,归并时去除
    while (i < ta && j < ln) {
        if (cmpt(rk, y[j], x[i])) sa[p++] = y[j++];
        else sa[p++] = x[i++];
    }
    while (i < ta) sa[p++] = x[i++];
    while (j < ln) sa[p++] = y[j++];
}
char s[N]; int n;
int main() {
    scanf("%s", s);
    for (n = 0; s[n]; n++)
        rk[n] = s[n]-'0'+1;
    //初始的rank是字符的ascii码
    //注意0是保留的,普通字符从1开始
    DC3(rk, sa, n, 75);
    for (int i = 0; i < n; i++)
        printf("%d ", sa[i]+1);
    puts("");
}

其他的问题

因为存在递归,时空复杂度 f(n)=\Theta(n)+f(\frac{2}{3}n)=\Theta(n)。(不知道怎么推的可以看一下这个)

又因为 n+(\frac{2}{3})n+(\frac{2}{3})^2n+(\frac{2}{3})^3n+...=3n,所以数组开三倍大小。