后缀数组

后缀数组

         后缀数组就是把一个文本串的所有后缀按字典序从小到大排放的数组。详细介绍见刘汝佳《算法竞赛训练指南》。

         AC自动机可以处理多模板的文本匹配问题，而后缀数组也可以处理多模板的文本匹配问题。那么它们有什么区别呢？

         AC自动机需要事先知道所有的模板，然后对于一个(在线输入的)文本串进行多模板匹配，也就是说模板一定要全部事先知道，需要匹配的文本可以动态的输入。

         后缀数组需要事先知道整个文本串，模板可以一个一个的动态输入。在实际应用中，你很多时候是无法事先知道要查询的模板的(如搜索引擎)。假设你要查找一篇(或多篇)文章里面有没有出现一个词组(模板)，你可以先预处理该文本，计算出它的后缀数组，然后用你输入的这个词组(模板)，对该文本的后缀数组进行二分查找(因为所有后缀已经按字典序排好)，最终通过O(mlogn)(n为文本长度，m为模板长度,后面还会介绍O(m+logn)时间复杂度的算法)的时间复杂度你可以知道该词组(模板)出现的所有位置。(如果此时用KMP去找匹配点的话，复杂度为O(n+m)，在文本串长度n远大于模板串长度m时代价太高)

         下面给出后缀数组的代码：

后缀数组(注释版)

#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;

const int maxn=20000+1000;
struct SuffixArray
{
    //保存原始字符串+‘\0’后形成的字符串
    //即原始字符串在s中表示[0,n-2]范围，
    //然后s[n-1]其实是人为添加的'\0'字符
    char s[maxn];

    //排名(后缀)数组,sa[i]==j表示字典序为i的后缀是后缀j
    //其中i从0到n-1,j从0到n-1范围
    int sa[maxn];

    //名次数组，rank[i]==j表示后缀i的字典序排名为j
    int rank[maxn];
    int height[maxn];

    //辅助数组用于x和y数组
    int t1[maxn],t2[maxn];

    //c[i]==j表示关键字<=i的关键字有j个
    int c[maxn];

    //s原始字符串+‘\0'字符后的长度
    //由于添加了尾0,所以n一般都>=2的
    int n;//n>=2,不能等于1,否则build_height()函数可能会出BUG

    //m大于s[]数组出现的任意字符的int值
    void build_sa(int m)
    {
        int i,*x=t1,*y=t2;

        //预处理每个后缀的长度为1的前缀,求出x数组和sa数组
        //此时x[i]==j表示第i个字符的绝对值(可以看成是名次数组)
        //但有可能x[1]=2,且x[3]=2,说明1字符和3字符完全一样。
        //此时算出的sa[i]==j表示当前长度为1的字符串的排名数组,
        //排名数组值不会一样
        //就算x[1]==x[3]==2,但是sa[1]=1,而sa[2]=3。
        //即就算1号字符和3号字符是完全一样的,
        //但是排名时第1名是1号字符,第2名才是3号字符
        for(i=0;i<m;i++) c[i]=0;
        for(i=0;i<n;i++) c[x[i]=s[i]]++;

        //此时c[i]表示关键字<=i的关键字一共有c[i]个
        for(i=1;i<m;i++) c[i]+=c[i-1];

        //计算当前长度(1)的排名数组
        for(i=n-1;i>=0;i--) sa[--c[x[i]]] = i;


        //每轮循环开始前我们通过之前的计算得到了x[]和sa[]：
        //每个后缀的长为k的前缀(即每个后缀的前[0,k-1]字符)的名次数组x[],
        //我们还知道每个后缀的长为k的前缀的排名数组sa[],
        //然后通过sa[]数组我们可以求得每个后缀的第[k,2*k-1]字符的排名数组y[],
        //然后通过k字符的x[]与k字符的y[],
        //我们可以求得每个后缀的长为2k的前缀字符串的sa[]排名数组
        //然后通过该sa[]排名数组，和k字符的x数组,我们可以求得2k字符的x[]数组
        //以上每轮的x[]名次数组都是可能有重复值大小的,但是sa[]值不会重复
        //比如表示k个字符的x[1]=2,x[4]=2时,
        //那么表示[1,k+1]字符串与[4,k+4]字符完全相同,且排名为2(最高排名为0)
        //当哪轮求出的x[]数组正好由n个值(即所有值都不重复时)
        //说明所有后缀已经排序完毕
        for(int k=1;k<=n;k<<=1)
        {
            //先计算每个后缀的前缀的[k,2*k-1]字符的排名数组y
            //即y是每个后缀的长为2k前缀的第二关键字
            int p=0;

            //y[p]==i表第二关键字为第p名的是后缀i
            //由于当前处理的是每个后缀的前缀的[k,2*k-1]字符
            //而后缀n-k到后缀n-1不存在第k个字符(想想是不是)
            //所以他们的第二关键字的名字自然优先
            for(i=n-k;i<n;i++) y[p++]=i;

            //除了上面那些后缀不存在第二关键字
            //x+k后缀的第1关键字排名-k 等于 x后缀的第2关键字排名
            for(i=0;i<n;i++)if(sa[i]>=k) y[p++]=sa[i]-k;

            //上面已经计算出了y[],(x[]数组上一轮已经算出)
            //下面通过第1关键字x[]名次数组和第2关键字y[]排名数组
            //计算综合后每个后缀的长2k前缀的sa[]数组
            for(i=0;i<m;i++) c[i]=0;
            for(i=0;i<n;i++) c[x[y[i]]]++;
            for(i=1;i<m;i++) c[i]+=c[i-1];
            for(i=n-1;i>=0;i--) sa[--c[x[y[i]]]] = y[i];

            //交换x和y,令y表示名次数组
            //计算综合后每个后缀的长2k前缀的x[]数组
            swap(x,y);
            //此时p用来记录x[]数组中不同值的个数
            p=1;x[sa[0]]=0;
            for(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-1:p++;
                //上面y[sa[i]+k]的sa[i]+k<=n-1,因为只要然后两个不同的后缀必然要分出大小
                //所以在他们y[sa[i]]==y[sa[i-1]],即这两个后缀的长k的第一关键字相同的情况下
                //他们必定还存在第二关键需要比较
            if(p>=n) break;
            m=p;
        }
    }

    //此函数详解见刘汝佳<<训练指南>>P222
    //height[i]表示sa[i-1]后缀与sa[i]后缀的最大公共前缀长度
    //即表示排名i-1和排名i的后缀的最大公共前缀LCP长度
    //所以height数组只有[1,n-1]是有效下标
    void build_height()//n不能等于1,否则出BUG
    {
        int i,j,k=0;
        for(i=0;i<n;i++)rank[sa[i]]=i;
        for(i=0;i<n;i++)
        {
            if(k)k--;
            j=sa[rank[i]-1];
            while(s[i+k]==s[j+k]) k++;
            height[rank[i]]=k;
        }
    }
}sa;

后缀数组RMQ版：

#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
const int maxn=1000000+100;
struct SuffixArray
{
    char s[maxn];
    int sa[maxn],rank[maxn],height[maxn];
    int t1[maxn],t2[maxn],c[maxn],n;
    int dmin[maxn][20];
    void build_sa(int m)
    {
        int i,*x=t1,*y=t2;
        for(i=0;i<m;i++) c[i]=0;
        for(i=0;i<n;i++) c[x[i]=s[i]]++;
        for(i=1;i<m;i++) c[i]+=c[i-1];
        for(i=n-1;i>=0;i--) sa[--c[x[i]]]=i;
        for(int k=1;k<=n;k<<=1)
        {
            int p=0;
            for(i=n-k;i<n;i++) y[p++]=i;
            for(i=0;i<n;i++)if(sa[i]>=k) y[p++]=sa[i]-k;
            for(i=0;i<m;i++) c[i]=0;
            for(i=0;i<n;i++) c[x[y[i]]]++;
            for(i=1;i<m;i++) c[i]+=c[i-1];
            for(i=n-1;i>=0;i--) sa[--c[x[y[i]]]] = y[i];
            swap(x,y);
            p=1,x[sa[0]]=0;
            for(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-1:p++;
            if(p>=n) break;
            m=p;
        }
    }
    void build_height()//n不能等于1,否则出BUG
    {
        int i,j,k=0;
        for(i=0;i<n;i++)rank[sa[i]]=i;
        for(i=0;i<n;i++)
        {
            if(k)k--;
            j=sa[rank[i]-1];
            while(s[i+k]==s[j+k])k++;
            height[rank[i]]=k;
        }
    }
    void initMin()
    {
        for(int i=1;i<=n;i++) dmin[i][0]=height[i];
        for(int j=1;(1<<j)<=n;j++)
            for(int i=1;i+(1<<j)-1<=n;i++)
                dmin[i][j]=min(dmin[i][j-1] , dmin[i+(1<<(j-1))][j-1]);
    }
    int RMQ(int L,int R)//取得范围最小值
    {
        int k=0;
        while((1<<(k+1))<=R-L+1)k++;
        return min(dmin[L][k] , dmin[R-(1<<k)+1][k]);
    }
    int LCP(int i,int j)//求后缀i和j的LCP最长公共前缀
    {
        int L=rank[i],R=rank[j];
        if(L>R) swap(L,R);
        L++;//注意这里
        return RMQ(L,R);
    }
}sa;

        下面介绍利用后缀数组的O(1)时间复杂度的LCP(i,j)操作来在O(m+logn)时间复杂度里面定位一个长m的模板串的算法。本算法的本质还是采用二分法，看mid这个后缀与模板T的相对大小，然后确定下次查询的范围是[L,mid-1]还是[mid+1,R]。不过比较mid与T的大小，我们不用从头开始比较了，我们利用ans保存与T最匹配的那个后缀的名次，用max_match保存ans与T的最大公共前缀长度。然后通过LCP(ans,mid)我们可以将每次比较mid与模板串T的次数大大减少。具体看下面：(引自许智磊《后缀数组》论文)

         多模式串的模式匹配问题

         给定一个固定待匹配串 S，长度为 n，然后每次输入一个模式串P，长度为m，要求返回 P 在 S 中的一个匹配或者返回匹配失败。所谓匹配指某个位置 i满足 1≤i≤n-m+1 使得 S[i..(i+m-1)]=P，也即 Suffix(i)的长m的前缀==P。

         我们知道，如果只有一个模式串，最好的算法就是 KMP 算法，时间复杂度为 O(n+m)，但是如果有多个模式串，我们就要考虑做适当的预处理使得对每个模式串进行匹配所花的时间小一些。最简单的预处理莫过于建立 S 的后缀数组（先在 S 的后面添加'$'），然后每次寻找匹配转化为用二分查找法在 SA 中找到和 P 的公共前缀最长的一个后缀，判断这个最长的公共前缀是否等于 m。

         这样，每次比较 P 和一个后缀的复杂度为O(m)，因为最坏情况下可能比较了 m 个字符。二分查找需要调用比较的次数为 O(logn)，因此总复杂度为O(mlogn)，于是每次匹配的复杂度从 O(n+m)变为 O(mlogn)，可以说改进了不少。

         可是这样仍然不能令我们满足。前面提到 LCP 可以增加后缀数组的威力，

         我们来试试用在这个问题上。

         我们分析原始的二分查找算法，大体有以下几步：

        Step1 令 left=1,right=n,max_match=0。

        Step2 令 mid=(left+right)/2（这里“/”表示取整除法）。

        Step 3 顺次比较 Suffix(SA[mid])和 P 的对应字符，找到两者的最长公共

前缀 r，并判断出它们的大小关系。若r>max_match 则令 max_match=r,ans=mid。

        Step4 若 Suffix(SA[mid])<P 则令 left=mid+1，若 Suffix(SA[mid])>P 则令right=mid-1，若 Suffix(SA[mid])=P 则转至 Step 6。

        Step5 若 left<right 则转至 Step 2，否则至 Step 6。

        Step6 若 max_match=m 则输出 ans，否则输出“无匹配”。

         注意力很快集中在Step 3，如果能够避免每次都从头开始比较Suffix(SA[mid])和 P 的对应字符，也许复杂度就可以进一步降低。

         类似于前面求 height 数组，我们考虑利用以前求得的最长公共前缀作为比较的“基础”，避免冗余的字符比较。

         在比较 Suffix(SA[mid])和 P 之前，我们先用常数时间计算LCP(mid,ans)，然后比较 LCP(mid,ans)和 max_match：

         情况一：LCP(mid,ans)<max_match，则说明 Suffix(SA[mid])和 P 的最长公共前缀就是 LCP(mid,ans)，即直接可以确定 Step 3 中的 r=LCP(mid,ans)，所以可以直接比较两者的第 r+1 个字符（结果一定不会是相等）就可以确定 Suffix(SA[mid])和 P 的大小。这种情况下，字符比较次数为 1 次。

         情况二： LCP(mid,ans)≥max_match， 则说明 Suffix(SA[mid])和 Suffix(SA[ans])的前 max_match 个字符一定是相同的， 于是 Suffix(SA[mid])和 P 的前 max_match个字符也是相同的，于是比较两者的对应字符可以从第 max_match+1 个开始，最后求出的 r 一定大于等于原先的 max_match，字符比较的次数为 rmax_match+1，不难看出 Step 3 执行过后 max_match 将等于 r。

         设每次 Step 3 执行之后 max_match 值增加的量为∆max。在情况一中，∆max=0，字符比较次数为 1=∆max+1；在情况二中，∆max=r-max_match，字符比较次数为 r-max_match+1，也是∆max+1。综上所述，每次 Step 3 进行字符比较的次数为∆max+1。

         总共的字符比较次数为所有的∆max 累加起来再加上 Step 3 执行的次数。所有∆max 累加的结果显然就是最后的 max_match 值，不会超过 len(P)=m，而 Step 3 执行的次数为O(logn)，因此总共的字符比较次数为 O(m+logn)。而整个算法的复杂度显然和字符比较次数同阶，为 O(m+logn)。

         至此，问题得到圆满解决，通过 O(nlogn)的时间进行预处理（构造后缀数组、名词数组，计算 height 数组，RMQ 预处理），之后就可以在 O(m+logn)的
时间内对一个长度为 m 的模式串 P

强烈建议看罗穗骞《后缀数组——处理字符串的有力工具》来理解后缀数组的应用场景。

后缀数组应用

POJ 1743Musical Theme(后缀数组)：找不重叠的两个相同子串。解题报告！

POJ 3261Milk Patterns(后缀数组)：找可重叠至少出现K次的子串。解题报告！

SPOJ 694Distinct Substrings(后缀数组)：找字符串的不同子串个数，字符串长<=1000。解题报告！

SPOJ 705New DistinctSubstrings (后缀数组): 找字符串的不同子串个数，字符串长<=50000。解题报告！

URAL1297 Palindrome(最长回文子串:后缀数组)：寻找最长回文子串。解题报告！

POJ 2406Power Strings(后缀数组)：后缀数组做法超时,直接用KMP做更快。解题报告！

POJ 2774Long Long Message(后缀数组:公共子串)：求两个串的最长公共连续子串的长度。解题报告！

URAL1517. Freedom of choice(后缀数组:最长公共连续子串) ：求两个串的最长公共连续子串的长度。解题报告！

POJ 3294Life Forms(后缀数组)：求n个字符串中超过一半字符串中的最长公共连续字串是什么,如果有多解,按字典序输出。解题报告！

SPOJ 220. Relevant Phrases of Annihilation(后缀数组)：给n个串，求在每个串中至少出现2次且不重叠的子串最大长度。解题报告！