zoukankan      html  css  js  c++  java
  • FFT求解字符串匹配

    一、最基础的字符串匹配

    给出1个长为n的串S和1个长为m的串T,询问T在S中出现的位置。

    这是kmp经典问题,但现在我们要用FFT解决

    [dis(S_i,T)=sum_{k=0}^{m-1}(S_{i+k}-T_k)^2 ]

    (dis(S_i,T)=0),则S中从i开始的m个字符和T匹配,把T串翻转,就变成了

    [dis(S_i,T)=sum_{k=0}^{m-1}(S_{i+k}-T_{m-k-1})^2 ]

    展开就是

    [dis(S_i,T)=sum_{k=0}^{m-1}(S_{i+k}*S_{i+k}+T_{m-k-1}*T_{m-k-1}-2*S_{i+k}*T_{m-k-1}) ]

    式子前2项可以用前缀和解决,最后一项是卷积形式,可以使用FFT

    二、带有通配符的字符串匹配

    给出1个长为n的串S和1个长为m的串T,两个串均含有通配符,询问T在S中出现的位置。

    借用上面的结论,通配符的位置看作0

    [dis(S_i,T)=sum_{k=0}^{m-1}(S_{i+k}-T_{m-k-1})^2*S_{i+k}*T_{m-k-1} ]

    展开就是

    [dis(S_i,T)=sum_{k=0}^{m-1}((S_{i+k})^3*T_{m-k-1}-2*(S_{i+k})^2*(T_{m-k-1})^2+S_{i+k}*(T_{m-k-1})^3) ]

    做3遍FFT
    (S_i)(T)匹配 (<=>) (F_{m+i-1}=0)

    例题:洛谷4147
    https://www.luogu.com.cn/problem/P4173

    #include<bits/stdc++.h>
    
    using namespace std;
    
    const int N=(1<<20)+10;
    const double eps=1e-9;
    const double pi=acos(-1);
    
    char s1[N],s2[N];
    int t1[N],t2[N];
    
    int rev[N];
    
    struct Complex
    {
    	double x,y;
    	Complex(double x_=0,double y_=0):x(x_),y(y_){}
    	Complex operator + (Complex P)
    	{
    		return Complex(x+P.x,y+P.y);
    	}
    	Complex operator - (Complex P)
    	{
    		return Complex(x-P.x,y-P.y);
    	}
    	Complex operator * (Complex P)
    	{
    		return Complex(x*P.x-y*P.y,x*P.y+y*P.x);
    	}
    };
    typedef Complex E;
    
    E a1[N],a2[N];
    
    double f[N];
    
    int pos[N];
    
    void fft(E *a,int len,int type)
    {
    	for(int i=0;i<len;++i)
    		if(i<rev[i]) swap(a[i],a[rev[i]]);
    	for(int i=1;i<len;i<<=1)
    	{
    		E wn(cos(pi/i),type*sin(pi/i));
    		for(int p=i<<1,j=0;j<len;j+=p)
    		{
    			E w(1,0);
    			for(int k=0;k<i;++k,w=w*wn)
    			{
    				E x=a[j+k],y=a[j+k+i]*w;
    				a[j+k]=x+y; a[j+k+i]=x-y;
    			}
    		}
    	}
    	if(type==-1)
    		for(int i=0;i<len;++i) a[i].x/=len;
    }
    
    void FFT(E *a,E *b,int len,int kk)
    {
    	fft(a,len,1);
    	fft(b,len,1);
    	for(int i=0;i<len;++i) a[i]=a[i]*b[i];
    	fft(a,len,-1);
    	for(int i=0;i<len;++i) f[i]+=kk*a[i].x;
    }
    
    int main()
    {
    	int m,n,ans=0;
    	scanf("%d%d",&m,&n);
    	scanf("%s%s",s1,s2);
    	reverse(s1,s1+m);
    	for(int i=0;i<m;++i) 
    		if(s1[i]!='*') t1[i]=s1[i]-'a'+1;
    	for(int i=0;i<n;++i)
    		if(s2[i]!='*') t2[i]=s2[i]-'a'+1;
    	int num=m+n-1,len=1,bit=0;
    	while(len<num) len<<=1,bit++;
    	for(int i=0;i<len;++i) rev[i]=(rev[i>>1]>>1)|((i&1)<<bit-1);
    	for(int i=0;i<n;++i) 
    	{
    		a1[i].x=t2[i]*t2[i]*t2[i];
    		a1[i].y=0;
    	}
    	for(int i=n;i<len;++i) a1[i].x=a1[i].y=0;
    	for(int i=0;i<m;++i)
    	{
    		a2[i].x=t1[i];
    		a2[i].y=0;
    	}
    	for(int i=m;i<len;++i) a2[i].x=a2[i].y=0;
    	FFT(a1,a2,len,1);
    	for(int i=0;i<n;++i) 
    	{
    		a1[i].x=t2[i];
    		a1[i].y=0;
    	}
    	for(int i=n;i<len;++i) a1[i].x=a1[i].y=0;
    	for(int i=0;i<m;++i)
    	{
    		a2[i].x=t1[i]*t1[i]*t1[i];
    		a2[i].y=0;
    	}
    	for(int i=m;i<len;++i) a2[i].x=a2[i].y=0;
    	FFT(a1,a2,len,1);
    	for(int i=0;i<n;++i) 
    	{
    		a1[i].x=t2[i]*t2[i];
    		a1[i].y=0;
    	}
    	for(int i=n;i<len;++i) a1[i].x=a1[i].y=0;
    	for(int i=0;i<m;++i)
    	{
    		a2[i].x=t1[i]*t1[i];
    		a2[i].y=0;
    	}
    	for(int i=m;i<len;++i) a2[i].x=a2[i].y=0;
    	FFT(a1,a2,len,-2);
    	for(int i=0;i+m-1<n;++i)
    		if(f[i+m-1]<0.5) pos[++ans]=i;
    	printf("%d
    ",ans);
    	for(int i=1;i<=ans;++i) printf("%d ",pos[i]+1); 
    }
    
    
    作者:xxy
    本文版权归作者和博客园共有,转载请用链接,请勿原文转载,Thanks♪(・ω・)ノ。
  • 相关阅读:
    ServletContext
    PS切图
    session实战案例
    Array Destruction
    Session详解
    No More Inversions 全网最详细 回文序列的逆序对
    Sum of Paths (DP、预处理)
    cookie详解
    web的状态管理
    简单最大流/最小割复习
  • 原文地址:https://www.cnblogs.com/TheRoadToTheGold/p/15128004.html
Copyright © 2011-2022 走看看