Beyond the Void
BYVoid
最長公共子串問題的後綴數組解法
本文正體字版由OpenCC轉換

[最長公共子串]

最長公共子串(Longest Common Substring ,簡稱LCS)問題,是指求給定的一組字符串長度最大的共有的子串的問題。例如字符串"abcb",“bca”,“acbc"的LCS就是"bc”。

求多串的LCS,顯然窮舉法是極端低效的算法。改進一些的算法是用一個串的每個後綴對其他所有串進行部分匹配,用KMP算法,時間複雜度爲O(NL^2),其中N爲字符串個數,L爲每個串的長度。更優秀的有廣義後綴樹的方法,時間可以達到 O(NL)。本文介紹一種基於後綴數組的LCS解法,利用二分查找技術,時間複雜度可以達到O(NLlogL)。

[最長公共子串問題的後綴數組解法]

關於後綴數組的構建方法以及Height數組的性質,本文不再具體介紹,可以參閱IOI國家集訓隊2004年論文《後綴數組》(許智磊)和IOI國家集訓隊2009年論文《後綴數組——處理字符串的有力工具》(羅穗騫)

回顧一下後綴數組,SA[i]表示排名第i的後綴的位置,Height[i]表示後綴SA[i]和SA[i-1]的最長公共前綴(Longest Common Prefix,LCP),簡記爲Height[i]=LCP(SA[i],SA[i-1])。連續的一段後綴SA[i..j]的最長公共前綴,就是H[i-1..j]的最小值,即LCP(SA[i..j])=Min(H[i-1..j])。

求N個串的最長公共子串,可以轉化爲求一些後綴的最長公共前綴的最大值,這些後綴應分屬於N個串。具體方法如下:

設N個串分別爲S1,S2,S3,…,SN,首先建立一個串S,把這N個串用不同的分隔符連接起來。S=S1[P1]S2[P2]S3…SN-1[PN-1]SN,P1,P2,…PN-1應爲不同的N-1個不在字符集中的字符,作爲分隔符(後面會解釋爲什麼)。

接下來,求出字符串S的後綴數組和Height數組,可以用倍增算法,或DC3算法。

然後二分枚舉答案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的長度爲NL+N-1,所以驗證的時間複雜度爲O(NL)。

到這裏,我們就可以理解爲什麼分隔符P1..PN-1必須是不同的N-1個不在字符集中的字符了,因爲這樣才能保證S的後綴的公共前綴不會跨出一個原有串的範圍。

後綴數組是一種處理字符串的強大的數據結構,配合LCP函數與Height數組的性質,後綴數組更是如虎添翼。利用後綴數組,容易地求出了多個串的LCS,而且時空複雜度也相當優秀了。雖然比起後綴樹的解法有所不如,但其簡明的思路和容易編程的特點卻在實際的應用中並不輸於後綴樹。

附:POI 2000 Repetitions 最長公共子串

/* 
 * Problem: POI2000 pow
 * Author: Guo Jiabao
 * Time: 2009.4.16 21:00
 * State: Solved
 * Memo: 多串最長公共子串 後綴數組 二分答案
*/
#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;
}

上次修改時間 2017-02-03

相關日誌