代码

#include <algorithm>
#include <math.h>
using namespace std;
const int N = 4e6 + 5;
const double PI = acos(-1.0);
struct Vec {
    double real, imag;
    Vec(double a = 0, double b = 0) : real(a), imag(b) {}
    Vec(const Vec &rhs) : real(rhs.real), imag(rhs.imag) {}
    Vec conj() const { return {real, -imag}; }
    Vec &operator=(const Vec &rhs) {
        real = rhs.real, imag = rhs.imag;
        return *this;
    }
    Vec &operator+=(const Vec &rhs) {
        real += rhs.real, imag += rhs.imag;
        return *this;
    }
    Vec &operator-=(const Vec &rhs) {
        real -= rhs.real, imag -= rhs.imag;
        return *this;
    }
    Vec &operator*=(const Vec &rhs) {
        double tmp1 = real * rhs.real - imag * rhs.imag, tmp2 = real * rhs.imag + imag * rhs.real;
        real = tmp1, imag = tmp2;
        return *this;
    }
    Vec &operator/=(const Vec &rhs) {
        double tmp1 = (real * rhs.real + imag * rhs.imag) /
                      (rhs.real * rhs.real + rhs.imag * rhs.imag),
               tmp2 = (imag * rhs.real - real * rhs.imag) /
                      (rhs.real * rhs.real + rhs.imag * rhs.imag);
        real = tmp1, imag = tmp2;
        return *this;
    }
    friend Vec operator+(const Vec &lhs, const Vec &rhs) {
        Vec res(lhs);
        return res += rhs;
    }
    friend Vec operator-(const Vec &lhs, const Vec &rhs) {
        Vec res(lhs);
        return res -= rhs;
    }
    friend Vec operator*(const Vec &lhs, const Vec &rhs) {
        Vec res(lhs);
        return res *= rhs;
    }
    friend Vec operator/(const Vec &lhs, const Vec &rhs) {
        Vec res(lhs);
        return res /= rhs;
    }
} eps[N], ans[N];
void init(int n) {
    int l = n >> 1;
    Vec w(cos(PI / l), sin(PI / l));
    eps[0] = {1, 0};
    for (int i = l + 1; i < n; ++i) eps[i] = eps[i - 1] * w;
    for (int i = l - 1; i > 0; --i) eps[i] = eps[i << 1];
}
void idft(int n, Vec x[]) { // dit
    for (int i = 0; i < n; ++i) x[i].imag = -x[i].imag;
    Vec u, v;
    for (int j = 0; j < n; j += 2) {
        u = x[j], v = x[j + 1];
        x[j] = u + v;
        x[j + 1] = u - v;
    }
    for (int i = 4; i <= n; i <<= 1) {
        for (int j = 0, l = i >> 1; j < n; j += i) {
            for (int k = 0; k < l; ++k) {
                u = x[k + j], v = eps[l + k] * x[k + j + l];
                x[k + j] = u + v;
                x[k + j + l] = u - v;
            }
        }
    }
    double c = 1.0 / n;
    for (int i = 0; i < n; ++i) x[i].real *= c, x[i].imag *= -c;
}
void dft(int n, Vec x[]) { // dif
    Vec u, v;
    for (int i = n; i >= 4; i >>= 1) {
        for (int j = 0, l = i >> 1; j < n; j += i) {
            for (int k = 0; k < l; ++k) {
                u = x[k + j], v = x[k + j + l];
                x[k + j] = u + v;
                x[k + j + l] = (u - v) * eps[l + k];
            }
        }
    }
    for (int j = 0; j < n; j += 2) {
        u = x[j], v = x[j + 1];
        x[j] = u + v;
        x[j + 1] = u - v;
    }
}
void poly_multiply(unsigned *a, int n, unsigned *b, int m, unsigned *c) {
    for (int i = 0; i <= n || i <= m; ++i) ans[i] = {a[i], b[i]};
    int len = 1;
    while (len < n + m + 2) len <<= 1;
    init(len);
    dft(len, ans);
    for (int i = 0; i < len; ++i) ans[i] *= ans[i];
    idft(len, ans);
    for (int i = 0; i <= n + m; ++i) c[i] = ans[i].imag / 2.0 + 0.5;
}

评测结果