BZOJ4259 : 殘缺的字串

Claris發表於2015-09-16
字串

假設字串是從第0位開始的，那麼對於兩個長度都為n的字串A,B，定義距離函式\[dis(A,B)=\sum_{i=0}^{n-1}(A[i]-B[i])^2[A[i]!='*'][B[i]!='*']\]

若把*號都設定為0，那麼有\[dis(A,B)=\sum_{i=0}^{n-1}(A[i]-B[i])^2A[i]B[i]\]

如果$dis(A,B)=0$，那麼A和B完全匹配。

對於這個問題，假設我們列舉B的末尾位置i，設$f[i]=dis(A,B[i-m+1,i])$，那麼B的這一個子串與A完全匹配，有

\[f[i]=\sum_{j=0}^{m-1}(A[j]-B[i-m+1+j])^2A[j]B[i-m+1+j]=0\]

如果把A串翻轉，並在後面不斷補0直至和B串等長的話，那麼有

\[\begin{eqnarray*}
f[i]&=&\sum_{j=0}^i(A[j]-B[i-j])^2A[j]B[i-j]\\
&=&\sum_{j=0}^i(A[j]^2-2A[j]B[i-j]+B[i-j]^2)A[j]B[i-j]\\
&=&\sum_{j=0}^iA[j]^3B[i-j]-2\sum_{j=0}^iA[j]^2B[i-j]^2+\sum_{j=0}^iA[j]B[i-j]^3
\end{eqnarray*}\]

顯然可以分成三段做FFT求出所有的f[i]，時間複雜度為$O(n\log n)$。

 

#include<cstdio>
#include<cmath>
#include<algorithm>
#define N 1048576
using namespace std;
char sa[N],sb[N];int n,m,i,j,k,a[N],b[N],ans,q[N];
struct comp{
  double r,i;comp(double _r=0,double _i=0){r=_r;i=_i;}
  comp operator+(const comp&x){return comp(r+x.r,i+x.i);}
  comp operator-(const comp&x){return comp(r-x.r,i-x.i);}
  comp operator*(const comp&x){return comp(r*x.r-i*x.i,r*x.i+i*x.r);}
}A[N],B[N],C[N];
const double pi=acos(-1.0);
void FFT(comp*a,int n,int t){
  for(int i=1,j=0;i<n-1;i++){
    for(int s=n;j^=s>>=1,~j&s;);
    if(i<j)swap(a[i],a[j]);
  }
  for(int d=0;(1<<d)<n;d++){
    int m=1<<d,m2=m<<1;
    double o=pi/m*t;comp _w(cos(o),sin(o));
    for(int i=0;i<n;i+=m2){
      comp w(1,0);
      for(int j=0;j<m;j++){
        comp &A=a[i+j+m],&B=a[i+j],t=w*A;
        A=B-t;B=B+t;w=w*_w;
      }
    }
  }
  if(t==-1)for(int i=0;i<n;i++)a[i].r/=n;
}
int main(){
  scanf("%d%d%s%s",&m,&n,sa,sb);
  for(i=0,j=m-1;i<j;i++,j--)swap(sa[i],sa[j]);
  for(i=0;i<m;i++)if(sa[i]!='*')a[i]=sa[i]-'a'+1;
  for(i=0;i<n;i++)if(sb[i]!='*')b[i]=sb[i]-'a'+1;
  for(k=1;k<n+m;k<<=1);
  for(i=0;i<k;i++)A[i]=comp(a[i]*a[i]*a[i],0),B[i]=comp(b[i],0);
  for(FFT(A,k,1),FFT(B,k,1),i=0;i<k;i++)C[i]=C[i]+A[i]*B[i];
  for(i=0;i<k;i++)A[i]=comp(a[i],0),B[i]=comp(b[i]*b[i]*b[i],0);
  for(FFT(A,k,1),FFT(B,k,1),i=0;i<k;i++)C[i]=C[i]+A[i]*B[i];
  for(i=0;i<k;i++)A[i]=comp(a[i]*a[i],0),B[i]=comp(b[i]*b[i],0);
  for(FFT(A,k,1),FFT(B,k,1),i=0;i<k;i++)C[i]=C[i]-A[i]*B[i]*comp(2,0);
  FFT(C,k,-1);
  for(i=m-1;i<n;i++)if(C[i].r<0.5)q[++ans]=i-m+2;
  for(printf("%d\n",ans),i=1;i<ans;i++)printf("%d ",q[i]);
  if(ans)printf("%d",q[ans]);
  return 0;
}

　　

相關文章