后缀数组解决经典问题

 

1.求最长回文子串

思路：穷举每一位，然后计算以这个字符为中心的最长回文子串。这里要分两种情况，1.回文子串的长度为奇数，2.长度为偶数。两种情况都可以转化为求一个后缀和一个反过来写的后缀的最长公共前缀。具体是将整个字符串反过来写在原字符串后面，中间用一个特殊的字符隔开。这样就成了求这个新的字符串的某两个后缀的最长公共前缀。然后在查找最长公共前缀lcp的时候利用到了RMQ，我们知道对于两个后缀j和k，设rank[j]<rank[k]，则易知后缀j，k的LCP长度等于height[rank[j]+1]，height[rank[j]+2]，...,height[rank[k]]中的最小值 ，就是RMQ(height, rank[j]+1, rank[k])，这是height[]数组的两两比较得到的
 

#include <algorithm>
#include <iostream>
#include <cstring>
#include <cstdio>
#include <queue>
using namespace std;
const int maxn = 2010;
int sa[maxn]; 
int t1[maxn], t2[maxn], c[maxn];
int rank[maxn], height[maxn];
int s[maxn], r[maxn];
char str[maxn];
int st[maxn][20];
void build_sa(int s[], int n, int m) {
    int i, j, p, *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 (j = 1; j <= n; j <<= 1) {
        p = 0;
        for (i = n-j; i < n; i++) y[p++] = i;
        for (i = 0; i < n; i++) 
            if (sa[i] >= j) 
                y[p++] = sa[i] - j;
        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-1]] == y[sa[i]] && y[sa[i-1]+j] == y[sa[i]+j] ? p-1 : p++;
 
        if (p >= n) break;
        m = p;
    }
}
void getHeight(int s[],int n) {
    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 ST(int n) {
    for (int i = 1; i <= n; i++)
        st[i][0] = i;
 
    for (int j = 1; (1<<j) <= n; j++) {
        for (int i = 1; i + (1<<j) <= n; i++) {
            int p = st[i][j-1];    
            int q = st[i+(1<<(j-1))][j-1];
            st[i][j] = height[p] > height[q] ? q : p;
        }
    }
}
 
int RMQ(int i, int j) {
    int k = 0;
    if (i > j)
        swap(i, j);
    i++;
 
    while ((1<<(k+1)) <= j-i+1) k++;
    i = st[i][k];
    j = st[j-(1<<k)+1][k];
 
    return min(height[i], height[j]);
}
 
int main() {
    while (scanf("%s", str) != EOF) {
        int len = strlen(str);
        int n = 2 * len + 1;
        for (int i = 0; i < len; i++)
            r[i] = str[i];
        r[len] = 1;
        for (int i = 0; i < len; i++)
            r[i+len+1] = str[len-1-i];
 
        r[n] = 0; //notice
        build_sa(r, n+1, 128);
        getHeight(r, n);
 
        ST(n);
        int mid = n >> 1;
        int ans = 0, cur = 0;
        for (int i = 0; i < mid; i++) {
            int j = RMQ(rank[i], rank[n-i-1]);  //奇对称
            if ((j<<1) - 1 > ans) { 
                ans = (j<<1) - 1;
                cur = i - j + 1;
            }
            if (i) {
                j = RMQ(rank[i], rank[n-i]); //偶对称
                if ((j << 1) > ans) {
                    ans = j << 1;
                    cur = i - j;
                }
            }
        }
 
        for (int i = cur; i < cur + ans; i++)
            printf("%c", r[i]);
        printf("\n");
    }
    return 0;
}

2.最长公共子串

先用倍增算法求出字符串S的后缀数组和Height数组。然后二分枚举答案A，假设N个串可以有长度为A的公共字串，并对A的可行性进行验证。如果验证A可行，A’(A’<A)也一定可行，尝试增大A，反之尝试缩小A。最终可以取得A的最大可行值，就是这N个串的最长公共子串的长度。可以证明，尝试次数是O(logL)的。于是问题就集中到了如何验证给定的长度A是否为可行解。方法是，找出在Height数组中找出连续的一段Height[i..j]，使得i<=k<=j均满足Height[k]>=A，并且i-1<=k<=j中，SA[k]分属于原有N个串S1..SN。如果能找到这样的一段，那么A就是可行解，否则A不是可行解。

具体查找i..j时，可以先从前到后枚举i的位置，如果发现Height[i]>=A，则开始从i向后枚举j的位置，直到找到了Height[j+1]<A，判断[i..j]这个区间内SA是否分属于S1..SN。如果满足，则A为可行解，然后直接返回，否则令i=j+1继续向后枚举。S中每个字符被访问了O(1)次，S的长度为N*L+N-1，所以时间复杂度为O(N*L)。

#include <iostream>
#include <cstdio>
#include <cstdlib>
#include <cmath>
#include <cstring>
const int MAXL=10011,MAXN=6;
using namespace std;
struct SuffixArray
{
	struct RadixElement
	{
		int id,k[2];
	}RE[MAXL],RT[MAXL];
	int N,A[MAXL],SA[MAXL],Rank[MAXL],Height[MAXL],C[MAXL];
	void RadixSort()
	{
		int i,y;
		for (y=1;y>=0;y--)
		{
			memset(C,0,sizeof(C));
			for (i=1;i<=N;i++) C[RE[i].k[y]]++;
			for (i=1;i<MAXL;i++) C[i]+=C[i-1];
			for (i=N;i>=1;i--) RT[C[RE[i].k[y]]--]=RE[i];
			for (i=1;i<=N;i++) RE[i]=RT[i];
		}
		for (i=1;i<=N;i++)
		{
			Rank[ RE[i].id ]=Rank[ RE[i-1].id ];
			if (RE[i].k[0]!=RE[i-1].k[0] || RE[i].k[1]!=RE[i-1].k[1])
				Rank[ RE[i].id ]++;
		}
	}
	void CalcSA()
	{
		int i,k;
		RE[0].k[0]=-1;
		for (i=1;i<=N;i++)
			RE[i].id=i,RE[i].k[0]=A[i],RE[i].k[1]=0;
		RadixSort();
		for (k=1;k+1<=N;k*=2)
		{
			for (i=1;i<=N;i++)
				RE[i].id=i,RE[i].k[0]=Rank[i],RE[i].k[1]=i+k<=N?Rank[i+k]:0;
			RadixSort();
		}
		for (i=1;i<=N;i++)
			SA[ Rank[i] ]=i;
	}
	void CalcHeight()
	{
		int i,k,h=0;
		for (i=1;i<=N;i++)
		{
			if (Rank[i]==1)
				h=0;
			else
			{
				k=SA[Rank[i]-1];
				if (--h<0) h=0;
				for (;A[i+h]==A[k+h];h++);
			}
			Height[Rank[i]]=h;
		}
	}
}SA;
int N,Ans,Bel[MAXL];
char S[MAXL];
void init()
{
	int i;
	freopen("pow.in","r",stdin);
	freopen("pow.out","w",stdout);
	scanf("%d",&N);
	SA.N=0;
	for (i=1;i<=N;i++)
	{
		scanf("%s",S);
		for (char *p=S;*p;p++)
		{
			SA.A[++SA.N]=*p-'a'+1;
			Bel[SA.N]=i;
		}
		if (i<N)
			SA.A[++SA.N]=30+i;
	}
}
bool check(int A)
{
	int i,j,k;
	bool ba[MAXN];
	for (i=1;i<=SA.N;i++)
	{
		if (SA.Height[i]>=A)
		{
			for (j=i;SA.Height[j]>=A && j<=SA.N;j++);
			j--;
			memset(ba,0,sizeof(ba));
			for (k=i-1;k<=j;k++)
				ba[Bel[SA.SA[k]]]=true;
			for (k=1;ba[k] && k<=N;k++);
			if (k==N+1)
				return true;
			i=j;
		}
	}
	return false;
}
void solve()
{
	int a,b,m;
	SA.CalcSA();
	SA.CalcHeight();
	a=0;b=SA.N;
	while (a+1<b)
	{
		m=(a+b)/2;
		if (check(m))
			a=m;
		else
			b=m-1;
	}
	if (check(b))
		a=b;
	Ans=a;
}
int main()
{
	init();
	solve();
	printf("%dn",Ans);
	return 0;
}