zoukankan      html  css  js  c++  java
  • BZOJ 3456: 城市规划 [多项式求逆元 组合数学 | 生成函数 多项式求ln]

    3456: 城市规划

    题意:n个点组成的无向连通图个数


    以前做过,今天复习一下


    (f[n])为n个点的无向连通图个数

    n个点的完全图个数为(2^{inom{n}{2}})

    和Bell数的推导很类似,枚举第一个cc的点的个数

    [2^{inom{n}{2}} = sum_{i=1}^n inom{n-1}{i-1} f(i) 2^{inom{n-i}{2}} ]

    整理后

    [frac{2^{inom{n}{2}}}{(n-1)!} = sum_{i=1}^n frac{f(i)}{(i-1)!}frac{2^{inom{n-i}{2}}}{(n-i)!} ]

    这是卷积的形式

    [C(x) = A(x)B(x) ightarrow A(x) = C(x)B^{-1}(x) ]

    多项式求逆就可做了

    注意(b_0=1)


    其实就是EGP求ln......

    简单的有标号集合计数

    [G(x) = sum_{ige 0} 2^{inom{i}{2}}frac{x^i}{i!} = sum_{i ge 0} frac{F(x)^i}{i!} = e^{F(x)}\ F(x) = ln G(x) = int frac{G'(x)}{G(x)}dx ]


    第一份是多项式求ln的代码,注意EGP要除以阶乘逆元
    #include <iostream>
    #include <cstdio>
    #include <cstring>
    #include <algorithm>
    #include <cmath>
    using namespace std;
    typedef long long ll;
    const int N = (1<<18) + 5;
    inline int read(){
    	char c=getchar();int x=0,f=1;
    	while(c<'0'||c>'9'){if(c=='-')f=-1;c=getchar();}
    	while(c>='0'&&c<='9'){x=x*10+c-'0';c=getchar();}
    	return x*f;
    }
    
    int P = 1004535809, inv2 = (P+1)/2;
    inline int Pow(ll a, int b) {
    	ll ans = 1;
    	for(; b; b>>=1, a=a*a%P)
    		if(b&1) ans=ans*a%P;
    	return ans;
    }
    
    ll inv[N];
    namespace ntt {
    	int g = 3, rev[N];
    	void dft(int *a, int n, int flag) {
    		int k = 0; while((1<<k) < n) k++;
    		for(int i=0; i<n; i++) {
    			rev[i] = (rev[i>>1]>>1) | ((i&1)<<(k-1));
    			if(i < rev[i]) swap(a[i], a[rev[i]]);
    		}
    		for(int l=2; l<=n; l<<=1) {
    			int m = l>>1, wn = Pow(g, flag == 1 ? (P-1)/l : P-1-(P-1)/l);
    			for(int *p = a; p != a+n; p += l)
    				for(int k=0, w=1; k<m; k++, w = (ll)w*wn %P) {
    					int t = (ll) w * p[k+m] %P, r = p[k];
    					p[k+m] = (r - t + P) %P;
    					p[k] = (r + t) %P;
    				}
    		}
    		if(flag == -1) {
    			ll inv = Pow(n, P-2);
    			for(int i=0; i<n; i++) a[i] = a[i] * inv %P;
    		}
    	}
    	
    	void inverse(int *a, int *b, int l) {
    		static int t[N];
    		if(l == 1) {b[0] = Pow(a[0], P-2); return;}
    		inverse(a, b, l>>1);
    		int n = l<<1;
    		for(int i=0; i<l; i++) t[i] = a[i], t[i+l] = 0;
    		dft(t, n, 1); dft(b, n, 1);
    		for(int i=0; i<n; i++) b[i] = (ll) b[i] * (2 - (ll) t[i] * b[i] %P + P) %P;
    		dft(b, n, -1); for(int i=l; i<n; i++) b[i] = 0;
    	}
    
    	void ln(int *a, int *b, int l) {
    		static int da[N], ia[N];
    		int n = l<<1;
    		for(int i=0; i<n; i++) da[i] = ia[i] = 0;
    		for(int i=0; i<l-1; i++) da[i] = (ll) (i+1) * a[i+1] %P;
    		inverse(a, ia, l);
    		dft(da, n, 1); dft(ia, n, 1);
    		for(int i=0; i<n; i++) b[i] = (ll) da[i] * ia[i] %P;
    		dft(b, n, -1);
    		for(int i=l-1; i>0; i--) b[i] = (ll) inv[i] * b[i-1] %P; b[0] = 0;
    		for(int i=l; i<n; i++) b[i] = 0;
    	}
    	
    }
    
    ll fac[N], facInv[N];
    int n, a[N], b[N];
    int main() {
    	freopen("in", "r", stdin);
    	n = read();
    	int len = 1;
    	while(len <= n) len <<= 1;
    	inv[1] = 1; fac[0] = facInv[0] = 1;
    	for(int i=1; i<=n; i++) {
    		if(i != 1) inv[i] = (P - P/i) * inv[P%i] %P;
    		fac[i] = fac[i-1] * i %P;
    		facInv[i] = facInv[i-1] * inv[i] %P;
    	}
    	a[0] = 1;
    	for(int i=1; i<=n; i++) a[i] = Pow(2, (ll) i * (i - 1) /2 %(P-1)) * facInv[i] %P;
    	ntt::ln(a, b, len);
    	//for(int i=0; i<=n; i++) printf("%d ", b[i]); puts("");
    	int ans = b[n] * fac[n] %P;
    	printf("%d", ans);
    }
    
    
    
    #include <iostream>
    #include <cstdio>
    #include <cstring>
    #include <algorithm>
    #include <cmath>
    using namespace std;
    typedef long long ll;
    const int N = (1<<18) + 5;
    inline int read(){
    	char c=getchar();int x=0,f=1;
    	while(c<'0'||c>'9'){if(c=='-')f=-1;c=getchar();}
    	while(c>='0'&&c<='9'){x=x*10+c-'0';c=getchar();}
    	return x*f;
    }
    
    ll P = 1004535809;
    inline ll Pow(ll a, int b) {
    	ll ans = 1;
    	for(; b; b>>=1, a=a*a%P)
    		if(b&1) ans=ans*a%P;
    	return ans;
    }
    inline void mod(int &x) {if(x>=P) x-=P; else if(x<0) x+=P;}
    namespace fnt {
    	int n, g=3, rev[N];
    	void dft(int *a, int n, int flag=1) {
    		for(int i=0; i<n; i++) if(i < rev[i]) swap(a[i], a[rev[i]]);
    
    		for(int l=2; l<=n; l<<=1) {
    			int m = l>>1;
    			ll wn = Pow(g, flag==1 ? (P-1)/l : P-1-(P-1)/l);
    			for(int *p=a; p!=a+n; p+=l) {
    				ll w = 1;
    				for(int k=0; k<m; k++) {
    					ll t = p[k+m] * w %P;
    					mod(p[k+m] = p[k] - t);
    					mod(p[k] = p[k] + t);
    					w = w * wn %P;
    				}
    			}
    		}
    		if(flag == -1) {
    			ll inv = Pow(n, P-2);
    			for(int i=0; i<n; i++) a[i] = a[i] * inv %P;
    		}
    	}
    	int t[N];
    	void inverse(int *a, int *b, int l) { // mod x^l
    		if(l == 1) {b[0] = Pow(a[0], P-2); return;}
    		inverse(a, b, (l+1)>>1);
    		int n = 1, k = 0; while(n < l<<1) n<<=1, k++;
    		for(int i=0; i<n; i++) rev[i] = (rev[i>>1]>>1) | ((i&1)<<(k-1));
    		for(int i=0; i<l; i++) t[i] = a[i]; for(int i=l; i<n; i++) t[i] = 0;
    		dft(t, n, 1); dft(b, n, 1);
    		for(int i=0; i<n; i++) b[i] = (ll) b[i] * (2 - (ll) t[i] * b[i] %P + P) %P;
    		dft(b, n, -1);
    		for(int i=l; i<n; i++) b[i] = 0;
    	}
    
    	void mul(int *a, int *b, int l) {
    		int n = 1; while(n < l<<1) n<<=1;
    		dft(a, n, 1); dft(b, n, 1);
    		for(int i=0; i<n; i++) a[i] = (ll) a[i] * b[i] %P;
    		dft(a, n, -1);
    	}
    }
    
    int n, a[N], b[N], bi[N];
    ll inv[N], fac[N], facInv[N];
    int main() {
    	freopen("in", "r", stdin);
    	n=read();
    	inv[1] = 1; fac[0] = facInv[0] = 1;
    	for(int i=1; i<=n; i++) {
    		if(i != 1) inv[i] = (P - P/i) * inv[P%i] %P;
    		fac[i] = fac[i-1] * i %P;
    		facInv[i] = facInv[i-1] * inv[i] %P;
    	}
    	for(int i=1; i<=n; i++) {
    		ll mi = Pow(2, (ll) i * (i-1) / 2 %(P-1));
    		a[i] = mi * facInv[i-1] %P;
    		b[i] = mi * facInv[i] %P;
    	}
    	b[0] = 1;
    
    	fnt::inverse(b, bi, n+1); 
    	fnt::mul(a, bi, n+1);
    	ll ans = a[n] * fac[n-1] %P;
    	printf("%lld", ans);
    }
    
    
  • 相关阅读:
    < java.util >-- Set接口
    Codeforces 627 A. XOR Equation (数学)
    Codeforces 161 B. Discounts (贪心)
    Codeforces 161 D. Distance in Tree (树dp)
    HDU 5534 Partial Tree (完全背包变形)
    HDU 5927 Auxiliary Set (dfs)
    Codeforces 27E. Number With The Given Amount Of Divisors (暴力)
    lght oj 1257
    Codeforces 219D. Choosing Capital for Treeland (树dp)
    Codeforces 479E. Riding in a Lift (dp + 前缀和优化)
  • 原文地址:https://www.cnblogs.com/candy99/p/6736081.html
Copyright © 2011-2022 走看看