后缀数组是后缀Trie的一个替代品。一个字符串的后缀Trie是把这个字符串所有的后缀给插入到一个Trie中。由于字符串的任意一个子串一定是这个字符串某个后缀的前缀,所以说可以直接在这个Trie里面进行查找就可以找到任意一个字符串是否在这个字符串中,但是最坏情况下这棵Trie的空间复杂度(或者说结点数)可以到达O(N^2)级别,因此需要优化,于是诞生了后缀树,而后缀数组便是后缀树的一个简单替代品。现在在这个字符串的后面加入一个字符$,并规定这个字符小于串中的所有其它字符,想象一下把一棵Trie的每个结点的所有儿子从小到大排序,那么不难发现所有的$都和这棵Trie的一个叶子一一对应。后缀数组(sa)储存的是Trie从左到右所有的的叶子在原串中对应的后缀的首字符的数组下标(规定排序的时候小的结点在左边,大的在右边)。
后缀数组本身就大致介绍完了,先贴出后缀数组的O(NlogN)构造算法:(感谢某写的书中代码很实用但是有时候充满了bug的Lrj)
(关于此算法以及其中用到的各种奇奇妙妙的技巧不做介绍,这只是小结)
1 int sa[maxn],c[maxn],t1[maxn<<1],t2[maxn<<1]; 2 void build_sa(int m) 3 { 4 int i,*x=t1,*y=t2; 5 for(i=0;i<m;i++) c[i]=0; 6 for(i=0;i<n;i++) c[x[i]=S[i]]++; 7 for(i=1;i<m;i++) c[i]+=c[i-1]; 8 for(i=n-1;i>=0;i--) sa[--c[x[i]]]=i; 9 for(int k=1;k<=n;k<<=1){ 10 int p=0; 11 for(i=n-k;i<n;i++) y[p++]=i; 12 for(i=0;i<n;i++) if(sa[i]>=k) y[p++]=sa[i]-k; 13 for(i=0;i<m;i++) c[i]=0; 14 for(i=0;i<n;i++) c[x[i]]++; 15 for(i=1;i<m;i++) c[i]+=c[i-1]; 16 for(i=n-1;i>=0;i--) sa[--c[x[y[i]]]]=y[i]; 17 swap(x,y); 18 p=1,x[sa[0]]=0; 19 for(i=1;i<n;i++) x[sa[i]]=y[sa[i-1]]==y[sa[i]]&&y[sa[i-1]+k]==y[sa[i]+k]?p-1:p++; 20 if(p==n) break; 21 m=p; 22 } 23 }
还有两个非常重要的东西——rank和height数组,rank[i]表示的是以原字符串i位置开头的后缀在sa中的数组下标,height[i]表示的是后缀sa[i-1]和后缀sa[i]的LCP(最长公共前缀),height的构造算法用到了递推来以复杂度O(N)的时间复杂度完成。具体证明请见lrj的大白书(ORZ实际上是我时间不多而且我比较懒)
给出height数组的构造代码:
1 void getHeight() 2 { 3 for(int i=0;i<n;i++) rank[sa[i]]=i; 4 int k=0; 5 for(int i=0;i<n;i++){ 6 if(!rank[i]){ height[0]=k=0; continue; } 7 if(k) k--; 8 int j=sa[rank[i]-1]; 9 while(S[j+k]==S[i+k]) k++; 10 height[rank[i]]=k; 11 } 12 }
于是我们就可以做很多事情辣!!!
Q1:计算任意两个后缀的最长公共前缀。
当你理解height数组实际上是后缀Trie上叶子i和叶子i-1的LCA的深度(规定后缀Trie的深度为0)之后,不难发现当你询问后缀x,y的最长公共前缀的时候你只需要在height中进行一次RMQ(rank[x]+1,rank[y])就可以了(假定rank[x]+1<=rank[y]),道理可以参考dfs序下的用dep数组求LCA的算法,或者用Trie来形象理解。
RMQ是静态的,用st表实现。给出st表的构造和查询以及解决问题的代码:
1 void getst() 2 { 3 for(int i=0;i<n;i++) d[i][0]=height[i]; 4 for(int i=n-1;i>=0;i--) 5 for(int j=1;(1<<j)<=n-i;j++) 6 d[i][j]=min(d[i][j-1],d[i+(1<<j-1)][j-1]); 7 } 8 int RMQ(int x,int y) 9 { 10 int k=0; 11 while((1<<k+1)<y-x+1) k++; 12 return min(d[x][k],d[y-(1<<k)+1][k]); 13 } 14 void work() 15 { 16 scanf("%d",&Q); 17 for(int i=1;i<=Q;i++){ 18 scanf("%d%d",&x,&y); 19 x=rank[x],y=rank[y]; 20 if(x>y) swap(x,y); 21 if(x!=y) printf("%d\n",RMQ(x+1,y));; 22 else printf("%d\n",x+1); 23 } 24 }
Q2:计算最长的可重叠重复子串和最长的不可重叠重复子串。
可重叠的最长重复子串你发现就是取height的最大值。不可重叠的采用二分,每一次二分一个长度mid,相当于所有height[i]<mid的都变成一块隔板,把sa数组变成了几段,求每一段里面sa[i]的最大值MAX和最小值MIN,如果存在一组MAX-MIN>=mid就说明mid是可行解,二分的区间左端点右移,否则二分的区间右端点左移(思考这个问题请注意height的意义)。
后一问代码:
1 void work() 2 { 3 int ans=0,A=0,B=n+1,mid; 4 while(A<B){ 5 mid=A+B>>1; bool ok=0; 6 int MAX=sa[0],MIN=sa[0]; 7 for(int i=1;i<n;i++){ 8 if(height[i]>=mid){ 9 MAX=max(MAX,sa[i]),MIN=min(MIN,sa[i]); 10 if(MAX-MIN>=mid) { ok=1; break; } 11 } 12 else MAX=MIN=sa[i]; 13 } 14 if(ok) A=mid+1,ans=mid; 15 else B=mid; 16 } 17 printf("%d\n",ans); 18 }
Q3:计算不同子串个数(。。。实际上就是看后缀Trie上有多少个结点,代码就不给了。)。
上面的问题全部都是关于某个字符串后缀的,怎么灵活应用呢?答案是:当你要把一堆串进行处理的时候,把它们接起来,每个字符串的后面丢一个分隔符(互不相同,没有在原字符串中出现过),比如下面这三个问题:
Q4:计算这两个字符串的最长公共子串,并输出字典序最小的那个。
对于给出的两个字符串(假设全部都是小写字母),将它们接起来,中间插入一个A,你可以发现任意跨越A位置的两个后缀的LCP正是这两个后缀的最长相同前缀(就是这两个位置开头的最长相同子串)。因为这两个后缀的LCP一定不包含A。可以发现当RMQ查询区间一段区间的时候,左端点不动,右端点右移,结果具有单调递减的性质。所以只需要扫一遍sa判断是不是这个和上一个对应的位置恰好在分隔符的两端,是的话对height取max就可以了。输出的时候再走一遍即可(利用sa的字典序性质)。
1 void work() 2 { 3 int ans=0; 4 for(int i=1;i<n;i++) 5 if(sa[i-1]<nn&&sa[i]>nn||sa[i-1]>nn&&sa[i]<nn) ans=max(ans,height[i]); 6 printf("%d\n",ans); 7 for(int i=1;i<n;i++) if((sa[i-1]<nn&&sa[i]>nn||sa[i-1]>nn&&sa[i]<nn)&&height[i]==ans1){ 8 for(int j=0;j<height[i];j++) tt[j]=S[sa[i]+j]; 9 break; 10 } 11 puts(tt); 12 }
Q5:计算长度不小于K的公共子串的个数(子串相同但位置不同算不同的子串)。
例如:T="xx",P="xx",K=1,则长度不小于K 的公共子串的个数是5。
再如:T="aababaa",P="abaabaa",K=2,则长度不小于K 的公共子串的个数是22。
可以发现这个问题是有单调性做法的。定义两个点(字符串两个下标)的种类不同当且仅当这两个点位于分隔点的两端。单独看两个点x,y的贡献就是LCP(x,y)-K+1。基本的思路是一路上计算每个点和前面的异种类点对答案产生的贡献。令cnt1表示当前遇到的分隔点前的点的数量,cnt2表示当前遇到的分隔点后的点的数量。可以发现当height一路单调递减的时候可以直接一路用当前的height来计算当前点和前面的异种类点对答案产生的共修贡献,然后更新cnt1和cnt2。但是有个问题,如果出现了height从某一段开始递减的情况怎么破?借助前面的思考,发现可以用单调栈来实现,当height[i]大于等于栈顶的时候直接丢到栈里面去,否则按照出栈的顺序依次合并栈顶元素的信息,计算贡献(相当于sa序列反着看原来递增的height数组变成了递减的),当height[i]<K的时候全部出栈,合并信息,然后把这个点丢进去垫底。(加入的分割符ascll码为1,好像终端里面打出来是个笑脸)
1 struct data{ int v,cnt1,cnt2; }stk[maxn]; int top; 2 void work() 3 { 4 long long ans=0; 5 for(int i=1;i<=n;i++){ 6 if(height[i]<K){ 7 while(top>1){ 8 ans+=(1ll*stk[top].cnt1*stk[top-1].cnt2+1ll*stk[top].cnt2*stk[top-1].cnt1)*(stk[top].v-K+1); 9 stk[top-1].cnt1+=stk[top].cnt1,stk[top-1].cnt2+=stk[top].cnt2; 10 top--; 11 } 12 stk[top=1]=(data){0,sa[i]<nn,sa[i]>nn}; 13 } 14 else{ 15 if(stk[top].v<=height[i]) stk[++top]=(data){height[i],sa[i]<nn,sa[i]>nn}; 16 else{ 17 while(stk[top].v>height[i]){ 18 ans+=(1ll*stk[top].cnt1*stk[top-1].cnt2+1ll*stk[top].cnt2*stk[top-1].cnt1)*(stk[top].v-K+1); 19 stk[top-1].cnt1+=stk[top].cnt1,stk[top-1].cnt2+=stk[top].cnt2; 20 top--; 21 } 22 stk[++top]=(data){height[i],sa[i]<nn,sa[i]>nn}; 23 } 24 } 25 } 26 cout<<ans<<‘\n‘; 27 }
Q6:输入N个小写字符串,你的任务是求出一个长度最大的字符串,使得它在超过一半的给出字符串中作为子串出现。如果有多解,按照字典序从小到大输出所有解,N<=100,每个字符串长度<=1000。(来自某大白书:the form of life!)
还是有单调做法。一样的思路,先用相同的分隔符隔开然后求一遍sa和height,原本是要让所有的分隔符互不相同才能保证正确性的,分隔符相同sa数组性质并没有太大实际影响,关键是height数组就有可能不对头了,所以需要手动调整一下height,记录len[i]表示第i个字符串加入之后当前整个字符串的长度(包括分隔符),belong[i]表示某个位置的字符属于哪个字符串(对于分隔符来说belong=0),于是可以通过计算sa[i]和sa[i-1]在原字符串中作为后缀的长度来和height[i]取min调整。之后的大致思路:当前[l,r]的最小值小于等于答案一定不会更优,当满足出现次数的限制的时候应该让l跳到[l,r]最小值最后一次出现的位置,查询两个后缀LCP长度的时候只有区间内最小的height起作用,利用这些不难脑补出一个利用单调队列完成的线性算法(主要利用对height查询的单调性)。
判掉N=1的情况后,注意一下细节就好。(原题多组数据)
1 void data_in() 2 { 3 n=0; 4 memset(S,0,sizeof(S)); 5 for(int i=1;i<=N;i++){ 6 gets(S+n); 7 len[i]=strlen(S),belong[len[i]]=0; 8 for(int j=n;j<strlen(S);j++) belong[j]=i; 9 n=strlen(S)+1,S[n-1]=‘A‘; 10 } 11 } 12 void work() 13 { 14 if(N==1){ 15 S[n-1]=‘\0‘; puts(S); putchar(‘\n‘); 16 return; 17 } 18 build_sa(200); 19 getHeight(); 20 for(int i=1;i<n;i++){ 21 int v1=max(0,len[belong[sa[i]]]-sa[i]),v2=max(0,len[belong[sa[i-1]]]-sa[i-1]); 22 height[i]=min(height[i],min(v1,v2)); 23 } 24 int l=0,r=0,cnt=0,ans=0; 25 memset(vis,0,sizeof(vis)); 26 while(r<n){ 27 if(height[r]<=ans){ 28 while(l<r) vis[belong[sa[l++]]]=0; 29 cnt=0; 30 if(belong[sa[r]]) vis[belong[sa[r]]]++,cnt++; 31 } 32 else{ 33 while(front!=rear&&mq[rear-1].v>height[r]) rear--; 34 if(front!=rear&&mq[rear-1].v==height[r]) mq[rear-1].pos=r; 35 else mq[rear++]=(data){r,height[r]}; 36 if(belong[sa[r]]&&vis[belong[sa[r]]]++==0) cnt++; 37 while(cnt*2>N){ 38 ans=mq[front].v; 39 while(l<mq[front].pos){ 40 if(belong[sa[r]]&&--vis[belong[sa[l]]]==0) cnt--; 41 l++; 42 } 43 front++; 44 } 45 } 46 r++; 47 } 48 if(!ans) { puts("?"); putchar(‘\n‘); return; } 49 for(int i=1;i<n;i++) if(height[i]>=ans){ 50 memset(vis,0,sizeof(vis)); cnt=0; 51 if(belong[sa[i-1]]&&vis[belong[sa[i-1]]]++==0) cnt++; 52 if(vis[belong[sa[i]]]++==0) cnt++; 53 while(i!=n-1&&height[i+1]>=ans) 54 if(vis[belong[sa[++i]]]++==0) cnt++; 55 if(cnt*2<=N) continue; 56 for(int j=sa[i];j<sa[i]+ans;j++) putchar(S[j]); 57 putchar(‘\n‘); 58 } 59 putchar(‘\n‘); 60 }