VEXOBEN
Vexoben
Feb 15, 2019
It takes 24 minutes to read this article.

本文介绍后缀数组的构造和简单应用。

大部分内容可以在罗穗骞2009年集训队中看到,本文进行了部分精简。

(刷新以获取数学公式)

定义

1、后缀(\(Suffix\)):对于字符串\(S[1..n]\),用\(Suffix(i)\)表示子串\(S[i..n]\)

2、后缀数组(\(Sa\)):后缀数组 \(Sa\) 是一个一维数组,它保存 \(1..n\) 的某个排列 \(SA[1], SA[2],……,SA[n],\)并且保证 \(Suffix(SA[i])\) < \(Suffix(SA[i+1])\),\(1≤i\)<\(n\)。

3、名次数组(\(Rank\)):名次数组 \(Rank[i]\) 保存的是 \(Suffix(i)\) 在所有后缀中从小到大排列的“名次”。

4、最长公共前缀(\(lcp\)):对于字符串\(A\),\(B\),定义\(lcp(A,B)\)为字符串\(A\),\(B\)的最长公共前缀。对于整数\(i\),\(j\),定义\(lcp(i,j)=lcp(Suffix(i),Suffix(j))\)

5、\(height\)数组:\(height[i]=lcp(sa[i-1],sa[i])\),即排名为\(i\)和\(i-1\)的两个后缀的最长公共前缀。

由\(height\)数组的定义可以得到它的一个性质:对于\(i\)和\(j\)不妨设\(rank[i]\)<\(rank[j]\),那么\(lcp(i,j)=min_{i+1≤k≤j}\{height[k]\}\)

倍增构造后缀数组

基数排序

对于\(n\)个双关键字数据,数据范围在\([1..m]\)之间,已知两个数组:\(x[i]\)表示第\(i\)个数第一关键字的值,\(y[i]\)表示第二关键字排名为\(i\)的数的位置(这里\(y[i]\)是两两不同的,如果两数两个关键字都相同,我们可以用序号比较出他们的第二关键字)我们可以用基数排序在\(O(n+m)\)的复杂度内完成对他们的排序。

对第一关键字开\(m\)个桶,然后按第二关键字从大到小的顺序将他们扔进桶中即可完成排序。

for (int i = 1; i <= m; ++i) cnt[i] = 0;
for (int i = 1; i <= n; ++i) ++cnt[x[i]];
for (int i = 2; i <= m; ++i) cnt[i] += cnt[i-1];
for (int i = n; i >= 1; --i) sa[cnt[x[y[i]]]--] = y[i];

倍增

我们进行\(O(logn)\)次排序.第\(x\)次排序以\(S[i..i+2^x)\)为第一关键字,以\(S[i+2^x..i+2^{x+1})\)为第二关键字.不过最开始要以第二关键字是序号(即\(y[i]=i\))做一次.

#include<iostream>
#include<cstdio>
#include<cstring>
#define rint register int
#define inv inline void
#define ini inline int
#define maxn 1000050
using namespace std;
char s[maxn];
int y[maxn],x[maxn],c[maxn],sa[maxn],rk[maxn],height[maxn],wt[30];
int n,m;
liti
inv putout(int x)
{
    if(!x) {putchar(48);return;}
    rint l=0;
    while(x) wt[++l]=x%10,x/=10;
    while(l) putchar(wt[l--]+48);
}

inv get_SA()
{
    for (rint i=1;i<=n;++i) ++c[x[i]=s[i]];
    //c数组是桶 
    //x[i]是第i个元素的第一关键字 
    for (rint i=2;i<=m;++i) c[i]+=c[i-1]; 
    //做c的前缀和,我们就可以得出每个关键字最多是在第几名 
    for (rint i=n;i>=1;--i) sa[c[x[i]]--]=i; 
    for (rint k=1;k<=n;k<<=1)
    {
        rint num=0;
        for (rint i=n-k+1;i<=n;++i) y[++num]=i;
        //y[i]表示第二关键字排名为i的数,第一关键字的位置 
        //第n-k+1到第n位是没有第二关键字的 所以排名在最前面 
        for (rint i=1;i<=n;++i) if (sa[i]>k) y[++num]=sa[i]-k;
        //排名为i的数 在数组中是否在第k位以后
        //如果满足(sa[i]>k) 那么它可以作为别人的第二关键字,就把它的第一关键字的位置添加进y就行了
        //所以i枚举的是第二关键字的排名,第二关键字靠前的先入队 
        for (rint i=1;i<=m;++i) c[i]=0;
        //初始化c桶 
        for (rint i=1;i<=n;++i) ++c[x[i]];
        //因为上一次循环已经算出了这次的第一关键字 所以直接加就行了 
        for (rint i=2;i<=m;++i) c[i]+=c[i-1];//第一关键字排名为1~i的数有多少个 
        for (rint i=n;i>=1;--i) sa[c[x[y[i]]]--]=y[i],y[i]=0;
        //因为y的顺序是按照第二关键字的顺序来排的 
        //第二关键字靠后的,在同一个第一关键字桶中排名越靠后 
        //基数排序 
        swap(x,y);
        //这里不用想太多,因为要生成新的x时要用到旧的,就把旧的复制下来,没别的意思 
        x[sa[1]]=1;num=1;
        for (rint 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;
        //因为sa[i]已经排好序了,所以可以按排名枚举,生成下一次的第一关键字 
        if (num==n) break;
        m=num;
        //这里就不用那个122了,因为都有新的编号了 
    }
    for (rint i=1;i<=n;++i) putout(sa[i]),putchar(' ');
}

int main()
{
    gets(s+1);
    n=strlen(s+1);m=122;
    //我们第一次读入字符直接不用转化,按原来的ascll码来就可以了 
    get_SA();
}

\(height\)数组的构造

性质

定义\(h[i]=height[rank[i]]\),也就是 \(suffix(i)\)和在它前一名的后缀的最长公共前缀。

\(h\) 数组有以下性质: \(h[i]≥h[i-1]-1\)

当\(h[i-1]=1\)时显然成立.当\(h[i-1]>1\)时,设 \(suffix(k)\) 是排在 \(suffix(i-1)\)前一名的后缀,则它们的最长公共前缀是\(h[i-1]\)。那么\(suffix(k+1)\)将排在\(suffix(i)\)的前面并且 \(suffix(k+1)\)和 \(suffix(i)\)的最长公共前缀是\(h[i-1]-1\),所以 \(suffix(i)\)和在它前一名的后缀的最长公共前缀至少是 \(h[i-1]-1\)。

构造

利用这个\(h\)数组的性质,按照 \(h[1],h[2],……,h[n]\)的顺序可以在线性时间内处理出\(height\)数组.

void calc_height() {
	for (int i = 1; i <= n; ++i) rnk[sa[i]] = i;
	for (int i = 1; i <= n; ++i) {
		int p = sa[rnk[i] - 1], k = max(0, hgt[rnk[i - 1]] - 1);
		while (s[p + k] == s[i + k]) ++k;
		hgt[rnk[i]] = k;
	}
}

例题

可重叠的 k 次最长重复子串

来源:POJ3261

题意

给定一个字符串,求至少出现 \(k\) 次的最长重复子串,这 \(k\) 个子串可以重叠。

保证\(s.length ≤ 100000\)

题解

后缀数组的一个技巧:用\(height\)数组将后缀分组。

题目等价于:求一个最长的子串,它是至少 \(k\) 个后缀的前缀。

先二分一个答案 \(mid\) ,如果 \(height[i] ≥ mid\), 将 \(sa[i]\) 和 \(sa[i - 1]\) 分在一组。只需记录有没有大小超过 \(k\) 的组即可。

#include<cstdio>
#include<iostream>
using namespace std;
const int N = 1e6 + 10;

int n, k;
int a[N], x[N], y[N], sa[N], rnk[N], cnt[N], hgt[N];

void calc_sa() {
	int m = 1e6 + 1;
	for (int i = 1; i <= n; ++i) ++cnt[x[i] = a[i]];
	for (int i = 2; i <= m; ++i) cnt[i] += cnt[i - 1];
	for (int i = n; i >= 1; --i) sa[cnt[x[i]]--] = i;
	for (int w = 1, p; w <= n; w <<= 1, m = p) {
		p = 0;
		for (int i = n - w + 1; i <= n; ++i) y[++p] = i;
		for (int i = 1; i <= n; ++i) if (sa[i] > w) y[++p] = sa[i] - w;
		for (int i = 1; i <= m; ++i) cnt[i] = 0;
		for (int i = 1; i <= n; ++i) ++cnt[x[i]];
		for (int i = 1; i <= m; ++i) cnt[i] += cnt[i - 1];
		for (int i = n; i >= 1; --i) sa[cnt[x[y[i]]]--] = y[i];
		swap(x, y);
		x[sa[1]] = 1; p = 1;
		for (int i = 2; i <= n; ++i) {
			int f = (y[sa[i]] == y[sa[i - 1]] && y[sa[i] + w] == y[sa[i - 1] + w]);
			x[sa[i]] = f ? p : ++p;
		}
		if (p == n) break;
	}
}

void calc_height() {
	for (int i = 1; i <= n; ++i) rnk[sa[i]] = i;
	for (int i = 1; i <= n; ++i) {
		int p = sa[rnk[i] - 1], k = max(0, hgt[rnk[i - 1]] - 1);
		while (a[p + k] == a[i + k]) ++k;
		hgt[rnk[i]] = k;
	}
}

int check(int mid) {
	int cnt = 0;
	for (int i = 2; i <= n; ++i) {
		if (hgt[i] < mid) cnt = 1;
		else if (++cnt >= k) return 1;
	}
	return 0;
}

int main() {
	cin >> n >> k;
	for (int i = 1; i <= n; ++i) {
		scanf("%d", &a[i]);
		++a[i];
	}
	calc_sa();
	calc_height();
	int l = 1, r = n, ans = 0;
	while (l <= r) {
		int mid = (l + r) >> 1;
		if (check(mid)) {
			l = mid + 1;
			ans = mid;
		}
		else r = mid - 1;
	}
	cout << ans << endl;
	return 0;
}

本质不同子串个数

来源:HDU4622

题意

给定字符串\(S\), \(q\)组询问,每次询问\(S[l..r]\)本质不同子串个数。

保证 \(s.length≤2000\), \(q ≤ 10000\)

题解

先考虑计算\(S[1..n]\)本质不同子串个数。

字符串\(suffix(i)\)的前缀有\(n-sa[i]+1\)个,其中\(height[i]\)个和\(suffix(i-1)\)重复,所以每次将 \(n-sa[i]+1-lcp(i, i - 1)\)累加进答案即可。由于保证了字典序,很显然这种做法是正确的。

而在询问一个区间时,每个后缀被删掉了一个后缀。如果\(suffix(sa[i]) < suffix(sa[i-1])\),一定有\(suffix(sa[i])\)变为了\(suffix(sa[i-1])\)的前缀,那么我们直接无视掉\(suffix(sa[i])\)。对于剩下的串,按顺序将\(n-sa[i]+1-lcp(i, lcp)\)累加进答案即可。

还有一种做法是对于每个前缀求一遍后缀数组,然后从而将问题转化为询问\(S[l..n]\)本质不同子串个数。


void calc_sa();

void calc_height();

void Solve() {
	scanf("%s", s + 1);
	n = strlen(s + 1);
	calc_sa();
	calc_height();
	cin >> q;
	while (q--) {
		int l, r;
		scanf("%d%d", &l, &r);
		int ans = 0, lcp = 0;
		for (int i = 1; i <= n; ++i) {
			lcp = min(lcp, hgt[i]);
			if (l <= sa[i] && sa[i] <= r) {
				if (r - sa[i] + 1 - lcp > 0) {
					ans += max(0, r - sa[i] + 1 - lcp);
					lcp = r - sa[i] + 1;
				}
			}
		}
		printf("%d\n", ans);
	}
}

最长回文子串

题意

给定字符串\(S\),求最长回文子串。

保证\(s.length≤100000\)

题解

manacher模板题,也可以用后缀数组做到线性复杂度。

将字符串翻转加在原来的字符串后面,中间用一个特殊符号连接。

比如:\(S="aabaaaab"\),\(S'="aabaaaab?baaaabaa"\)

原来的一个回文串就可以转化为\(S'\)中两个后缀的最长公共前缀。

比如以\(S[3]\)为中心的回文串,对应\("aaaab?baaaabaa"\)和\("aa"\)的公共前缀。

以\(S[3]和S[4]\)为中心的回文串,对应\("aaaab?baaaabaa"\)和\("baa"\)的公共前缀。

用DC3求后缀数组,\(O(N)\)预处理RMQ就可以做到线性复杂度了。

大概是没有人愿意打的。

最长公共子串

来源:POJ2774

题意

给定字符串\(A\),\(B\),求最长公共子串。

保证\(a,b.length≤100000\)

题解

将字符串拼成 \(A?B\),找在\('?'\)两边的最大\(height\)即可。


void calc_sa();

void calc_height();

int main() {
	scanf("%s", s + 1); n = strlen(s + 1);
	scanf("%s", t + 1); m = strlen(t + 1);
	s[n + 1] = 'A';
	for (int i = 1; i <= m; ++i) s[n + i + 1] = t[i];
	calc_sa(n + m + 1);
	calc_height(n + m + 1);
	int ans = 0;
	for (int i = 1; i <= n + m + 1; ++i) {
		if (sa[i] == n + 1) continue;
		if (sa[i - 1] == n + 1) continue;
		int flag = (sa[i] <= n) ^ (sa[i - 1] <= n);
		if (flag) ans = max(ans, hgt[i]);
	}
	cout << ans << endl;
	return 0;
}

长度不小于 k 的公共子串的个数

来源:POJ3415

题意

给定字符串\(A\),\(B\),求三元组\((i,j,k)\)的数目,使得\(A[i..i+k-1]=B[j..j+k-1],k>=K\)

保证\(a,b.length≤100000\)

题解

将字符串拼成 \(A?B\)。加入\(suffix(sa[i])\)时,对于\(j\)<\(i\),产生的贡献是\(([sa_i≤n] xor [sa_j≤n])*min_{j+1≤k≤i}\{height_k\}\)

对\(sa_i≤n\)和\(sa_i>n\)分别维护单调栈。每次将单调栈中所有大于\(height_{i+1}\)的元素更改为\(height_{i+1}\)。注意要把相同的元素压在一起,这样每弹出一次栈中元素就少一个,每次只加入一个元素,就保证了这部分复杂度线性。


void calc_sa();

void calc_height();

int main() {
	while (scanf("%d", &k) && k) {
		scanf("%s", s + 1); n = strlen(s + 1);
		scanf("%s", t + 1); m = strlen(t + 1);
		s[n + 1] = '!';
		for (int i = 1; i <= m; ++i) s[n + i + 1] = t[i];
		calc_sa(n + m + 1);
		calc_height(n + m + 1);
		stack<pair<int, int> > st[2];
		long long ans = 0, val[2] = {0, 0};
		for (int i = 1; i <= n + m + 1; ++i) {
			if (sa[i] == n + 1) continue;
			int p = sa[i] <= n, v = hgt[i + 1], cnt[2] = {0, 0};
			ans += val[p];			
			++cnt[p ^ 1];
			val[p ^ 1] += max(0, v - k + 1);
			for (int j = 0; j < 2; ++j) {
				while (!st[j].empty() && st[j].top().first > v) {
					pair<int, int> tmp = st[j].top();
					st[j].pop();
					val[j] -= 1LL * tmp.second * max(0, tmp.first - k + 1);
					val[j] += 1LL * tmp.second * max(0, v - k + 1);
					cnt[j] += tmp.second;
				}
				if (cnt[j]) st[j].push(make_pair(v, cnt[j]));
			}
 		}
		printf("%lld\n", ans);
	}
	return 0;
}