algorithm/or_convolution.cpp

View this file on GitHub
Last update: 2022-01-28 13:07:07+09:00

Depends on

Code

#pragma once
#include <cassert>
#include <vector>
#include "../utility/rep.cpp"
#include "subset_mobius_transform.cpp"
#include "subset_zeta_transform.cpp"

template <class R>
std::vector<typename R::Type> or_convolution(std::vector<typename R::Type> a, std::vector<typename R::Type> b) {
    assert(a.size() == b.size());
    subset_zeta_transform<typename R::Sum>(a);
    subset_zeta_transform<typename R::Sum>(b);
    for (const int i : rep(a.size())) a[i] = R::Product::operation(a[i], b[i]);
    subset_mobius_transform<typename R::Sum>(a);
    return a;
}

#line 2 "algorithm/or_convolution.cpp"
#include <cassert>
#include <vector>
#line 2 "utility/rep.cpp"
#include <algorithm>

class Range {
    struct Iter {
        int itr;
        constexpr Iter(const int pos) noexcept : itr(pos) {}
        constexpr void operator++() noexcept { ++itr; }
        constexpr bool operator!=(const Iter& other) const noexcept { return itr != other.itr; }
        constexpr int operator*() const noexcept { return itr; }
    };
    const Iter first, last;

  public:
    explicit constexpr Range(const int first, const int last) noexcept : first(first), last(std::max(first, last)) {}
    constexpr Iter begin() const noexcept { return first; }
    constexpr Iter end() const noexcept { return last; }
};

constexpr Range rep(const int l, const int r) noexcept { return Range(l, r); }
constexpr Range rep(const int n) noexcept { return Range(0, n); }
#line 5 "algorithm/subset_mobius_transform.cpp"

template <class G> void subset_mobius_transform(std::vector<typename G::Type>& f) {
    const int n = f.size();
    assert((n & (n - 1)) == 0);
    for (int i = n; i >>= 1;)
        for (const int j : rep(n))
            if (j & i) f[j] = G::operation(G::inverse(f[j & ~i]), f[j]);
}
#line 5 "algorithm/subset_zeta_transform.cpp"

template <class S> void subset_zeta_transform(std::vector<typename S::Type>& f) {
    const int n = f.size();
    assert((n & (n - 1)) == 0);
    for (int i = 1; i < n; i <<= 1)
        for (const int j : rep(n))
            if (j & i) f[j] = S::operation(f[j & ~i], f[j]);
}
#line 7 "algorithm/or_convolution.cpp"

template <class R>
std::vector<typename R::Type> or_convolution(std::vector<typename R::Type> a, std::vector<typename R::Type> b) {
    assert(a.size() == b.size());
    subset_zeta_transform<typename R::Sum>(a);
    subset_zeta_transform<typename R::Sum>(b);
    for (const int i : rep(a.size())) a[i] = R::Product::operation(a[i], b[i]);
    subset_mobius_transform<typename R::Sum>(a);
    return a;
}