結果
| 問題 |
No.2507 Yet Another Subgraph Counting
|
| コンテスト | |
| ユーザー |
|
| 提出日時 | 2023-08-31 21:26:30 |
| 言語 | C++17(gcc12) (gcc 12.3.0 + boost 1.87.0) |
| 結果 |
AC
|
| 実行時間 | 88 ms / 2,000 ms |
| コード長 | 9,336 bytes |
| コンパイル時間 | 11,023 ms |
| コンパイル使用メモリ | 250,448 KB |
| 最終ジャッジ日時 | 2025-02-16 15:51:14 |
|
ジャッジサーバーID (参考情報) |
judge4 / judge6 |
(要ログイン)
| ファイルパターン | 結果 |
|---|---|
| other | AC * 52 |
ソースコード
// TLE O(4^n n^2)
#pragma GCC target("avx2")
#pragma GCC optimize("O3")
#pragma GCC optimize("unroll-loops")
#include <algorithm>
#include <array>
#include <cassert>
#include <iostream>
#include <limits>
#include <utility>
#include <vector>
#ifdef _MSC_VER
# include <intrin.h>
#else
# include <x86intrin.h>
#endif
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 <size_t N>
void push_front(std::array<uint64_t, N> &a, uint64_t v) {
for (size_t i = N; i --> 1;) a[i] = a[i - 1];
a[0] = v;
}
template <size_t N>
void muleq(std::array<uint64_t, N>& p, const std::array<uint64_t, N>& q, size_t siz) {
for (size_t i = siz; i --> 0;) {
uint64_t val = 0;
for (size_t j = 0; j <= i; ++j) val += p[i - j] * q[j];
p[i] = val;
}
}
template <size_t N, typename InputIterator>
std::vector<std::array<uint64_t, N + 1>> add_rank(InputIterator first, InputIterator last) {
const size_t n = last - first;
std::vector<std::array<uint64_t, N + 1>> fs(n);
for (size_t i = 0; i < n; ++i) {
fs[i][bits::popcount(i)] = first[i];
}
return fs;
}
template <size_t N>
std::vector<std::array<uint64_t, N + 1>> add_rank(const std::vector<uint64_t> &a) {
return add_rank<N>(a.begin(), a.end());
}
template <size_t N>
std::vector<uint64_t> remove_rank(const std::vector<std::array<uint64_t, N + 1>>& polys) {
const size_t n = polys.size();
std::vector<uint64_t> a(n);
for (size_t i = 0; i < n; ++i) a[i] = polys[i][bits::popcount(i)];
return a;
}
}
template <size_t N>
std::vector<uint64_t> subset_exp(const std::vector<uint64_t>& f) {
assert(f[0] == 0);
const size_t n = bits::bit_length(f.size()) - 1;
auto rf = details::add_rank<N>({ 1 });
rf.reserve(size_t(1) << n);
for (size_t k = 0; k < n; ++k) {
auto rg = details::add_rank<N>(f.begin() + (1 << k), f.begin() + (1 << (k + 1)));
{
const size_t n = rg.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) {
const size_t s = i, t = i + b;
for (size_t q = 0; q <= k; ++q) {
rg[t][q] += rg[s][q];
}
}
}
for (size_t j = 0; j < size_t(1) << k; ++j) {
details::push_front(rg[j], 1);
details::muleq(rg[j], rf[j], k + 2);
rf.push_back(rg[j]);
}
}
{
const size_t k = n;
const size_t n = rf.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) {
const size_t s = i, t = i + b;
for (size_t q = 0; q <= k; ++q) {
rf[t][q] -= rf[s][q];
}
}
}
return details::remove_rank<N>(rf);
}
}
}
constexpr size_t N_MAX = 13;
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) if (c[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<N_MAX>(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<N_MAX>(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';
}