【算法与数据结构】—— 后缀数组

后缀数组



—— 摘要 ——


后缀数组是处理字符串的有力工具。后缀数组是后缀树的一个非常精巧的替代品,它比后缀树容易编程实现,能够实现后缀树的很多功能且时间复杂度也并不逊色,而且它比后缀树所占用的内存空间小很多。可以说,在信息学竞赛中后缀数组比后缀树要更为实用。




—— 基本概念 ——


在对后缀数组进行定义之前,需要先了解以下几个基础名词:

  1. 子串:字符串S的子串 str[i , j] ( i ≤ j )表示从字符串S中提取出S[i],S[i+1],……,S[j],并将这些单字符依次序连接形成的字符串。
  2. 后缀:后缀是指从某个位置i开始到整个字符串末尾结束的一个特殊子串。如果我们将某个字符串的后缀用Suffix(i)来表示,那么Suffix(i)=str[i,strlen(S)],其中strlen(S)表示字符串S的长度。(在叙述部分,我们假设字符串的索引从1开始,下同)。
  3. 字符串大小比较:字符串的大小比较是一种基于字典序进行的比较。比如对于两个字符串a,b,我们可以令i从1开始顺次比较a[i]和b[i]的大小,如果a[i]=b[i],则令i加1,继续往后比较;否则,一旦出现a[i]<b[i]则认为a<b,反之则认为a>b,结束比较。如果i>strlen(a)或i>strlen(b)都仍未比出结果,则认为字符串长度更短的那个字符串更小。
    从上面字符串的大小比较来看,如果两个字符串的长度不一致,那么将这两个字符串进行比较的结果一定是不等的,因为不满足a=b的必要条件:strlen(a)=stalen(b)。因此对于某个字符串S而言,S的两个不同开头位置的后缀必不相等。
  4. 后缀数组:后缀数组sa是一个一维数组,它保存的是1-n的某个全排列,并且保证Suffix(sa[i])<Suffix(sa[i+1]),1≤i≤n。也就是将S的n个后缀从小到大进行排序后,将排好序的后缀的开头位置顺次放进sa中。
  5. 名次数组:名次数组rank[i]保存的是Suffix(i)在所有后缀中按从小到大排列的“名次”。

简单说来,后缀数组是“排第几的是谁?”,名次数组是“你排第几?”。容易看出,后缀数组和名次数组为互逆运算。如下图所示:
Alt
设字符串的长度为n。为了方便比较大小,可以在字符串的最后面再添加一个字符,要求这个字符没有在前面的字符串中出现过,而且比前面的字符都要小。
任意两个后缀如果直接比较大小,最多需要比较字符n次,也就是说最迟在比较第n个字符时一定能分出“胜负”。而在求出名次数组后,可以仅用0(1)的时间比较出任意两个后缀的大小。但我们的程序在进行求解时并未直接求解rank数组,而是在对后缀数组sa进行求解。实际上,由于后缀数组和rank数组之间存在着互逆性,所以我们在求出后缀数组或名次数组中的其中一个以后,便可以用0(n)的时间求出另外一个。




—— 具体实现 ——


对于后缀数组的实现,主要有两种算法:

  1. 倍增算法
  2. DC3算法

DC3算法的代码量稍微有点大,在竞赛中遇到的时候几乎不会用这算法(当然,更多的原因是这玩意儿有点难,我确实是搞不懂),因此下面仅对倍增算法进行介绍。
在对倍增算法进行学习前,要求读者要懂得基数排序的基本思想,否则会很难理解这个算法。如果有不了解这个算法的同学,请先移步这篇博客——基数排序(后缀数组基础)进行学习后再来。
倍增算法的主要思路是:用倍增的方法对每个字符开始的长度为2k的子字符串进行排序,求出排名,即rank值。k从0开始,每次加1,当2k大于n以后,每个字符开始的长度为2k的子字符串便相当于所有的后缀。并且这些子字符串一定都已比较出大小,即rank值中没有相同的值,那么此时的rank值就是最后的结果。每一次排序都利用上次长度为2k-1的字符串的rank值,那么长度为2k的字符串就可以用两个长度为2k-1的字符串的排名作为关键字表示,然后进行基数排序,最后便可得出长度为2k的字符串的rank值。以字符串“aabaaaab”为例,整个过程如下图所示。其中x、y是表示长度为2k的字符串的两个关键字。
Alt
其具体的解题过程是:
(1)首先计算S[0],S[1],…,S[n-1]的排名(注意这个单个字符的排序)。比如上面,对于aabaaaab,排序后为:1,1,2,1,1,1,1,2;
(2)计算子串S[0,1],S[1,2],S[2,3],…,S[n-2,n-1],S[n-1,null] 的排名(注意最后一个的第二个字符为空),由于我们知道了单个字符的排名, 那么每个子串可以用一个二元组来表示,比如S[0,1]={1,1},S[1,2]={1,2},S[2,3]={2,1},等等,也就是aa,ab,ba,aa,aa,aa,ab,bε(ε表示空)的排名,排序后为:1,2,4,1,1,1,2,3
(3)计算子串S[0,1,2,3],S[1,2,3,4],S[2,3,4,5],……,S[n-4,n-3,n-2,n-1],S[n-3,n-2,n-1],S[n−2,n−1,],S[n−1,]的排名,方法与上面相同,也是用一个二元组来表示。
接下来的执行流程和(2)(3)类似:每次使用两个2x-1长度的子串来计算2x长度的子串的排名,直到某一次排序后n个数字各不相同,最后就能将rank数组得到。
那在代码中我们怎么通过两个2x-1长度的子串来计算2x长度的子串的排名呢?还是基数排序!我们通过基数排序来对具体关键字所在的索引进行排序,这里面有个关键点是:对关键字所在索引进行排序。仔细想,这样的操作最终不就会得到我们的后缀数组么?确实如此。虽然上面图示的流程给你的感觉是在对rank数组进行求解,但实际上rank数组在倍增法中只是作为一个求后缀数组sa的跳板,我们最终的目的是为了得到后缀数组。
纸上谈兵也许会很空洞,下面我们结合具体的代码进行分析:

const int N=15010;
class SuffixArray{
	private:
		static const int MAX=N;
		int wa[MAX],wb[MAX],wd[MAX],r[MAX];
		bool isSame(int *r,int a,int b,int len)
		{ return r[a]==r[b] && r[a+len]==r[b+len]; }
		void da(int n,int m)
		{
			//第一次桶排序,求出字符长度为1时的sa数组 
			int *x=wa,*y=wb,*t;
			for(int i=0;i<m;i++) wd[i]=0;
			for(int i=0;i<n;i++) wd[x[i]=r[i]]++;
			for(int i=1;i<m;i++) wd[i]+=wd[i-1];
			for(int i=n-1;i>=0;i--) sa[--wd[x[i]]]=i;
			for(int j=1,p=1;p<n;j<<1,m=p){
				//对第二关键字排序
				p=0;
				for(int i=n-j;i<n;i++) y[p++]=i;
				for(int i=0;i<n;i++) if(sa[i]>=j) y[p++]=sa[i]-j;
				//对第一关键字排序
				for(int i=0;i<m;i++) wd[i]=0;
				for(int i=0;i<n;i++) wd[x[i]]++;
				for(int i=1;i<m;i++) wd[i]+=wd[i-1];
				for(int i=n-1;i>=0;i--) sa[--wd[x[y[i]]]]=y[i];
				//利用指针操作调换两数组的内容 
				t=x,x=y,y=t;
				//更新x数组 
				p=1,x[sa[0]]=0;
				for(int i=1;i<n;i++) x[sa[i]]=isSame(y,sa[i-1],sa[i],j)?p-1:p++;
			}
		}
	public:
		int sa[MAX],n;					//sa表示后缀数组,n表示传递进来的字符串的长度 
		void calSuffixArray(char *s)	//计算后缀数组sa 
		{
			n=0;						//初始化n 
			while(*s){					//遍历字符串 
				r[n++]=*s-'a'+1;		//初始化r数组:将字符型数组转化为数值型数组 
				s++;
			} 
			r[n]=0;						//最后在r数组的末尾再添加一个0 
			da(n+1,27);					//开始求sa数组 
		} 
};

我们以该类中的函数调用顺序来进行说明。
首先是public成员中的calSuffixArray函数,该函数作为其类的唯一入口函数,实际上也就是用于求解后缀数组的唯一方法。该函数接受的参数是一个字符串,这个字符串即是需要被计算后缀数组的目标字符串。在calSuffixArray函数里,首先是通过一个while循环来对字符串s进行遍历,遍历的目的有两个:
① 得到字符串长度,并将其存进变量n中;
② 将字符串的每个字符单独提取出,用相应的数字进行代替,并存进数组r中(比如上面代码的意思是用1-26分别表示字母a-z,当然,本程序的处理对象仅仅是小写字符串,如果含大写字母你可以用*s-‘A’+1来替代,如果含特殊字符(如!、{、%、?、$等)可以用*s-‘A’+1来替代。具体的替代方案需根据你的处理对象并结合ASCII码表得出)。
通过上面的步骤,我们便得到了字符串s的等效替代数组r,这个数组在后面的da函数中才会发挥作用。紧接着,我们就直接调用da函数。这个函数有两个参数:
① n表示你传递的等效替代数组r的长度。由于我们为了方便后面的处理,在数组r的末尾添加了一个0,因此数组r的长度实际上是比字符串s多1的。所以我们在调用da函数时传递进的参数为n+1。
② m指示了你在用基数排序时所用到的“桶”的数量(基数排序是基于桶排序的)。由于我给出的示例是处理仅含小写字母的字符串,而小写英文字母又只有26个,因此在calSuffixArray函数中传递进的参数为27。如果你的处理对象含大小写字母,那么你可以将m改为58;如果你的处理对象含大小写字母和特殊字符,那么你可以将m改为94(具体数值可以根据你的处理对象并结合ASCII码表得出)。

接下来到了最关键的部分:da函数
进入函数首先是进行一个基数排序,代码如下:

for(i=0;i<m;i++) wd[i]=0;
for(i=0;i<n;i++) wd[x[i]=r[i]]++;
for(i=1;i<m;i++) wd[i]+=wd[i-1];
for(i=n-1;i>=0;i--) sa[--wd[x[i]]]=i;

其作用是对长度为1的字符串进行排序,比如对字符串”aabaaaab”进行排序,并将排序后的结果放进后缀数组sa中。而上面的x数组保存的值相当于就算rank值。这里的工作可被视为初始化工作。
接下来是若干次的基数排序,目的是对每个长度为2k的字符串所形成的二元组进行排序。这里的长度是以2的倍数形式递增,所以我们可以用一个循环来实现,然后在每次循环中都对当前长度下的二元组进行排序处理。因此在这部分的排序可分为两次:第一次是对第二关键字排序,第二次是对第一关键字排序。对第二关键字排序的结果实际上可由上一次求得的sa直接算出,没有必要再算一次。代码如下:

for(p=0,i=n-j;i<n;i++) y[p++]=i;
for(i=0;i<n;i++)if(sa[i]>=j) y[p++]=sa[i]-j;

其中变量j表示当前字符串的长度,数组y保存的是对第二关键字排序的结果。然后要对第一关键字进行排序,代码如下:

for(int i=0;i<m;i++) wd[i]=0;
for(int i=0;i<n;i++) wd[x[i]]++;
for(int i=1;i<m;i++) wd[i]+=wd[i-1];
for(int i=n-1;i>=0;i--) sa[--wd[x[y[i]]]]=y[i];

这样便求出了新的sa值。在求出sa后,下一步是更新rank值。这里要注意的是,可能有多个字符串的rank值是相同的,所以必须比较两个字符串是否完全相同,y数组的值已经没有必要保存,为了节省空间,这里用y数组保存rank值。这里又有一个小优化,将x和y定义为指针类型,复制整个数组的操作可以用交换指针的值代替,不必将数组中值一个一个的复制。代码如下:

t=x,x=y,y=t;
p=1,x[sa[0]]=0;
for(int i=1;i<n;i++) x[sa[i]]=isSame(y,sa[i-1],sa[i],j)?p-1:p++;

其中isSame函数的代码是:

bool isSame(int *r,int a,int b,int len)
{ return r[a]==r[b] && r[a+len]==r[b+len]; }

这里可以看到规定r[n-1]=0的好处,如果r[a]=r[b],说明以r[a]或r[b],开头的长度为1的字符串肯定不包括字符r[n-1],所以调用变量r[a+1]和r[b+1]不会导致数组下标越界,这样就不需要做特殊判断。执行完上面的代码后,rank值保存在x数组中,而变量p的结果实际上就是不同的字符串的个数。因此,如果p等于n,那么函数可以结束。因为在当前长度的字符串中,已经没有相同的字符串,接下来的排序不会再改变rank数组,从而也就不会再改变后缀数组sa。
下面给出执行calSuffixArray函数的算法流程图:
Alt
下面给出一个求解后缀数组sa的演示程序:

#include<iostream>
using namespace std;

const int N=15010;
class SuffixArray{
	private:
		static const int MAX=N;
		int wa[MAX],wb[MAX],wd[MAX],r[MAX];
		bool isSame(int *r,int a,int b,int len)
		{ return r[a]==r[b] && r[a+len]==r[b+len]; }
		void da(int n,int m)
		{
			//第一次桶排序,求出字符长度为1时的sa数组 
			int *x=wa,*y=wb,*t;
			for(int i=0;i<m;i++) wd[i]=0;
			for(int i=0;i<n;i++) wd[x[i]=r[i]]++;
			for(int i=1;i<m;i++) wd[i]+=wd[i-1];
			for(int i=n-1;i>=0;i--) sa[--wd[x[i]]]=i;
			for(int j=1,p=1;p<n;j<<1,m=p){
				//对第二关键字排序
				p=0;
				for(int i=n-j;i<n;i++) y[p++]=i;
				for(int i=0;i<n;i++) if(sa[i]>=j) y[p++]=sa[i]-j;
				//对第一关键字排序
				for(int i=0;i<m;i++) wd[i]=0;
				for(int i=0;i<n;i++) wd[x[i]]++;
				for(int i=1;i<m;i++) wd[i]+=wd[i-1];
				for(int i=n-1;i>=0;i--) sa[--wd[x[y[i]]]]=y[i];
				//利用指针操作调换两数组的内容 
				t=x,x=y,y=t;
				//更新x数组 
				p=1,x[sa[0]]=0;
				for(int i=1;i<n;i++) x[sa[i]]=isSame(y,sa[i-1],sa[i],j)?p-1:p++;
			}
		}
	public:
		int sa[MAX],n;					//sa表示后缀数组,n表示传递进来的字符串的长度 
		void calSuffixArray(char *s)	//计算后缀数组sa 
		{
			n=0;						//初始化n 
			while(*s){					//遍历字符串 
				r[n++]=*s-'a'+1;		//初始化r数组:将字符型数组转化为数值型数组 
				s++;
			} 
			r[n]=0;						//最后在r数组的末尾再添加一个0 
			da(n+1,27);					//开始求sa数组 
		} 
};

int main()
{
	char chs[N];
	cin>>chs;
	SuffixArray suffixArray;
	suffixArray.calSuffixArray(chs);
	for(int i=1;i<=suffixArray.n;i++) cout<<suffixArray.sa[i]+1<<" ";
	cout<<endl; 
	return 0;
}

该程序的作用是对某个仅含小写字母的字符串,求出其对应的后缀数组。比如输入字符串”lljjllj”,效果如下:
Alt
数据结构的存在是为了给算法奠基,而算法的意义则是为解决一些实际问题。
下面我们来看一下后缀数组能用于处理些什么问题:



—— 应用一:不同子串个数 ——


给定一个仅含小写字母的字符串,问其有多少个不同的子串。(如:对于字符串abad而言,其不同子串有:a、b、d、ab、ba、ad、aba、bad、abad共9个)。
如果给的字符串长度在1000以内,我们是可以用暴力枚举的方法来进行求解的,但是长度范围一旦过大,暴力法就不奏效了。此时,我们的后缀数组就派上用场了。
在解释用后缀数组求不同字串数量之前,我们需要先了解一个东西:公共前缀。
什么是公共前缀?比如对于字符串”abcdef”和”abcxyz”而言,其公共前缀显然是子串”abc”;而对于字符串”abcdef”和”uvwxyz”而言,其公共前缀则为空。
我们定义height[i]=suffix(sa[i-1])和suffix(sa[i])的公共前缀长度。即:height数组为排名相邻的两个后缀的公共前缀长度。比如对于字符串”aabaaaab”而言,我们不难得出其height数组的内容如下:
Alt
那么对于任意的j和k(rank[j]<rank[k]),suffix(j)和suffix(k)的公共前缀长度为:
min( { height[rank[j]+1],height[rank[j]+2], height[rank[j]+3],…,height[rank[k] } )
比如对于字符串”aabaaaab”而言,其后缀”abaaaab”和”aaab”的公共前缀长度为:
min( { height[3],height[4],height[5],height[6] } ) = min( { 2,3,1,2 } ) = 1

到这里有同学可能有点懵。我们不是来求不同子串个数么?现在咋在整这什么height数组哦。
实际上,对于任何一个子串而言,其都能被表达为某个后缀的前缀。那么求解某个字符串的不同字串个数实际上就等价于求所有后缀之间,不同前缀的个数。如果所有的后缀按照suffix(sa[1]),suffix(sa[2]),suffix(sa[3]),… , suffix(sa[n])的顺序计算,不难发现,对于每一次新加进来的后缀suffix(sa[k]), 它将产生n-sa[k]+1个新的前缀。但是其中有height[k]个是和前面的字符串的前缀是相同的。所以suffix(sa[k])将“贡献”出n-sa[k]+1-height[k]个不同的子串。累加后便是原问题的答案。这个做法的时间复杂度为0(n)。

于是,现在我们的任务是先求出height数组。
如果按height[2],height[3], … height [n]的顺序计算,最坏情况下时间复杂度为0(n2)。这样做并没有利用字符串的性质。如果我们定义h[i]=height[rank[i]],也就是suffix(i)和在它前一名的后缀的公共前缀。那么h数组具有性质:h[i]≥h[i-1]。
此时如果我们按照h[1], h[2], … h[n]的顺序计算,并利用h数组的性质,则时间复杂度可降为0(n)。
在具体实现时,其实没有必要保存h数组,只须按照h[1], h[2], … h[n]的顺序计算即可。代码如下:

void calRank()						//计算名次数组rank 
{  for(int i=1;i<=n;i++) rank[sa[i]]=i;  }
void calHeight()					//计算相邻的两个后缀的最长公共前缀长度数组height 
{
	int j,k=0;
	calRank();
	for(int i=0;i<n;i++){
		if(k) k--;
		int j=sa[rank[i]-1];
		while(r[i+k]==r[j+k]) k++;
		height[rank[i]]=k;
	}
}

这样就可以通过调用calHeight()函数将height数组求出。

在得到了height数组后,我们便可根据前面的分析,直接写出求解某个字符串的不同子串的算法如下:

long long calSubstringNum(char *s)	//对于某个字符串str,计算其不同子串的个数 
{
	calSuffixArray(s);
	calHeight();
	long long ans=0;
	for(int i=1;i<=n;i++)
		ans += n - sa[i] -height[i];
	return ans;
}

下面给出整个过程的完整代码:
#include<iostream>
using namespace std;

const int N=15010;
class SuffixArray{
	private:
		static const int MAX=N;
		int wa[MAX],wb[MAX],wd[MAX],r[MAX];
		bool isSame(int *r,int a,int b,int len)
		{ return r[a]==r[b] && r[a+len]==r[b+len]; }
		void da(int n,int m)						//n表示传递进的字符串的长度,m表示最大桶的数量 
		{
			int *x=wa,*y=wb,*t;
			for(int i=0;i<m;i++) wd[i]=0;
			for(int i=0;i<n;i++) wd[x[i]=r[i]]++;
			for(int i=1;i<m;i++) wd[i]+=wd[i-1];
			for(int i=n-1;i>=0;i--) sa[--wd[x[i]]]=i;
			for(int j=1,p=1;p<n;j<<1,m=p){
				//对第二关键字排序
				p=0;
				for(int i=n-j;i<n;i++) y[p++]=i;
				for(int i=0;i<n;i++) if(sa[i]>=j) y[p++]=sa[i]-j;
				//对第一关键字排序
				for(int i=0;i<m;i++) wd[i]=0;
				for(int i=0;i<n;i++) wd[x[i]]++;
				for(int i=1;i<m;i++) wd[i]+=wd[i-1];
				for(int i=n-1;i>=0;i--) sa[--wd[x[y[i]]]]=y[i];
				//利用指针操作调换两数组的内容
				t=x,x=y,y=t;
				//更新x数组 
				p=1,x[sa[0]]=0;
				for(int i=1;i<n;i++) x[sa[i]]=isSame(y,sa[i-1],sa[i],j)?p-1:p++;
			}
		}
	public:
		int sa[MAX],rank[MAX],height[MAX],n;//分别表示后缀数组、名次数组、排名相邻的两个后缀的公共前缀长度数组、传递进来的初始字符串长度 
		void calSuffixArray(char *s)		//计算后缀数组sa 
		{
			n=0;							//初始化n 
			while(*s){						//遍历字符串
				r[n++]=*s-'!'+1;			//初始化r数组:将字符型数组转化为数值型数组
				s++;
			}
			r[n]=0; 						//对于r数组而言,其还需要在最后加一个0方便处理 
			da(n+1,100);					//由于上面多加了个0,因此传递进da函数中的r数组长度需再加1 
		}
		void calRank()						//计算名次数组rank 
		{  
			for(int i=1;i<=n;i++) 
			rank[sa[i]]=i;  
		}
		void calHeight()					//计算相邻的两个后缀的最长公共前缀长度数组height 
		{
			calRank();
			int k=0;
			for(int i=0;i<n;i++){
				if(k) k--;
				int j=sa[rank[i]-1];
				while(r[i+k]==r[j+k]) k++;
				height[rank[i]]=k;
			}
		}
		long long calSubStringNum(char *s)	//对于某个字符串str,计算其不同子串的个数 
		{
			calSuffixArray(s);
			calHeight();
			long long ans=0;
			for(int i=1;i<=n;i++)
				ans += n - sa[i] - height[i];
			return ans;
		}
};

int main()
{
	char chs[N];cin>>chs;
	SuffixArray suffixArray;
	cout<<"字符串:"<<chs<<"的不同字串个数为:"<<suffixArray.calSubstringNum(chs)<<endl;
	for(int i=1;i<=suffixArray.n;i++) cout<<suffixArray.sa[i]<<" ";		//sa数组的索引从1开始,但是值从0开始 
	cout<<endl;
	for(int i=0;i<suffixArray.n;i++) cout<<suffixArray.rank[i]<<" ";	//rank数组的索引从0开始,但是值从1开始 
	cout<<endl;
	for(int i=1;i<=suffixArray.n;i++) cout<<suffixArray.height[i]<<" ";	//height数组的索引从1开始 
	cout<<endl; 
	return 0;
}

趁热打铁:蓝桥杯 算法提高 着急的WYF(不同子串个数)



—— 应用二:最长公共子串 ——


我们先来了解下什么是公共子串:如果字符串L同时出现在字符串A和字符串B中,则称字符串L是字符串A和字符串B的公共子串。
现在如果给出两个字符串A和B,让你输出其最长公共子串(Longest Common Substring ,简称LCS)。你要怎么做?
当然,你是可以用暴力法的,但是那样的时间复杂度达到了0(n3),在字符长度超过1000时就难以胜任了。这时,我们的后缀数组再次出场!

我们知道字符串的任何一个子串都可表达为这个字符串的某个后缀的前缀,因此求A和B的最长公共子串等价于求A的后缀和B的后缀的最长公共前缀的最大值。如果枚举A和B的所有的后缀,那么这样做显然效率低下。由于要计算A的后缀和B的后缀的最长公共前缀,所以可先将第二个字符串写在第一个字符串后面,中间用一个没有出现过的字符隔开,再求这个新的字符串的后缀数组。观察一下,看看能不能从这个新的字符串的后缀数组中找到一些规律。以A=”aaabaa”,B=”abaa”为例,首先我们用一个特殊字符”$”来将这两个字串连接起,得到新的字符串”aaaba$abaa”,然后我们可以得到其height数组的内容如下(其中属于字符串A的后缀用蓝色作为背景、属于字符串B的后缀用黄色作为背景):

Alt
在上表中,容易看出height数组里的最大值为3,并且对于字符串A=”aaabaa”,B=”abaa”而言,其最长公共字串确实是”aba”,长度为3。那是不是所有height值中的最大值就是答案呢?从上表的构建逻辑来看确实是成立的,但是这样的成立建立在一个基础之上:这两个后缀不在同一个字符串中的!所以实际上,只有当suffix(sa[i-1])和suffix(sa[i])不是同一个字符串中的两个后缀时,height[i]才满足要求。而这个条件很容易满足,用一个变量k存放特殊字符”$”所在位置,然后将每次欲更新最大值的那个height[i]所在的后缀用于和k进行一个比较即可。
若记字符串A和字符串B的长度分别为len(A)和len(B)。则求新的字符串的后缀数组和height数组的时间是O(len(A)+ len(B)),然后求排名相邻但原来不在同一个字符串中的两个后缀的height值的最大值,时间也是O(len(A)+ len(B)),所以整个做法的时间复杂度为O(len(A)+ len(B))。时间复杂度已经取到下限,由此看出,这是一个非常优秀的算法。

下面给出“求最长公共子串”的代码:

bool check(int i,int k)							//判断height[i]的两个比较后缀是否不在同一个字符串中(即是否在k的两边)
{
	if( (sa[i-1]-k)*(sa[i]-k)<0 ) return true;	//结果小于0说明异号,异号说明sa[i-1]和sa[i]在k的两边(即要么sa[i-1]<k && sa[i]>k;要么sa[i]<k && sa[i-1]>k) 
	return false;
}

int calLongestCommonSubString(char *a,char *b,char *s)	//计算某两个字符串a,b的最长公共子串s,并返回其长度
{
	int max=0;											//用于保存最长公共子串的长度 
	int pos=strlen(a);									//得到字符串a的长度(之后将用它标记最大height的位置)
	int flagPos=pos;									//记录两个字符中间插入的特殊字符的位置 
	a[pos++]='$';										//在字符串a的末尾插入特殊字符 
	strcpy(a+pos,b);									//将字符串b追加到字符串a的后面
	calSuffixArray(a);									//计算字符串a的后缀数组sa 
	calHeight();										//计算字符串a的rank数组和height数组 
	for(int i=2;i<=n;i++)								//寻找最大的height
		if(height[i]>max && check(i,flagPos)){
			pos=i;										//更新pos 
			max=height[i];								//更新max 
		}
	int x=0,maxPos=max+sa[pos];
	for(int i=sa[pos];i<maxPos;i++)
		s[x++] = a[i];
	return max;
}

最后给出整个过程的完整代码(对于给定的某两个字符串(仅含小写字母),输出这两个字符串的最长公共子串):
#include<iostream>
#include<cstring>
using namespace std;

const int N=15010;
class SuffixArray{
	private:
		static const int MAX=15010;
		int wa[MAX],wb[MAX],wd[MAX],r[MAX];
		bool isSame(int *r,int a,int b,int len)
		{ return r[a]==r[b] && r[a+len]==r[b+len]; }
		bool check(int i,int k)
		{
			if( (sa[i-1]-k)*(sa[i]-k)<0 ) return true;	//结果小于0说明异号,异号说明sa[i-1]和sa[i]在k的两边 
			return false;
		} 
		void da(int n,int m)
		{
			int *x=wa,*y=wb,*t;
			for(int i=0;i<m;i++) wd[i]=0;
			for(int i=0;i<n;i++) wd[x[i]=r[i]]++;
			for(int i=1;i<m;i++) wd[i]+=wd[i-1];
			for(int i=n-1;i>=0;i--) sa[--wd[x[i]]]=i;
			for(int j=1,p=1;p<n;j<<2,m=p){
				p=0;
				for(int i=n-j;i<n;i++) y[p++]=i;
				for(int i=0;i<n;i++) if(sa[i]>=j) y[p++]=sa[i]-j;
				for(int i=0;i<m;i++) wd[i]=0;
				for(int i=0;i<n;i++) wd[x[i]]++;
				for(int i=1;i<m;i++) wd[i]+=wd[i-1];
				for(int i=n-1;i>=0;i--) sa[--wd[x[y[i]]]]=y[i];
				t=x,x=y,y=t;
				p=1,x[sa[0]]=0;
				for(int i=1;i<n;i++) x[sa[i]]=isSame(y,sa[i-1],sa[i],j)?p-1:p++;
			}
		}
	public:
		int sa[MAX],rank[MAX],height[MAX],n;
		void calSuffixArray(char *s)
		{
			n=0;
			while(*s){
				r[n++]=*s-'a'+1;
				s++;
			}
			r[n]=0;
			da(n+1,100);
		}
		void calRank()
		{
			for(int i=1;i<=n;i++)
				rank[sa[i]]=i;
		}
		void calHeight()
		{
			calRank();
			int k=0;
			for(int i=0;i<n;i++){
				if(k) k--;
				int j=sa[rank[i]-1];
				while(r[i+k]==r[j+k]) k++;
				height[rank[i]]=k;
			}
		}
		long long calSubStringNum(char *s)
		{
			calSuffixArray(s);
			calHeight();
			long long ans=0;
			for(int i=1;i<=n;i++)
				ans += n-sa[i]-height[i];
			return ans;
		}
		int calLongestCommonSubString(char *a,char *b,char *s)	//计算某两个字符串a,b的最长公共子串s,并返回其长度 
		{
			int max=0;											//用于保存最长公共子串的长度 
			int pos=strlen(a);									//得到字符串a的长度(之后将用它标记最大height的位置)
			int flagPos=pos;									//记录两个字符中间插入的特殊字符的位置 
			a[pos++]='$';										//在字符串a的末尾插入特殊字符 
			strcpy(a+pos,b);									//将字符串b追加到字符串a的后面
			calSuffixArray(a);									//计算字符串a的后缀数组sa 
			calHeight();										//计算字符串a的rank数组和height数组 
			for(int i=2;i<=n;i++)								//寻找最大的height
				if(height[i]>max && check(i,flagPos)){
					pos=i;										//更新pos 
					max=height[i];								//更新max 
				}
			int x=0,maxPos=max+sa[pos];
			for(int i=sa[pos];i<maxPos;i++)
				s[x++] = a[i];
			return max;
		}
};

int main()
{
	char chs[N],chs_[N],ans[N];
	cin>>chs>>chs_;
	SuffixArray suffixArray;
	cout<<suffixArray.calLongestCommonSubString(chs,chs_,ans)<<endl;
	cout<<ans<<endl;
	return 0;
}


本文参考自:
《国家集训队论文》罗穗骞
《后缀数组:倍增法和DC3的简单理解》Only the Strong Survive 的博客
感谢他们在文章中精彩的推论和分析,使我能在巨人的肩膀上学习。


  • 11
    点赞
  • 23
    收藏
    觉得还不错? 一键收藏
  • 打赏
    打赏
  • 6
    评论

“相关推荐”对你有帮助么?

  • 非常没帮助
  • 没帮助
  • 一般
  • 有帮助
  • 非常有帮助
提交
评论 6
添加红包

请填写红包祝福语或标题

红包个数最小为10个

红包金额最低5元

当前余额3.43前往充值 >
需支付:10.00
成就一亿技术人!
领取后你会自动成为博主和红包主的粉丝 规则
hope_wisdom
发出的红包

打赏作者

theSerein

你的鼓励将是我创作的最大动力

¥1 ¥2 ¥4 ¥6 ¥10 ¥20
扫码支付:¥1
获取中
扫码支付

您的余额不足,请更换扫码支付或充值

打赏作者

实付
使用余额支付
点击重新获取
扫码支付
钱包余额 0

抵扣说明:

1.余额是钱包充值的虚拟货币,按照1:1的比例进行支付金额的抵扣。
2.余额无法直接购买下载,可以购买VIP、付费专栏及课程。

余额充值