洛谷 P4173 残缺的字符串 fft

题目描述

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

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

输入输出格式

输入格式：
第一行包含两个正整数 $m$ , $n$ ，分别表示 AA 串和 BB 串的长度。

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

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

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

输出格式：
第一行包含一个整数 $k$ ，表示 $B$ 串中可以完全匹配 $A$ 串的位置个数。
若 $k>0$ ，则第二行输出 $k$ 个正整数，从小到大依次输出每个可以匹配的开头位置（下标从 $1$ 开始）。

输入输出样例

输入样例#1：
3 7
a*b
aebr*ob
输出样例#1：
2
1 5
说明

$100％$的数据满足$1≤m≤n≤300000$。

分析：
一般字符串的距离可以设为

$$dis(A,B)=\sum_{i=1}^{len}(A_i-B_i)^2$$
显然只有当相同时贡献为$0$，平方为了防止正负抵消，$dis(A,B)=0$时，两串相同。
上面这种式子是可以用$fft$做的，把模板串取反。
对于本题，由于串中有$*$，把这个东西在多项式中设为$0$，那么我们可以重新设置
$$dis'(A,B)=\sum_{i=1}^{len}(A_i-B_i)^2*A_i*B_i$$
同样取反模板串，得
$$dis'(A,B)=\sum_{i=1}^{len}(A_i-B_{len-i+1})^2*A_i*B_{len-i+1}$$
随便化简一下，
$$dis'(A,B)=\sum_{i=1}^{len}{A_i}^3*B_{len-i+1}-2*\sum_{i=1}^{len}{A_i}^2*(B_{len-i+1})^2+\sum_{i=1}^{len}{A_i}*(B_{len-i+1})^3$$
然后就是计算$A^3*B$，$A^2*B^2$，$A*B^3$，跑三次$fft$跑出来。
由于常数较大，我们考虑减小常数，我们可以把复数部设为第二个多项式，这个和任意模差不多，然后就可以一次$dft$跑出两个多项式了。我的代码要做$6$次变换。其实算一个多项式可以先保留点值多项式，最后一次逆变换即可。这样只会做$4$次逆变换。

代码：

// luogu-judger-enable-o2
#include <iostream>
#include <cstdio>
#include <cmath>
#include <algorithm>
#define LL long long

const int maxn=3e5+7;
const int maxp=1048580;
const double pi=acos(-1);

using namespace std;

int n,m,len;
char s[maxn],t[maxn];
LL f[maxn],g[maxn],a[maxn],b[maxn],c[maxn],F[maxn],G[maxn];
int r[maxp],ans[maxp];

struct rec{
    double x,y;
}w[maxp],dft[maxp],idft[maxp];

rec operator +(rec a,rec b)
{
    return (rec){a.x+b.x,a.y+b.y};
}

rec operator -(rec a,rec b)
{
    return (rec){a.x-b.x,a.y-b.y};
}

rec operator *(rec a,rec b)
{
    return (rec){a.x*b.x-a.y*b.y,a.x*b.y+a.y*b.x};
}

rec operator !(rec a)
{
    return (rec){a.x,-a.y};
}

void fft(rec *a,int f)
{
    for (int i=0;i<len;i++)
    {
        if (i<r[i]) swap(a[i],a[r[i]]);
    }
    w[0]=(rec){1,0};
    for (int i=2;i<=len;i*=2)
    {
        rec wn=(rec){cos(2*pi/i),f*sin(2*pi/i)};
        for (int j=i/2;j>=0;j-=2) w[j]=w[j/2];
        for (int j=1;j<i/2;j+=2) w[j]=w[j-1]*wn;
        for (int j=0;j<len;j+=i)
        {
            for (int k=0;k<i/2;k++)
            {
                rec u=a[j+k],v=a[j+k+i/2]*w[k];
                a[j+k]=u+v;
                a[j+k+i/2]=u-v;
            }
        }
    }
}

void FFT(LL *a,LL *b,LL *c,int n,int m)
{
    len=1;
    int k=0;
    while (len<(n+m)) len*=2,k++;
    for (int i=0;i<len;i++)
    {
        r[i]=(r[i>>1]>>1)|((i&1)<<(k-1));
    }
    for (int i=0;i<len;i++)
    {
        if (i<n) dft[i].x=a[i]; else dft[i].x=0;
        if (i<m) dft[i].y=b[i]; else dft[i].y=0;
    }
    fft(dft,1); 
    for (int i=0;i<len;i++)
    {
        int j=(len-1)&(len-i);
        rec da,db;
        da=(dft[i]+(!dft[j]))*(rec){0.5,0};
        db=(dft[i]-(!dft[j]))*(rec){0,-0.5};
        idft[i]=da*db;
    }
    fft(idft,-1);
    for (int i=0;i<m;i++) c[i]=trunc(idft[i].x/len+0.5);
}

int main()
{
    scanf("%d%d",&n,&m);
    scanf("%s",s);
    scanf("%s",t);
    for (int i=0;i<n;i++)
    {
        if (s[i]=='*') f[i]=0;
                  else f[i]=s[i]-'a'+1;
        F[i]=f[i]*f[i]*f[i];
    }
    for (int i=0;i<m;i++)
    {
        if (t[i]=='*') g[i]=0;
                  else g[i]=t[i]-'a'+1;
        G[i]=g[i];
    }   
    reverse(f,f+n);
    reverse(F,F+n);
    FFT(F,G,a,n,m);
    for (int i=0;i<n;i++)
    {
        if (f[i]!=0) F[i]=F[i]/f[i];
    }
    for (int i=0;i<m;i++) G[i]=G[i]*g[i];
    FFT(F,G,b,n,m);
    for (int i=0;i<n;i++)
    {
        if (f[i]!=0) F[i]=F[i]/f[i];
    }
    for (int i=0;i<m;i++) G[i]=G[i]*g[i];
    FFT(F,G,c,n,m);
    for (int i=n-1;i<=m-1;i++)
    {
        LL d=a[i]-2*b[i]+c[i];
        if (!d)
        {
            ans[++ans[0]]=i-(n-2);
        }
    }
    printf("%d\n",ans[0]);
    for (int i=1;i<=ans[0];i++) printf("%d ",ans[i]);
}