后缀数组详解
后缀数组定义:让人懵逼的有力工具 ,指某一字符串后缀按照字典序的一个排列。sa[i] = j 的含义为所有后缀按照字典序排列,排在第 i 个的是后缀 j 。
1.各变量的含义
此算法涉及到很多坨数组及变量,故在这里做一个罗列:
sa[ i ]:表示排名为 i 的后缀的起始位置下标;
rak[ i ]:本应写作rank数组,但因与编译器关键字重名,这里写作rak[ i ] ; 表示起始位置为 i 的后缀的排名;
tp[ i ]:基数排序的第二关键字,表示第二关键字排名为 i 的后缀的起始下标;
tax[ i ]:i 号元素出现了多少次,用于辅助基数排序;
s :字符串,s[i]表示字符串中第 i 个字符串;
lcp(x,y):字符串x与字符串y的最长公共前缀,这里指排名为x与排名为y的最长公共前缀;
height[ i ]:lcp(sa[i],sa[i−1]),即排名为 i 的后缀的字符串与排名为 i−1 的后缀的字符串的最长公共前缀;
H[ i ]:height[ rak[ i ] ],即 i 号前缀与前一名(不一定是i-1号)的最长公共前缀。
在说一下 ran 数组和 ra 数组便于深刻理解:他们有下面的关系等式:ran[sa[i]]=i
和sa[ran[i]]=i
,仔细揣摩一下就能理解其含义~。
2.具体思想
求 sa
数组其实关键在于排序,如果直接sort快排的话,每一次的时间复杂度是 O(nlogn) ,再加上比较的时间复杂度,其总时间复杂度会达到 O(n*n*logn) 。
我们考虑到是按字典序排序,所以我们用基数排序,每次对单字符排序,之后用单字符对多字符排序,以此类推,直至所有后缀顺序都不相同为止。这里在基数排序的过程中,不是一个一个的字符相加,而是用倍增的思想,因为如果一个字符一个字符相加,会出现某些字符重复排序的情况,所以直接倍增即可。考虑到裸基数排序的时间复杂度是 O(n),这里在基数排序基础上加个二分(倍增),其时间复杂度进一步降低为O(logn),在加上比较的时间复杂度为O(n),所以总的时间复杂度为O(nlogn),相当可观。
下面是倍增基数排序的图解:
3.代码实现及讲解
1. 基数排序
const int maxn=1e6+10;
char s[maxn];
int rak[maxn],sa[maxn],tax[maxn],tp[maxn];
void Qsort()
{
//M为桶的个数,及字符集的个数,len为字符串的长度+1。
for (int i = 0; i <= M; i++) tax[i] = 0; //把桶清零
for (int i = 1; i <= len; i++) tax[rak[i]]++; //统计每个名词出现的个数
for (int i = 1; i <= M; i++) tax[i] += tax[i - 1]; //做前缀和
//可以快速定位每个位置应有的排名,可以统计比当前名次小的后缀有多少个
for (int i = len; i >= 1; i--) sa[ tax[rak[tp[i]]]-- ] = tp[i]; //&@#$%......
//i从大到小依次枚举,那么sa[tax[rak[tp[i]]]−−]的意思就是说:用rak[i]定位到第一关键字的大小;那么tax[rak[tp[i]]]就表示当第一关键字相同时,第二关键字较大的这个后缀的排名是啥到了排名,我们也就能更新sa了,--表示减去自身,其他应有多少个排名。
}
2.倍增
void SuffixSort()
{
M = 1010; //字符集的大小,一共需要多少个桶
for (int i=1;i<=len;i++)
rak[i] = s[i]-'0'+1,tp[i]=i; //初始化rak和tp,注意rak的加1
Qsort();
for (int w=1;w<=len;w<<=1)
{
//w:当前倍增的长度,w = x表示已经求出了长度为x的后缀的排名,现在要更新长度为2x的后缀的排名
//p表示不同的后缀的个数,很显然原字符串的后缀都是不同的,因此p = N时可以退出循环
int p = 0; //这里的p仅仅是一个计数器
for (int i=len-w+1;i<=len;i++)
tp[++p]=i;
for (int i=1;i<=len;i++){
if (sa[i]>w)
tp[++p] = sa[i]-w; //这两个for是后缀数组的核心部分,是对第二关键字排序
}
/*假设我们现在需要得到的长度为w,那么sa[i]表示的实际是长度为w/2的后缀中排名为i的位置(也就是上一轮的结果)
我们需要得到的tp[i]表示的是:长度为w的后缀中,第二关键字排名为i的位置。
之所以能这样更新,是因为i号后缀的前w/2个字符形成的字符串是i−w/2号后缀的后w/2个字符形成的字符串*/
Qsort(); //此时我们已经更新出了第二关键字,利用上一轮的rak更新本轮的sa
memcpy(tp,rak,sizeof(rak)); //这里原本tp已经没有用了
rak[sa[1]]=p=1;
for (int i=2;i<=len;i++)
rak[sa[i]]= (tp[sa[i-1]]==tp[sa[i]]&&tp[sa[i-1]+w]==tp[sa[i]+w])?p:++p;
//这里当两个后缀上一轮排名相同时本轮也相同,至于为什么大家可以思考一下
if(p==len) break;
M=p;
}
}
4.精髓:height数组
height数组在变量含义中也提过:它表示lcp( sa[i],sa[i-1] );
H数组也同上说过:他表示height[ rak[ i ] ],即从i开始的后缀与排名前一名的后缀的最长公共前缀。
性质:\(H[i] >= H[i-1]-1\)
证明?记了也不懂系列,干脆不计了,如果确实想了解的话,可以参考博客:https://www.cnblogs.com/zwfymqz/p/8413523.html#_label4
求height数组代码
之所以能够线性求出height数组,还要依靠上面那条重要性质,show code:
void GetHeight()
{
int j,k=0;
for(int i=1;i<=len;i++)
{
if(k) k--;
int j=sa[rak[i]-1];
while(s[i+k]==s[j+k]) k++;
Height[rak[i]]=k;
//printf("%d\n", k);
}
}
5.模板(无注释版)
上面看懂了吗?没看懂?没关系,会用就行了(虽然蒟蒻我也没看懂,待我再研究研究),下面上无注释版模板:
基排:
void Qsort()
{
for (int i=0;i<=M;i++) tax[i] = 0;
for (int i=1;i<=len;i++) tax[rak[i]]++;
for (int i=1;i<=M;i++) tax[i]+=tax[i - 1];
for (int i=len;i>= 1;i--) sa[ tax[rak[tp[i]]]-- ]=tp[i];
}
倍增:
void suffixsort()
{
M=1010;
for(int i=1;i<=len;++i){
rak[i]=s[i];
tp[i]=i;
}
Qsort();
for(int w=1;w<=len;w<<=1)
{
int p=0;
for(int i=len-w+1;i<=len;++i)
tp[++p] = i;
for(int i=1;i<=len;++i)
if(sa[i]>w)
tp[++p]=sa[i]-w;
Qsort();
memcpy(tp,rak,sizeof(rak));
rak[sa[1]]=p=1;
for(int i=2;i<=len;++i)
rak[sa[i]]=( (tp[sa[i-1]]==tp[sa[i]]) && (tp[sa[i-1]+w]==tp[sa[i]+w]) )?p:++p;
if(p==len) break;
M=p;
}
}
线性求height数组:
void GetHeight()
{
int j,k=0;
for(int i=1;i<=len;i++)
{
if(k) k--;
int j=sa[rak[i]-1];
while(s[i+k]==s[j+k]) k++;
Height[rak[i]]=k;
//printf("%d\n", k);
}
}
AC代码(求LCP,连接取max height值即可):
#include<iostream>
#include<cstdlib>
#include<iomanip>
#include<algorithm>
#include<cstring>
using namespace std;
const int maxn=1e6+10;
char s[maxn],str[maxn];
int rak[maxn],tax[maxn],tp[maxn],sa[maxn],Height[maxn];
int len1,len,M;
void Qsort()
{
for(int i=0;i<=M;++i) tax[i]=0;
for(int i=1;i<=len;++i) tax[rak[i]]++;
for(int i=1;i<=M;++i) tax[i]+=tax[i-1];
for(int i=len;i>=1;i--) sa[ tax[rak[tp[i]]]-- ]=tp[i];
}
void suffixsort()
{
M=1010; //桶的个数
for(int i=1;i<=len;++i){
rak[i]=s[i];
tp[i]=i;
}
Qsort();
for(int w=1;w<=len;w<<=1)
{
int p=0;
for(int i=len-w+1;i<=len;++i)
tp[++p] = i;
for(int i=1;i<=len;++i)
if(sa[i]>w)
tp[++p]=sa[i]-w;
Qsort();
memcpy(tp,rak,sizeof(rak));
rak[sa[1]]=p=1;
for(int i=2;i<=len;++i)
rak[sa[i]]=( (tp[sa[i-1]]==tp[sa[i]]) && (tp[sa[i-1]+w]==tp[sa[i]+w]) )?p:++p;
if(p==len) break;
M=p;
}
}
void getheight()
{
int k=0;
for(int i=1;i<=len;++i)
{
if(k) k--;
int j=sa[rak[i]-1];
while(s[i+k]==s[j+k]) k++;
Height[rak[i]]=k;
}
}
//yeshowmuchiloveyoumydearmotherreallyicannotbelieveit#yeaphowmuchiloveyoumydearmother
int main()
{
ios::sync_with_stdio(false);
cin>>str+1;
len1=strlen(str+1);
str[len1+1]='#';
cin>>str+len1+2;
len=strlen(str+1);
for(int i=1;i<=len;++i)
s[i]=str[i];
suffixsort();
getheight();
int res=-1;
for(int i=1;i<=len;++i)
{
if((sa[i]<=len1&&sa[i-1]>len1+1)||(sa[i]>len1+1&&sa[i-1]<=len1))
res=max(res,Height[i]);
}
cout<<res<<endl;
system("pause");
return 0;
}
poj 1743:求不可重叠最长公共子串(后缀数组+二分):
#include<iostream>
#include<cstdlib>
#include<algorithm>
#include<cstring>
#include<stdio.h>
using namespace std;
const int maxn=2e4+10;
int rak[maxn],sa[maxn],tax[maxn],tp[maxn],Height[maxn];
int a[maxn]; //这里其实每个数组元素就相当于一个字符,连起来就相当于一个字符串,不要死板的以为只能输char数组,有很多转换变形
int M,len;
void Qsort()
{
for(int i=0;i<=M;++i) tax[i]=0;
for(int i=1;i<=len;++i) tax[rak[i]]++;
for(int i=1;i<=M;++i) tax[i]+=tax[i-1];
for(int i=len;i>=1;--i) sa[tax[rak[tp[i]]]--] = tp[i];
}
void suffixsort()
{
M=210;
for(int i=1;i<=len;++i)
{
rak[i]=a[i];
tp[i]=i;
}
Qsort();
for(int w=1;w<=len;w<<=1)
{
int p=0;
for(int i=len-w+1;i<=len;++i)
tp[++p]=i;
for(int i=1;i<=len;++i)
if(sa[i]>w)
tp[++p]=sa[i]-w;
Qsort();
swap(rak,tp);
rak[sa[1]]=p=1;
for(int i=2;i<=len;++i){
rak[sa[i]]=((tp[sa[i-1]]==tp[sa[i]])&&(tp[sa[i-1]+w]==tp[sa[i]+w]))?p:++p;
if(p==len) break;
M=p;
}
}
}
void getHeight()
{
int j,k=0;
for(int i=1;i<=len;++i){
if(k) k--;
int j=sa[rak[i]-1];
while(a[i+k]==a[j+k]) k++;
Height[rak[i]] = k;
}
}
int check(int x)
{
int mx=sa[1],mi=sa[1]; //mx为sa最大值,mi为sa最小值
for(int i=2;i<=len;++i){
if(Height[i]<x) mx=mi=sa[i];
else{
if(sa[i]<mi) mi=sa[i];
if(sa[i]>mx) mx=sa[i];
if(mx-mi>x) return 1; //sa最大值与最小值差x个(不能等于x,否则首尾有一个重叠)说明不重叠,满足条件
}
}
return 0;
}
int main()
{
//ios::sync_with_stdio(false);
while(scanf("%d",&len)!=EOF)
{
if(len==0) break;
for(int i=1;i<=len;++i)
scanf("%d",&a[i]);
for(int i=1;i<len;++i)
a[i]=a[i+1]-a[i]+90; //防止出现负数
len--; //所有元素差分后长度减1
suffixsort();
getHeight();
//下面为二分代码
int res=0;
int l=1,r=len,mid;
while(l<r)
{
mid=(l+r)>>1;
if(check(mid)){ //满足则向右更新l值,求最大res
res=mid;
l=mid+1;
}
else{
r=mid;
}
}
if(res<4) printf("0\n");
else{
printf("%d\n",res+1); //因为差分过,所以只需判断是否长度为4即可,最后满足的答案也要加1
}
}
//system("pause");
return 0;
}
6.应用
1:给定一个字符串,求它们的两个后缀的最长公共前缀
解:假设它们是某一个字符串的后缀,其排名分别为 j 和 k ,那么最长公共前缀就是min(Height[rak[j]+1],Height[rak[j]+1],......Height[rak[k]])
2:最长可重复子串
解:求height数组,取其中最大值即可。因为任意两个后缀的最长公共前缀都是Height数组里面某一段的最小值,这个值一定不大于height数组里面的最大值。
3:最长不可重叠重复子串
解:先二分答案,将题目变为判定性问题:判断是否存在两个长度为 k 的子串是相同且不重叠的。解决这个问题还得用 height 数组。把排序后的后缀分成若干组,其中每组的后缀之间的 height 值都不小于 k 。然后判断每组后缀,其最大 sa 值与最小 sa 值之差是否大于等于k。如有一组满足,则存在长度为 K 的子串重复且不重叠。
4:可重叠的最少出现k次的最长重复子串
解:和上例类似,先二分答案,然后将后缀数组分成若干组。不同的是,这里要判断的是有没有一个组的后缀个数大于等于k。如果有,那么存在k个相同的子串满足条件,否则不存在。
5:不相同的子串的个数
解:即求所有后缀之间不相同的前缀的个数。我们如果按 suffix(sa[1])
、suffix(sa[2])
......suffix(sa[n])
的顺序计算,可以发现,每次新加进来的后缀 suffix(sa[k])
,都将产生 n-sa[k]+1 个新的后缀,而其中 height[k] 个适合前面字符串相同的。所以suffix([sa[k]])
的真是贡献是n-sa[k]+1-height[k]
个不同的子串,累加一下即可。
6:最长回文子串
解:将一个不可能出现的符号加在这个字符串后面,再讲该字符串反着来加到这个字符串后面。这样问题就变成了求这个新的字符串的最长公共前缀。这样做的时间复杂度为 O(nlogn),其中用 RMQ 算法做预处理可将时间复杂度变为O(n)。(其实也可以用 Manacher 来求。)
7:连续重复子串
解:问题为一个字符串是由某个字符串s重复R次得到的,求R的最大值。我们穷举k即可,判断条件 k能整除字符串长度L整除以及suffix(1)与suffix(k+1)的最长公共前缀是否等于n-k。(用kmp也可以求)
8:重复次数最多的连续重复子串
解:先穷举长度 L,然后求长度为 L 的子串最多能连续出现几次。首先连续出现1 次是肯定可以的,所以这里只考虑至少 2 次的情况。假设在原字符串中连续出现 2 次,记这个子字符串为 S,那么 S 肯定包括了字符 r[0], r[L], r[L*2],
r[L*3], ……中的某相邻的两个。所以只须看字符 r[L*i]和 r[L*(i+1)]往前和往后各能匹配到多远,记这个总长度为 K,那么这里连续出现了 K/L+1 次。最后
看最大值是多少。时间复杂度为O(nlogn)。
9:两个字符串的最长公共子串
用一个不可能出现的字符将两个字符串接起来,求一下height数组最大且不是同一个字符串中的即可(判断sa下标和两个字符串之间长度即可)。O(lena+lenb)。