P4173 残缺的字符串 字符串匹配 快速傅里叶变换fft

若字符串中不含通配符如何进行匹配

设$A$是模式串，长度$m$，$B$是原串，长度$n$。

定义匹配函数 $P(x)={\sum }_{i=0}^{m-1}[A(i)-B(x-m+i+1)]$

当$P(x)=0$时，字符串$B$的以第$x$为结束的连续$m$位，与$A$完全匹配

因为上述匹配函数是加和的形式进行判断，无法判断顺序

考虑将匹配函数进行平方，消除顺序的影响

$P(x)={\sum }_{i=0}^{m-1}[A(i)-B(x-m+i+1){]}^2$

上述函数无法使用$fft$进行求解，将$A$逆置

$P(x)={\sum }_{i=0}^{m-1}[S(m-i-1)-B(x-m+i+1){]}^2$

拆项

$P(x)={\sum }_{i=0}^{m-1}S(m-i-1{)}^2+{\sum }_{i=0}^{m-1}B(x-m+i+1{)}^2-2{\sum }_{i=0}^{m-1}S(m-i-1)B(x-m+i+1)$

观察最后一项：$m-i-1+x-m+i+1=x$

$P(x)={\sum }_{i=0}^{m-1}S(i{)}^2+{\sum }_{i=0}^{m-1}B(i{)}^2-2{\sum }_{i+j=x}^{m-1}S(i)B(j)$

题意

给你模式串和原串，串中存在通配符*，输出原串中可以和模式串完全匹配的每个开头位置。

通过上述对基本字符串的匹配，应该对如何用$fft$有了简单了解。

带有通配符？?

定义匹配函数
$P(x)={\sum }_{i=0}^{m-1}[A(i)-B(x-m+i+1){]}^{2}A(i)B(x-m+i+1)$

$S$是$A$的翻转：$S(i)=A(m-i-1)$

$P(x)={\sum }_{i=0}^{m-1}[S(m-i-1)-B(x-m+i+1){]}^{2}S(m-i-1)B(x-m+i+1)$

展开化简得到

$f_i={\sum }_{i+j=x}S(i{)}^{3}B(j)-2{\sum }_{i+j=x}S(i{)}^{2}B(j{)}^2+{\sum }_{i+j=x}S(i)B(j{)}^3$

使用$fft$计算上述柿子，得到答案。

代码

#include <iostream>
#include <string>
#include <cstring>
#include <cmath>
#include <algorithm>
#include <queue>
#include <iomanip>
#include <map>
#include <cstdio>
#include <stack>
#include <set>
using namespace std;
typedef long long ll;
typedef unsigned long long ull;
typedef pair<int ,int > pii;
#define endl '\n'
ll gcd(ll a, ll b){
    return b == 0 ? a : gcd(b, a % b);
}
void input(){
    freopen("in.txt", "r", stdin);
    freopen("out.txt", "w", stdout);
}
inline int read(){
	int x=0,f=1;char c=getchar();
	while(c<'0'||c>'9') {if(c=='-') f=-1;c=getchar();}
	while (c>='0'&&c<='9') x=(x<<3)+(x<<1)+(c^48),c=getchar();
	return x*f;
}
const int N = 2e6+10, M = N * 2, inf = 1e8;
const double PI = acos(-1);
int rev[N], bit, tot;
struct Complex{
    double x, y;
    Complex operator+ (const Complex& t) const { return {x + t.x, y + t.y}; }
    Complex operator- (const Complex& t) const{ return {x - t.x, y - t.y}; }
    Complex operator* (const Complex& t) const{ return {x * t.x - y * t.y, x * t.y + y * t.x}; }
}a[N], b[N], p[N];
void fft(Complex a[], int inv){
    for(int i = 0; i < tot; i++)
        if(i < rev[i]) swap(a[i], a[rev[i]]);
    for(int mid = 1; mid < tot; mid <<= 1){
        Complex w1 = Complex({cos(PI/mid), inv * sin(PI/mid)});
        for(int i = 0; i < tot; i += mid * 2){
            Complex wk = Complex({1, 0});
            for(int j = 0; j < mid; j ++, wk = wk * w1){
                Complex x = a[i + j], y = wk * a[i + j + mid];
                a[i + j] = x + y, a[i + j + mid] = x - y;
            }
        }
    }
}
int n, m;
double A[N], B[N];
char s1[N], s2[N];
int main(){
    ios::sync_with_stdio(0); cin.tie(0); cout.tie(0);
    cin>>m>>n>>s1>>s2;
    while((1 << bit) < n + m + 1) bit++;
    tot = 1 << bit;
    for(int i = 0; i < tot; i++) rev[i] = (rev[i >> 1] >> 1) | ((i & 1) << (bit - 1));
    for(int i = 0; i < m; i++) if(s1[i] != '*') A[i] = s1[i] - 'a' + 1;
    for(int i = 0; i < n; i++) if(s2[i] != '*') B[i] = s2[i] - 'a' + 1;
    reverse(A, A + m);

    for(int i = 0; i < m; i++) a[i] = {A[i]*A[i]*A[i], 0};
    for(int i = 0; i < n; i++) b[i] = {B[i], 0};
    fft(a, 1), fft(b, 1);
    for(int i = 0; i < tot; i++) p[i] = p[i] + a[i] * b[i];
    for(int i = 0; i < tot; i++) a[i] = {0, 0}, b[i] = {0, 0};

    for(int i = 0; i < m; i++) a[i] = {A[i]*A[i], 0};
    for(int i = 0; i < n; i++) b[i] = {B[i]*B[i], 0};
    fft(a, 1), fft(b, 1);
    Complex f2 = {2, 0};
    for(int i = 0; i < tot; i++) p[i] = p[i] - a[i] * b[i] * f2;
    for(int i = 0; i < tot; i++) a[i] = {0, 0}, b[i] = {0, 0};

    for(int i = 0; i < m; i++) a[i] = {A[i], 0};
    for(int i = 0; i < n; i++) b[i] = {B[i]*B[i]*B[i], 0};
    fft(a, 1), fft(b, 1);
    for(int i = 0; i < tot; i++) p[i] = p[i] + a[i] * b[i];
    fft(p, -1);

    vector<int> v;
    for(int i = m - 1; i < n; i++)
        if(!(int)(p[i].x/(double)tot+0.5)) v.push_back(i-m+2);
    cout<<v.size()<<endl;
    for(auto i : v) cout<<i<<" ";
    return 0;
}