您的位置:首页 > 其它

洛谷3763:[TJOI2017]DNA——题解

2018-06-19 16:37 405 查看

https://www.luogu.org/problemnew/show/P3763

加里敦大学的生物研究所,发现了决定人喜不喜欢吃藕的基因序列S,有这个序列的碱基序列就会表现出喜欢吃藕的性状,但是研究人员发现对碱基序列S,任意修改其中不超过3个碱基,依然能够表现出吃藕的性状。现在研究人员想知道这个基因在DNA链S0上的位置。所以你需要统计在一个表现出吃藕性状的人的DNA序列S0上,有多少个连续子串可能是该基因,即有多少个S0的连续子串修改小于等于三个字母能够变成S。

是的这篇代码过不了BZOJ(因为懒得卡了/不想写SAM/不会DC3……众多原因)

其实容错三次匹配并不吓人,我们可以先跳跃匹配到匹配不上的地方,然后cnt++,继续跳跃……直到匹配完全或者cnt>3为止。

这个跳跃完全可以枚举起点,然后用SA来求lcp进而实现跳跃匹配以此变成$O(n)$的。

所以总复杂度是$O(Tnlogn)$的……只要卡卡就能过洛谷。

当然为了过BZOJ,要么常数优秀(写SAM,然后遍历,每次选择一个节点往下走,如果和当前节点匹配不上则cnt++,匹配复杂度不变但是常数小),要么就学DC3,要么……其实后缀数组卡卡也能过。

#include<map>
#include<cmath>
#include<stack>
#include<queue>
#include<cstdio>
#include<cctype>
#include<vector>
#include<cstdlib>
#include<cstring>
#include<iostream>
#include<algorithm>
using namespace std;
typedef long long ll;
const int N=2e5+5;
inline int turn(char ch){
if(ch==0)return 0;
if(ch=='A')return 1;
if(ch=='G')return 2;
if(ch=='C')return 3;
if(ch=='T')return 4;
return 5;
}
char s
,p
;
int n,rk
,height
,w
,sa
;
inline bool pan(int *x,int i,int j,int k){
int ti=i+k<n?x[i+k]:-1;
int tj=j+k<n?x[j+k]:-1;
return ti==tj&&x[i]==x[j];
}
void SA_init(){
int *x=rk,*y=height,r=6;
for(int i=0;i<r;i++)w[i]=0;
for(int i=0;i<n;i++)w[turn(s[i])]++;
for(int i=1;i<r;i++)w[i]+=w[i-1];
for(int i=n-1;i>=0;i--)sa[--w[turn(s[i])]]=i;
r=1;x[sa[0]]=0;
for(int i=1;i<n;i++)
x[sa[i]]=s[sa[i]]==s[sa[i-1]]?r-1:r++;
for(int k=1;r<n;k<<=1){
int yn=0;
for(int i=n-k;i<n;i++)y[yn++]=i;
for(int i=0;i<n;i++)
if(sa[i]>=k)y[yn++]=sa[i]-k;
for(int i=0;i<r;i++)w[i]=0;
for(int i=0;i<n;i++)w[x[y[i]]]++;
for(int i=1;i<r;i++)w[i]+=w[i-1];
for(int i=n-1;i>=0;i--)sa[--w[x[y[i]]]]=y[i];
swap(x,y);r=1;x[sa[0]]=0;
for(int i=1;i<n;i++)
x[sa[i]]=pan(y,sa[i],sa[i-1],k)?r-1:r++;
}
}
void height_init(){
int i,j,k=0;
for(int i=1;i<=n;i++)rk[sa[i]]=i;
for(int i=0;i<n;i++){
if(k)k--;
int j=sa[rk[i]-1];
while(s[i+k]==s[j+k])k++;
height[rk[i]]=k;
}
}
int f
[20],lg
;
inline int qpow(int a){return 1<<a;}
void st_init(){
for(int i=1;i<=n;i++){
f[i-1][0]=height[i];
lg[i]=lg[i-1];
if((1<<lg[i]+1)==i)lg[i]++;
}
for(int j=1;j<=lg
;j++){
for(int i=0;i<n;i++){
if(i+qpow(j)-1>=n)break;
f[i][j]=min(f[i][j-1],f[i+qpow(j-1)][j-1]);
}
}
}
int lcp(int i,int j){
int l=rk[i],r=rk[j];if(l>r)swap(l,r);
l--;r--;if(r<0)return 0;l++;
int len=r-l+1,k=lg[len],h=qpow(k);
return min(f[l][k],f[r-h+1][k]);
}
int main(){
int t;scanf("%d",&t);
while(t--){
scanf("%s%s",s,p);
n=strlen(s);int m=strlen(p);
s[n++]='#';
for(int i=0;i<m;i++)s[n++]=p[i];
s[n++]=0;SA_init();n--;height_init();st_init();
int ans=0;
for(int i=0;i<n-2*m;i++){
int cnt=0;
for(int j=0;j<m&&cnt<=3;){
if(s[i+j]!=s[n-m+j])cnt++,j++;
else j+=lcp(i+j,n-m+j);
}
if(cnt<=3)ans++;
}
printf("%d\n",ans);
}
return 0;
}

+++++++++++++++++++++++++++++++++++++++++++

 +本文作者:luyouqi233。               +

 +欢迎访问我的博客:http://www.cnblogs.com/luyouqi233/+

+++++++++++++++++++++++++++++++++++++++++++

内容来自用户分享和网络整理,不保证内容的准确性,如有侵权内容,可联系管理员处理 点击这里给我发消息
标签: