結果
問題 |
No.2507 Yet Another Subgraph Counting
|
ユーザー |
|
提出日時 | 2023-08-22 16:27:55 |
言語 | C++17(gcc12) (gcc 12.3.0 + boost 1.87.0) |
結果 |
AC
|
実行時間 | 234 ms / 2,000 ms |
コード長 | 8,738 bytes |
コンパイル時間 | 6,490 ms |
コンパイル使用メモリ | 221,568 KB |
最終ジャッジ日時 | 2025-02-16 12:18:54 |
ジャッジサーバーID (参考情報) |
judge3 / judge4 |
(要ログイン)
ファイルパターン | 結果 |
---|---|
other | AC * 52 |
ソースコード
#include <cassert> #include <iostream> #include <limits> #include <utility> #include <vector> #include <immintrin.h> namespace library { namespace bits { template <typename T, std::enable_if_t<std::is_integral_v<T>, std::nullptr_t> = nullptr> T kth_bit(T v, size_t k) { return (v >> k) & 1; } template <typename T, std::enable_if_t<std::is_integral_v<T>, std::nullptr_t> = nullptr> size_t bit_length(const T v) { if constexpr (std::numeric_limits<std::make_unsigned_t<T>>::digits <= 32) { return 32 - __builtin_clz(v); } else { return 64 - __builtin_clzll(v); } } template <typename T, std::enable_if_t<std::is_integral_v<T>, std::nullptr_t> = nullptr> size_t popcount(const T v) { if constexpr (std::numeric_limits<std::make_unsigned_t<T>>::digits <= 32) { return __builtin_popcount(v); } else { return __builtin_popcountll(v); } } template <typename T, std::enable_if_t<std::is_integral_v<T>, std::nullptr_t> = nullptr> size_t count_tz(const T v) { if constexpr (std::numeric_limits<std::make_unsigned_t<T>>::digits <= 32) { return __builtin_ctz(v); } else { return __builtin_ctzll(v); } } // See https://www.intel.com/content/www/us/en/docs/intrinsics-guide/index.html#text=_pdep&ig_expand=4939 template <typename T, std::enable_if_t<std::is_integral_v<T>, std::nullptr_t> = nullptr> __attribute__((target("bmi2"))) T pdep(const T src, const T mask) { /* T dst = 0; for (size_t i = 0, j = 0; i < BIT_NUM; ++i) { if (kth_bit(mask, i)) dst |= kth_bit(src, j) << i, ++j; } return dst; */ if constexpr (std::numeric_limits<std::make_unsigned_t<T>>::digits <= 32) { return _pdep_u32(src, mask); } else { return _pdep_u64(src, mask); } } } namespace subset_transform { template <typename T> void zeta(std::vector<T>& x) { const size_t n = x.size(); for (size_t b = 1; b < n; b *= 2) for (size_t l = 0; l < n; l += 2 * b) for (size_t i = l; i < l + b; ++i) { x[i + 1 * b] += x[i + 0 * b]; } } template <typename T> void mobius(std::vector<T>& x) { const size_t n = x.size(); for (size_t b = 1; b < n; b *= 2) for (size_t l = 0; l < n; l += 2 * b) for (size_t i = l; i < l + b; ++i) { x[i + 1 * b] -= x[i + 0 * b]; } } } namespace set_power_series { namespace details { template <typename T> struct polynomial : std::vector<T> { using std::vector<T>::vector; polynomial& operator+=(const polynomial& q) { for (size_t i = 0; i < q.size(); ++i) (*this)[i] += q[i]; return *this; } polynomial& operator-=(const polynomial& q) { for (size_t i = 0; i < q.size(); ++i) (*this)[i] -= q[i]; return *this; } polynomial& operator*=(const polynomial& q) { const size_t n = this->size(); polynomial r(n); for (size_t i = 0; i < n; ++i) for (size_t j = 0; i + j < n; ++j) r[i + j] += (*this)[i] * q[j]; return *this = std::move(r); } polynomial operator+(const polynomial& q) { polynomial p = *this; p += q; return p; } polynomial operator-(const polynomial& q) { polynomial p = *this; p -= q; return p; } polynomial operator*(const polynomial& q) { polynomial p = *this; p *= q; return p; } }; template <typename T> std::vector<polynomial<T>> add_rank(const std::vector<T>& a) { const size_t n = a.size(); std::vector fs(n, polynomial<T>(bits::count_tz(n) + 1, T{ 0 })); for (size_t i = 0; i < n; ++i) fs[i][bits::popcount(i)] = a[i]; return fs; } template <typename T> std::vector<T> remove_rank(const std::vector<polynomial<T>>& polys) { const size_t n = polys.size(); std::vector<T> a(n); for (size_t i = 0; i < n; ++i) a[i] = polys[i][bits::popcount(i)]; return a; } } template <typename T> std::vector<T> subset_convolution(const std::vector<T>& a, const std::vector<T>& b) { const size_t n = a.size(); auto ra = details::add_rank(a); auto rb = details::add_rank(b); subset_transform::zeta(ra); subset_transform::zeta(rb); for (size_t i = 0; i < n; ++i) ra[i] *= rb[i]; subset_transform::mobius(ra); return details::remove_rank(ra); } template <typename T> std::vector<T> subset_exp(const std::vector<T>& f) { assert(f[0] == 0); const size_t n = bits::bit_length(f.size()) - 1; std::vector<T> g{ 1 }; for (size_t i = 0; i < n; ++i) { std::vector<T> h = subset_convolution(g, std::vector<T>(f.begin() + (1 << i), f.begin() + (1 << (i + 1)))); std::move(h.begin(), h.end(), std::back_inserter(g)); } return g; } } } using vertex = size_t; using vertex_set = size_t; using edge = std::pair<vertex, vertex>; using Int = uint64_t; using set_power_series = std::vector<Int>; std::vector<Int> count_cycles(const size_t n, const size_t, const std::vector<edge>& edges) { // adjacency list std::vector adj(n, std::vector<vertex>{}); for (const auto& [u, v] : edges) adj[u].push_back(v), adj[v].push_back(u); // "c" mentioned in the editorial std::vector<Int> c(1u << n); // dp[S: vertex set][v: vertex] := # simple paths from min S to v passing vertices in S (but not passing vertices not in S) std::vector dp(1u << n, std::vector<Int>(n)); // base cases for (vertex v = 0; v < n; ++v) { dp[1u << v][v] = 1; } for (vertex_set S = 1; S < 1u << n; ++S) { // min S const vertex start = library::bits::count_tz(S); for (vertex cur = 0; cur < n; ++cur) for (const vertex nxt : adj[cur]) { if (start == nxt) { c[S] += dp[S][cur]; } else if (start < nxt and not library::bits::kth_bit(S, nxt)) { const vertex_set T = S | (1u << nxt); dp[T][nxt] += dp[S][cur]; } } } for (vertex_set S = 1; S < 1u << n; ++S) { const size_t card = library::bits::popcount(S); if (card == 1) c[S] = 1; if (card == 2) c[S] = 0; if (card >= 3) c[S] /= 2; } return c; } Int solve(const size_t n, const size_t m, const std::vector<edge>& edges) { // E[S: vertex set] := # edges connecting vertices in S. std::vector<Int> E(1u << n); for (const auto& [u, v] : edges) ++E[(1u << u) | (1u << v)]; library::subset_transform::zeta(E); // "c" mentioned in the editorial const set_power_series c = count_cycles(n, m, edges); // "f" mentioned in the editorial set_power_series f(1u << n); for (vertex_set C = 1; C < 1u << n; ++C) { // max C const vertex t = library::bits::bit_length(C) - 1; // {0, ..., t} - C const vertex_set S = ((1u << (t + 1)) - 1) ^ C; const size_t k = library::bits::popcount(S); // "g_C" mentioned in the editorial set_power_series g(1u << k); for (vertex_set A = 0; A < 1u << k; ++A) { // For more information about pdep, see https://www.intel.com/content/www/us/en/docs/intrinsics-guide/index.html#text=_pdep&ig_expand=4939. const vertex_set T = library::bits::pdep(A, S); g[A] = f[T] * (E[T | C] - E[T] - E[C]); } // "h_C" mentioned in the editorial const set_power_series h = library::set_power_series::subset_exp(g); for (vertex_set A = 0; A < 1u << k; ++A) { const vertex_set X = library::bits::pdep(A, S) | C; f[X] += c[C] * h[A]; } } return library::set_power_series::subset_exp(f).back(); } int main() { size_t n, m; std::cin >> n >> m; std::vector<edge> edges(m); for (auto& [u, v] : edges) { std::cin >> u >> v; --u, --v; } std::cout << solve(n, m, edges) << '\n'; }