Submission #17532188

---

Source Code Expand

#ifndef SET_POWER_SERIES
#define SET_POWER_SERIES

// The class SubsetFunction<T> represents a function on the subsets of
// {0, 1, ..., N-1} with values in T. A subset is represented by its bitmask.
// The following operations are supported:
//  xor_convolution
//  and_convolution
//  or_convolution
//  subset_convolution
//  complement
//  exp
//  log
//
// Only the high-level API (that is SubsetFunction) shall be used.
// With respect to copy and assignment, the class SubsetFunction behaves like
// a vector.
// 
// Internally, all operations are implemented directly on arrays (and, whenever
// possible, in place). There are no memory leaks.
// Operating directly on arrays simplifies the implementation (since often one
// has to treat a subarray as an array and with vectors this is not possible).

// TODO: Maybe use some global arrays as temporary memory for subset_convolution and
//       inverse_subset_convolution (to avoid many allocations).

// TODO: Remove iostream.

// TODO: Remove main() and test.
// TODO: Add comments.

#include <iostream>

#include <vector>
#include <cassert>
#include <cstring>
#include <algorithm>
using namespace std;

// All functions in internal operate on pointers and are, therefore, a bit
// uncomfortable to use for the end-user.
// It is much better to use the fancier APIs offered by the class SubsetFunction.
namespace internal {

template<typename T>
void clear(int N, T* A) { fill(A, A+(1<<N), 0); }

// In place xor-transform. See xor_convolution to understand its importance.
// A'[x] = TODO.
// Complexity: O(N*2^N).
template<typename T>
void xor_transform(int N, T* A, bool inverse = false) {
    for (int len = 1; 2 * len <= (1<<N); len *= 2) {
        for (int x = 0; x < (1<<N); x += 2 * len) {
            for (int y = 0; y < len; y++) {
                T u = A[x + y];
                T v = A[x + y + len];
                A[x + y] = u + v;
                A[x + y + len] = u - v;
            }
        }
    }
    if (inverse) for (int i = 0; i < (1<<N); i++) A[i] /= (1<<N);
}

// In place and-transform. See and_convolution to understand its importance.
//   A'[x] = \sum_{x\subseteq y} A[y].
// Complexity: O(N*2^N).
template<typename T>
void and_transform(int N, T* A, bool inverse = false) {
    for (int len = 1; 2 * len <= (1<<N); len *= 2) {
        for (int x = 0; x < (1<<N); x += 2 * len) {
            for (int y = 0; y < len; y++) {
                if (inverse) A[x+y] -= A[x+y+len];
                else A[x+y] += A[x+y+len];
            }
        }
    }
}

// In place or-transform. See or_convolution to understand its importance.
//   A'[x] = \sum_{y\subseteq x} A[y].
// Complexity: O(N*2^N).
template<typename T>
void or_transform(int N, T* A, bool inverse = false) {
    for (int len = 1; 2 * len <= (1<<N); len *= 2) {
        for (int x = 0; x < (1<<N); x += 2 * len) {
            for (int y = 0; y < len; y++) {
                if (inverse) A[x+y+len] -= A[x+y];
                else A[x+y+len] += A[x+y];
            }
        }
    }
}

// C[x] = \sum_{y&z = 0, y|z = x} A[y]*B[z].
// Complexity: O(N²*2^N).
template<typename T>
void subset_convolution(int N, const T* A, const T* B, T* C) {
    T** A_popcnt = new T*[N+1];
    T** B_popcnt = new T*[N+1];
    T* tmp = new T[1<<N];
    for (int i = 0; i <= N; i++) {
        A_popcnt[i] = new T[1<<N];
        B_popcnt[i] = new T[1<<N];
        clear(N, A_popcnt[i]);
        clear(N, B_popcnt[i]);
    }
    for (int x = 0; x < (1<<N); x++) {
        int q = __builtin_popcount(x);
        A_popcnt[q][x] = A[x];
        B_popcnt[q][x] = B[x];
    }
    for (int i = 0; i <= N; i++) {
        or_transform(N, A_popcnt[i]);
        or_transform(N, B_popcnt[i]);
    }
    for (int l = 0; l <= N; l++) {
        clear(N, tmp);
        for (int i = 0; i <= l; i++) {
            int j = l-i;
            for (int x = 0; x < (1<<N); x++)
                tmp[x] += A_popcnt[i][x]*B_popcnt[j][x];
        }
        or_transform(N, tmp, true);
        for (int x = 0; x < (1<<N); x++)
            if (__builtin_popcount(x) == l) C[x] = tmp[x];
    }
    for (int i = 0; i <= N; i++) {
        delete[] A_popcnt[i];
        delete[] B_popcnt[i];
    }
    delete[] A_popcnt;
    delete[] B_popcnt;
    delete[] tmp;
}

// B = exp(A).
//
// For x != 0, it holds
// B[x] = \sum_{0 < y_1 < y_2 < ... < y_k: y_i disjoint, y_1|y_2|...|y_k = x}
//                    A[y_1]*A[y_2]*...*A[y_k].
// The value B[0] is (arbitrarily) set to 1.
//
// The value of A[0] is disregarded.
template<typename T>
void exp(int N, const T* A, T* B) {
    B[0] = 1;
    for (int n = 0; n < N; n++) 
        subset_convolution(n, A + (1<<n), B, B + (1<<n));
}

// subset_convolution(A, C) = B.
//
// ACHTUNG: The value A[0] must be invertible.
template<typename T>
void inverse_subset_convolution(int N, const T* A, const T* B, T* C) {
    T inv = 1/A[0];
    
    T** A_popcnt = new T*[N+1];
    T** C_popcnt = new T*[N+1];
    for (int i = 0; i <= N; i++) {
        A_popcnt[i] = new T[1<<N];
        C_popcnt[i] = new T[1<<N];
        clear(N, A_popcnt[i]);
        clear(N, C_popcnt[i]);
    }
    for (int x = 0; x < (1<<N); x++) {
        int q = __builtin_popcount(x);
        A_popcnt[q][x] = A[x];
    }
    for (int i = 0; i <= N; i++) or_transform(N, A_popcnt[i]);
    
    for (int l = 0; l <= N; l++) {
        for (int i = 0; i <= l; i++) {
            int j = l-i;
            for (int x = 0; x < (1<<N); x++)
                C_popcnt[l][x] += A_popcnt[i][x]*C_popcnt[j][x];
        }
        or_transform(N, C_popcnt[l], true);
        for (int x = 0; x < (1<<N); x++) {
            int q = __builtin_popcount(x);
            if (q != l) C_popcnt[l][x] = 0;
            else {
                C_popcnt[l][x] = (B[x] - C_popcnt[l][x])*inv;
                C[x] = C_popcnt[l][x];
            }
        }
        or_transform(N, C_popcnt[l]);
    }

    for (int i = 0; i <= N; i++) {
        delete[] A_popcnt[i];
        delete[] C_popcnt[i];
    }
    delete[] A_popcnt;
    delete[] C_popcnt;
}

// exp(B) = A.
// The value of B[0] is (arbitrarily) set to 0.
//
// ACHTUNG: It must hold A[0] = 1.
template<typename T>
void log(int N, const T* A, T* B) {
    assert(A[0] == 1);
    B[0] = 0;
    for (int n = 0; n < N; n++)
        inverse_subset_convolution(n, A, A+(1<<n), B+(1<<n));
}

// A'[x] = A[complement of x]
template<typename T>
void complement(int N, T* A) {
    int pot = 1<<(N-1);
    for (int x = 0; x < pot; x++) swap(A[x], A[x|pot]);
}


// C[x] = \sum_{y^z = x} A[y]B[z]
// Complexity: O(N2^N).
template<typename T>
void xor_convolution(int N, const T* A, const T* B, T* C) {
    T* B2 = new T[1<<N];
    for (int x = 0; x < (1<<N); x++) C[x] = A[x], B2[x] = B[x];
    walsh_hadamard_transform(N, C);
    walsh_hadamard_transform(N, B2);
    for (int x = 0; x < (1<<N); x++) C[x] *= B2[x];
    walsh_hadamard_transform(N, C, true);
    delete[] B2;
}

// C[x] = \sum_{y&z = x} A[y]B[z]
// Complexity: O(N2^N).
template<typename T>
void and_convolution(int N, const T* A, const T* B, T* C) {
    T* B2 = new T[1<<N];
    for (int x = 0; x < (1<<N); x++) C[x] = A[x], B2[x] = B[x];
    and_transform(N, C);
    and_transform(N, B2);
    for (int x = 0; x < (1<<N); x++) C[x] *= B2[x];
    and_transform(N, C, true);
    delete[] B2;
}

// C[x] = \sum_{y|z = x} A[y]B[z]
// Complexity: O(N2^N).
template<typename T>
void or_convolution(int N, const T* A, const T* B, T* C) {
    T* B2 = new T[1<<N];
    for (int x = 0; x < (1<<N); x++) C[x] = A[x], B2[x] = B[x];
    or_transform(N, C);
    or_transform(N, B2);
    for (int x = 0; x < (1<<N); x++) C[x] *= B2[x];
    or_transform(N, C, true);
    delete[] B2;
}

}; // namespace internal

template<typename T>
struct SubsetFunction {
    int N;
    T* arr;
    SubsetFunction(int N): N(N) {
        arr = new T[1<<N];
        clear();
    }
    ~SubsetFunction() { delete[] arr; }
    SubsetFunction(const SubsetFunction& other) {
        N = other.N;
        arr = new T[1<<N];
        for (int x = 0; x < (1<<N); x++) arr[x] = other[x];
    }
    SubsetFunction& operator=(const SubsetFunction& other) {
        if (this != &other) {
            if (N != other.N) {
                N = other.N;
                delete[] arr;
                arr = new T[1<<N];
            }
            for (int x = 0; x < (1<<N); x++) arr[x] = other[x];
        }
        return *this;
    }
    T& operator[](int index) { return arr[index]; }
    const T& operator[](int index) const { return arr[index]; }
    void clear() { internal::clear(N, arr); }
};

#define CONVOLUTION_DEF(OP)                                         \
template<typename T>                                                \
SubsetFunction<T> OP##_convolution(const SubsetFunction<T>& A,      \
                                   const SubsetFunction<T>& B) {    \
    assert(A.N == B.N);                                             \
    SubsetFunction<T> C(A.N);                                       \
    internal::OP##_convolution(A.N, A.arr, B.arr, C.arr);           \
    return C;                                                       \
}

CONVOLUTION_DEF(xor)              // Complexity O(N2^N).
CONVOLUTION_DEF(and)              // Complexity O(N2^N).
CONVOLUTION_DEF(or)               // Complexity O(N2^N).
CONVOLUTION_DEF(subset)           // Complexity O(N²2^N).
CONVOLUTION_DEF(inverse_subset)   // Complexity O(N²2^N).

#define UNARY_OPERATOR_DEF(OP)                          \
template<typename T>                                    \
SubsetFunction<T> OP(const SubsetFunction<T>& A) {      \
    SubsetFunction<T> B(A.N);                           \
    internal::OP(A.N, A.arr, B.arr);                    \
    return B;                                           \
}

UNARY_OPERATOR_DEF(complement)    // Complexity O(N).
UNARY_OPERATOR_DEF(exp)           // Complexity O(N²2^N).
UNARY_OPERATOR_DEF(log)           // Complexity O(N2^N).

// #include "../arithmetic_lib.cpp"
typedef unsigned long long int ULL;
typedef long long int LL;


// Computes the inverse of n modulo m.
// Precondition: n >= 0, m > 0 and gcd(n, m) == 1.
// The returned value satisfies 0 <= x < m (Inverse(0, 1) = 0).
// ACHTUNG: It must hold max(m, n) < 2**31 to avoid integer overflow.
LL Inverse(LL n, LL m) {
    n %= m;
    if (n <= 1) return n; // Handles properly (n = 0, m = 1).
    return m - ((m * Inverse(m, n) - 1) / n);
}

// Fast exponentiation modulo mod. Returns x**e modulo mod.
// Assumptions: 0 <= x < mod
//              mod < 2**31.
//              0 <= e < 2**63
LL pow(LL x, LL e, LL mod) {
    LL res = 1;
    for (; e >= 1; e >>= 1) {
        if (e & 1) res = res * x % mod;
        x = x * x % mod;
    }
    return res;
}

// Struct that computes x % mod faster than usual, if mod is always the same.
// It gives a x1.8 speed up over the % operator (with mod ~ 1e9 and x large).
// It is an implementation of the Barrett reduction, see
//    https://en.wikipedia.org/wiki/Barrett_reduction .
// If fast_mod is an instance of the class, then fast_mod(x) will return
// x % mod. There are no restrictions on the values of mod and x, provided
// that they fit in an unsigned long long (and mod != 0).
//
// ACHTUNG: The integer type __uint128_t must be available.
struct FastMod {
    ULL mod;
    ULL inv;
    FastMod(ULL mod) : mod(mod), inv(-1ULL / mod) {}
    ULL operator()(ULL x) {
        ULL q = (ULL)((__uint128_t(inv) * x) >> 64);
        ULL r = x - q * mod;
        return r - mod * (r >= mod);
    }
};

// Class for integers modulo mod.
// It supports all expected operations: +, -, *, /, pow, ==, < and >.
// It is as fast as it can be.
// The modulo mod shall be set through set_mod().
//
// Assumptions: mod < (1<<30).
// ACHTUNG: The integer type __uint128_t must be available.
//
// Remark: To handle larger moduli (up to 1<<62), one has to:
//          1. replace int in the definitions of mod, n.
//          2. The parameter of fast_mod must be __uint128_t, so it must be
//             changed in the definition of fast_mod and in the definition of
//             the operators * and *=.
//          3. fast_mod must be a naive modulo operation, no barrett reduction.
//          4. In Inverse, __int128_t shall be used.
struct ModularInteger {
    static int mod;
    static __uint128_t inv_mod; // Necessary for fast_mod.
    int n; // 0 <= n < mod at all times
    static void set_mod(int _mod) {
        mod = _mod;
        inv_mod = -1ULL / mod;
    }
    ModularInteger(): n(0) {}
    ModularInteger(LL _n): n(_n % mod) {
        n += (n<0)*mod;
    }
    LL get() const { return n; }
    static int fast_mod(ULL x) { // Barrett reduction.
        ULL q = (inv_mod * x) >> 64;
        int m = x - q * mod;
        m -= mod * (m >= mod);
        return m;
    }

    ModularInteger operator-() const { return n?mod-n:0; }
};
int ModularInteger::mod;
__uint128_t ModularInteger::inv_mod;

ModularInteger operator +(const ModularInteger& A, const ModularInteger& B) {
    ModularInteger C;
    C.n = A.n + B.n;
    C.n -= (C.n >= ModularInteger::mod)*ModularInteger::mod;
    return C;
}

void operator +=(ModularInteger& A, const ModularInteger& B) {
    A.n += B.n;
    A.n -= (A.n >= ModularInteger::mod)*ModularInteger::mod;
}

ModularInteger operator -(const ModularInteger& A, const ModularInteger& B) {
    ModularInteger C;
    C.n = A.n - B.n;
    C.n += (C.n < 0)*ModularInteger::mod;
    return C;
}

void operator -=(ModularInteger& A, const ModularInteger& B) {
    A.n -= B.n;
    A.n += (A.n < 0)*ModularInteger::mod;
}

ModularInteger operator *(const ModularInteger& A, const ModularInteger& B) {
    ModularInteger C;
    C.n = ModularInteger::fast_mod(((ULL)A.n) * B.n);
    return C;
}

void operator *=(ModularInteger& A, const ModularInteger& B) {
    A.n = ModularInteger::fast_mod(((ULL)A.n) * B.n);
}

// Assumption: gcd(B, mod) = 1.
ModularInteger operator /(const ModularInteger& A, const ModularInteger& B) {
    return A * Inverse(B.n, ModularInteger::mod);
}

// Assumption: gcd(B, mod) = 1.
void operator/=(ModularInteger& A, const ModularInteger& B) {
    A *= Inverse(B.n, ModularInteger::mod);
}

ModularInteger power(ModularInteger A, ULL e) {
    ModularInteger res = 1;
    for (; e >= 1; e >>= 1) {
        if (e & 1) res *= A;
        A *= A;
    }
    return res;
}

bool operator ==(const ModularInteger& A, const ModularInteger& B) {
    return A.n == B.n;
}
bool operator !=(const ModularInteger& A, const ModularInteger& B) {
    return A.n != B.n;
}
bool operator <(const ModularInteger& A, const ModularInteger& B) {
    return A.n < B.n;
}
bool operator >(const ModularInteger& A, const ModularInteger& B) {
    return A.n > B.n;
}
bool operator <=(const ModularInteger& A, const ModularInteger& B) {
    return A.n <= B.n;
}
bool operator >=(const ModularInteger& A, const ModularInteger& B) {
    return A.n >= B.n;
}

ostream& operator <<(ostream& out, const ModularInteger& A) {
    out << A.n;
    return out;
}

istream& operator >>(istream& in, ModularInteger& A) {
  LL a;
  in >> a;
  A = ModularInteger(a);
  return in;
}

typedef ModularInteger mint;


template<typename T>
SubsetFunction<T> slow_subset_convolution(const SubsetFunction<T>& A,
                                          const SubsetFunction<T>& B) {
    int N = A.N;
    SubsetFunction<T> C(N);
    for (int x = 0; x < (1<<N); x++) for (int y = 0; y < (1<<N); y++) {
        if (x&y) continue;
        C[x|y] += A[x] * B[y];
    }
    return C;
}


template<typename T>
SubsetFunction<T> slow_or_convolution(const SubsetFunction<T>& A,
                                      const SubsetFunction<T>& B) {
    int N = A.N;
    SubsetFunction<T> C(N);
    for (int x = 0; x < (1<<N); x++) for (int y = 0; y < (1<<N); y++) {
        C[x|y] += A[x] * B[y];
    }
    return C;
}

int main() {
    mint::set_mod(998244353);

    // ARC105F
    int N, M;
    cin >> N >> M;
    vector<vector<bool>> aa(N, vector<bool>(N, false));
    for (int i = 0; i < M; i++) {
        int a, b;
        cin >> a >> b;
        a--, b--;
        aa[a][b] = true;
        aa[b][a] = true;
    }
    mint inv2 = mint(1)/2;
    SubsetFunction<mint> A(N);
    for (int x = 0; x < (1<<N); x++) {
        A[x] = 1;
        for (int i = 0; i < N; i++) {
            if ((x&(1<<i)) == 0) continue;
            for (int j = i+1; j < N; j++) {
                if ((x&(1<<j)) == 0 or !aa[i][j]) continue;
                A[x] *= inv2;
            }
        }
    }

    auto B = subset_convolution(A, A);
    for (int x = 0; x < (1<<N); x++) B[x] /= A[x];
    auto C = log(B);
    cout << C[(1<<N)-1]/2 << endl;


    // GENERAL TESTING
    
    // int N = 2;
    // SubsetFunction<mint> A(N);
    // for (int x = 0; x < (1<<N); x++) A[x] = rand() % 4;

    // SubsetFunction<mint> B(N);
    // for (int x = 0; x < (1<<N); x++) B[x] = rand() % 4;

    // SubsetFunction<mint> C = subset_convolution(A, B);

    // auto C2 = slow_subset_convolution(A, B);

    // auto B2 = inverse_subset_convolution(A, C);

    // for (int x = 0; x < (1<<N); x++) assert(C[x] == C2[x]);
    // for (int x = 0; x < (1<<N); x++) assert(B[x] == B2[x]);

    // auto E = exp(A);
    // auto L = log(E);
    // for (int x = 1; x < (1<<N); x++) assert(L[x] == A[x]);
    // assert(L[0] == 0);

    // auto D = or_convolution(A, B);
    // auto D2 = slow_or_convolution(A, B);
    // for (int x = 0; x < (1<<N); x++) cout << A[x] << " ";
    // cout << endl;
    // for (int x = 0; x < (1<<N); x++) cout << B[x] << " ";
    // cout << endl;
    // for (int x = 0; x < (1<<N); x++) cout << D[x] << " ";
    // cout << endl;
    // for (int x = 0; x < (1<<N); x++) cout << D2[x] << " ";
    // cout << endl;
    // for (int x = 0; x < (1<<N); x++) assert(D[x] == D2[x]);

    // internal::or_transform(N, D2.arr);
    // internal::or_transform(N, A.arr);
    // internal::or_transform(N, B.arr);
    // cout << endl;
    // for (int x = 0; x < (1<<N); x++) cout << A[x] << " ";
    // cout << endl;
    // for (int x = 0; x < (1<<N); x++) cout << B[x] << " ";
    // cout << endl;
    // for (int x = 0; x < (1<<N); x++) cout << D2[x] << " ";
    // cout << endl;

    // YOSUPO subset convolution
    // int N;
    // cin >> N;
    // SubsetFunction<mint> A(N), B(N);
    // for (int i = 0; i < (1<<N); i++) cin >> A[i];
    // for (int i = 0; i < (1<<N); i++) cin >> B[i];
    // auto C = subset_convolution(A, B);
    // for (int i = 0; i < (1<<N); i++) cout << C[i] << " ";
    // cout << "\n";
}


#endif // SET_POWER_SERIES

Submission Info

Judge Result