代码

#include <algorithm>
#include <math.h>
using namespace std;
const int N = 4e6 + 5;
const double PI = acos(-1.0);
struct Complex {
    double real, imag;
    Complex(double a = 0, double b = 0) : real(a), imag(b) {}
    Complex(const Complex &rhs) : real(rhs.real), imag(rhs.imag) {}
    Complex conj() const { return {real, -imag}; }
    Complex &operator=(const Complex &rhs) {
        real = rhs.real, imag = rhs.imag;
        return *this;
    }
    Complex &operator+=(const Complex &rhs) {
        real += rhs.real, imag += rhs.imag;
        return *this;
    }
    Complex &operator-=(const Complex &rhs) {
        real -= rhs.real, imag -= rhs.imag;
        return *this;
    }
    Complex &operator*=(const Complex &rhs) {
        double tmp1 = real * rhs.real - imag * rhs.imag, tmp2 = real * rhs.imag + imag * rhs.real;
        real = tmp1, imag = tmp2;
        return *this;
    }
    Complex &operator/=(const Complex &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 Complex operator+(const Complex &lhs, const Complex &rhs) {
        Complex res(lhs);
        return res += rhs;
    }
    friend Complex operator-(const Complex &lhs, const Complex &rhs) {
        Complex res(lhs);
        return res -= rhs;
    }
    friend Complex operator*(const Complex &lhs, const Complex &rhs) {
        Complex res(lhs);
        return res *= rhs;
    }
    friend Complex operator/(const Complex &lhs, const Complex &rhs) {
        Complex res(lhs);
        return res /= rhs;
    }
} eps[N], ans[N];
void init(int n) {
    int l = n >> 1;
    Complex w = {cos(PI / l), sin(PI / l)};
    eps[l] = {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, Complex x[]) { // dit
    for (int i = 0; i < n; ++i) x[i].imag = -x[i].imag;
    Complex 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, Complex x[]) { // dif
    Complex 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;
}

评测结果