結果
問題 | No.310 2文字しりとり |
ユーザー |
👑 ![]() |
提出日時 | 2024-02-20 22:16:17 |
言語 | C++17 (gcc 13.3.0 + boost 1.87.0) |
結果 |
AC
|
実行時間 | 958 ms / 6,000 ms |
コード長 | 11,473 bytes |
コンパイル時間 | 2,106 ms |
コンパイル使用メモリ | 129,964 KB |
最終ジャッジ日時 | 2025-02-19 18:10:29 |
ジャッジサーバーID (参考情報) |
judge3 / judge2 |
(要ログイン)
ファイルパターン | 結果 |
---|---|
other | AC * 28 |
ソースコード
#include<iostream>#include<vector>#include<random>namespace po167{std::random_device po_seed_gen;std::mt19937 po_random_gen(po_seed_gen());long long randLong(long long l = 0, long long r = (1ll << 62)){return std::uniform_int_distribution<long long>(l, r - 1)(po_random_gen);}template<class T>std::vector<T> rand_vec(int sz, long long l =0, long long r = (1ll << 62)){std::vector<T> res(sz);for(auto &x: res) x = T(randLong(l, r));return res;}//https://ikatakos.com/pot/programming_algorithm/number_theory/barlekamp_masseytemplate<class T>std::vector<T> Berlekamp_Massey(std::vector<T> s){std::vector<T> C={1}, B={1};int L = 0;int m = 1;T b = 1;for (int n = 0; n < int(s.size()); n++){T d = 0;for (int j = 0; j < L + 1; j++){d += C[j] * s[n - j];}if (d == 0){m++;continue;}T iv = d / b;if (2 * L <= n){auto tC = C;while (int(C.size()) <= n + 1 - L) C.push_back(0);for (int i = 0; i < int(B.size()); i++) C[i + m] -= iv * B[i];B = tC;L = n + 1 - L;b = d;m = 1;}else{for(int i = 0; i < int(B.size()); i++){C[i + m] -= iv * B[i];}m++;}}return C;}// https://nyaannyaan.github.io/library/matrix/black-box-linear-algebra.hpptemplate<class T>T det_sparse_matrix(std::vector<std::vector<T>> A, T mj = 0){int n = A.size();std::vector<T> D;struct pos_mat{int x;int y;T val;};std::vector<pos_mat> p;for (int i = 0; i < n; i++){for(int j = 0; j < n; j++) if (A[i][j] != mj) {p.push_back({i, j, A[i][j] - mj});}}while (true){while (true){D = rand_vec<T>(n);bool ok = 1;for (auto x: D) if (x == 0) ok = 0;if (ok) break;}std::vector<pos_mat> AD = p;for (int i = 0; i < int(AD.size()); i++) AD[i].val *= D[AD[i].y];std::vector<T> u = rand_vec<T>(n), v = rand_vec<T>(n);std::vector<T> b(n);std::vector<T> a(2 * n + 1);b = u;for (int i = 0; i < 2 * n + 1; i++){T sum = 0;for (int j = 0; j < n; j++){sum += u[j] * D[j];a[i] += u[j] * v[j];}sum *= mj;for (int j = 0; j < n; j++){u[j] = sum;}for (pos_mat tmp: AD){u[tmp.x] += tmp.val * b[tmp.y];}b = u;}auto mp = Berlekamp_Massey(a);if (mp.back() == 0) return 0;if (int(mp.size()) != n + 1) continue;T res = mp.back();if (n & 1) res *= -1;T tmp = 1;for (auto d: D) tmp *= d;return res / tmp;}}}namespace po167{template<class S>std::vector<S> count_outdeg(std::vector<std::vector<S>> &G){int N = G.size();std::vector<S> outdeg(N);for (int i = 0; i < N; i++) for (int j = 0; j < N; j++) outdeg[i] += G[i][j];return outdeg;}template<class S>std::vector<S> count_indeg(std::vector<std::vector<S>> &G){int N = G.size();std::vector<S> indeg(N);for (int i = 0; i < N; i++) for (int j = 0; j < N; j++) indeg[j] += G[i][j];return indeg;}// O(|M|^{3})template<class T>T det_matrix(std::vector<std::vector<T>> M){int N = M.size();if (N == 0) return 1;if (N == 1) return M[0][0];if (N == 2) return M[0][0] * M[1][1] - M[0][1] * M[1][0];T res = 1;for (int i = 0; i < N; i++){for (int j = i; j < N; j++){if (M[j][i] != 0){if (j != i){swap(M[i], M[j]);res *= -1;}break;}}if (M[i][i] == 0) return 0;res *= M[i][i];if (i + 1 == N) break;T v = 1 / M[i][i];for (int j = i + 1; j < N; j++){T t = M[j][i] * v;for (int k = i; k < N; k++){M[j][k] -= M[i][k] * t;}}}return res;}template<class T,class S>std::vector<std::vector<T>> Directed_Matrix_tree_Theorem_sub(std::vector<std::vector<S>> &G, int u = 0){int N = G.size();std::vector L(N - 1, std::vector<T>(N - 1, 0));std::vector<S> outdeg(N);for (int i = 0; i < N; i++) for (int j = 0; j < N; j++) if (i != j) outdeg[i] += G[i][j];for (int i = 0; i < N; i++) for (int j = 0; j < N; j++){if (i == u || j == u) continue;int a = i, b = j;if (u < i) a--;if (u < j) b--;if (i == j) L[a][b] = outdeg[i];else L[a][b] -= G[i][j];}return L;}template<class T,class S>T Directed_Matrix_tree_Theorem(std::vector<std::vector<S>> &G, int u = 0){if (int(G.size()) == 1) return T(1);return det_matrix(Directed_Matrix_tree_Theorem_sub<T>(G, u));}// s.t// forall i,j// 0 <= G[i][j]// if not connected// use remove isotemplate<class T,class S>std::pair<T, std::vector<std::vector<T>>> Count_Euler_Circuit_sub(std::vector<std::vector<S>> G, std::vector<T> &fact_base, bool fact_in = true){int N = G.size();std::vector<S> outdeg = count_outdeg(G);std::vector<S> indeg = count_indeg(G);// indeg == outdeg ?for (int i = 0; i < N; i++) if (indeg[i] != outdeg[i]) return {0, {{1}}};// connected ?std::vector<bool> seen(N);std::vector<int> order={0};seen[0] = 1;for (int i = 0; i < N; i++){if (i == int(order.size())) return {0, {{1}}};for (int j = 0; j < N; j++){if (G[order[i]][j] != 0 && !seen[j]){seen[j] = 1;order.push_back(j);}}}T res = 1;if (fact_in) for (int i = 0; i < N; i++) res *= fact_base[outdeg[i] - 1];return {res ,Directed_Matrix_tree_Theorem_sub<T, S>(G)};}template<class T,class S>T Count_Euler_Circuit(std::vector<std::vector<S>> G, std::vector<T> &fact_base, bool fact_in = true){auto tmp = Count_Euler_Circuit_sub(G, fact_base, fact_in);if (tmp.first == 0) return T(0);return tmp.first * det_matrix(tmp.second);}template<class S>std::vector<std::vector<S>> remove_isolated_vertex(std::vector<std::vector<S>> G){int N = G.size();std::vector<int> seen(N, -1);for (int i = 0; i < N; i++) for (int j = 0; j < N; j++){if (G[i][j] != 0) seen[i] = 1, seen[j] = 1;}int ind = 0;for (int i = 0; i < N; i++) if (seen[i] == 1) seen[i] = ind++;std::vector res(ind, std::vector<S>(ind));for (int i = 0; i < N; i++) if (seen[i] != -1){for (int j = 0; j < N; j++) if (seen[j] != -1){res[seen[i]][seen[j]] = G[i][j];}}return res;}// s.t// forall i,j// 0 <= G[i][j]// if not connected// use remove isotemplate<class T,class S>std::pair<T, std::vector<std::vector<T>>> Count_Eulerian_Trail_sub(std::vector<std::vector<S>> G,std::vector<T> &fact_base, bool fact_in = true){bool ok = 1;int N = G.size();std::vector<S> outdeg = count_outdeg(G), indeg = count_indeg(G);S sum = 0;for (int i = 0; i < N; i++){sum += outdeg[i];}int st = -1, ed = -1;for (int i = 0; i < N; i++) if (outdeg[i] != indeg[i]){if (std::abs(outdeg[i] - indeg[i]) > 1){ok = 0;break;}if (outdeg[i] > indeg[i]){if (ed != -1) ok = 0;ed = i;} else {if (st != -1) ok = 0;st = i;}}if (!ok) return {0, {{1}}};if ((st == -1) ^ (ed == -1)) return {0, {{1}}};if (st == -1){auto tmp = Count_Euler_Circuit_sub(G, fact_base, fact_in);return {T(sum) * tmp.first, tmp.second};}G[st][ed]++;return Count_Euler_Circuit_sub(G, fact_base, fact_in);}template<class T,class S>T Count_Eulerian_Trail(std::vector<std::vector<S>> G,std::vector<T> &fact_base, bool fact_in = true){auto tmp = Count_Eulerian_Trail_sub<T, S>(G, fact_base, fact_in);if (tmp.first == 0) return 0;return tmp.first * det_matrix(tmp.second);}}namespace po167{template<class T>struct Binomial{std::vector<T> fact_vec, fact_inv_vec;void extend(int m = -1){int n = fact_vec.size();if (m == -1) m = n * 2;if (n >= m) return;fact_vec.resize(m);fact_inv_vec.resize(m);for (int i = n; i < m; i++){fact_vec[i] = fact_vec[i - 1] * T(i);}fact_inv_vec[m - 1] = T(1) / fact_vec[m - 1];for (int i = m - 1; i > n; i--){fact_inv_vec[i - 1] = fact_inv_vec[i] * T(i);}}Binomial(int MAX = 2){fact_vec.resize(1, T(1));fact_inv_vec.resize(1, T(1));extend(MAX + 1);}T fact(int i){if (i < 0) return 0;while (int(fact_vec.size()) <= i) extend();return fact_vec[i];}T invfact(int i){if (i < 0) return 0;while (int(fact_inv_vec.size()) <= i) extend();return fact_inv_vec[i];}T C(int a, int b){if (a < b || b < 0) return 0;return fact(a) * invfact(b) * invfact(a - b);}T invC(int a, int b){if (a < b || b < 0) return 0;return fact(b) * fact(a - b) *invfact(a);}T P(int a, int b){if (a < b || b < 0) return 0;return fact(a) * invfact(a - b);}T inv(int a){if (a < 0) return inv(-a) * T(-1);if (a == 0) return 1;return fact(a - 1) * invfact(a);}};}struct mint{static constexpr int m = 1000000007;int x;mint() : x(0){}mint(long long x_):x(x_ % m){if (x < 0) x += m;}int val(){return x;}mint &operator+=(mint b){if ((x += b.x) >= m) x -= m; return *this;}mint &operator-=(mint b){if ((x -= b.x) < 0) x += m; return *this;}mint &operator*=(mint b){x= (long long)(x) * b.x % m; return *this;}mint pow(long long e) const {mint r = 1,b =*this;while (e){if (e & 1) r *= b;b *= b;e >>= 1;}return r;}mint inv(){return pow(m - 2);}mint &operator/=(mint b){return *this *= b.pow(m - 2);}friend mint operator+(mint a, mint b){return a += b;}friend mint operator-(mint a, mint b){return a -= b;}friend mint operator/(mint a, mint b){return a /= b;}friend mint operator*(mint a, mint b){return a *= b;}friend bool operator==(mint a, mint b){return a.x == b.x;}friend bool operator!=(mint a, mint b){return a.x != b.x;}};int solve_yuki_0310(int N, int M, std::vector<int> A, std::vector<int> B){if (N * N == M) return 1;std::vector G(N ,std::vector<int>(N, 1));for (int i = 0; i < M; i++){G[A[i] - 1][B[i] - 1] = 0;}po167::Binomial<mint> bin(N + 100);auto tmp = po167::Count_Eulerian_Trail_sub(po167::remove_isolated_vertex(G), bin.fact_vec);return (tmp.first * po167::det_sparse_matrix<mint>(tmp.second , -1)).val();}int main(){int N, M;std::cin >> N >> M;std::vector<int> A(M), B(M);for (int i = 0; i < M; i++) std::cin >> A[i] >> B[i];std::cout << solve_yuki_0310(N, M, A, B) << "\n";}