【后缀数组】【字符串】【数据结构】后缀数组小结
后缀数组小结
后缀数组是啥？
简而言之，后缀数组的基础内容就是对于给定的一个字符串，将其后缀按字典序排序后所形成的一个编号序列。
例如：
对于字符串ababa:
其后缀数组为SA[]={5,3,1,4,2} ，因为a < aba < ababa <ba <baba
接下来介绍一下如何求suffix array及其衍生数组(rank,height)等。
后缀数组如何求？
最朴素的想法是直接提取所有的后缀并排序，使用快速排序的复杂度为\(O(n\log n * T (cmp))\ =\ O (n^2 log n)\).
考虑如何优化这一过程。
我们观察字符串的比较过程，都是从首位开始进行比较的。
考虑一种比较方式，每次利用每个后缀的第1位先比较，然后确定大致位置，再分组在组内比较第2位，这显然是一种可取的思想，如果使用基数排序是O(n)的，否则是\(O(n \log n)\)的。考虑后缀的特性，发现了什么？
我们每一次比较都在做与之前相同的事，通俗的解释一下就是，假设第i个后缀与第j个后缀的首字母相同，事实上我接下来就是在比较第i+1个后缀与第j+1个后缀的首字母，而这个相对顺序，在上一轮排序中，我应该已经知道了才对。利用这一想法，发现第二轮比较的事物可以结合上一轮进行，然后我们又发现了什么？第二轮比较之后，事实上我获得的是每个后缀的前2个字母的相对顺序，那么由此，我就可以得到前4个，前8个……
这是一种典型的倍增思想，也是十分适合于解决求后缀数组，下面介绍一下倍增法求后缀数组的具体实现。
首先，在上述的讨论中，我们已经发现了，首先预处理出每个后缀首字母的rank，然后开始倍增:第i轮，拥有的信息是所有后缀的前\(2^{i-1}\)首字母的相对顺序，要求的是前\(2^i\)的相对顺序，对于一个后缀j，我们考虑要利用的是其上一轮本身的rank，以及后缀\(j+2^{i-1}\)上一轮的rank，显然是一个双关键字排序的过程，第一关键字为本身的rank，第二关键字为后缀\(j+2^{i-1}\)上一轮的rank，然后去重，直至每个后缀都拥有其独有的rank位置，显然总复杂度为\(O(n \log n T(cmp))\)当使用快速排序时，为\(O(n \log^2 n)\)，尽管这样以及足够优秀，但是仍不足以满足对效率的要求，事实上，考虑每轮rank数不超过n，使用基数排序就可以优化对于排序的效率要求，复杂度转变为\(O(n\log n)\).
由于网络上的代码大多简单难懂，缺乏注释，先来看一个比较直观好理解的代码(by 学弟 Pink_Rabbit)
#include<cstdio>
#include<cstring>
int len, sig=127;
char str[100005];
int rank[100005], rk2[100005];
int buk[100005], buk2[100005], nxt[100005], to[100005], tot;
inline void ins(int x,int y){nxt[++tot]=buk[x];to[tot]=y;buk[x]=tot;}
inline void ins2(int x,int y){nxt[++tot]=buk2[x];to[tot]=y;buk2[x]=tot;}
int SA[100005];
void getSA(){
for(int i=1;i<=len;++i) rank[i]=str[i];
for(int j=1,p;j<=len;j<<=1){
tot=p=0;
memset(buk,-1,sizeof buk);
for(int i=1;i<=len;++i)
rk2[i] = i+j<=len ? rank[i+j] : 0,
ins(rk2[i],i);
memset(buk2,-1,sizeof buk2);
for(int i=sig;i>=0;--i){
for(int k=buk[i];~k;k=nxt[k]){
ins2(rank[to[k]],to[k]);
}
}
for(int i=0;i<=sig;++i){
int lst=-1;
for(int k=buk2[i];~k;k=nxt[k]){
if(rk2[to[k]]==lst) rank[to[k]]=p;
else rank[to[k]]=++p, lst=rk2[to[k]];
}
}
sig=p;
for(int i=1;i<=len;++i) printf("%d ",rank[i]); puts("");
}
for(int i=1;i<=len;++i) SA[rank[i]]=i;
for(int i=1;i<=len;++i) printf("%d ",SA[i]); puts("");
}
int main(){
scanf("%s",str+1);
len=strlen(str+1);
getSA();
return 0;
}
/*使用链表是最直观的书写方式，如果对于网上的绝大多数模板不甚理解，建议这样写，虽然常数大，但这是最直观的写法。*/
接下来给出一个常见的模板
/*
poj 2774 SA
copyright (c) Melacau.
All rights reserved.
*/
/*
task:
Given 2 string(1<=|S|<=10^5)
Ask about the max length of the same substrng of the two string.
solution:
1.Get a new string s = a+" "+b
2.Get height of s
3.Get the max height[i] in case of i belongs to a and i+1 belongs to b or ...
*/
#include <stdio.h>
#include <string.h>
#define R register
#define MN 200005
inline int max(int a,int b){return a>b?a:b;}
int sa[2][MN],rk[2][MN],num[MN],vl[MN],n,g,la,lb,ans,height[MN];char a[MN],b[MN];
inline void getsa(int *sa,int *rk,int *SA,int *RK,R int k){
R int i;for (i=1; i<=n; ++i) num[rk[sa[i]]]=i;//base sort 1st key
for (i=n; i; --i) if (sa[i]>k) SA[num[rk[sa[i]-k]]--]=sa[i]-k;//get sa in [1..n-k] based on 2nd key(loop) and 1st key(num[])
for (i=n-k; ++i<=n; ) SA[num[rk[i]]--]=i;//get sa in (n-k,n]
for (i=1; i<=n; ++i) RK[SA[i]]=RK[SA[i-1]]+(rk[SA[i]]!=rk[SA[i-1]]||rk[SA[i]+k]!=rk[SA[i-1]+k]);//get all new rk
}
void build(){
R int i;for (i=1; i<=n; ++i) ++num[vl[i]];
for (i=1; i<=27; ++i) num[i]+=num[i-1];//get base sort R1
for (i=1; i<=n; ++i) sa[0][num[vl[i]]--]=i;//get sa R1
for (i=1; i<=n; ++i) rk[0][sa[0][i]]=rk[0][sa[0][i-1]]+(vl[sa[0][i]]!=vl[sa[0][i-1]]);//get all rk R1
for (i=1; rk[g][sa[g][n]]<n; i<<=1,g^=1) getsa(sa[g],rk[g],sa[g^1],rk[g^1],i);//multiply get sa
}
void get_height(){
R int i,j,k=0;for (i=1;i<=n;++i)
if (rk[g][i]!=1){
j=sa[g][rk[g][i]-1];//as h[i]>=h[i-1]-1 get pos
while(vl[i+k]==vl[j+k]) ++k;//get the length of height
height[rk[g][i]]=k;if (k) --k;//for h[i]>=h[i-1]-1 , if k > 0 then --k
}
}
int main(){
R int i;
scanf("%s",a+1);la=strlen(a+1);
scanf("%s",b+1);lb=strlen(b+1);
for (i=1; i<=la; ++i) vl[i]=a[i]-'a'+1;
vl[la+1]=27;n=la+lb+1;
for (i=la+2; i<=n; ++i) vl[i]=b[i-la-1]-'a'+1;
build();get_height();
for (R int i=2; i<=n; ++i)
if (sa[g][i]<=la&&sa[g][i-1]>la+1||sa[g][i]>la+1&&sa[g][i-1]<=la)
ans=max(ans,height[i]);
printf("%d\n",ans);
return 0;
}
接下来对上述代码分段解释：
for (i=1; i<=n; ++i) ++num[vl[i]]; for (i=1; i<=27; ++i) num[i]+=num[i-1]; for (i=1; i<=n; ++i) sa[0][num[vl[i]]--]=i; for (i=1; i<=n; ++i) rk[0][sa[0][i]]=rk[0][sa[0][i-1]]+(vl[sa[0][i]]!=vl[sa[0][i-1]]);首先第1行，典型的计数统计；
其次第2行，是为了计算每个字母的排名范围；
接下来第3行是比较难理解的，首先明确一点，这里的sa意味着，排名对应的后缀的起始位置，因此，这句话的实质就是基数排序。简单的解释一下就是，我现在给你个数字，这个数字对应的排名区间你已经知道了，所以这个区间内的某个排名就对应了这个数的位置，现在枚举数，求排名，手算或者观察还是很好理解的。
第4行是算每个位置的排名，不多解释。
接下来看倍增部分，外层循环不说，滚动数组优化，使用传指针的方式简化代码，提高效率只是技巧，观察实现本身：
inline void getsa(int *sa,int *rk,int *SA,int *RK,R int k){ R int i;for (i=1; i<=n; ++i) num[rk[sa[i]]]=i; for (i=n; i; --i) if (sa[i]>k) SA[num[rk[sa[i]-k]]--]=sa[i]-k; for (i=n-k; ++i<=n; ) SA[num[rk[i]]--]=i; for (i=1; i<=n; ++i) RK[SA[i]]=RK[SA[i-1]]+(rk[SA[i]]!=rk[SA[i-1]]||rk[SA[i]+k]!=rk[SA[i-1]+k]); }首先解释第1行，这样基数排序的目的是，对上一次的rank，计算对应的区间(为什么？)这层循环的意思就是找到每个rank对应的终止排名。
接下来看第2行，考虑原本的第i名后缀j，如果它可以作为新一轮的第2关键字的话，由于循环是倒序的，显然相同第一关键字越早出现的被分配到的排名越后，这是符合实现的，因此按照这个顺序，对存在第2关键字的后缀分配排名。
什么，还有没有第2关键字的？显然是存在的，即后缀本身长度已经不大于当前的倍增长度的，我们认为他们的第2关键字为0，这也就是第3行的作用，对这部分后缀分配新排名。
注意：第2、3两行是无法颠倒顺序的，第2行循环的顺序也很重要，原因是这是基于第2关键字的实现。
最后第4行，计算排名，浅显易懂。
至此，我们已经完成了SA的计算与rank的计算。
辅助数组height
接下来，讲一下一个辅助数组，height的作用：
它的定义是,height[i]表示排名为i的后缀与排名为i-1的后缀的最长公共前缀。
显然暴力求的话是\(O(n^2)\)十分睿智。
介绍两个性质:
\(height[i] >= height[rk[sa[i]-1]]-1\)
\(|lcp(suffix[i],suffix[j])|=min\{ height[k] \}(rk[i] \leq k \leq rk[j])\)
性质2比较显然。
下面给出性质1证明：
如果height[rk[sa[i]-1]]=1/0,这是显然成立的。
分析大于1的情况：
首先设sa[rk[sa[i]-1]-1]=k
显然suffix[k]与suffix[sa[i]-1]在去掉首字母之后仍有公共前缀，此时他们分别是suffix[k+1],suffix[sa[i]]，那么这两个后缀的最长公共前缀的长度显然等于height[rk[sa[i]-1]]-1。又根据性质2，则有性质1。
于是，根据性质2，我们就比较容易写出代码：
void get_height(){ R int i,j,k=0;for (i=1;i<=n;++i) if (rk[g][i]!=1){ j=sa[g][rk[g][i]-1]; while(vl[i+k]==vl[j+k]) ++k; height[rk[g][i]]=k;if (k) --k; } }
总结
利用这些数组，我们就可以比较轻松的解决一些字符串问题，说白了，后缀数组及其辅助数组height不过是一个工具，具体如何使用还需大家自行琢磨。本文只是简单的解释了一下后缀数组模板的首先及代码书写的缘由，并不准备介绍例题，给出的模板是poj 2774 查询两个字符串的lcs。可以结合代码自行理解使用技巧。
