bzoj4259: 残缺的字符串

传送门：http://www.lydsy.com/JudgeOnline/problem.php?id=4259

思路：FFT是一种高效的字符串匹配算法

只要想到FFT就好办了

令为‘*’的位置值为0

那么两个字符串匹配当且仅当

Σ(a[i]-b[i])^2*a[i]*b[i]=0

式子的含义很好理解，前面就是字符相同，后面就是该位是否为通配符

为了FFT能做，我们要把A串反过来才有一个卷积的形式

设F[i]就是B串以i结尾的子串B[i-m+1,i]与A匹配得到的值

f[i]=Σ(a[j]-b[i-j])^2*a[j]*b[i-j]

    =Σ(a[j]^3*b[i-j]-2*(a[j])^2*(b[i-j])^2)+a[j]*b[i-j]^3)

分成三个部分，分别做FFT，加起来就是答案了

（公式写的比较挫，可以去看claris的博客http://www.cnblogs.com/clrs97/p/4814499.html）

#include<cmath>
#include<cstdio>
#include<cstring>
#include<iostream>
#include<algorithm>
const int maxn=300010,maxt=(1<<19)+10;
using namespace std;
int n,m,ans[maxn],a[maxn],b[maxn];char s1[maxn],s2[maxn];
struct com{double r,i;}A[maxt],B[maxt],C[maxt];
com operator +(com a,com b){return (com){a.r+b.r,a.i+b.i};}
com operator -(com a,com b){return (com){a.r-b.r,a.i-b.i};}
com operator *(com a,com b){return (com){a.r*b.r-a.i*b.i,a.r*b.i+a.i*b.r};}
com operator /(com a,int b){return (com){a.r/b,a.i/b};}
com operator *(com a,int b){return (com){a.r*b,a.i*b};}
double sqr(int x){return x*x;}
double cube(int x){return x*x*x;}

struct FFT{
	int rev[maxt],n;
	void clear(com a[]){memset(a,0,sizeof(*a)*(n+1));}
	void init(int nn){
		for (n=1;n<nn;n<<=1);
		for (int i=1;i<n;i++) rev[i]=rev[i-(i&-i)]+(n>>1)/(i&-i);
	}
	void fft(com a[],int op){
		for (int i=1;i<n;i++) if (rev[i]<i) swap(a[rev[i]],a[i]);
		for (int sz=2;sz<=n;sz<<=1){
			com w0=(com){cos(2*M_PI/sz),sin(2*M_PI*op/sz)};
			for (int bg=0;bg<n;bg+=sz){
				com w=(com){1,0};
				for (int pos=bg;pos<bg+(sz>>1);pos++){
					com x=a[pos],y=w*a[pos+(sz>>1)];
					a[pos]=x+y,a[pos+(sz>>1)]=x-y;
					w=w*w0;
				}
			}
		}
		if (op==-1) for (int i=0;i<n;i++) a[i]=a[i]/n;
	}
}T;

void init(){
	scanf("%d%d",&m,&n);
	scanf("%s%s",s1,s2);
	for (int i=0;i<m;i++) a[i]=s1[i]=='*'?0:s1[i]-'a'+1;
	for (int i=0;i<n;i++) b[i]=s2[i]=='*'?0:s2[i]-'a'+1;
	for (int i=0;i<(m>>1);i++) swap(a[i],a[m-i-1]);
	//for (int i=0;i<m;i++) printf("%d ",(int)a[i]);puts("");
	//for (int i=0;i<n;i++) printf("%d ",(int)b[i]);puts("");
}

void work(){
	T.init(n);
	
	T.clear(A),T.clear(B);
	for (int i=0;i<m;i++) A[i]=(com){a[i],0};
	for (int i=0;i<n;i++) B[i]=(com){cube(b[i]),0};
	T.fft(A,1),T.fft(B,1);
	for (int i=0;i<T.n;i++) C[i]=C[i]+A[i]*B[i];
	//for (int i=0;i<T.n;++i) printf("r=%.3lf i=%.3lf\n",C[i].r,C[i].i);puts("");
	
	T.clear(A),T.clear(B);
	for (int i=0;i<m;i++) A[i]=(com){sqr(a[i]),0};
	for (int i=0;i<n;i++) B[i]=(com){sqr(b[i]),0};
	T.fft(A,1),T.fft(B,1);
	for (int i=0;i<T.n;i++) C[i]=C[i]-A[i]*B[i]*2;
	//for (int i=0;i<T.n;++i) printf("r=%.3lf i=%.3lf\n",C[i].r,C[i].i);puts("");
	
	T.clear(A),T.clear(B);
	for (int i=0;i<m;i++) A[i]=(com){cube(a[i]),0};
	for (int i=0;i<n;i++) B[i]=(com){b[i],0};
	T.fft(A,1),T.fft(B,1);
	for (int i=0;i<T.n;i++) C[i]=C[i]+A[i]*B[i];
	//for (int i=0;i<T.n;++i) printf("r=%.3lf i=%.3lf\n",C[i].r,C[i].i);puts("");
	
	T.fft(C,-1);
	int cnt=0;
	for (int i=0;i<n-m+1;i++)
		if (fabs(C[m-1+i].r)<1e-1)
			ans[++cnt]=i+1;
	printf("%d\n",cnt);
	for (int i=1;i<cnt;i++) printf("%d ",ans[i]);
	printf("%d\n",ans[cnt]);
}

int main(){
	init(),work();
	return 0;
}