后缀数组

后缀数组

后缀数组定义： 让人懵逼的有力工具 ，指某一字符串后缀按照字典序的一个排列。sa[i] = j 的含义为所有后缀按照字典序排列，排在第 i 个的是后缀 j 。

1.各变量的含义

此算法涉及到很多坨数组及变量，故在这里做一个罗列：

sa[ i ]：表示排名为 i 的后缀的起始位置下标；

rak[ i ]：本应写作rank数组，但因与编译器关键字重名，这里写作rak[ i ] ; 表示起始位置为 i 的后缀的排名；

tp[ i ]：基数排序的第二关键字，表示第二关键字排名为 i 的后缀的起始下标；

tax[ i ]：i 号元素出现了多少次，用于辅助基数排序；

s ：字符串，s[i]表示字符串中第 i 个字符串；

lcp(x,y)：字符串x与字符串y的最长公共前缀，这里指排名为x与排名为y的最长公共前缀；

height[ i ]：lcp(sa[i],sa[i−1])，即排名为 i 的后缀的字符串与排名为 i−1 的后缀的字符串的最长公共前缀；

H[ i ]：height[ rak[ i ] ]，即 i 号前缀与前一名（不一定是i-1号）的最长公共前缀。

在说一下 ran 数组和 ra 数组便于深刻理解：他们有下面的关系等式：ran[sa[i]]=i和sa[ran[i]]=i，仔细揣摩一下就能理解其含义~。

2.具体思想

求 sa 数组其实关键在于排序，如果直接sort快排的话，每一次的时间复杂度是 O(nlogn) ，再加上比较的时间复杂度，其总时间复杂度会达到 O(n*n*logn) 。

我们考虑到是按字典序排序，所以我们用基数排序，每次对单字符排序，之后用单字符对多字符排序，以此类推，直至所有后缀顺序都不相同为止。这里在基数排序的过程中，不是一个一个的字符相加，而是用倍增的思想，因为如果一个字符一个字符相加，会出现某些字符重复排序的情况，所以直接倍增即可。考虑到裸基数排序的时间复杂度是 O(n)，这里在基数排序基础上加个二分（倍增），其时间复杂度进一步降低为O(logn)，在加上比较的时间复杂度为O(n)，所以总的时间复杂度为O(nlogn)，相当可观。

下面是倍增基数排序的图解：

3.代码实现及讲解

1. 基数排序

const int maxn=1e6+10;
char s[maxn];
int rak[maxn],sa[maxn],tax[maxn],tp[maxn];
void Qsort() 
{
    //M为桶的个数，及字符集的个数，len为字符串的长度+1。
    for (int i = 0; i <= M; i++) tax[i] = 0;                //把桶清零
    for (int i = 1; i <= len; i++) tax[rak[i]]++;           //统计每个名词出现的个数
    for (int i = 1; i <= M; i++) tax[i] += tax[i - 1];      //做前缀和
    			//可以快速定位每个位置应有的排名，可以统计比当前名次小的后缀有多少个
    for (int i = len; i >= 1; i--) sa[ tax[rak[tp[i]]]-- ] = tp[i];		//&@#$%......
    			
    //i从大到小依次枚举，那么sa[tax[rak[tp[i]]]−−]的意思就是说：用rak[i]定位到第一关键字的大小；那么tax[rak[tp[i]]]就表示当第一关键字相同时，第二关键字较大的这个后缀的排名是啥到了排名，我们也就能更新sa了，--表示减去自身，其他应有多少个排名。
}

2.倍增

void SuffixSort() 
{
    M = 1010;                                 //字符集的大小，一共需要多少个桶
    for (int i=1;i<=len;i++)
        rak[i] = s[i]-'0'+1,tp[i]=i;		//初始化rak和tp,注意rak的加1
    Qsort();
    for (int w=1;w<=len;w<<=1)
    {
        //w:当前倍增的长度，w = x表示已经求出了长度为x的后缀的排名，现在要更新长度为2x的后缀的排名
        //p表示不同的后缀的个数，很显然原字符串的后缀都是不同的，因此p = N时可以退出循环
        int p = 0;								//这里的p仅仅是一个计数器
        for (int i=len-w+1;i<=len;i++) 
            	tp[++p]=i;
        for (int i=1;i<=len;i++){ 
            	if (sa[i]>w) 
                    tp[++p] = sa[i]-w; 		//这两个for是后缀数组的核心部分，是对第二关键字排序
        }
        
/*假设我们现在需要得到的长度为w，那么sa[i]表示的实际是长度为w/2的后缀中排名为i的位置(也就是上一轮的结果)
我们需要得到的tp[i]表示的是：长度为w的后缀中，第二关键字排名为i的位置。
之所以能这样更新，是因为i号后缀的前w/2个字符形成的字符串是i−w/2号后缀的后w/2个字符形成的字符串*/
            
        Qsort();					//此时我们已经更新出了第二关键字，利用上一轮的rak更新本轮的sa
        memcpy(tp,rak,sizeof(rak));	//这里原本tp已经没有用了
        rak[sa[1]]=p=1;
        for (int i=2;i<=len;i++)
            rak[sa[i]]= (tp[sa[i-1]]==tp[sa[i]]&&tp[sa[i-1]+w]==tp[sa[i]+w])?p:++p;
        //这里当两个后缀上一轮排名相同时本轮也相同，至于为什么大家可以思考一下
        if(p==len)	break;
        M=p;
    }
}

4.精髓：height数组

height数组在变量含义中也提过：它表示lcp( sa[i],sa[i-1] )；

H数组也同上说过：他表示height[ rak[ i ] ]，即从i开始的后缀与排名前一名的后缀的最长公共前缀。

性质：$H[i]>=H[i-1]-1$

证明？记了也不懂系列，干脆不计了，如果确实想了解的话，可以参考博客：https://www.cnblogs.com/zwfymqz/p/8413523.html#_label4

求height数组代码

之所以能够线性求出height数组，还要依靠上面那条重要性质，show code:

void GetHeight() 
{
    int j,k=0;
    for(int i=1;i<=len;i++) 
    {
        if(k) 	k--;
        int j=sa[rak[i]-1];
        while(s[i+k]==s[j+k]) k++;
        Height[rak[i]]=k;
        //printf("%d\n", k);
    }
}

5.模板（无注释版）

上面看懂了吗？没看懂？没关系，会用就行了（虽然蒟蒻我也没看懂，待我再研究研究），下面上无注释版模板：

基排：

void Qsort() 
{
    for (int i=0;i<=M;i++) 		tax[i] = 0;   
    for (int i=1;i<=len;i++)	tax[rak[i]]++;           
    for (int i=1;i<=M;i++) 		tax[i]+=tax[i - 1];    
    for (int i=len;i>= 1;i--) 	sa[ tax[rak[tp[i]]]-- ]=tp[i];
}

倍增：

void suffixsort()
{
    M=1010;       
    for(int i=1;i<=len;++i){
        rak[i]=s[i];
        tp[i]=i;
    }
    Qsort();
    for(int w=1;w<=len;w<<=1)
    {
        int p=0;
        for(int i=len-w+1;i<=len;++i)
            tp[++p] = i;
        for(int i=1;i<=len;++i)
            if(sa[i]>w)
                tp[++p]=sa[i]-w;
        Qsort();
        memcpy(tp,rak,sizeof(rak));
        rak[sa[1]]=p=1;
        for(int i=2;i<=len;++i)
            rak[sa[i]]=( (tp[sa[i-1]]==tp[sa[i]]) && (tp[sa[i-1]+w]==tp[sa[i]+w]) )?p:++p;
        if(p==len)  break;
        M=p;
    }
}

线性求height数组：

void GetHeight() 
{
    int j,k=0;
    for(int i=1;i<=len;i++) 
    {
        if(k) 	k--;
        int j=sa[rak[i]-1];
        while(s[i+k]==s[j+k]) k++;
        Height[rak[i]]=k;
        //printf("%d\n", k);
    }
}

AC代码（求LCP,连接取max height值即可）：

#include<iostream>
#include<cstdlib>
#include<iomanip>
#include<algorithm>
#include<cstring>

using namespace std;
const int maxn=1e6+10;
char s[maxn],str[maxn];
int rak[maxn],tax[maxn],tp[maxn],sa[maxn],Height[maxn];
int len1,len,M;
void Qsort()
{
    for(int i=0;i<=M;++i)  tax[i]=0;
    for(int i=1;i<=len;++i)     tax[rak[i]]++;
    for(int i=1;i<=M;++i)   tax[i]+=tax[i-1];
    for(int i=len;i>=1;i--) sa[ tax[rak[tp[i]]]-- ]=tp[i];
}
void suffixsort()
{
    M=1010;       //桶的个数
    for(int i=1;i<=len;++i){
        rak[i]=s[i];
        tp[i]=i;
    }
    Qsort();
    for(int w=1;w<=len;w<<=1)
    {
        int p=0;
        for(int i=len-w+1;i<=len;++i)
            tp[++p] = i;
        for(int i=1;i<=len;++i)
            if(sa[i]>w)
                tp[++p]=sa[i]-w;
        Qsort();
        memcpy(tp,rak,sizeof(rak));
        rak[sa[1]]=p=1;
        for(int i=2;i<=len;++i)
            rak[sa[i]]=( (tp[sa[i-1]]==tp[sa[i]]) && (tp[sa[i-1]+w]==tp[sa[i]+w]) )?p:++p;
        if(p==len)  break;
        M=p;
    }
}
void getheight()
{
    int k=0;
    for(int i=1;i<=len;++i)
    {
        if(k)   k--;
        int j=sa[rak[i]-1];
        while(s[i+k]==s[j+k])   k++;
        Height[rak[i]]=k;
    }
}
//yeshowmuchiloveyoumydearmotherreallyicannotbelieveit#yeaphowmuchiloveyoumydearmother

int main()
{
    ios::sync_with_stdio(false);

    cin>>str+1;
    len1=strlen(str+1);
    str[len1+1]='#';
    cin>>str+len1+2;
    len=strlen(str+1);
    for(int i=1;i<=len;++i)
        s[i]=str[i];
    suffixsort();
    getheight();
    int res=-1;
    for(int i=1;i<=len;++i)
    {   
        if((sa[i]<=len1&&sa[i-1]>len1+1)||(sa[i]>len1+1&&sa[i-1]<=len1))
            res=max(res,Height[i]);
    }
    cout<<res<<endl;

    system("pause");
    return 0;
}

poj 1743:求不可重叠最长公共子串（后缀数组+二分）：

#include<iostream>
#include<cstdlib>
#include<algorithm>
#include<cstring>
#include<stdio.h>

using namespace std;
const int maxn=2e4+10;
int rak[maxn],sa[maxn],tax[maxn],tp[maxn],Height[maxn];
int a[maxn];           //这里其实每个数组元素就相当于一个字符，连起来就相当于一个字符串，不要死板的以为只能输char数组,有很多转换变形 
int M,len;
void Qsort()
{
    for(int i=0;i<=M;++i)   tax[i]=0;
    for(int i=1;i<=len;++i)    tax[rak[i]]++;
    for(int i=1;i<=M;++i)   tax[i]+=tax[i-1];
    for(int i=len;i>=1;--i)     sa[tax[rak[tp[i]]]--] = tp[i];
}
void suffixsort()
{
    M=210;
    for(int i=1;i<=len;++i)
    {
        rak[i]=a[i];
        tp[i]=i;
    }
    Qsort();
    for(int w=1;w<=len;w<<=1)
    {
        int p=0;
        for(int i=len-w+1;i<=len;++i)
            tp[++p]=i;
        for(int i=1;i<=len;++i)
            if(sa[i]>w)
                tp[++p]=sa[i]-w;
        Qsort();
        swap(rak,tp);
        rak[sa[1]]=p=1;
        for(int i=2;i<=len;++i){
            rak[sa[i]]=((tp[sa[i-1]]==tp[sa[i]])&&(tp[sa[i-1]+w]==tp[sa[i]+w]))?p:++p;
            if(p==len)  break;
            M=p;
        }
    }
}
void getHeight()
{
    int j,k=0;
    for(int i=1;i<=len;++i){
        if(k)   k--;
        int j=sa[rak[i]-1];
        while(a[i+k]==a[j+k])   k++;
        Height[rak[i]] = k;
    }
}
int check(int x)
{
    int mx=sa[1],mi=sa[1];              //mx为sa最大值，mi为sa最小值
    for(int i=2;i<=len;++i){
        if(Height[i]<x) mx=mi=sa[i];
        else{
            if(sa[i]<mi)    mi=sa[i];
            if(sa[i]>mx)    mx=sa[i];
            if(mx-mi>x)     return 1;   //sa最大值与最小值差x个（不能等于x，否则首尾有一个重叠）说明不重叠，满足条件
        }
    }
    return 0;
}

int main()
{
    //ios::sync_with_stdio(false);

    while(scanf("%d",&len)!=EOF)
    {
        if(len==0)      break;
        for(int i=1;i<=len;++i)
            scanf("%d",&a[i]);
        for(int i=1;i<len;++i)
            a[i]=a[i+1]-a[i]+90;        //防止出现负数
        len--;                          //所有元素差分后长度减1
        suffixsort();
        getHeight();
        //下面为二分代码
        int res=0;
        int l=1,r=len,mid;
        while(l<r)
        {
            mid=(l+r)>>1;
            if(check(mid)){             //满足则向右更新l值，求最大res
                res=mid;
                l=mid+1;
            }
            else{
                r=mid;
            }
        }
        if(res<4)   printf("0\n");
            else{
                printf("%d\n",res+1);           //因为差分过，所以只需判断是否长度为4即可，最后满足的答案也要加1
            }
    }

    //system("pause");
    return 0;
}

6.应用

1：给定一个字符串，求它们的两个后缀的最长公共前缀

​ 解：假设它们是某一个字符串的后缀，其排名分别为 j 和 k ，那么最长公共前缀就是min(Height[rak[j]+1],Height[rak[j]+1],......Height[rak[k]])

2：最长可重复子串

​ 解：求height数组，取其中最大值即可。因为任意两个后缀的最长公共前缀都是Height数组里面某一段的最小值，这个值一定不大于height数组里面的最大值。

3：最长不可重叠重复子串

​ 解：先二分答案，将题目变为判定性问题：判断是否存在两个长度为 k 的子串是相同且不重叠的。解决这个问题还得用 height 数组。把排序后的后缀分成若干组，其中每组的后缀之间的 height 值都不小于 k 。然后判断每组后缀，其最大 sa 值与最小 sa 值之差是否大于等于k。如有一组满足，则存在长度为 K 的子串重复且不重叠。

4:可重叠的最少出现k次的最长重复子串

​ 解：和上例类似，先二分答案，然后将后缀数组分成若干组。不同的是，这里要判断的是有没有一个组的后缀个数大于等于k。如果有，那么存在k个相同的子串满足条件，否则不存在。

5：不相同的子串的个数

​ 解：即求所有后缀之间不相同的前缀的个数。我们如果按 suffix(sa[1])、suffix(sa[2])…suffix(sa[n])的顺序计算，可以发现，每次新加进来的后缀 suffix(sa[k])，都将产生 n-sa[k]+1 个新的后缀，而其中 height[k] 个适合前面字符串相同的。所以suffix([sa[k]])的真是贡献是n-sa[k]+1-height[k]个不同的子串，累加一下即可。

6：最长回文子串

​ 解：将一个不可能出现的符号加在这个字符串后面，再讲该字符串反着来加到这个字符串后面。这样问题就变成了求这个新的字符串的最长公共前缀。这样做的时间复杂度为 O(nlogn)，其中用 RMQ 算法做预处理可将时间复杂度变为O(n)。（其实也可以用 Manacher 来求。）

7：连续重复子串

​ 解：问题为一个字符串是由某个字符串s重复R次得到的，求R的最大值。我们穷举k即可，判断条件 k能整除字符串长度L整除以及suffix(1)与suffix(k+1)的最长公共前缀是否等于n-k。（用kmp也可以求）

8：重复次数最多的连续重复子串

​ 解：先穷举长度 L，然后求长度为 L 的子串最多能连续出现几次。首先连续出现1 次是肯定可以的，所以这里只考虑至少 2 次的情况。假设在原字符串中连续出现 2 次，记这个子字符串为 S，那么 S 肯定包括了字符 r[0], r[L], r[L*2],
r[L*3], ……中的某相邻的两个。所以只须看字符 r[L*i]和 r[L*(i+1)]往前和往后各能匹配到多远，记这个总长度为 K，那么这里连续出现了 K/L+1 次。最后
看最大值是多少。时间复杂度为O(nlogn)。

9：两个字符串的最长公共子串

​ 用一个不可能出现的字符将两个字符串接起来，求一下height数组最大且不是同一个字符串中的即可(判断sa下标和两个字符串之间长度即可)。O(lena+lenb)。# 后缀数组