Convolution (Index multiplication mod P)
(Mylib/Convolution/convolution_multiply.cpp)

• View this file on GitHub
• Last update: 2021-04-23 23:44:44+09:00
• Link: View error logs on GitHub Actions

Operations

• convolution_multiply(f, g, P)
  • Pは素数
  • $h_k = \sum_{k = i * j \pmod{P}} f_i * g_j$

Requirements

Notes

Problems

References

• https://drken1215.hatenablog.com/entry/2020/11/04/172300

Depends on

• Mod pow (Mylib/Number/Mod/mod_pow.cpp)
• Primitive root (Mylib/Number/Prime/primitive_root.cpp)

Verified with

• test/yukicoder/931/main.test.cpp

Code

#pragma once
#include <vector>
#include "Mylib/Number/Prime/primitive_root.cpp"

namespace haar_lib {
  template <typename T, const auto &convolve>
  std::vector<T> convolution_multiply(const std::vector<T> &A, const std::vector<T> &B, int P) {
    const int p_root = primitive_root(P);

    std::vector<int> index(P);

    {
      int64_t s = 1;

      for (int i = 0; i < P; ++i) {
        index[s] = i;
        s *= p_root;
        if (s >= P) s %= P;
      }
    }

    std::vector<T> a(P), b(P);
    for (int i = 0; i < (int) A.size(); ++i) a[index[i % P]] = A[i];
    for (int i = 0; i < (int) B.size(); ++i) b[index[i % P]] = B[i];

    auto c = convolve(a, b);
    std::vector<T> ret(P);

    {
      int64_t s = 1;

      for (auto x : c) {
        ret[s] += x;
        s *= p_root;
        if (s >= P) s %= P;
      }
    }

    return ret;
  }
};  // namespace haar_lib

#line 2 "Mylib/Convolution/convolution_multiply.cpp"
#include <vector>
#line 2 "Mylib/Number/Mod/mod_pow.cpp"
#include <cstdint>

namespace haar_lib {
  constexpr int64_t mod_pow(int64_t n, int64_t p, int64_t m) {
    int64_t ret = 1;
    while (p > 0) {
      if (p & 1) (ret *= n) %= m;
      (n *= n) %= m;
      p >>= 1;
    }
    return ret;
  }
}  // namespace haar_lib
#line 3 "Mylib/Number/Prime/primitive_root.cpp"

namespace haar_lib {
  constexpr int primitive_root(int p) {
    int pf[30] = {};
    int k      = 0;
    {
      int n = p - 1;
      for (int64_t i = 2; i * i <= p; ++i) {
        if (n % i == 0) {
          pf[k++] = i;
          while (n % i == 0) n /= i;
        }
      }
      if (n != 1)
        pf[k++] = n;
    }

    for (int g = 2; g <= p; ++g) {
      bool ok = true;
      for (int i = 0; i < k; ++i) {
        if (mod_pow(g, (p - 1) / pf[i], p) == 1) {
          ok = false;
          break;
        }
      }

      if (not ok) continue;

      return g;
    }
    return -1;
  }
}  // namespace haar_lib
#line 4 "Mylib/Convolution/convolution_multiply.cpp"

namespace haar_lib {
  template <typename T, const auto &convolve>
  std::vector<T> convolution_multiply(const std::vector<T> &A, const std::vector<T> &B, int P) {
    const int p_root = primitive_root(P);

    std::vector<int> index(P);

    {
      int64_t s = 1;

      for (int i = 0; i < P; ++i) {
        index[s] = i;
        s *= p_root;
        if (s >= P) s %= P;
      }
    }

    std::vector<T> a(P), b(P);
    for (int i = 0; i < (int) A.size(); ++i) a[index[i % P]] = A[i];
    for (int i = 0; i < (int) B.size(); ++i) b[index[i % P]] = B[i];

    auto c = convolve(a, b);
    std::vector<T> ret(P);

    {
      int64_t s = 1;

      for (auto x : c) {
        ret[s] += x;
        s *= p_root;
        if (s >= P) s %= P;
      }
    }

    return ret;
  }
};  // namespace haar_lib

Back to top page