【快速傅里叶变换应用(字符串匹配)】残缺的字符串

【题目描述】
很久很久以前，在你刚刚学习字符串匹配的时候，有两个仅包含小写字母的字符串A和B，其中A串长度为m，B串长度为n。可当你现在再次碰到这两个串时，这两个串已经老化了，每个串都有不同程度的残缺。

你想对这两个串重新进行匹配，其中A为模板串，那么现在问题来了，请回答，对于B的每一个位置i，从这个位置开始连续m个字符形成的子串是否可能与A串完全匹配?

【输入】
第一行包含两个正整数m,n(1<=m<=n<=300000)，分别表示A串和B串的长度。

第二行为一个长度为m的字符串A。

第三行为一个长度为n的字符串B。

两个串均仅由小写字母和号组成，其中号表示相应位置已经残缺。

【输出】
第一行包含一个整数k，表示B串中可以完全匹配A串的位置个数。

若k>0，则第二行输出k个正整数，从小到大依次输出每个可以匹配的开头位置（下标从1开始）。

【样例输入】
3 7
ab
aebrob

【样例输出】
2
1 5

思路

这道题很有意思。字符串本身是不确定的，所以各个字符串的算法都显得十分无力。
1.一般情况下字符串匹配
那么此时，我们考虑最暴力的做法。对于每一个字符我们都给它赋予一个值。定义c表示从b的第i个位置起与a的匹配情况。那么在没有通配符的条件下，我们可以得到这样一个式子。
$$c[i]=\sum _{j=0}^{m-1}(b[i+j]-a[j])$$
但是这样依然存在问题，比如ba和ab会被判定为匹配，这都是因为万恶的负数。所以我们考虑把每一个值平方，这样当且仅当$c[i]==0$时字符串匹配。即：
$$c[i]=\sum _{j=0}^{m-1}(b[i+j]-a[j]{)}^2$$
可是这样将会达到O($n^2$)的复杂度，显然是不够优秀的。这时候，我们注意到a和b的下标的差为定值，这足以让我们联想到多项式乘法中两个多项式的和为定值。所以我们考虑把a数组翻转。这样上面这个式子就变成了：
$$c[i]=\sum _{j=0}^{m-1}(b[i+j]-a[len-j-1]{)}^2$$
显然，此时对于每一个枚举的i,a和b下标的和就是定值了。
那么我们把上面这个式子拆开：
$$c[i]=\sum _{j=0}^{m-1}(b[i+j{]}^2+a[len-j-1{]}^2-2\ast b[i+j]\ast a[len-j-1])$$
显然，对于前两项，我们可以直接简单处理得到，而对于第三项，它正是一个多项式乘法，所以就可以用快速傅里叶变换完成了。

2.含通配符的字符串匹配
仿照上面的思路，我们同样定义一个c数组，但是考虑当某一位置为通配符时，怎么样让它在上式中乘法的贡献为0。简单而粗暴地，我们把通配符的值赋为0，那么：
$$c[i]=\sum _{j=0}^{m-1}(b[i+j]-a[len-j-1]{)}^2\ast b[i+j]\ast a[len-j-1]$$
于是问题又巧妙地解决了。
我们把式子拆开：
$$c[i]=\sum _{j=0}^{m-1}(b[i+j{]}^3\ast a[len-j-1]+a[len-j-1{]}^3\ast b[i+j]-2\ast b[i+j{]}^2\ast a[len-j-1{]}^2)$$
显然，虽然有高次，但是式子的值依然只与下标有关，所以使用快速傅里叶变换做三次多项式乘法。

#include<bits/stdc++.h>
#define re register
using namespace std;
const int N=1.1e6+5;
int n,m;
int a[N],b[N]; 
inline int red()
{
    int data=0;int w=1; char ch=0;
    ch=getchar();
    while(ch!='-' && (ch<'0' || ch>'9')) ch=getchar();
    if(ch=='-') w=-1,ch=getchar();
    while(ch>='0' && ch<='9') data=(data<<3)+(data<<1)+ch-'0',ch=getchar();
    return data*w;
}
struct cp{
	double x,y;
	cp(double _x=0,double _y=0){x=_x,y=_y;}
	friend cp operator +(cp a,cp b){return cp(a.x+b.x,a.y+b.y);}
	friend cp operator -(cp a,cp b){return cp(a.x-b.x,a.y-b.y);}
	friend cp operator *(cp a,cp b){return cp(a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x);}
}p1[N],p2[N],p[N];
int lim=1,l=0,rev[N];
char s2[N>>1|1],s1[N>>1|1];
double pi=3.14159265358979323846;
inline void pre()
{
	while(lim<=(n+m))lim<<=1,l++;
	for(int re i=0;i<lim;i++)rev[i]=(rev[i>>1]>>1)|((i&1)<<(l-1));
}
inline void fft(cp *f,int op)
{
	for(int re i=0;i<lim;i++)if(i>rev[i])swap(f[i],f[rev[i]]);
	for(int re mid=1;mid<lim;mid<<=1)
	{
		cp tmp(cos(pi/mid),op*sin(pi/mid));
		for(int re j=0;j<lim;j+=(mid<<1))
		{
			cp w(1,0);
			for(int re k=0;k<mid;k++,w=w*tmp)
			{
				cp x=f[j+k],y=f[j+k+mid]*w;
				f[j+k]=x+y;f[j+k+mid]=x-y;
			}
		}
	}
	if(op==-1)for(int re i=0;i<lim;i++)f[i].x/=lim;
}
int ans[N],cnt=0;
int main()
{
	m=red();n=red();scanf("%s%s",s1,s2);pre();
    reverse(s1, s1 + m);
	for(int re i=0;i<m;i++)a[i]=(s1[i]!='*')*(s1[i]-'a'+1);
	for(int re i=0;i<n;i++)b[i]=(s2[i]!='*')*(s2[i]-'a'+1);
	
	for(int re i=0;i<m;i++)p1[i].x=a[i]*a[i]*a[i];
	for(int re i=0;i<n;i++)p2[i].x=b[i];
	fft(p1,1);fft(p2,1);
	for(int re i=0;i<lim;i++)p[i]=p[i]+p1[i]*p2[i];
	
	memset(p1,0,sizeof(p1));
	memset(p2,0,sizeof(p2));
	for(int re i=0;i<m;i++)p1[i].x=a[i];
	for(int re i=0;i<n;i++)p2[i].x=b[i]*b[i]*b[i];
	fft(p1,1);fft(p2,1);
	for(int re i=0;i<lim;i++)p[i]=p[i]+p1[i]*p2[i];
	
	memset(p1,0,sizeof(p1));
	memset(p2,0,sizeof(p2));
	for(int re i=0;i<m;i++)p1[i].x=a[i]*a[i]<<1;
	for(int re i=0;i<n;i++)p2[i].x=b[i]*b[i];
	fft(p1,1);fft(p2,1);
	for(int re i=0;i<lim;i++)p[i]=p[i]-p1[i]*p2[i];
	
	fft(p,-1);
	for(int re i=m-1;i<n;i++)
		if((long long)(p[i].x+0.5)==0)
			ans[++cnt]=i+2-m;
	printf("%d\n",cnt);
	for(int re i=1;i<=cnt;i++)
		printf("%d ",ans[i]);
}