P4173 残缺的字符串(FFT字符串匹配)

P4173 残缺的字符串(FFT字符串匹配)
P4173

解题思路：

经典套路将模式串翻转,将*设为0，设以目标串的x位置匹配结束的匹配函数为(P(x)=sum^{m-1}_{i=0}[A(m-1-i)-B(x-(m-1-i))]^2A(m-1-i)B(x-(m-1-i))]),展开之后化简为(P(x)=sum_{i+j=x}A^3(i)B(j)-2sum_{i+j=x}A^2(i)B^2(j)+sum_{i+j=x}A(i)B^3(j))
做三次FFT即可,然后交题就出了一堆玄学错误

#include <bits/stdc++.h>
using namespace std;
/*    freopen("k.in", "r", stdin);
    freopen("k.out", "w", stdout); */
//clock_t c1 = clock();
//std::cerr << "Time:" << clock() - c1 <<"ms" << std::endl;
//#pragma comment(linker, "/STACK:1024000000,1024000000")
#define de(a) cout << #a << " = " << a << endl
#define rep(i, a, n) for (int i = a; i <= n; i++)
#define per(i, a, n) for (int i = n; i >= a; i--)
typedef long long ll;
typedef unsigned long long ull;
typedef pair<int, int> PII;
typedef pair<double, double> PDD;
typedef vector<int, int> VII;
#define inf 0x3f3f3f3f
const ll INF = 0x3f3f3f3f3f3f3f3f;
const ll MAXN = 2e6 + 10;
const ll MAXM = 5e6 + 7;
const ll MOD = 1e9 + 7;
const double eps = 1e-6;
const double pi = acos(-1.0);
struct Complex
{
    double x, y;
    Complex(double xx = 0, double yy = 0) { x = xx, y = yy; }
} a[MAXN], b[MAXN], c[MAXN];
Complex operator+(Complex a, Complex b) { return Complex(a.x + b.x, a.y + b.y); }
Complex operator-(Complex a, Complex b) { return Complex(a.x - b.x, a.y - b.y); }
Complex operator*(Complex a, Complex b) { return Complex(a.x * b.x - a.y * b.y, a.x * b.y + a.y * b.x); } //不懂的看复数的运算那部分
int l = 0, r[MAXN];
int limit = 1;
void FFT(Complex *A, int type)
{
    for (int i = 0; i < limit; i++)
        if (i < r[i])
            swap(A[i], A[r[i]]); //求出要迭代的序列
    for (int mid = 1; mid < limit; mid <<= 1)
    {                                                    //待合并区间的长度的一半
        Complex Wn(cos(pi / mid), type * sin(pi / mid)); //单位根
        for (int R = mid << 1, j = 0; j < limit; j += R)
        {                    //R是区间的长度，j表示前已经到哪个位置了
            Complex w(1, 0); //幂
            for (int k = 0; k < mid; k++, w = w * Wn)
            {                                                 //枚举左半部分
                Complex x = A[j + k], y = w * A[j + mid + k]; //蝴蝶效应
                A[j + k] = x + y;
                A[j + mid + k] = x - y;
            }
        }
    }
   /*  if (type == -1)
        for (int i = 0; i < limit; i++)
            a[i].x /= limit; */
}
char s[MAXN], t[MAXN];
int ta[MAXN] = {0}, tb[MAXN] = {0};
int ans[MAXN] = {0};
int main()
{
    int n, m;
    scanf("%d%d%s%s", &m, &n, t, s);
    for (int i = 0; i < m; i++)
        ta[m - i - 1] = t[i] == '*' ? 0 : (t[i] - 'a' + 1);
    for (int i = 0; i < n; i++)
        tb[i] = s[i] == '*' ? 0 : (s[i] - 'a' + 1);
    while (limit <= m + n)
        limit <<= 1, l++;
    for (int i = 0; i < limit; i++)
        r[i] = (r[i >> 1] >> 1) | ((i & 1) << (l - 1));
    for (int i = 0; i < limit; i++)
    {
        a[i] = Complex(ta[i], 0);
        b[i] = Complex(tb[i] * tb[i] * tb[i], 0);
    }
    FFT(a, 1), FFT(b, 1);
    for (int i = 0; i < limit; i++)
        c[i] = c[i] + a[i] * b[i];
    for (int i = 0; i < limit; i++)
    {
        a[i] = Complex(ta[i] * ta[i], 0);
        b[i] = Complex(tb[i] * tb[i], 0);
    }
    FFT(a, 1), FFT(b, 1);
    for (int i = 0; i < limit; i++)
        c[i] = c[i] - a[i] * b[i] * Complex(2.0, 0);
    for (int i = 0; i < limit; i++)
    {
        a[i] = Complex(ta[i] * ta[i] * ta[i], 0);
        b[i] = Complex(tb[i], 0);
    }
    FFT(a, 1), FFT(b, 1);
    for (int i = 0; i < limit; i++)
        c[i] = c[i] + a[i] * b[i];
    FFT(c, -1);
    int cnt = 0;
    for (int i = m - 1; i < n; i++)
    {
        if (int(c[i].x / limit + 0.5) == 0)
            ans[cnt++] = i - m + 2;
    }
    printf("%d
", cnt);
    for (int i = 0; i < cnt; i++)
        printf("%d ", ans[i]);
    printf("
");
    return 0;
}