[BZOJ 3456]城市规划
转载请注明出处:http://danihao123.is-programmer.com/
aji多项式求逆毁我青春,,,
设\(f_i\)表示\(i\)个点的有标号无向联通图,考虑所有可能的图(记\(F_i\)为\(i\)个点的有标号无向图,显然\(F_i = 2^{\binom{i}{2}}\))和\(f\)的关系(使用图计数的经典套路:枚举1所在的联通块大小):
\[F_n = \sum_{i = 1}^n \binom{n-1}{i-1} f_i F_{n - i}\]
看起来事卷积?但是这个卷积没有办法直接用FFT/NTT求(当然分离一下项啊,移下项就可以分治NTT力)。
考虑进一步化简柿子。完全展开后会发现右边有个非常碍眼的\((n-1)!\),所以两边除一下:
\[\frac{F_n}{(n-1)!} =\sum_{i = 1}^n \frac{f_i}{(i-1)!}\cdot \frac{F_{n - i}}{(n - i)!}\]
然后这个问题就很毒瘤了:我们要求的答案多项式(除上那个阶乘)和一个已知多项式做卷积,可以得到另一个已知多项式……这样就需要多项式除法了,于是乎多项式逆元派上了用场。
代码:
#include <cstdio> #include <cstring> #include <cstdlib> #include <cctype> #include <algorithm> #include <utility> typedef long long ll; const int maxn = 1020010; const ll ha = 1004535809LL; const ll bs = 3LL; ll pow_mod(ll a, ll b) { ll ans = 1LL, res = a % ha; while(b) { if(1LL & b) ans = (ans * res) % ha; res = (res * res) % ha; b >>= 1; } return ans; } ll inv(ll x) { return pow_mod(x, ha - 2LL); } int flip(int bi, int x) { int ans = 0; for(int i = 0; i < bi; i ++) { if((1 << i) & x) { ans += (1 << (bi - i - 1)); } } return ans; } void NTT(ll *A, int bi, bool flag = false) { int n = 1 << bi; for(int i = 0; i < n; i ++) { int v = flip(bi, i); if(v < i) std::swap(A[v], A[i]); } for(int L = 1; L < n; L <<= 1) { ll xi = pow_mod(3LL, (ha - 1LL) / (ll(L << 1))); if(flag) xi = inv(xi); for(int i = 0; i < n; i += (L << 1)) { ll w = 1LL; for(int j = i; j < i + L; j ++) { ll v1 = A[j], v2 = A[j + L]; A[j] = (v1 + (w * v2) % ha) % ha; A[j + L] = (v1 - (w * v2) % ha + ha) % ha; w = (w * xi) % ha; } } } } void poly_mul(ll *A, ll *B, int bi, ll *C) { static ll T1[maxn], T2[maxn]; int n = (1 << bi); std::copy(A, A + n, T1); std::copy(B, B + n, T2); NTT(T1, bi); NTT(T2, bi); #ifdef LOCAL puts("poly_mul :"); for(int i = 0; i < (n); i ++) { printf("%lld ", A[i]); } puts(""); for(int i = 0; i < (n); i ++) { printf("%lld ", B[i]); } puts(""); for(int i = 0; i < (n); i ++) { printf("%lld ", T1[i]); } puts(""); for(int i = 0; i < (n); i ++) { printf("%lld ", T2[i]); } puts(""); #endif for(int i = 0; i < n; i ++) { T1[i] = (T1[i] * T2[i]) % ha; } #ifdef LOCAL for(int i = 0; i < (n); i ++) { printf("%lld ", T1[i]); } puts(""); #endif NTT(T1, bi, true); ll inv_n = inv(n); for(int i = 0; i < n; i ++) { T1[i] = (T1[i] * inv_n) % ha; } std::copy(T1, T1 + n, C); #ifdef LOCAL for(int i = 0; i < (n); i ++) { printf("%lld ", C[i]); } puts(""); #endif } void poly_inv(int mod, ll *B, ll *BB) { if(mod == 1) { BB[0] = inv(B[0]); } else { poly_inv((mod + 1) >> 1, B, BB); int bi = 0, sz = 1; while(sz <= ((mod * 2) + 1)) { bi ++; sz <<= 1; } ll inv_sz = inv(sz); static ll tmp[maxn]; std::copy(B, B + mod, tmp); std::fill(tmp + mod, tmp + sz, 0LL); NTT(tmp, bi); NTT(BB, bi); for(int i = 0; i < sz; i ++) { tmp[i] = (tmp[i] * BB[i]) % ha; tmp[i] = (tmp[i] * (ha - 1LL)) % ha; tmp[i] = (tmp[i] + 2LL) % ha; tmp[i] = (tmp[i] * BB[i]) % ha; } NTT(tmp, bi, true); for(int i = 0; i < sz; i ++) { tmp[i] = (tmp[i] * inv_sz) % ha; } std::copy(tmp, tmp + mod, BB); std::fill(BB + mod, BB + sz, 0LL); } } int main() { static ll fac[maxn], A[maxn], B[maxn], BB[maxn]; int n; scanf("%d", &n); int bi = 0, sz = 1; while(sz <= n + 1) { bi ++; sz <<= 1; } fac[0] = 1LL; for(int i = 1; i <= n; i ++) { fac[i] = (fac[i - 1] * (ll(i))) % ha; } B[0] = 1LL; for(int i = 1; i <= n; i ++) { B[i] = pow_mod(2LL, (ll(i)) * (ll(i - 1)) / 2LL); B[i] = (B[i] * inv(fac[i])) % ha; } for(int i = 1; i <= n; i ++) { A[i] = pow_mod(2LL, (ll(i)) * (ll(i - 1)) / 2LL); A[i] = (A[i] * inv(fac[i - 1])) % ha; } poly_inv(n + 1, B, BB); poly_mul(A, BB, bi, A); printf("%lld\n", (A[n] * fac[n - 1]) % ha); return 0; }