[学习笔记]后缀数组

参考:

后缀数组 最详细讲解

上面一篇是转载这一篇的:

后缀数组 学习笔记

 

一、

后缀:suff(i),后缀要排序

sa[i],排名为i的后缀开始位置

rk[i],i开始位置的后缀的排名。

rk[sa[i]]=i,sa[rk[i]]=i

 

后缀要排序。

暴力排序:O(n^2logn)

可以利用倍增+基数排序优化到O(nlogn)

以下图片转载自:后缀数组 最详细讲解

 

为了避免常数过大,不能用vector当桶来用。

以及各种优化具体代码实现在下面。

 

二、

不能光排序。没有用。

设LCP(i,j)为,排名为i和排名为j的后缀的最长公共前缀长度。(注意i,j是排名)

即:LCP(suff[sa[i]],suff[sa[j]])(这个LCP是最长公共前缀的意思)

 

性质:

1.

LCP(i,j)=LCP(j,i)

LCP(i,i)=len(sa[i])=n-sa[i]+1

这个太平凡了。

 

 

 2.LCP(i,k)=min(LCP(i,j),LCP(j,k)) 对于任意1<=i<=j<=k<=n,都满足

证明:

设p=min(LCP(i,j),LCP(j,k))

那么有:LCP(i,j)>=p,LCP(j,k)>=p

那么必然有LCP(i,k)>=p

假如LCP(i,k)=q>p

那么有:s[i+p+1]=s[k+p+1],并且s[i+p+1]!=s[j+p+1],s[j+p+1]!=s[k+p+1]

那么,i,k的字典序显然更接近。不满足:1<=i<=j<=k<=n

所以LCP(i,k)<=p

所以,LCP(i,k)=p

证毕

 

3.LCP(i,k)=min(LCP(j,j-1)) 对于所有的1<i<=j<=k<=n取一个min

有:LCP(i,k)=min(LCP(i,i+1),LCP(i+1,k)),然后不断拆下去,即可证明。

 

三、LCP求法

根据以上性质,

人们设计了:height[i]表示,LCP(suff[sa[i]],suff[sa[i-1]])

设h[i]表示,height[rk[i]]

关键性质:

h[i]>=h[i-1]-1;

证明:

设sa[rk[i-1]-1]=k

h[i-1]=p

那么有:suff[k],suff[i-1]的前p位相同。

1.s[k]!=s[i-1]

即p=0,显然成立。

2.s[k]=s[i-1]

suff[k],suff[i-1]的前p位相同。

所以,suff[k+1],suff[i]的前p-1位相同,

LCP(suff[i],suff[k+1])=p-1

由于k+1对于i位置来讲,suff[k+1]比较平凡,不如suff[sa[rk[i]-1]]和suff[i]相似度更大。

所以,一定有:h[i]>=p-1=h[i-1]-1

证毕

 

有了这个性质,就可以算出height来了。或者h

然后利用LCP(i,k)=min(LCP(j,j-1)) 对于所有的1<i<=j<=k<=n取一个min

这一条就是LCP(i,k)=min(height[j])对于所有的1<i<j<=k<=n取一个min

经常求出height之后,用RMQ等可以直接找到LCP(i,k)

 

代码:

#include<bits/stdc++.h>
#define reg register int
#define il inline
#define numb (ch^'0')
using namespace std;
typedef long long ll;
il void rd(int &x){
    char ch;x=0;bool fl=false;
    while(!isdigit(ch=getchar()))(ch=='-')&&(fl=true);
    for(x=numb;isdigit(ch=getchar());x=x*10+numb);
    (fl==true)&&(x=-x);
}
namespace Miracle{
const int N=1e6+6;
char s[N];
int n,m;
int x[2*N],y[2*N],sa[N];//不开两倍,下面求x的时候会RE(当然求x的时候,对sa[i]+k和n+1取min也可以,总之要是0 
int c[N],num;
int rk[N],hei[N];
void get_SA(){
    m=122;//m表示剩下的排名种类个数 
    for(reg i=1;i<=n;++i) ++c[x[i]=s[i]];//c是桶,x是第一关键字 
    for(reg i=2;i<=m;++i) c[i]+=c[i-1];//c做前缀和,得到第一关键字为x的后缀排名最少是多少。 
    for(reg i=1;i<=n;++i) sa[c[x[i]]--]=i;//处理当前sa 
    for(reg k=1;k<=n;k<<=1){
        num=0;
        for(reg i=n-k+1;i<=n;++i) y[++num]=i;//y[i]表示第二关键字排在 第i位的后缀,起始位置是多少 
        //这里直接把第二关键字是0的放进去,因为它们肯定在最前面 
        for(reg i=1;i<=n;++i) if(sa[i]>k) y[++num]=sa[i]-k;//如果满足sa[i]>k,就有i会作为sa[i]-k的第二关键字。加进去 
        for(reg i=1;i<=m;++i) c[i]=0;//清空桶 
        for(reg i=1;i<=n;++i) ++c[x[i]];//累计第一关键字个数 
        for(reg i=2;i<=m;++i) c[i]+=c[i-1];//前缀和 
        for(reg i=n;i>=1;--i) sa[c[x[y[i]]]--]=y[i],y[i]=0;//处理当前sa,注意倒序这样第一关键字相同的,第二关键字一定升序排列 
        swap(x,y);//复制一下,便于后面处理x[i] 
        num=1;
        x[sa[1]]=1;
        for(reg 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_HEI(){
    for(reg i=1;i<=n;++i) rk[sa[i]]=i;//刷rank 
    int k=0;
    for(reg i=1;i<=n;++i){
        if(rk[i]==1) continue;//h[sa[1]]=0 
        if(k) --k;//h[i]>=h[i-1]-1
        int j=sa[rk[i]-1];
        while(j+k<=n&&i+k<=n&&s[j+k]==s[i+k]) ++k;//暴力匹配,因为每次k最多-1,k最多是n,所以复杂度O(n) 
        hei[rk[i]]=k;
    }
}
int main(){
    scanf("%s",s+1);
    n=strlen(s+1);
    get_SA();
    for(reg i=1;i<=n;++i){
        printf("%d ",sa[i]);
    }return 0;
}

}
int main(){
    Miracle::main();
    return 0;
}

/*
   Author: *Miracle*
   Date: 2018/11/14 21:24:52
*/

 

 

四、例题

P2408 不同子串个数

子串是所有后缀的所有前缀

对于出现两次的子串,必然起始位置的lcp大于等于长度。

某个子串每多出现一次,就会在某个hei[i]中出现。

所以,个数就是:n*(n-1)/2+n-∑hei[i]

每个出现多次的子串,都一定会被减掉多余的次数次。

 

更好的考虑方法是:

每个后缀贡献的前缀:就是n-sa[i]+1-height[i]

总共的前缀个数,减去和上一个排名的相同的前缀(即子串)个数。

剩下的一定都是之前没有出现过的。

所以所有之前出现出现过的都减去了。

当然,这个式子∑一下就是上面的那个。

 

BZOJ 4310: 跳蚤

求出本质不同的子串个数。

二分排名mid

从小到大枚举sa,求出第mid大的本质不同的子串出现位置。

具体来讲,每个位置的sa贡献的不同子串个数是:n-sa[i]+1-height[i]

而且新出现的一定是没有出现的最小的那几个。

所以,如果贡献的子串小于k个,就k-=贡献,继续找。否则就是剩下的第k个位置。

(当然,可以对贡献做一个前缀和,然后二分即可。O(n)->O(logn))

之后贪心砍

从后往前,由于是子串的最大子串,所以,每当得到一个子串比第mid个子串大,必须砍掉。

不能不砍,而如果在之前砍的话,留下的部分字典序一定更大。不优。

所以贪心正确。

比较两个串字典序大小,利用后缀数组直接RMQ求LCP即可。

注意,二分的mid是long long级别。re也要开long long

#include<bits/stdc++.h>
#define reg register int
#define il inline
#define numb (ch^'0')
using namespace std;
typedef long long ll;
il void rd(int &x){
    char ch;x=0;bool fl=false;
    while(!isdigit(ch=getchar()))(ch=='-')&&(fl=true);
    for(x=numb;isdigit(ch=getchar());x=x*10+numb);
    (fl==true)&&(x=-x);
}
namespace Miracle{
const int N=1e5+5;
int n,m,p;
char s[N];
int sa[N],rk[N],hei[N];
int x[2*N],y[2*N];
int c[N],num;
void SA(){
    m=230;
    for(reg i=1;i<=n;++i) ++c[x[i]=s[i]];
    for(reg i=2;i<=m;++i) c[i]+=c[i-1];
    for(reg i=n;i>=1;--i) sa[c[x[i]]--]=i;
    for(reg k=1;k<=n;k<<=1){
        num=0;
        for(reg i=n-k+1;i<=n;++i) y[++num]=i;
        for(reg i=1;i<=n;++i) if(sa[i]>k) y[++num]=sa[i]-k;
        for(reg i=1;i<=m;++i) c[i]=0;
        for(reg i=1;i<=n;++i) ++c[x[i]];
        for(reg i=2;i<=m;++i) c[i]+=c[i-1];
        for(reg i=n;i>=1;--i) sa[c[x[y[i]]]--]=y[i],y[i]=0;
        swap(x,y);
        num=1;x[sa[1]]=1;
        for(reg 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 HEI(){
    for(reg i=1;i<=n;++i) rk[sa[i]]=i;
    int k=0;
    for(reg i=1;i<=n;++i){
        if(rk[i]==1) continue;
        if(k) --k;
        int j=sa[rk[i]-1];
        while(j+k<=n&&i+k<=n&&s[j+k]==s[i+k]) ++k;
        hei[rk[i]]=k;
    }
}
int f[N][20];
int lg[N];
int LCP(int i,int j){
    if(i>j) swap(i,j);
    if(i==j) return n-sa[i]+1;
    int len=j-i;
    int ret=min(f[i+1][lg[len]],f[j-(1<<lg[len])+1][lg[len]]);
    return ret;
}
bool cmp(int s1,int d1,int s2,int d2){//s1<=s2 return 1
    int lcp=LCP(rk[s1],rk[s2]);
    //cout<<" lcp "<<lcp<<endl;
    int l1=d1-s1+1,l2=d2-s2+1;
    if(l1<=l2){
        if(lcp>=l1) return true;
        if(lcp<l1) return s[s1+lcp]<=s[s2+lcp];
    } 
    else{
        if(lcp>=l2) return false;
        if(lcp<l2) return s[s1+lcp]<=s[s2+lcp];
    }
}
ll l,r;
int st,nd;
bool che(ll mid){
    //cout<<" mid ------------------------"<<mid<<endl;
    int tot=0;
       st=1,nd=0;
    ll re=mid;
    while(re>n-sa[st]+1-hei[st]){
        re-=n-sa[st]+1-hei[st];
        ++st;
    }
    nd=sa[st]+hei[st]-1+re;
    st=sa[st];
    //cout<<" st "<<st<<" nd "<<nd<<endl;
    int rr=n;
    for(reg i=n;i>=1;--i){
        if(tot>p) return false;
        while(i>=1&&cmp(i,rr,st,nd)) --i;
        ++tot;
        rr=i;
        //cout<<" tot "<<tot<<" : "<<i<<" "<<rr<<endl;
        ++i;
    }
    return tot<=p;
}
int main(){
    scanf("%d",&p);
    scanf("%s",s+1);
    n=strlen(s+1);
    SA();HEI();
    for(reg i=1;i<=n;++i) {
        lg[i]=(i>>(lg[i-1]+1))?lg[i-1]+1:lg[i-1];
        f[i][0]=hei[i];
        r+=n-sa[i]+1-hei[i];
    }
    for(reg j=1;(1<<j)<=n;++j){
        for(reg i=1;i+(1<<j)-1<=n;++i){
            f[i][j]=min(f[i][j-1],f[i+(1<<(j-1))][j-1]);
        }
    }
    l=1;
    ll ans=0;
//    for(reg i=1;i<=n;++i){
//        cout<<i<<" : "<<rk[i]<<" "<<sa[i]<<" "<<hei[i]<<endl;
//    }
//    cout<<" judge "<<cmp(4,5,1,4)<<endl;
    while(l<=r){
        ll mid=(l+r)>>1;
        if(che(mid)){
            ans=mid;r=mid-1;            
        }
        else l=mid+1;
    }
    ans=che(ans);
    for(reg i=st;i<=nd;++i) putchar(s[i]);
    return 0;
}

}
int main(){
    Miracle::main();
    return 0;
}

/*
   Author: *Miracle*
   Date: 2018/11/15 21:46:33
*/
跳蚤

 

 

P2852 [USACO06DEC]牛奶模式Milk Patterns

对于出现k次的子串,考虑sa序列,枚举最长子串的出现位置中排名最靠后的位置的起始位置,为了找到最长的子串,必然找前面连续k个位置所代表的后缀,

不断求LCP,也就是求LCP(i,i-k+1),就是这个子串长度

单调队列优化即可。

#include<bits/stdc++.h>
#define reg register int
#define il inline
#define numb (ch^'0')
using namespace std;
typedef long long ll;
il void rd(int &x){
    char ch;x=0;bool fl=false;
    while(!isdigit(ch=getchar()))(ch=='-')&&(fl=true);
    for(x=numb;isdigit(ch=getchar());x=x*10+numb);
    (fl==true)&&(x=-x);
}
namespace Miracle{
const int N=1000000+5;
int s[N];
int n,m,p;
int sa[N],rk[N],hei[N],c[N];
int x[2*N],y[2*N];
int num;
void SA(){
    num=0;
    m=1000000;
    for(reg i=1;i<=n;++i) ++c[x[i]=s[i]];
    for(reg i=2;i<=m;++i) c[i]+=c[i-1];
    for(reg i=1;i<=n;++i) sa[c[x[i]]--]=i;
    for(reg k=1;k<=n;k<<=1){
        num=0;
        for(reg i=n-k+1;i<=n;++i) y[++num]=i;
        for(reg i=1;i<=n;++i) if(sa[i]>k) y[++num]=sa[i]-k;
        for(reg i=1;i<=m;++i) c[i]=0;
        for(reg i=1;i<=n;++i) ++c[x[i]];
        for(reg i=2;i<=m;++i) c[i]+=c[i-1];
        for(reg i=n;i>=1;--i) sa[c[x[y[i]]]--]=y[i],y[i]=0;
        swap(x,y);
        num=1;
        x[sa[1]]=1;
        for(reg 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;
    }
//    for(reg i=1;i<=n;++i){
//        cout<<i<<" : "<<sa[i]<<endl;
//    }
}
void HEI(){
    for(reg i=1;i<=n;++i) rk[sa[i]]=i;
    int k=0;
    for(reg i=1;i<=n;++i){
        if(rk[i]==1) continue;
        if(k) --k;
        int j=sa[rk[i]-1];
        while(j+k<=n&&i+k<=n&&s[j+k]==s[i+k]) ++k;
        hei[rk[i]]=k;
    }
}
int q[N],l,r;
int wrk(){
    l=1;
    int ret=0;
    for(reg i=2;i<=p-1;++i){
        while(l<=r&&hei[q[r]]>=hei[i]) --r;
        q[++r]=i;
    }
    //if(l<=r) ret=max(ret,hei[q[l]]);
    for(reg i=p;i<=n;++i){
        while(l<=r&&hei[q[r]]>=hei[i]) --r;
        q[++r]=i;
        while(l<=r&&q[l]<i-p+2)++l;
        if(l<=r) ret=max(ret,hei[q[l]]);
    }
    return ret;
}
int main(){
    rd(n);rd(p);
    for(reg i=1;i<=n;++i) rd(s[i]);
    SA();HEI();
    int ans=wrk();
    printf("%d\n",ans);
    return 0;
}

}
int main(){
    Miracle::main();
    return 0;
}

/*
   Author: *Miracle*
   Date: 2018/11/15 8:52:42
*/
牛奶模式

 

P3181 [HAOI2016]找相同字符

另一篇:[HAOI2016]找相同字符

 

P2178 [NOI2015]品酒大会

求出height,显然要先求出LCP=r的个数以及最大值。然后前缀处理一下。

求LCP=r的个数及最大值:

1.分治。类似于“找相同字符”那个题,注意,mid左边选择一个最值,mid右边选择一个最值。维护两个指针扫过的位置的最值。

最后合并即可。

#include<bits/stdc++.h>
#define reg register int
#define int long long 
#define il inline
#define numb (ch^'0')
using namespace std;
typedef long long ll;
il void rd(int &x){
    char ch;x=0;bool fl=false;
    while(!isdigit(ch=getchar()))(ch=='-')&&(fl=true);
    for(x=numb;isdigit(ch=getchar());x=x*10+numb);
    (fl==true)&&(x=-x);
}
namespace Miracle{
const int N=300000+5;
const int inf=0x3f3f3f3f3f3f3f3f;
int n,m;
char s[N];
ll a[N];
struct node{
    ll ans;
    ll mx;
    void init(){
        ans=0;mx=-inf;
    }
    node operator +(const node &b){
        node c;c.init();
        c.ans=ans+b.ans;
        c.mx=max(mx,b.mx);
        return c;
    }
}b[N];
int sa[N],rk[N],hei[N];
int y[2*N],x[2*N];
int c[N],num;
void SA(){
    m=122;
    for(reg i=1;i<=n;++i) ++c[x[i]=s[i]];
    for(reg i=2;i<=m;++i) c[i]+=c[i-1];
    for(reg i=n;i>=1;--i) sa[c[x[i]]--]=i;
    for(reg k=1;k<=n;k<<=1){
        num=0;
        for(reg i=n-k+1;i<=n;++i) y[++num]=i;
        for(reg i=1;i<=n;++i) if(sa[i]>k) y[++num]=sa[i]-k;
        for(reg i=1;i<=m;++i) c[i]=0;
        for(reg i=1;i<=n;++i) ++c[x[i]];
        for(reg i=2;i<=m;++i) c[i]+=c[i-1];
        for(reg i=n;i>=1;--i) sa[c[x[y[i]]]--]=y[i],y[i]=0;
        swap(x,y);
        num=1;x[sa[1]]=1;
        for(reg 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 HEI(){
    for(reg i=1;i<=n;++i) rk[sa[i]]=i;
    int k=0;
    for(reg i=1;i<=n;++i){
        if(rk[i]==1) continue;
        if(k) --k;
        int j=sa[rk[i]-1];
        while(j+k<=n&&i+k<=n&&s[j+k]==s[i+k]) ++k;
        hei[rk[i]]=k;
    }
}
void divi(int l,int r){
//    cout<<" divi "<<l<<" and "<<r<<" ---------------------------"<<endl;
    if(l==r) return;
    int mid=(l+r)/2;
//    cout<<" mid "<<mid<<endl;
    ll mxl=-inf,mil=inf;
    ll mxr=-inf,mir=inf;
    
    int mi=0x3f3f3f3f;
    int ptr=mid+1;
    for(reg i=mid+1;i<=r;++i){
        mi=min(mi,hei[i]);
        while(ptr>=l+1&&hei[ptr]>=mi){
            //tmp.upda(a[sa[ptr-1]]);
            mil=min(mil,a[sa[ptr-1]]);
            mxl=max(mxl,a[sa[ptr-1]]);
            ptr--;
        }
        b[mi].ans+=(mid-ptr+1);
        if(mil!=inf)b[mi].mx=max(b[mi].mx,mil*a[sa[i]]);
        if(mxl!=-inf)b[mi].mx=max(b[mi].mx,mxl*a[sa[i]]);
        
//        cout<<" ii "<<i<<" ptr "<<ptr<<endl;
//        for(reg i=0;i<=2;++i){
//            cout<<i<<" : "<<b[i].val[0]<<" "<<b[i].val[1]<<" || "<<b[i].val[2]<<" "<<b[i].val[3]<<endl;
//        }
    }
//    cout<<" half "<<endl;
    
    mi=inf;
    ptr=mid;
    for(reg i=mid;i>=l;--i){
        mi=min(mi,hei[i+1]);
        while(ptr+1<=r&&hei[ptr+1]>mi){
            //tmp.upda(a[sa[ptr+1]]);
            mir=min(mir,a[sa[ptr+1]]);
            mxr=max(mxr,a[sa[ptr+1]]);
            ++ptr;
        }
        b[mi].ans+=(ptr-mid);
        if(mir!=inf)b[mi].mx=max(b[mi].mx,mir*a[sa[i]]);
        if(mxr!=-inf)b[mi].mx=max(b[mi].mx,mxr*a[sa[i]]);
    }
//    cout<<" end "<<endl;
//    for(reg i=0;i<=2;++i){
//        cout<<i<<" : "<<b[i].ans<<" "<<b[i].mx<<endl;
//    }
    divi(l,mid);
    divi(mid+1,r);
}
int main(){
    scanf("%lld",&n);
    scanf("%s",s+1);
    for(reg i=1;i<=n;++i)scanf("%lld",&a[i]);
    SA();HEI();
//    for(reg i=1;i<=n;++i){
//        cout<<i<<" "<<sa[i]<<" "<<hei[i]<<endl;
//    }
//    cout<<"--------------------------------"<<endl;
//    
//    cout<<"hei : ";
//    for(reg i=1;i<=n;++i){
//        cout<<hei[i]<<" ";
//    }cout<<endl;
//    cout<<"a[i]: ";
//    for(reg i=1;i<=n;++i){
//        cout<<a[i]<<" ";
//    }cout<<endl;
//    cout<<" sa : ";
//    for(reg i=1;i<=n;++i){
//        cout<<sa[i]<<" ";
//    }cout<<endl;
//    
    for(reg i=0;i<=n;++i)b[i].init();
    divi(1,n);
    for(reg i=n-1;i>=0;--i){
        b[i]=b[i]+b[i+1];
    }
    for(reg i=0;i<=n-1;++i){
        if(b[i].mx==-inf) b[i].mx=0;
        printf("%lld %lld\n",b[i].ans,b[i].mx);
    }
    return 0;
}

}
signed main(){
    Miracle::main();
    return 0;
}

/*
   Author: *Miracle*
   Date: 2018/11/15 15:54:36
*/
品酒大会

2.并查集。按照height排序,每次合并i,i+1两个排名所代表的字符串集合。设lcp长度为k,两边的集合两两组合,lcp一定都为k。

把两棵树的size相乘,加入ans[k][0]。每个并查集维护集合内的a[sa[i]]的最大值和最小值,然后ans[k][1]=max(ans[k][1],mi1*mi2,mx1*mx2)

本质是考虑两个要求LCP的i,j,取得min的位置在哪里。在取得min的位置的那里,就会合并这两个集合并且统计答案。

3.SAM做法

利用Parent树DP,类似“差异”一题。对反串建SAM,在LCA处统计LCP个数以及最大值,最小值。然后按照len从后往前处理一下即可。

 

最后都要前缀处理一下。

 

五、总结

后缀数组可以对一个串求出height,

就可以快速支持查找i,j位置开始的lcp长度。即LCP(rk[i],rk[j]),然后用RMQ即可优化。

和区间问题比较像。但是不是直观地在原始字符串上处理区间,而是对sa序列进行区间处理。

由于是区间,又是找一个最小值,所以可以利用到的方法就比较多了。

分治,并查集,单调栈,单调队列都可以尝试。

并且相邻两个排名LCP长度是最长的,这个性质有时也会用到。

 

posted @ 2018-11-15 15:26  *Miracle*  阅读(273)  评论(3编辑  收藏  举报