详解后缀数组(SA)数组

算法能解决的问题

• 求解最长公共前缀$lcp$
• 在短时间内将$n$个后缀按照字典序排序

算法原理

倍增法$O(n\log n)$

对于一个字符串有很多后缀, 第$i$个后缀

• $sa$数组, $sa(i)$表示排名第$i$位的后缀
• $rk$数组, $rk(i)$表示第$i$个后缀的排名
• *$height(i)$表示$sa(i)$与$sa(i-1)$两个后缀的最长公共前缀

假设当前字符串$aababb$, 后缀如下

1. $aababb$
2. $ababb$
3. $babb$
4. $abb$
5. $bb$
6. $b$

---

第一次只关注每个后缀的第一个字符按照字典序排序, 如果没字符用$null$表示, 比任意字符都小

第一次排序的结果, 第一个字符相同相对顺序不变

1. $aababb$
2. $ababb$
3. $abb$
4. $babb$
5. $bb$
6. $b$

可以用基数排序算法时间复杂度$O(n)$

---

假设执行某一轮后已经按照前$k$个字符排好序了, 考虑前$2k$个字符

• 将每个后缀的前$k$个字符当作第一关键字, 将后$k$个字符当作第二关键字进行排序, 双关键字排序

• 每一轮结束之后, 将第一关键字离散化为一个整数方便排序, 也是按照基数排序来排

---

实现细节

基数排序范围$1\sim n$

统计每个数出现次数

因为求的是排名, 求的是$\le w$的数字个数, 就是对出现次数求前缀和

然后倒序枚举每个数字"(为了保证相同值的先后顺序不发生变化)

---

小于等于$3$的数字个数是$5$, 然后对应位置$-1$

小于等于$1$的数字个数$2$

小于等于$1$的数字个数$1$

小于等于$2$的数字个数$4$

小于等于$2$的数字个数$3$

小于等于$5$的数字个数$6$

---

这样就得到了小于等于某个数的排名

双关键字类似, 为了保证排序是稳定的

1. 现将所有元素按照第二关键字排序
2. 从后往前再按照第一关键字排, 保证第二关键字相对顺序不变

注意到

第$i$个后缀的第二关键字就是第$i+k$后缀的第一关键字

基数排序部分

将$y$清空, 将$x$信息存在$y$中

然后对前$2k$的字符进行离散化, 离散化范围是$1\sim n$

$sa(i)$个后缀的第二关键字是$sa(i+k)$后缀的第一关键字

逻辑很简单, 如果当前离散化的值和前一个字符串离散化后的值相等, 直接复制过来

如果比前面的值大则$+1$

这样就可以在$O(n\log n)$时间复杂度内求出$sa$数组

---

现在问题变成如何求$height$数组?

假设已经将后缀排好序, $height(i)$表示
排名第$i$位的后缀和排名第$i-1$位的后缀的最长公共前缀的长度, 也就是

$$height(i)=lcp(i,i-1)$$

$i,j$是排名第$i,j$的后缀
$lcp$的一些性质

• $lcp(i,j)=lcp(j,i)$
• $lcp(i,i)=len(i)$
• $l\le k\le j$有
  $$lcp(i,j)=\min (lcp(i,k),lcp(k,j))$$

证明:

(1) 首先证明$lcp(i,j)\ge \min (lcp(i,k),lcp(k,j))$

黑色部分是$k,j$前缀的最长公共前缀, 红色部分是$i,k$部分的最长公共前缀

因为$lcp$是最长的, 因此$lcp(i,j)\ge \min (lcp(i,k),lcp(k,j))$, (1)证毕

(2) 再证明$lcp(i,j)\le \min (lcp(i,k),lcp(k,j))$

假设排名$i,j$的后缀的最长公共前缀是$A,C$, 那么有$A=C$, 再取出排名$k$的后缀的前缀是$B$

因为$i,k,j$是按照字典序排名的, 因此$A\le B\le C$

因为$A=C$, 夹逼定理得到$A=B=C$

得到$lcp(i,k)\ge A$, 并且$lcp(k,j)\ge C$, 也就是$\min (lcp(i,k),lcp(k,j))\ge lcp(i,j)$, (2)证毕

由(1)和(2)得到$lcp(i,j)=\min (lcp(i,k),lcp(k,j))$, 证毕

暴力的计算
$$lcp(i,j)=\min (lcp(i,i+1),lcp(i,i+2),...,lcp(i,j-1))$$
算法时间复杂度$O(n^2)$非常慢

优化

$h[i]=height[rk[i]]$, 第$i$个后缀和排名在第$i$个后缀前一名的后缀的最长公共前缀的长度

猜想(1): $h[i]\ge h[i-1]-1$

如果这个猜想是正确的, 那么枚举的时候直接从$i-1$开始枚举

对于每个$i$来说$k$最多减$1$, $k\le n$, $k$移动的范围是$[-n,n]$, 移动次数最多是$2n$, 因此算法时间复杂度降低到$O(n)$

设字典序第$i$个后缀在第$i-1$个后缀的下方

分为两种情况

• 如果$k$和$i-1$没有公共前缀, 那么$h[i-1]=0$, $h[i]\ge -1$显然成立

• 假设第一个字母相同

第$k+1$个后缀在第$i$个后缀的上面

第$k+1$个后缀和第$i$个后缀的最长公共前缀等于, 第$i-1$个后缀和第$k$个后缀的最长公共前缀去掉第一个字符

因此第$k+1$个后缀和第$i$个后缀的最长公共前缀的长度就是$h[i-1]-1$

但是中间的后缀数量未知, 我们设中间某个后缀是$t$, 那么由$lcp(i,j)=\min (lcp(i,k),lcp(k,j))$可知$lcp(k+1,i)=\min (lcp(t,i))$

因为$lcp(k+1,i)$是取的最小值因此
$lcp(t,i)\ge lcp(k+1,i)$, 也就是有
$$h[i]\ge h[i-1]-1$$
证毕

核心代码

void get_sa() {
   
    for (int i = 1; i <= n; ++i) c[x[i] = s[i]]++;
    for (int i = 2; i <= m; ++i) c[i] += c[i - 1];
    for (int i = n; i; --i) sa[c[x[i]]--] = i;

    for (int k = 1; k <= n; k <<= 1) {
   
        int num = 0;
        for (int i = n - k + 1; i <= n; ++i) y[++num] = i;

        // 希望求出排名为t的第二关键字, 实际求是求排名t + k的第一关键字, 枚举排名为t + k的后缀
        for (int i = 1; i <= n; ++i) {
   
            if (sa[i] > k) y[++num] = sa[i] - k;
        }
        // 清空基数排序数组, 再按照第一关键字排序
        for (int i = 1; i <= m; ++i) c[i] = 0;
        for (int i = 1; i <= n; ++i) c[x[i]]++;
        for (int i = 2; i <= m; ++i) c[i] += c[i - 1];
        // y存放t + k的第一关键字, 实际就是t的第二关键字, c是基数排序的桶, 获取到排名后数量减少
        for (int i = n; i; --i) sa[c[x[y[i]]]--] = y[i], y[i] = 0;
        // 将x信息存在y中
        swap(x, y);
        x[sa[1]] = 1, num = 1;
        // 现在已经将前2k个字符排好序了, 然后将每个后缀的前2k个字符构成的子串离散化在下一轮排序的时候可以作为关键字参与排名
        for (int i = 2; i <= n; ++i) {
   
            x[sa[i]] = y[sa[i]] == y[sa[i - 1]] && y[sa[i] + k] == y[sa[i - 1] + k] ? num : ++num;
        }
        if (num == n) break;
        m = num;
    }
}

求$height$数组的代码

void get_height() {
   
    for (int i = 1; i <= n; ++i) rk[sa[i]] = i;
    for (int i = 1, k = 0; i <= n; ++i) {
   
        if (rk[i] == 1) continue;

        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$, $k$每次倍增

• 对于当前后缀$t$需要先按照第二关键字排序, $t$的第二关键字就是$t+k$的第一关键字, 枚举$t+k$

• 再按照第二关键字逆序排第一关键字, 得到的相对顺序不会发生改变

• 这样就处理好了前$2k$个字符, 将前$2k$个字符的后缀离散化

• 重复该过程直到所有的后缀的排好序

示例代码

#include <bits/stdc++.h>

using namespace std;

const int N = 1e6 + 10;

int n, m;
string s;
int sa[N], x[N], y[N], c[N], rk[N], height[N];

void get_sa() {
   
    for (int i = 1; i <= n; ++i) c[x[i] = s[i]]++;
    for (int i = 2; i <= m; ++i) c[i] += c[i - 1];
    for (int i = n; i; --i) sa[c[x[i]]--] = i;

    for (int k = 1; k <= n; k <<= 1) {
   
        int num = 0;
        for (int i = n - k + 1; i <= n; ++i) y[++num] = i;

        // 希望求出排名为t的第二关键字, 实际求是求排名t + k的第一关键字, 枚举排名为t + k的后缀
        for (int i = 1; i <= n; ++i) {
   
            if (sa[i] > k) y[++num] = sa[i] - k;
        }
        // 清空基数排序数组, 再按照第一关键字排序
        for (int i = 1; i <= m; ++i) c[i] = 0;
        for (int i = 1; i <= n; ++i) c[x[i]]++;
        for (int i = 2; i <= m; ++i) c[i] += c[i - 1];
        // y存放t + k的第一关键字, 实际就是t的第二关键字, c是基数排序的桶, 获取到排名后数量减少
        for (int i = n; i; --i) sa[c[x[y[i]]]--] = y[i], y[i] = 0;
        // 将x信息存在y中
        swap(x, y);
        x[sa[1]] = 1, num = 1;
        // 现在已经将前2k个字符排好序了, 然后将每个后缀的前2k个字符构成的子串离散化在下一轮排序的时候可以作为关键字参与排名
        for (int i = 2; i <= n; ++i) {
   
            x[sa[i]] = y[sa[i]] == y[sa[i - 1]] && y[sa[i] + k] == y[sa[i - 1] + k] ? num : ++num;
        }
        if (num == n) break;
        m = num;
    }
}

void get_height() {
   
    for (int i = 1; i <= n; ++i) rk[sa[i]] = i;
    for (int i = 1, k = 0; i <= n; ++i) {
   
        if (rk[i] == 1) continue;
        
        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;
    }
}

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

    cin >> s;
    n = s.size();
    s = ' ' + s;
    m = 122;

    get_sa();
    get_height();

    for (int i = 1; i <= n; ++i) cout << sa[i] << ' ';
    cout << '\n';
    for (int i = 1; i <= n; ++i) cout << height[i] << ' ';
    cout << '\n';

    return 0;
}

*NOIP2015品酒大会

两个问, 第一个是有多少对$r$相似的酒对, 第二问是在所有$r$相似的酒对中的美味度乘积的最大值

$n\le 3\times 10^5$

对$height$分段

对于$height[i]<r$的情况作为分界线, 将真个数组进行分段

因为$lcp(i,j)=\min (lcp(i,k),lcp(k,j))$, $r$相似的酒只能发生在段内, 如果跨段就会导致最小值$<r$

假设分段之后每段的数量是$s$, 那么答案就是$\sum C_s^2$

乘积最大化三种情况

1. 两个正数相乘, 最大$\times$次大
2. 负数相乘, 最小$\times$次小
3. 一正一负只有两个数, 最大$\times$次大

维护每一个段内的最大值, 次小值, 最小值, 次小值

因为$r$如果从$0$开始枚举, 最开始所有$height$都可以考虑, 但是随着$r$增大需要分裂, 不好操作

$r$需要枚举所有情况, 按照从大到小枚举所有$r$, 每次变小需要合并一些段, 这样就好维护最大值, 次小值, 最小值, 次小值

动态维护集合的合并可以用并查集维护

算法时间复杂度$O(n\log n)$

#include <bits/stdc++.h>

#define x first
#define y second

using namespace std;

typedef long long LL;
typedef pair<LL, LL> PLL;

const int N = 3e5 + 10;
const LL INF = 2e18;

int n, m;
string s;
int sa[N], rk[N], x[N], y[N], c[N], height[N];
int w[N], p[N], sz[N];
LL max1[N], max2[N], min1[N], min2[N];
vector<int> hs[N];
PLL ans[N];


int find(int x) {
   
    if (p[x] != x) p[x] = find(p[x]);
    return p[x];
}

void get_sa() {
   
    for (int i = 1; i <= n; ++i) c[x[i] = s[i]]++;
    for (int i = 2; i <= m; ++i) c[i] += c[i - 1];
    for (int i = n; i; --i) sa[c[x[i]]--] = i;

    for (int k = 1; k <= n; k <<= 1) {
   
        int num = 0;
        for (int i = n - k + 1; i <= n; ++i) y[++num] = i;
        for (int i = 1; i <= n; ++i) {
   
            if (sa[i] > k) y[++num] = sa[i] - k;
        }
        for (int i = 1; i <= m; ++i) c[i] = 0;
        for (int i = 1; i <= n; ++i) c[x[i]]++;
        for (int i = 2; i <= m; ++i) c[i] += c[i - 1];
        for (int i = n; i; --i) sa[c[x[y[i]]]--] = y[i], y[i] = 0;
        swap(x, y);
        x[sa[1]] = 1, num = 1;
        for (int i = 2; i <= n; ++i) {
   
            x[sa[i]] = y[sa[i]] == y[sa[i - 1]] && y[sa[i] + k] == y[sa[i - 1] + k] ? num : ++num;
        }
        if (num == n) break;
        m = num;
    }
}

void get_height() {
   
    for (int i = 1; i <= n; ++i) rk[sa[i]] = i;
    for (int i = 1, k = 0; i <= n; ++i) {
   
        if (rk[i] == 1) continue;
        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;
    }
}

LL C(LL a) {
   
    return (a - 1ll) * a / 2ll;
}

PLL solve(int r) {
   
    static LL cnt = 0, maxv = -INF;

    for (auto x : hs[r]) {
   
        int a = find(x - 1), b = find(x);
        cnt -= C(sz[a]) + C(sz[b]);
        p[a] = b, sz[b] += sz[a];
        cnt += C(sz[b]);

        vector<LL> tmp = {
   max1[a], max2[a], max1[b], max2[b]};
        sort(tmp.begin(), tmp.end(), greater<LL>());
        max1[b] = tmp[0], max2[b] = tmp[1];
        tmp.clear();
        tmp = {
   min1[a], min2[a], min1[b], min2[b]};
        sort(tmp.begin(), tmp.end());
        min1[b] = tmp[0], min2[b] = tmp[1];

        maxv = max({
   maxv, max1[b] * max2[b], min1[b] * min2[b]});
    }
    if (maxv == -INF) return {
   cnt, 0};
    return {
   cnt, maxv};
}


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

    cin >> n, m = 122;
    cin >> s;
    s = ' ' + s;
    for (int i = 1; i <= n; ++i) cin >> w[i];

    get_sa();
    get_height();

    for (int i = 2; i <= n; ++i) hs[height[i]].push_back(i);

    for (int i = 1; i <= n; ++i) {
   
        p[i] = i, sz[i] = 1;
        max1[i] = min1[i] = w[sa[i]];
        max2[i] = -INF, min2[i] = INF;
    }

    for (int r = n - 1; r >= 0; --r) ans[r] = solve(r);
    for (int i = 0; i < n; ++i) cout << ans[i].x << ' ' << ans[i].y << '\n';

    return 0;
}

*生成魔咒

动态的过程

• 每次对字符串结尾添加一个字符
• 查询当前字符串子串的数量

首先考虑如何枚举子串, 首先枚举起点, 再枚举终点, 其实就是每个后缀$i$的所有前缀的集合就是子串的集合, 因此可以使用后缀数组

首先考虑简单的情况, 对于固定的字符串, 有多少个不同的子串

因此对于固定的字符串, 结果就是
$$\sum len(i)-height(i)$$

但是当前问题讨论的是有修改的情况, 如果在后面加一个字符, 会影响所有的后缀

但是可以反向考虑, 将字符串反转因为原串的不同子串的数量等价于反转之后的不同子串的数量

这样的好处是, 每次在前面添加一个字符, 最多添加一个后缀, 对已有的后缀是没有影响的

然后反向做这个问题, 先将所有串找出, 每次删除前面一个字符, 也就是删除一个前缀, 需要边删除边维护$height$数组

---

问题变为如何动态维护$height(i)$?

每次删除一个后缀相当于在排好序的后缀中删除一个序列

$sa$数组可以用一个双链表维护, 支持在$O(1)$时间内删除

假设希望删除的序列是虚线部分$k$, 对于$height(j)=lcp(j,k)$, $height(k)=lcp(k,i)$

现在要删除$k$, $height(j)=lcp(j,i)$, 由性质$lcp(i,j)=\min (lcp(i,k),lcp(k,j))$可知, 新的$lcp(j,i)=\min (lcp(k,i),lcp(j,k))$, 也可以在$O(1)$时间复杂度内计算出

示例代码

#include <bits/stdc++.h>

using namespace std;

typedef long long LL;
const int N = 1e5 + 10;

int n, m;
int s[N];
int sa[N], x[N], y[N], c[N], rk[N], height[N];
int u[N], d[N];
LL ans[N];

int get(int x) {
   
    static unordered_map<int, int> hash;
    if (!hash.count(x)) hash[x] = ++m;
    return hash[x];
}

void get_sa() {
   
    for (int i = 1; i <= n; ++i) c[x[i] = s[i]]++;
    for (int i = 2; i <= m; ++i) c[i] += c[i - 1];
    for (int i = n; i; --i) sa[c[x[i]]--] = i;

    for (int k = 1; k <= n; k <<= 1) {
   
        int num = 0;
        for (int i = n - k + 1; i <= n; ++i) y[++num] = i;
        for (int i = 1; i <= n; ++i) {
   
            if (sa[i] > k) y[++num] = sa[i] - k;
        }
        for (int i = 1; i <= m; ++i) c[i] = 0;
        for (int i = 1; i <= n; ++i) c[x[i]]++;
        for (int i = 2; i <= m; ++i) c[i] += c[i - 1];
        for (int i = n; i; --i) sa[c[x[y[i]]]--] = y[i], y[i] = 0;
        swap(x, y);
        x[sa[1]] = 1, num = 1;
        for (int i = 2; i <= n; ++i) {
   
            x[sa[i]] = y[sa[i]] == y[sa[i - 1]] && y[sa[i] + k] == y[sa[i - 1] + k] ? num : ++num;
        }
        if (num == n) break;
        m = num;
    }
}

void get_height() {
   
    for (int i = 1; i <= n; ++i) rk[sa[i]] = i;
    for (int i = 1, k = 0; i <= n; ++i) {
   
        if (rk[i] == 1) continue;
        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;
    }
}

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

    cin >> n;
    for (int i = n; i; --i) cin >> s[i], s[i] = get(s[i]);

    get_sa();
    get_height();

    LL res = 0;
    for (int i = 1; i <= n; ++i) {
   
        res += n - sa[i] + 1 - height[i];
        u[i] = i - 1, d[i] = i + 1;
    }
    d[0] = 1, u[n + 1] = n;

    for (int i = 1; i <= n; ++i) {
   
        ans[i] = res;
        int k = rk[i], j = d[k];

        res -= n - sa[k] + 1 - height[k];
        res -= n - sa[j] + 1 - height[j];
        height[j] = min(height[j], height[k]);
        res += n - sa[j] + 1 - height[j];
        d[u[k]] = d[k], u[d[k]] = u[k];
    }

    for (int i = n; i; --i) cout << ans[i] << '\n';

    return 0;
}