[BZOJ 3145]Str(SAM+启发式合并)

[BZOJ 3145]Str(SAM+启发式合并)

题面

Arcueid得知了两者的DNA片段，想寻求一个DNA片段，使得其在两者的DNA中都出现过。我们知道公主的脑袋有点不太灵活，如果两个DNA片段只有一个位置不同，她也会将其认为是相同的。所以请您找出这样的最长的DNA片段吧。

分析

考虑两个串中那个不同的位置,不妨设在第一个串中的位置为(i),第二个串中的位置为(j).那么最长长度就是(operatorname{lcs}(i-1,j-1)+1+operatorname{lcp}(i+1,j+1)).其中(operatorname{lcs(x,y)})分别表示前缀(x)和前缀(y)的最长公共后缀,和后缀(x)和后缀(y)的最长公共前缀(operatorname{lcp(x,y)})

那么我们可以把两个串接在一起,中间加一个分割符，得到一个新串(all)。然后对新串建后缀自动机。那么对于(operatorname{lcs}),就是该后缀自动机上(x)和(y)对应节点,在fail树上的LCA对应的len.因为沿着fail树向上跳相当于从长到短遍历后缀。(operatorname{lcp}),在新串的反串构成的SAM上同理。

那么问题就转化成,给出两棵树,且树的节点间有对应关系。求一些节点的lca深度和的最大值。不妨记这两棵树为(T_1,T_2). 我们遍历树(T_1),对于当前子树中的每个节点(x),我们找出子树里每个节点,它对应(all)中的一个位置(pos)的前缀。用set维护其子树内所有的合法前缀在(T_2)上的dfs序。具体地，我们找到(all)的反串中(pos)对应的位置,把这个后缀节点插入(x)对应的set. 也就是说. 按照dfs序是因为动态维护树上一堆点的lca时,插入一个点,只需要用和这个点dfs序相邻的点去更新答案。这样对每个子树启发式合并,就得到了(x)子树节点在(T_2)中的LCA.而(x)子树中的节点在(T_1)中的LCA就是(x). 因此答案就是(T_2)中的LCA的len+1+x的len.

当然也有可能只有前缀或后缀,这种情况直接枚举一下分割点即可。

感觉语言很难描述,可以看代码。注意常数优化。

代码

//sdfsdfs
#include<iostream>
#include<cstdio>
#include<cstring>
#include<set>
#include<vector>
#define maxn 400000
#define maxlogn 18
#define maxc 27
using namespace std;
int n,m,tot;
char s[maxn+5],t[maxn+5];
char all[maxn+5];;//s和t连起来 
struct SAM{
#define link(x) t[x].link 
#define len(x) t[x].len
	struct node{
		int ch[maxc];
		int link;
		int len;
		int right;
		bool flag;//表示这一位是否是结尾,如果是结尾就不能模糊下一位 
	}t[maxn*2+5];
	const int root=1;
	int last=root;
	int ptr=1;
	int extend(int c,int pos,int f){
		int p=last,cur=++ptr;
		len(cur)=len(p)+1;
		t[cur].right=pos;
		t[cur].flag=f;
		while(p&&t[p].ch[c]==0){
			t[p].ch[c]=cur;
			p=link(p);
		}
		if(p==0) link(cur)=root;
		else{
			int q=t[p].ch[c];
			if(len(p)+1==len(q)) link(cur)=q;
			else{
				int clo=++ptr;
				len(clo)=len(p)+1;
				for(int i=0;i<maxc;i++) t[clo].ch[i]=t[q].ch[i];
				link(clo)=link(q);
				link(q)=link(cur)=clo;
				while(p&&t[p].ch[c]==q){
					t[p].ch[c]=clo;
					p=link(p);
				}
			}
		} 
		last=cur;
		return cur;
	}
	
	int tim=0;
	vector<int>fail[maxn+5];
	int dfn[maxn+5];
	int deep[maxn+5];
//	int anc[maxn+5][maxlogn+5];
	int sz[maxn+5];
	int top[maxn+5],fa[maxn+5],son[maxn+5];
	void dfs1(int x,int f){
		deep[x]=deep[f]+1;
		sz[x]=1;
		fa[x]=f;
//		anc[x][0]=fa;
//		for(int i=1;i<=maxlogn;i++) anc[x][i]=anc[anc[x][i-1]][i-1];
		for(int y : fail[x]){
			dfs1(y,x);
			sz[x]+=sz[y];
			t[x].flag|=t[y].flag;
			if(sz[y]>sz[son[x]]) son[x]=y;
		}
	}
	void dfs2(int x,int t){
		dfn[x]=++tim;
		top[x]=t;
		if(son[x]) dfs2(son[x],t);
		for(int y : fail[x]){
			if(y!=son[x]) dfs2(y,y); 
		} 
	}
	void build(){
		for(int i=2;i<=ptr;i++) fail[link(i)].push_back(i);
		dfs1(1,0);
		dfs2(1,1);
	}
	int lca(int x,int y){
		while(top[x]!=top[y]){
			if(deep[top[x]]<deep[top[y]]) swap(x,y);
			x=fa[top[x]];
		}
		if(deep[x]<deep[y]) return x;
		else return y;
	}
	int lcp(int x,int y){//fail树上LCA的len就是LCP长度 
		return len(lca(x,y));
	}
	
}T1,T2,T3;
//T1是正序插入all,T2是倒序插入all,T3是正序插入S 

struct cmp{ 
	bool operator () (int a,int b){
		return T2.dfn[a]<T2.dfn[b];
	} 
}; 
int res;
set<int,cmp>p1[maxn+5],p2[maxn+5];//s串和t串在总的自动机上对应的节点 
int ans[maxn+5];
int id1[maxn+5],id2[maxn+5];//正串和反串的前i位在自动机上的编号 
void merge(int x,int y){
	if(p1[x].size()+p2[x].size()<p1[y].size()+p2[y].size()){
		p1[x].swap(p1[y]);
		p2[x].swap(p2[y]);
	}
	for(int t : p1[y]){
		auto it=p2[x].lower_bound(t);
		if(it!=p2[x].end()) ans[x]=max(ans[x],T2.lcp(*it,t));
		if(it!=p2[x].begin()) ans[x]=max(ans[x],T2.lcp(*(--it),t));//和dfs序相邻的更新lca 
	}
	for(int t : p2[y]){
		auto it=p1[x].lower_bound(t);
		if(it!=p1[x].end()) ans[x]=max(ans[x],T2.lcp(*it,t));
		if(it!=p1[x].begin()) ans[x]=max(ans[x],T2.lcp(*(--it),t));//和dfs序相邻的更新lca 
	}
	for(int t : p1[y]) p1[x].insert(t);
	for(int t : p2[y]) p2[x].insert(t);
//	p1[y].clear();
//	p2[y].clear(); 
}
int solve1(){
	int x=T3.root;
	int matlen=0;
	int ans=0;
	for(int i=n+2;i<=tot;i++){
		int c=all[i]-'a';
		if(T3.t[x].ch[c]){
			x=T3.t[x].ch[c];
			matlen++;
		}else{
			while(x&&T3.t[x].ch[c]==0) x=T3.link(x);
			if(x==0){
				x=T3.root;
				matlen=0;
			}else{
				matlen=T3.len(x)+1;
				x=T3.t[x].ch[c];
			}
		}
		ans=max(ans,matlen+T3.t[x].flag);
	}
	return ans;
}
void solve2(int x){
	int pos=T1.t[x].right;//pos之后是一个模糊点 
	if(pos>0&&pos+2<=n) p1[x].insert(id2[pos+2]);//pos属于s 
	if(pos>=n+2&&pos<=tot-2) p2[x].insert(id2[pos+2]);//pos属于t
	for(int y : T1.fail[x]){//S1上x的子树中的节点,x是它们的公共祖先,len(x)一定是lcp长度 
		solve2(y);
		merge(x,y); 
	} 
	if(ans[x]) res=max(res,ans[x]+T1.len(x)+1);//对反串S2合并得到lcs长度 
}

int main(){
#ifdef LOCAL
	freopen("4.in","r",stdin);
#endif
	scanf("%s",s+1);
	n=strlen(s+1);
	scanf("%s",t+1);
	m=strlen(t+1);
	for(int i=1;i<=n;i++) all[++tot]=s[i];
	all[++tot]='z'+1;
	for(int i=1;i<=m;i++) all[++tot]=t[i];
	for(int i=1;i<=tot;i++) id1[i]=T1.extend(all[i]-'a',i,1);
	for(int i=tot;i>=1;i--) id2[i]=T2.extend(all[i]-'a',tot-i+1,1);
	for(int i=1;i<=n;i++) T3.extend(all[i]-'a',i,i!=n);
	T1.build();
	T2.build();
	T3.build();
	res=max(res,solve1());
	solve2(1);
	printf("%d
",res);
}
/*
abac
acac
*/