代码

#include <algorithm>
#include <math.h>
using namespace std;
const int N = 4e6 + 5;
const double PI = acos(-1.0);
struct Complex {
  typedef double value_t;
  value_t real, imag;
  Complex() = default;
  Complex(value_t a, value_t b) : real(a), imag(b) {}
  Complex(const Complex &rhs) : real(rhs.real), imag(rhs.imag) {}
  ~Complex() = default;
  Complex conj() const { return {real, -imag}; }
  Complex &operator=(const Complex &rhs) { return real = rhs.real, imag = rhs.imag, *this; }
  Complex &operator+=(const Complex &rhs) { return real += rhs.real, imag += rhs.imag, *this; }
  Complex &operator-=(const Complex &rhs) { return real -= rhs.real, imag -= rhs.imag, *this; }
  Complex &operator*=(const Complex &rhs) {
    value_t tmp1 = real * rhs.real - imag * rhs.imag, tmp2 = real * rhs.imag + imag * rhs.real;
    return real = tmp1, imag = tmp2, *this;
  }
  Complex &operator/=(const Complex &rhs) {
    value_t 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);
    return real = tmp1, imag = tmp2, *this;
  }
  friend Complex operator+(const Complex &lhs, const Complex &rhs) { return Complex(lhs) += rhs; }
  friend Complex operator-(const Complex &lhs, const Complex &rhs) { return Complex(lhs) -= rhs; }
  friend Complex operator*(const Complex &lhs, const Complex &rhs) { return Complex(lhs) *= rhs; }
  friend Complex operator/(const Complex &lhs, const Complex &rhs) { return Complex(lhs) /= 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] = {(double)a[i], (double)b[i]};
  int len = 1;
  while (len < n + m + 1) 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;
}

评测结果