BZOJ4259 : 残缺的字符串

假设字符串是从第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;
}