Test/Math/Convolution_and.test.cpp

• View this file on GitHub
• Last update: 2025-12-20 02:32:16+09:00
• Problem: https://judge.yosupo.jp/problem/bitwise_and_convolution

Depends on

• Library/Math/Bit.hpp
• Library/Math/Convolution.hpp
• Library/Math/EuclideanAlgorithm.hpp
• Library/Math/Math.hpp
• Library/Math/Mobius.hpp
• Library/Math/ModInt.hpp
• Library/Math/Zeta.hpp
• Library/Range/istream.hpp
• Library/Utility/Tuple.hpp
• Library/Utility/io.hpp

Code

#define PROBLEM "https://judge.yosupo.jp/problem/bitwise_and_convolution"
#include <iostream>

// begin:tag includes
#include "../../Library/Math/Convolution.hpp"
#include "../../Library/Math/ModInt.hpp"
#include "../../Library/Range/istream.hpp"
// end:tag includes

int main() {
  std::cin.tie(0);
  std::ios::sync_with_stdio(0);

  using mint = mtd::ModInt<998244353>;

  auto [n] = mtd::io::in<int>();
  auto a = mtd::io::_input<std::vector<mint>>(1LL << n);
  auto b = mtd::io::_input<std::vector<mint>>(1LL << n);

  auto ans = mtd::convolution::bitwise_and(a, b);
  for (auto x : ans) { std::cout << x << " "; }
  std::cout << std::endl;
}

#line 1 "Test/Math/Convolution_and.test.cpp"
#define PROBLEM "https://judge.yosupo.jp/problem/bitwise_and_convolution"
#include <iostream>

// begin:tag includes
#line 2 "Library/Math/Convolution.hpp"

#include <ranges>
#include <vector>

#line 2 "Library/Math/Mobius.hpp"

#line 5 "Library/Math/Mobius.hpp"

#line 2 "Library/Math/Bit.hpp"

namespace mtd {

  constexpr unsigned clz(unsigned int x) {
    int c = 0;
    if (x == 0) { return 0; }
    if (x & 0xffff0000) {
      x &= 0xffff0000;
      c |= 0x10;
    }
    if (x & 0xff00ff00) {
      x &= 0xff00ff00;
      c |= 0x08;
    }
    if (x & 0xf0f0f0f0) {
      x &= 0xf0f0f0f0;
      c |= 0x04;
    }
    if (x & 0xcccccccc) {
      x &= 0xcccccccc;
      c |= 0x02;
    }
    if (x & 0xaaaaaaaa) { c |= 0x01; }
    return c + 1;
  }

  constexpr unsigned ctz(unsigned int n) {
    if (!n) return -1;
    unsigned int c = 32;
    n &= -static_cast<signed int>(n);
    if (n) c--;
    if (n & 0x0000FFFF) c -= 16;
    if (n & 0x00FF00FF) c -= 8;
    if (n & 0x0F0F0F0F) c -= 4;
    if (n & 0x33333333) c -= 2;
    if (n & 0x55555555) c -= 1;
    return c;
  }

  constexpr unsigned long long popcount(unsigned long long x) {
    x = x - ((x >> 1) & 0x5555555555555555);
    x = (x & 0x3333333333333333) + ((x >> 2) & 0x3333333333333333);
    x = (x + (x >> 4)) & 0x0f0f0f0f0f0f0f0f;
    x = x + (x >> 8);
    x = x + (x >> 16);
    x = x + (x >> 32);
    return x & 0x0000007f;
  }

}  // namespace mtd
#line 7 "Library/Math/Mobius.hpp"

namespace mtd::mobius {

  template <class T>
  auto n(const std::vector<T>& a) {
    auto ret = a;
    for (auto i : std::views::iota(static_cast<size_t>(1), a.size())) {
      ret[i] = a[i] - a[i - 1];
    }
    return ret;
  }

  template <class T>
  auto bit_subset(const std::vector<T>& a) {
    auto ret = a;
    int size = clz(a.size());
    for (auto b : std::views::iota(0, size)) {
      for (auto bit : std::views::iota(0, 1LL << size)) {
        if (((bit >> b) & 1) && bit < a.size()) {
          ret[bit] -= ret[bit ^ (1LL << b)];
        }
      }
    }
    return ret;
  }

  template <class T>
  auto bit_supset(const std::vector<T>& a) {
    auto ret = a;
    int size = clz(a.size());
    for (auto b : std::views::iota(0, size)) {
      for (auto bit : std::views::iota(0, 1LL << size)) {
        if (((bit >> b) & 1) && bit < a.size()) {
          ret[bit ^ (1LL << b)] -= ret[bit];
        }
      }
    }
    return ret;
  }

}  // namespace mtd::mobius
#line 2 "Library/Math/Zeta.hpp"

#line 5 "Library/Math/Zeta.hpp"

#line 7 "Library/Math/Zeta.hpp"

namespace mtd::zeta {

  template <class T>
  auto n(const std::vector<T>& a) {
    auto ret = a;
    for (auto i : std::views::iota(static_cast<size_t>(1), a.size())) {
      ret[i] += ret[i - 1];
    }
    return ret;
  }

  template <class T>
  auto bit_subset(const std::vector<T>& a) {
    auto ret = a;
    int size = clz(a.size());
    for (auto b : std::views::iota(0, size)) {
      for (auto bit : std::views::iota(0, 1LL << size)) {
        if (((bit >> b) & 1) && bit < a.size()) {
          ret[bit] += ret[bit ^ (1LL << b)];
        }
      }
    }
    return ret;
  }

  template <class T>
  auto bit_supset(const std::vector<T>& a) {
    auto ret = a;
    int size = clz(a.size());
    for (auto b : std::views::iota(0, size)) {
      for (auto bit : std::views::iota(0, 1LL << size)) {
        if (((bit >> b) & 1) && bit < a.size()) {
          ret[bit ^ (1LL << b)] += ret[bit];
        }
      }
    }
    return ret;
  }
}  // namespace mtd::zeta
#line 8 "Library/Math/Convolution.hpp"

namespace mtd::convolution {

  template <class T>
  auto bitwise_and(const std::vector<T>& a, const std::vector<T>& b) {
    auto za = mtd::zeta::bit_supset(a);
    auto zb = mtd::zeta::bit_supset(b);
    auto zab = std::vector<T>();
    for (auto i : std::views::iota(static_cast<size_t>(0), a.size())) {
      zab.emplace_back(za[i] * zb[i]);
    }
    auto zma = mtd::mobius::bit_supset(za);
    return mtd::mobius::bit_supset(zab);
  }

}  // namespace mtd::convolution
#line 2 "Library/Math/ModInt.hpp"

#line 4 "Library/Math/ModInt.hpp"
#include <iterator>

#line 2 "Library/Math/Math.hpp"

#include <cmath>

#include <numeric>

#include <optional>

#line 7 "Library/Math/Math.hpp"
#include <unordered_map>

#line 9 "Library/Math/Math.hpp"

#line 2 "Library/Math/EuclideanAlgorithm.hpp"

#line 5 "Library/Math/EuclideanAlgorithm.hpp"
#include <tuple>


namespace mtd {

  class EuclideanAlgorithm {
    using T = long long;

    // 大きすぎるとオーバーフローしてしまう

    const static inline T m_mx = 1e9;

    const T m_a;
    const T m_b;
    const T m_c;

    T m_gcd;
    T m_x;
    T m_y;

    auto excludedEuclidAlgorithm(T a, T b) -> std::tuple<T, T, T> {
      if (a < 0) {
        auto [g, x, y] = excludedEuclidAlgorithm(-a, -b);
        return {g, -x, -y};
      }
      if (b == 0) { return {a, 1, 0}; }
      auto [g, y, x] = excludedEuclidAlgorithm(b, a % b);
      y -= a / b * x;
      return {g, x, y};
    }

    auto kRange(T x, T b, T l) const -> std::pair<T, T> {
      // x + b * k >= l を満たす k の範囲を求める

      T xd = (l - x);
      if (b == 0 && x >= l) { return {-m_mx, m_mx}; }
      if (b == 0 && x < l) { return {m_mx, -m_mx}; }
      if (b > 0 && xd < 0) { return {xd / b, m_mx}; }
      if (b > 0 && xd >= 0) { return {(xd + b - 1) / b, m_mx}; }
      if (b < 0 && xd < 0) { return {-m_mx, (-xd) / (-b)}; }
      if (b < 0 && xd >= 0) { return {-m_mx, -(xd - b - 1) / (-b)}; }
      return {m_mx, -m_mx};
    }

  public:
    auto debug() const {
      std::cout << m_a << " * " << m_x << " + " << m_b << " * " << m_y << " = "
                << m_c << std::endl;
      std::cout << "calc: " << m_a * m_x + m_b * m_y << " = " << m_c
                << std::endl;
    }

    EuclideanAlgorithm(T a, T b, T c) : m_a(a), m_b(b), m_c(c) {
      if (a == 0 && b == 0) { throw std::runtime_error(""); }
      auto [g, x, y] = excludedEuclidAlgorithm(a, b);
      if (c % g > 0) {
        throw std::runtime_error(
            "There is no solution to the equation. c must be divisible by "
            "gcd(a,b).");
      }
      m_gcd = g;
      m_x = c / g * x;
      m_y = c / g * y;
    }
    EuclideanAlgorithm(T a, T b) : EuclideanAlgorithm(a, b, std::gcd(a, b)) {}

    auto gcd() const { return m_gcd; }
    auto get(T x, T y) const { return m_a * x + m_b * y; }
    auto get(T k) const -> std::pair<T, T> {
      if (m_b == 0) { return {m_x, m_y - k}; }
      if (m_a == 0) { return {m_x + k, m_y}; }
      return {m_x + m_b * k, m_y - m_a * k};
    }
    // x>=x_lとなるようなkの範囲

    auto getMinX(T x_l = 0) const -> std::pair<T, T> {
      return kRange(m_x, m_b, x_l);
    }
    // y>=y_lとなるようなkの範囲

    auto getMinY(T y_l = 0) const -> std::pair<T, T> {
      return kRange(m_y, -1 * m_a, y_l);
    }
    // x>=x_l, y>=y_lとなるようなkの範囲

    auto getMin(T x_l = 0, T y_l = 0) const -> std::pair<T, T> {
      auto [xl, xr] = getMinX(x_l);
      auto [yl, yr] = getMinY(y_l);
      return {std::max(xl, yl), std::min(xr, yr)};
    }
  };

}  // namespace mtd

#line 11 "Library/Math/Math.hpp"

namespace mtd {
  template <class T>
  class Math {
    const std::vector<T> m_fac;
    const std::vector<T> m_finv;

    auto constructFac(int s) {
      std::vector<T> fac(s);
      fac[0] = fac[1] = 1;
      for (long long i = 2; i < s; ++i) { fac[i] = fac[i - 1] * i; }
      return fac;
    }
    auto constructInv(int s) {
      std::vector<T> finv(s);
      finv[s - 1] = 1 / m_fac[s - 1];
      for (long long i = s - 2; i >= 0; --i) {
        finv[i] = finv[i + 1] * (i + 1);
      }
      return finv;
    }

  public:
    constexpr Math(int size = 3 * static_cast<int>(1e6))
        : m_fac(constructFac(size)), m_finv(constructInv(size)) {}

    /* O(log b) */
    static constexpr T pow(T a, long long b) {
      T ans = 1;
      while (b > 0) {
        if (b & 1) { ans *= a; }
        b >>= 1;
        a *= a;
      }
      return ans;
    }

    /* O(log mod) */
    template <class S>
    static constexpr std::optional<long long> log(S x, S y, S mod) {
      x %= mod;
      y %= mod;

      if (mod == 1) { return 0; }
      if (x == 0 && y == 0) { return 1; }
      if (x == 0 && y == 1) { return 0; }
      if (x == 0) { return std::nullopt; }
      if (y == 1) { return 0; }

      if (auto g = std::gcd(x, mod); g > 1) {
        if (y % g) { return std::nullopt; }
        auto nx = x / g;
        auto nmod = mod / g;
        auto ea = mtd::EuclideanAlgorithm(nx, -nmod, 1);
        auto [t, _] = ea.getMinX();
        auto [nx_inv, __] = ea.get(t);
        nx_inv %= nmod;
        if (auto ans = log(x, y / g * nx_inv, nmod); ans) {
          return ans.value() + 1;
        } else {
          return ans;
        }
      }

      auto s = static_cast<S>(std::sqrt(mod));
      S xe = y;
      std::unordered_map<S, S> map;
      map.reserve(s);
      for (auto i : std::views::iota(0, s)) {
        (xe *= x) %= mod;
        map[xe] = i + 1;
      }

      S xs = 1;
      for ([[maybe_unused]] auto _ : std::views::iota(0, s)) {
        (xs *= x) %= mod;
      }
      S xse = 1;
      for (auto i : std::views::iota(0, mod / s + 5)) {
        (xse *= xs) %= mod;
        if (map.contains(xse)) { return s * (i + 1) - map[xse]; }
      }
      return std::nullopt;
    }
    constexpr std::optional<long long> log(long long x,
                                           long long y) requires requires {
      typename T::value_type;
      T::mod();
    }
    { return log(x, y, T::mod()); }

    constexpr auto fact(int n) const { return (n < 0) ? 0 : m_fac[n]; }
    constexpr auto factInv(int n) const { return (n < 0 ? 0 : m_finv[n]); }
    constexpr auto comb(int n, int r) const {
      return fact(n) * factInv(r) * factInv(n - r);
    }
    constexpr auto perm(int n, int r) const { return fact(n) * factInv(n - r); }
  };
}  // namespace mtd

#line 7 "Library/Math/ModInt.hpp"

namespace mtd {

  template <int MOD, class T = long long>
  class ModInt {
  public:
    using value_type = T;
    T x;

    constexpr ModInt(T _x) : x(_x >= 0 ? _x % MOD : MOD + (_x % MOD)) {}
    constexpr ModInt() : ModInt(0) {}

    // 四則演算
    constexpr auto& operator+=(const ModInt<MOD, T>& m) {
      x += m.x;
      if (x >= MOD) { x -= MOD; }
      return *this;
    }
    constexpr auto& operator-=(const ModInt<MOD, T>& m) {
      x -= m.x;
      if (x < 0) { x += MOD; }
      return *this;
    }
    constexpr auto& operator*=(const ModInt<MOD, T>& m) {
      x *= m.x;
      if (x >= MOD) { x %= MOD; }
      return *this;
    }
    constexpr auto& operator/=(const ModInt<MOD, T>& m) {
      x *= mtd::Math<ModInt<MOD, T>>::pow(m.x, MOD - 2).x;
      if (x >= MOD) { x %= MOD; }
      return *this;
    }

    constexpr auto operator+(const ModInt<MOD, T>& m) const {
      auto t = *this;
      t += m;
      return t;
    }
    constexpr auto operator-(const ModInt<MOD, T>& m) const {
      auto t = *this;
      t -= m;
      return t;
    }
    constexpr auto operator*(const ModInt<MOD, T>& m) const {
      auto t = *this;
      t *= m;
      return t;
    }
    constexpr auto operator/(const ModInt<MOD, T>& m) const {
      auto t = *this;
      t /= m;
      return t;
    }

    constexpr auto& operator+=(const T& t) {
      return *this += ModInt<MOD, T>(t);
    }
    constexpr auto& operator-=(const T& t) {
      return *this -= ModInt<MOD, T>(t);
    }
    constexpr auto& operator*=(const T& n) {
      return *this *= ModInt<MOD, T>(n);
    }
    constexpr auto& operator/=(const T& n) {
      return *this /= ModInt<MOD, T>(n);
    }
    constexpr auto operator+(const T& t) const {
      return *this + ModInt<MOD, T>(t);
    }
    constexpr auto operator-(const T& t) const {
      return *this - ModInt<MOD, T>(t);
    }
    constexpr auto operator*(const T& t) const {
      return *this * ModInt<MOD, T>(t);
    }
    constexpr auto operator/(const T& t) const {
      return *this / ModInt<MOD, T>(t);
    }
    constexpr friend auto operator+(const T& t, const ModInt<MOD, T>& m) {
      return m + t;
    }
    constexpr friend auto operator-(const T& t, const ModInt<MOD, T>& m) {
      return -m + t;
    }
    constexpr friend auto operator*(const T& t, const ModInt<MOD, T>& m) {
      return m * t;
    }
    constexpr friend auto operator/(const T& t, const ModInt<MOD, T>& m) {
      return ModInt<MOD, T>(1) / m * t;
    }

    // 単項演算
    constexpr auto operator-() const { return ModInt<MOD, T>(0 - x); }

    // 比較演算
    constexpr auto operator!=(const ModInt<MOD, T>& m) const {
      return x != m.x;
    }
    constexpr auto operator==(const ModInt<MOD, T>& m) const {
      return !(x != m.x);
    }

    // 入出力
    constexpr friend std::ostream& operator<<(std::ostream& os,
                                              const ModInt<MOD, T>& m) {
      return os << m.x;
    }
    constexpr friend std::istream& operator>>(std::istream& is,
                                              ModInt<MOD, T>& m) {
      return is >> m.x;
    }

    constexpr auto val() const { return x; }
    static constexpr auto mod() { return MOD; }
  };

}  // namespace mtd
#line 2 "Library/Range/istream.hpp"

#line 4 "Library/Range/istream.hpp"

#line 2 "Library/Utility/io.hpp"

#line 5 "Library/Utility/io.hpp"
#include <type_traits>
#line 7 "Library/Utility/io.hpp"

#line 2 "Library/Utility/Tuple.hpp"

#include <functional>

namespace mtd {
  namespace util {
    template <class F, class T>
    constexpr auto tuple_transform(F&& f, T&& t) {
      return std::apply(
          [&]<class... Ts>(Ts&&... elems) {
            return std::tuple<std::invoke_result_t<F&, Ts>...>(
                std::invoke(f, std::forward<Ts>(elems))...);
          },
          std::forward<T>(t));
    }
    template <class F, class T>
    constexpr auto tuple_for_each(F&& f, T&& t) {
      std::apply(
          [&]<class... Ts>(Ts&&... elems) {
            (std::invoke(f, std::forward<Ts>(elems)), ...);
          },
          std::forward<T>(t));
    }
  }  // namespace util
}  // namespace mtd
#line 9 "Library/Utility/io.hpp"

namespace mtd {
  namespace io {

    namespace __details {
      template <typename T>
      concept is_vec = std::same_as<
          T, std::vector<typename T::value_type, typename T::allocator_type>>;
      template <typename T>
      concept is_mat = is_vec<T> && is_vec<typename T::value_type>;

    }  // namespace __details

    template <class T>
    constexpr auto _input() {
      T x;
      std::cin >> x;
      return x;
    }
    template <typename T>
    requires requires { typename std::tuple_size<T>::type; }
    constexpr auto _input() {
      T x;
      util::tuple_for_each([](auto&& i) { std::cin >> i; }, x);
      return x;
    }
    template <__details::is_vec T>
    constexpr auto _input(int n) {
      std::vector<typename T::value_type> v;
      v.reserve(n);
      for (auto i : std::views::iota(0, n)) {
        v.emplace_back(_input<typename T::value_type>());
      }
      return v;
    }
    template <__details::is_mat T>
    constexpr auto _input(int h, int w) {
      T mat;
      mat.reserve(h);
      for (auto i : std::views::iota(0, h)) {
        mat.emplace_back(_input<typename T::value_type>(w));
      }
      return mat;
    }

    template <int N, class Tuple, class T, class... Args, class... Sizes>
    constexpr auto _tuple_input(Tuple& t, Sizes... sizes);
    template <int N, class Tuple, __details::is_vec T, class... Args,
              class Size, class... Sizes>
    constexpr auto _tuple_input(Tuple& t, Size size, Sizes... sizes);
    template <int N, class Tuple, __details::is_mat T, class... Args,
              class Size, class... Sizes>
    constexpr auto _tuple_input(Tuple& t, Size size_h, Size size_w,
                                Sizes... sizes);

    template <int N, class Tuple, class T, class... Args, class... Sizes>
    constexpr auto _tuple_input(Tuple& t, Sizes... sizes) {
      std::get<N>(t) = _input<T>();
      if constexpr (sizeof...(Args) > 0) {
        _tuple_input<N + 1, Tuple, Args...>(t, sizes...);
      }
    }
    template <int N, class Tuple, __details::is_vec T, class... Args,
              class Size, class... Sizes>
    constexpr auto _tuple_input(Tuple& t, Size size, Sizes... sizes) {
      std::get<N>(t) = _input<T>(size);
      if constexpr (sizeof...(Args) > 0) {
        _tuple_input<N + 1, Tuple, Args...>(t, sizes...);
      }
    }
    template <int N, class Tuple, __details::is_mat T, class... Args,
              class Size, class... Sizes>
    constexpr auto _tuple_input(Tuple& t, Size size_h, Size size_w,
                                Sizes... sizes) {
      std::get<N>(t) = _input<T>(size_h, size_w);
      if constexpr (sizeof...(Args) > 0) {
        _tuple_input<N + 1, Tuple, Args...>(t, sizes...);
      }
    }

    template <class... Args, class... Sizes>
    requires(std::convertible_to<Sizes, size_t>&&...) constexpr auto in(
        Sizes... sizes) {
      auto base = std::tuple<Args...>();
      _tuple_input<0, decltype(base), Args...>(base, sizes...);
      return base;
    }

  }  // namespace io

}  // namespace mtd
#line 6 "Library/Range/istream.hpp"

namespace mtd {
  namespace ranges {

    constexpr int _inf = 1e9;

    template <class... Args>
    struct istream_view
        : public std::ranges::view_interface<istream_view<Args...>> {
      class iterator {
        int count;
        std::tuple<Args...> val;

      public:
        using difference_type = int;
        using value_type = decltype(val);
        using iterator_concept = std::input_iterator_tag;

        constexpr iterator() = default;
        constexpr explicit iterator(int _count) : count(_count) {
          operator++();
        }

        constexpr auto operator*() const { return val; }
        constexpr auto& operator++() {
          --count;
          if (count >= 0) { val = io::in<Args...>(); }
          return *this;
        }
        constexpr auto operator++(int) { return ++*this; }

        constexpr auto operator==(const iterator& s) const {
          return count == s.count;
        }
        constexpr auto operator==(std::default_sentinel_t) const {
          return count < 0 || std::cin.eof() || std::cin.fail() ||
                 std::cin.bad();
        }
        constexpr friend auto operator==(std::default_sentinel_t s,
                                         const iterator& li) {
          return li == s;
        }
      };

      int count;

    public:
      constexpr explicit istream_view(int _count) : count(_count) {}
      constexpr explicit istream_view() : istream_view(_inf) {}
      constexpr auto begin() const { return iterator(count); }
      constexpr auto end() const { return std::default_sentinel; }
    };
  }  // namespace ranges

  namespace views {
    namespace __detail {
      template <typename... _Args>
      concept __can_istream_view = requires {
        ranges::istream_view(std::declval<_Args>()...);
      };
    }  // namespace __detail

    template <class... Args>
    struct _Istream {
      template <class... _Tp>
      requires __detail::__can_istream_view<_Tp...>
      constexpr auto operator() [[nodiscard]] (_Tp&&... __e) const {
        return ranges::istream_view<Args...>(std::forward<_Tp>(__e)...);
      }
    };

    template <class... Args>
    inline constexpr _Istream<Args...> istream{};
  }  // namespace views

}  // namespace mtd
#line 8 "Test/Math/Convolution_and.test.cpp"
// end:tag includes

int main() {
  std::cin.tie(0);
  std::ios::sync_with_stdio(0);

  using mint = mtd::ModInt<998244353>;

  auto [n] = mtd::io::in<int>();
  auto a = mtd::io::_input<std::vector<mint>>(1LL << n);
  auto b = mtd::io::_input<std::vector<mint>>(1LL << n);

  auto ans = mtd::convolution::bitwise_and(a, b);
  for (auto x : ans) { std::cout << x << " "; }
  std::cout << std::endl;
}

Back to top page