后缀数组 学习笔记
定义
定义字符串从第 \(i\) 个字符开始的后缀为 \(suf_i\)。
后缀数组一种是处理字符串的有力工具。将字符串的所有后缀按字典序排序,则 \(sa_i\) 表示排名为 \(i\) 的后缀的起始位置,\(rk_i\) 为第 \(i\) 个后缀的排名。
倍增算法
直接排序,由于要进行字符串的比较,复杂度 \(O(n^2\log n)\)。
然而,所有后缀中有大量重合的地方。考虑一个倍增:枚举 \(2^k\),求所有长 \(2^k\) 的子串的排名,超出的部分视为无穷小的字符,当 \(2^k\geq n\) 时就得到了后缀数组。而长 \(2^k\) 的子串可以被拆成两个长 \(2^{k-1}\) 的子串,因此刚开始对长 \(1\) 的子串排序,后面只需要对旧的 \(rk_i,rk_{i+2^{k-1}}\) 做双关键字排序就得到了新的 \(sa_i,rk_i\)。
字符串的一个特点是字符集不会很大,因此可以使用不基于比较的排序,如计数排序、基数排序。复杂度 \(O(n\log n)\)。
void getsa(int n,int m,char s[],int sa[])
{
int x[n+5],y[n+5],b[max(n,m)+5];
memset(x,0,sizeof(x)),memset(y,0,sizeof(y)),memset(b,0,sizeof(b));
for(int i=1;i<=n;i++)b[x[i]=s[i]]++;
for(int i=1;i<=m;i++)b[i]+=b[i-1];
for(int i=n;i>=1;i--)sa[b[x[i]]--]=i;
for(int k=1;k<=n;k<<=1)
{
int p=0;
for(int i=n-k+1;i<=n;i++)y[++p]=i;
for(int i=1;i<=n;i++)if(sa[i]>k)y[++p]=sa[i]-k;
memset(b,0,sizeof(b));
for(int i=1;i<=n;i++)b[x[i]]++;
for(int i=1;i<=m;i++)b[i]+=b[i-1];
for(int i=n;i>=1;i--)sa[b[x[y[i]]]--]=y[i];
for(int i=1;i<=n;i++)swap(x[i],y[i]);
p=0;
for(int i=1;i<=n;i++)x[sa[i]]=y[sa[i]]==y[sa[i-1]]&&y[sa[i]+k]==y[sa[i-1]+k]?p:++p;
if(p==n)break;
m=p;
}
}
下面逐句解释代码:
int x[n+5],y[n+5],b[max(n,m)+5];
memset(x,0,sizeof(x)),memset(y,0,sizeof(y)),memset(b,0,sizeof(b));
for(int i=1;i<=n;i++)b[x[i]=s[i]]++;
for(int i=1;i<=m;i++)b[i]+=b[i-1];
for(int i=n;i>=1;i--)sa[b[x[i]]--]=i;
x,y,b
分别表示第 \(i\) 个后缀的第一关键字、第二关键字排名为 \(i\) 的后缀的起始位置、桶。桶同时用于第一轮的计数排序和后面的基数排序,而下面的值域可能来到 \(n\),因此开 \(\max(n,m)\)。
下面的是计数排序,对桶做前缀和之后,逆序扫一遍,当值相同时给大的编号分配靠后的位置,因此稳定。
int p=0;
for(int i=n-k+1;i<=n;i++)y[++p]=i;
for(int i=1;i<=n;i++)if(sa[i]>k)y[++p]=sa[i]-k;
这一部分对第二关键字排序。其中 \(suf_{n-2^{k}+1}\sim suf_n\) 在最前面,因为第二关键字为空。否则按排名枚举第二关键字,如果这个子串能作为第二关键字,就加入它的第一关键字。
memset(b,0,sizeof(b));
for(int i=1;i<=n;i++)b[x[i]]++;
for(int i=1;i<=m;i++)b[i]+=b[i-1];
for(int i=n;i>=1;i--)sa[b[x[y[i]]]--]=y[i];
这一部分对第一关键字排序,倒序枚举第二关键字。
for(int i=1;i<=n;i++)swap(x[i],y[i]);
p=0;
for(int i=1;i<=n;i++)x[sa[i]]=y[sa[i]]==y[sa[i-1]]&&y[sa[i]+k]==y[sa[i-1]+k]?p:++p;
if(p==n)break;
m=p;
这一部分计算新的第一关键字,也就是 \(rk\)。按当前的 \(sa\) 枚举,如果当前的串和上一个串不同,则排名增加。p
就是新的值域。
为什么 \(y\) 不用开两倍?因为虽然 \(sa_i+2^k\) 或 \(sa_{i-1}+2^k\) 可能超 \(n\),但当进入这一部分时有 \(s[sa_i,sa_i+2^k-1]=s[sa_{i-1},sa_{i-1}+2^k-1]\)。而如果 \(sa_i+2^k\) 或 \(sa_{i-1}+2^k\),则超出的部分都为空白,且空白的长度不同,所以此时不会进入后面。
height 数组
定义 \(\operatorname{lcp}\) 表示最长公共前缀的长度,\(ht_i=\operatorname{lcp}(suf_{sa_{i-1}},suf_{sa_i})\),也就是排名 \(i-1,i\) 的后缀的最长公共前缀。\(ht_1=0\)。
结论一:设 \(rk_i<rk_j\),则 \(\operatorname{lcp}(suf_i,suf_j)=\min_{k=rk_i+1}^{rk_j}ht_k\)。
考虑感性证明一下:在 \(sa\) 中,\(i,j\) 之间的后缀的字典序在 \(suf_i,suf_j\) 之间,则这些后缀的前 \(\operatorname{lcp}(suf_i,suf_j)\) 个字符相同。而在从 \(i\) 变化到 \(j\) 的过程中,至少有一步改变了 \(\operatorname{lcp}(suf_i,suf_j)+1\)。
结论二:\(ht_{rk_i}\geq ht_{rk_{i-1}}-1\)。
当 \(ht_{rk_{i-1}}\leq 1\),显然满足。
当 \(ht_{rk_{i-1}}>1\),\(sa_{rk_{i-1}-1}+1\leq n\)。令 \(suf_{i-1}=aAB,suf_{sa_{rk_{i-1}-1}}=aAC,suf_i=AB,suf_{sa_{rk_{i-1}-1}+1}=AC\)。其中 \(a\) 是一个字符,\(A,B,C\) 都是串,\(B<C\),则 \(suf_{sa_{rk_{i-1}-1}+1}\leq suf_{sa_{rk_i-1}}<suf_i\),\(suf_{sa_{rk_i}-1},suf_i\) 至少有公共前缀 \(A\)。
运用性质二,可以暴力匹配求 \(ht\)。因为 \(ht\) 不超过 \(n\),而每次只少一,因此复杂度 \(O(n)\)。
for(int i=1;i<=n;i++)rk[sa[i]]=i;
for(int i=1,k=0;i<=n;i++){
if(k)k--;
while(s[i+k]==s[sa[rk[i]-1]+k])k++;
ht[rk[i]]=k;
}
例题
经典结论:本质不同子串数为 \(\frac{n(n-1)}{2}-\sum_{i=1}^n ht_i\)。
证明:用全部子串减去重复。子串等于后缀的前缀,按字典序枚举所有后缀,每个后缀产生的重复子串数为这个后缀与之前的后缀的公共前缀,也就是 \(\max_{j=1}^{i-1}\operatorname{lcp}(suf_{sa_i},suf_{sa_j})\)。而根据结论一可知这个值就是 \(ht_i\)。
枚举 \(S_0\) 的每个左端点和 \(S\) 匹配。直接找到下一个失配的位置,需要快速求 \(S_0\) 的一段后缀和 \(S\) 的一段后缀的最长公共前缀。这可以将 \(S_0\) 和 \(S\) 拼接后算出 \(ht\),然后 ST 表维护。
#include<bits/stdc++.h>
using namespace std;
int c,n,m,sa[200005],rk[200005],ht[200005];
char s[200005],temp[100005],t[100005];
template<typename T,int maxn1,int maxn2>struct ST{
T st[maxn1][maxn2];
int lg[maxn1];
void build(int n,T a[]){
for(int i=2;i<=n;i++)lg[i]=lg[i>>1]+1;
for(int i=1;i<=n;i++)st[i][0]=a[i];
for(int j=1;1<<j<=n;j++)for(int i=1;i+(1<<j)-1<=n;i++)st[i][j]=min(st[i][j-1],st[i+(1<<(j-1))][j-1]);
}
T query(int l,int r){
return min(st[l][lg[r-l+1]],st[r-(1<<lg[r-l+1])+1][lg[r-l+1]]);
}
};
ST<int,200005,20>a;
void getsa(int n,int m,char s[],int sa[],int rk[],int ht[])
{
int x[n+5],y[n+5],b[max(m,n)+5];
memset(x,0,sizeof(x)),memset(y,0,sizeof(y)),memset(b,0,sizeof(b));
for(int i=1;i<=n;i++)b[x[i]=s[i]]++;
for(int i=1;i<=m;i++)b[i]+=b[i-1];
for(int i=n;i>=1;i--)sa[b[x[i]]--]=i;
for(int k=1;k<=n;k<<=1)
{
int p=0;
for(int i=n-k+1;i<=n;i++)y[++p]=i;
for(int i=1;i<=n;i++)if(sa[i]>k)y[++p]=sa[i]-k;
memset(b,0,sizeof(b));
for(int i=1;i<=n;i++)b[x[i]]++;
for(int i=1;i<=m;i++)b[i]+=b[i-1];
for(int i=n;i>=1;i--)sa[b[x[y[i]]]--]=y[i];
for(int i=1;i<=n;i++)swap(x[i],y[i]);
p=0;
for(int i=1;i<=n;i++)x[sa[i]]=y[sa[i]]==y[sa[i-1]]&&y[sa[i]+k]==y[sa[i-1]+k]?p:++p;
if(p==n)break;
m=p;
}
for(int i=1;i<=n;i++)rk[sa[i]]=i;
for(int i=1,k=0;i<=n;i++){
if(k)k--;
while(s[i+k]==s[sa[rk[i]-1]+k])k++;
ht[rk[i]]=k;
}
}
int main(){
cin>>c;
while(c--){
cin>>temp+1>>t+1,n=strlen(temp+1),m=strlen(t+1);
for(int i=1;i<=n;i++)s[i]=temp[i];
for(int i=1;i<=m;i++)s[n+i]=t[i];
getsa(n+m,84,s,sa,rk,ht),a.build(n+m,ht);
int ans=0;
for(int i=1;i<=n-m+1;i++){
int p=1;
for(int j=1;j<=4;j++){
if(rk[i+p-1]<rk[n+p])p+=a.query(rk[i+p-1]+1,rk[n+p])+(j!=4);
else p+=a.query(rk[n+p]+1,rk[i+p-1])+(j!=4);
if(p>m)break;
}
if(p>m)ans++;
}
cout<<ans<<'\n';
}
return 0;
}
问题即为求最长的重复出现至少 \(k\) 次的子串。先求出 \(ht_i\)。发现这 \(k\) 个子串所在的后缀一定是相邻的排名,否则不优。因此对每个长 \(k-1\) 的区间的最小 \(ht\) 取 \(\max\) 即可,可以单调队列。
#include<bits/stdc++.h>
using namespace std;
int n,k,s[20005],sa[20005],rk[20005],ht[20005],fp,bp=1,q[20005],ans;
void rebuild(int n,int a[],int cnt=0){
int temp[n+5];
for(int i=1;i<=n;i++)temp[i]=a[i];
sort(temp+1,temp+n+1),cnt=unique(temp+1,temp+n+1)-temp-1;
for(int i=1;i<=n;i++)a[i]=lower_bound(temp+1,temp+cnt+1,a[i])-temp;
}
void getsa(int n,int m,int s[],int sa[],int rk[],int ht[])
{
int x[n+5],y[n+5],b[max(m,n)+5];
memset(x,0,sizeof(x)),memset(y,0,sizeof(y)),memset(b,0,sizeof(b));
for(int i=1;i<=n;i++)b[x[i]=s[i]]++;
for(int i=1;i<=m;i++)b[i]+=b[i-1];
for(int i=n;i>=1;i--)sa[b[x[i]]--]=i;
for(int k=1;k<=n;k<<=1)
{
int p=0;
for(int i=n-k+1;i<=n;i++)y[++p]=i;
for(int i=1;i<=n;i++)if(sa[i]>k)y[++p]=sa[i]-k;
memset(b,0,sizeof(b));
for(int i=1;i<=n;i++)b[x[i]]++;
for(int i=1;i<=m;i++)b[i]+=b[i-1];
for(int i=n;i>=1;i--)sa[b[x[y[i]]]--]=y[i];
for(int i=1;i<=n;i++)swap(x[i],y[i]);
p=0;
for(int i=1;i<=n;i++)x[sa[i]]=y[sa[i]]==y[sa[i-1]]&&y[sa[i]+k]==y[sa[i-1]+k]?p:++p;
if(p==n)break;
m=p;
}
for(int i=1;i<=n;i++)rk[sa[i]]=i;
for(int i=1,k=0;i<=n;i++){
if(k)k--;
while(s[i+k]==s[sa[rk[i]-1]+k])k++;
ht[rk[i]]=k;
}
}
int main(){
cin>>n>>k,k--;
for(int i=1;i<=n;i++)cin>>s[i];
rebuild(n,s),getsa(n,n,s,sa,rk,ht);
for(int i=1;i<=n;i++){
while(fp>=bp&&ht[q[fp]]>ht[i])fp--;
q[++fp]=i;
while(q[bp]<=i-k)bp++;
if(i>=k)ans=max(ans,ht[q[bp]]);
}
return cout<<ans<<'\n',0;
}
转调后相同也就是差分数组相同。
二分答案,此时 \(ht\geq mid\) 的排名形成一段段区间,每个区间有长度为 \(mid\) 的公共子串。对于每个区间,判断最远的两个后缀的距离是否大于 \(mid\) 即可。
#include<bits/stdc++.h>
using namespace std;
int n,a[5005],sa[5005],rk[5005],ht[5005];
void getsa(int n,int m,int s[],int sa[],int rk[],int ht[])
{
int x[n+5],y[n+5],b[max(m,n)+5];
memset(x,0,sizeof(x)),memset(y,0,sizeof(y)),memset(b,0,sizeof(b));
for(int i=1;i<=n;i++)b[x[i]=s[i]]++;
for(int i=1;i<=m;i++)b[i]+=b[i-1];
for(int i=n;i>=1;i--)sa[b[x[i]]--]=i;
for(int k=1;k<=n;k<<=1)
{
int p=0;
for(int i=n-k+1;i<=n;i++)y[++p]=i;
for(int i=1;i<=n;i++)if(sa[i]>k)y[++p]=sa[i]-k;
memset(b,0,sizeof(b));
for(int i=1;i<=n;i++)b[x[i]]++;
for(int i=1;i<=m;i++)b[i]+=b[i-1];
for(int i=n;i>=1;i--)sa[b[x[y[i]]]--]=y[i];
for(int i=1;i<=n;i++)swap(x[i],y[i]);
p=0;
for(int i=1;i<=n;i++)x[sa[i]]=y[sa[i]]==y[sa[i-1]]&&y[sa[i]+k]==y[sa[i-1]+k]?p:++p;
if(p==n)break;
m=p;
}
for(int i=1;i<=n;i++)rk[sa[i]]=i;
for(int i=1,k=0;i<=n;i++){
if(k)k--;
while(s[i+k]==s[sa[rk[i]-1]+k])k++;
ht[rk[i]]=k;
}
}
bool check(int mid){
for(int i=1,minn,maxn;i<=n;i++){
if(ht[i]>=mid){
minn=min(minn,sa[i]),maxn=max(maxn,sa[i]);
if(maxn-minn>mid)return 1;
}
else minn=maxn=sa[i];
}
return 0;
}
int main(){
cin>>n;
for(int i=1;i<=n;i++)cin>>a[i];
n--;
for(int i=1;i<=n;i++)a[i]=a[i]-a[i+1]+88;
getsa(n,175,a,sa,rk,ht);
int l=1,r=n;
while(l<r){
int mid=(l+r+1)>>1;
if(check(mid))l=mid;
else r=mid-1;
}
return cout<<(l<4?0:l+1)<<'\n',0;
}
差分后转化为求多个串的最长公共子串。
将这些串用不在字符集中且互不相同的特殊字符连接。二分答案,对于每个极大的 \(ht\geq mid\) 的区间,如果这一段区间的后缀来自全部 \(n\) 个串,就说明有长 \(mid\) 的最长公共子串。
#include<bits/stdc++.h>
using namespace std;
int n,m,b[10010],sa[10010],rk[10010],ht[10010],vis[10];
char a[10010],temp[2005];
void getsa(int n,int m,char s[],int sa[],int rk[],int ht[])
{
int x[n+5],y[n+5],b[max(m,n)+5];
memset(x,0,sizeof(x)),memset(y,0,sizeof(y)),memset(b,0,sizeof(b));
for(int i=1;i<=n;i++)b[x[i]=s[i]]++;
for(int i=1;i<=m;i++)b[i]+=b[i-1];
for(int i=n;i>=1;i--)sa[b[x[i]]--]=i;
for(int k=1;k<=n;k<<=1)
{
int p=0;
for(int i=n-k+1;i<=n;i++)y[++p]=i;
for(int i=1;i<=n;i++)if(sa[i]>k)y[++p]=sa[i]-k;
memset(b,0,sizeof(b));
for(int i=1;i<=n;i++)b[x[i]]++;
for(int i=1;i<=m;i++)b[i]+=b[i-1];
for(int i=n;i>=1;i--)sa[b[x[y[i]]]--]=y[i];
for(int i=1;i<=n;i++)swap(x[i],y[i]);
p=0;
for(int i=1;i<=n;i++)x[sa[i]]=y[sa[i]]==y[sa[i-1]]&&y[sa[i]+k]==y[sa[i-1]+k]?p:++p;
if(p==n)break;
m=p;
}
for(int i=1;i<=n;i++)rk[sa[i]]=i;
for(int i=1,k=0;i<=n;i++){
if(k)k--;
while(s[i+k]==s[sa[rk[i]-1]+k])k++;
ht[rk[i]]=k;
}
}
bool check(int mid){
for(int i=1;i<=n;i++)vis[i]=0;
for(int i=1,t=1,num=0;i<=m;i++){
if(ht[i]<mid)t++,num=0;
if(vis[b[sa[i]]]!=t)num++,vis[b[sa[i]]]=t;
if(num==n)return 1;
}
return 0;
}
int main(){
cin>>n;
for(int i=1,t;i<=n;i++){
cin>>temp+1,t=strlen(temp+1);
for(int j=1;j<=t;j++)a[++m]=temp[j],b[m]=i;
a[++m]=i+'A';
}
getsa(m,122,a,sa,rk,ht);
int l=0,r=2000;
while(l<r){
int mid=(l+r+1)>>1;
if(check(mid))l=mid;
else r=mid-1;
}
return cout<<l<<'\n',0;
}
[[字符串]]