【BZOJ】3456: 城市规划 动态规划+多项式求逆
【题意】求n个点的带标号无向连通图个数 mod 1004535809。n<=130000。
【算法】动态规划+多项式求逆
【题解】设$g_n$表示n个点的无向图个数,那么显然
$$g_n=2^{\frac{n(n-1)}{2}}$$
设$f_n$表示n个点的无向连通图个数,通过枚举1号点所属连通块大小很容易得到$g_n$的等式:
$$g_n=\sum_{i=1}^{n}\binom{n-1}{i-1}*f_i*g_{n-i}$$
特别的,$g_0=1$。
将组合数拆分一下,即可得到:
$$\frac{g_n}{(n-1)!}=\sum_{i=1}^{n}\frac{f_i}{(i-1)!}*\frac{g_{n-i}}{(n-i)!}$$
多项式求逆即可,注意下标从1开始需要强制F(0)=0,以及由于$g_0=1$所以右边额外+1。(这个式子恰好是多项式乘法的形式)
复杂度O(n log n)。
#include<cstdio> #include<cstring> #include<cctype> #include<cmath> #include<queue> #include<stack> #include<set> #include<vector> #include<algorithm> #define ll long long #define lowbit(x) x&-x using namespace std; int read(){ char c;int s=0,t=1; while(!isdigit(c=getchar()))if(c=='-')t=-1; do{s=s*10+c-'0';}while(isdigit(c=getchar())); return s*t; } int min(int a,int b){return a<b?a:b;} int max(int a,int b){return a<b?b:a;} int ab(int x){return x>0?x:-x;} //int MO(int x){return x>=MOD?x-MOD:x;} //void insert(int u,int v){tot++;e[tot].v=v;e[tot].from=first[u];first[u]=tot;} /*------------------------------------------------------------*/ const int inf=0x3f3f3f3f,maxn=270010,MOD=1004535809; int n,F[maxn],G[maxn],H[maxn]; int power(int x,int k){int ans=1;while(k){if(k&1)ans=1ll*ans*x%MOD;x=1ll*x*x%MOD;k>>=1;}return ans;} int inv(int x){return power(x,MOD-2);} int M(int x){return x>=MOD?x-MOD:x;} void ntt(int *a,int n,int f){ int k=0; for(int i=0;i<n;i++){ if(i<k)swap(a[i],a[k]); for(int j=n>>1;(k^=j)<j;j>>=1); } for(int l=2;l<=n;l<<=1){ int m=l>>1,wn=(~f?power(3,(MOD-1)/l):power(3,MOD-1-(MOD-1)/l)); for(int *p=a;p!=a+n;p+=l){ int w=1; for(int i=0;i<m;i++){ int t=1ll*w*p[i+m]%MOD; p[i+m]=M(p[i]-t+MOD); p[i]=M(p[i]+t); w=1ll*w*wn%MOD; } } } if(f==-1){ int o=inv(n); for(int i=0;i<n;i++)a[i]=1ll*a[i]*o%MOD; } } int h[maxn]; void pinv(int *f,int *g,int n){ if(n==1){g[0]=inv(f[0]);return;} pinv(f,g,n>>1);n<<=1; for(int i=0;i<n/2;i++)h[i]=f[i]; //for(int i=n/2;i<n;i++)h[i]=g[i]=0; ntt(h,n,1);ntt(g,n,1); for(int i=0;i<n;i++)g[i]=1ll*g[i]*(2-1ll*h[i]*g[i]%MOD+MOD)%MOD; ntt(g,n,-1); for(int i=n/2;i<n;i++)g[i]=0;//! } int fac[maxn],fav[maxn]; int main(){ scanf("%d",&n);n++; fac[0]=1;for(int i=1;i<=n;i++)fac[i]=1ll*fac[i-1]*i%MOD; fav[n]=inv(fac[n]);for(int i=n;i>=1;i--)fav[i-1]=1ll*fav[i]*i%MOD; for(int i=1;i<n;i++)H[i]=1ll*power(2,1ll*i*(i-1)/2%(MOD-1))*fav[i-1]%MOD; for(int i=0;i<n;i++)G[i]=1ll*power(2,1ll*i*(i-1)/2%(MOD-1))*fav[i]%MOD;//n change int N=1;while(N<n+n)N<<=1; pinv(G,F,N>>1);//n>>1 ntt(H,N,1);ntt(F,N,1); for(int i=0;i<N;i++)F[i]=1ll*H[i]*F[i]%MOD; ntt(F,N,-1); printf("%lld",1ll*F[n-1]*fac[n-2]%MOD); return 0; }
注意:求逆元的时候传进去一倍即可,里面会再翻倍。逆元每次最后g数组后半部分要清零。n++之后对应的n要改变。
NTT最好预处理omega[],不然NTT太多次会变得很慢,即:
#include<cstdio> #include<cstring> #include<cctype> #include<cmath> #include<queue> #include<stack> #include<set> #include<vector> #include<algorithm> #define ll long long #define lowbit(x) x&-x using namespace std; int read(){ char c;int s=0,t=1; while(!isdigit(c=getchar()))if(c=='-')t=-1; do{s=s*10+c-'0';}while(isdigit(c=getchar())); return s*t; } int min(int a,int b){return a<b?a:b;} int max(int a,int b){return a<b?b:a;} int ab(int x){return x>0?x:-x;} //int MO(int x){return x>=MOD?x-MOD:x;} //void insert(int u,int v){tot++;e[tot].v=v;e[tot].from=first[u];first[u]=tot;} /*------------------------------------------------------------*/ const int inf=0x3f3f3f3f,maxn=270010,MOD=1004535809; int n,F[maxn],G[maxn],H[maxn]; int power(int x,int k){int ans=1;while(k){if(k&1)ans=1ll*ans*x%MOD;x=1ll*x*x%MOD;k>>=1;}return ans;} int inv(int x){return power(x,MOD-2);} int M(int x){return x>=MOD?x-MOD:x;} int o[maxn],oi[maxn]; void init(int n){ int x=1,y=power(3,(MOD-1)/n); for(int i=0;i<=n;i++){ o[i]=oi[n-i]=x; x=1ll*x*y%MOD; } } void ntt(int *a,int n,int *o,int f){ int k=0; for(int i=0;i<n;i++){ if(i<k)swap(a[i],a[k]); for(int j=n>>1;(k^=j)<j;j>>=1); } for(int l=2;l<=n;l<<=1){ int m=l>>1; for(int *p=a;p!=a+n;p+=l){ for(int i=0;i<m;i++){ int t=1ll*o[n/l*i]*p[i+m]%MOD; p[i+m]=M(p[i]-t+MOD); p[i]=M(p[i]+t); } } } if(f==-1){ int o=inv(n); for(int i=0;i<n;i++)a[i]=1ll*a[i]*o%MOD; } } int h[maxn]; void pinv(int *f,int *g,int n){ if(n==1){g[0]=inv(f[0]);return;} pinv(f,g,n>>1);n<<=1; init(n); for(int i=0;i<n/2;i++)h[i]=f[i]; //for(int i=n/2;i<n;i++)h[i]=g[i]=0;//? ntt(h,n,o,1);ntt(g,n,o,1); for(int i=0;i<n;i++)g[i]=1ll*g[i]*(2-1ll*h[i]*g[i]%MOD+MOD)%MOD; ntt(g,n,oi,-1); for(int i=n/2;i<n;i++)g[i]=0;//? } int fac[maxn],fav[maxn]; int main(){ scanf("%d",&n);n++; fac[0]=1;for(int i=1;i<=n;i++)fac[i]=1ll*fac[i-1]*i%MOD; fav[n]=inv(fac[n]);for(int i=n;i>=1;i--)fav[i-1]=1ll*fav[i]*i%MOD; for(int i=1;i<n;i++)H[i]=1ll*power(2,1ll*i*(i-1)/2%(MOD-1))*fav[i-1]%MOD; for(int i=0;i<n;i++)G[i]=1ll*power(2,1ll*i*(i-1)/2%(MOD-1))*fav[i]%MOD;//! int N=1;while(N<n+n)N<<=1; pinv(G,F,N>>1); init(N); ntt(H,N,o,1);ntt(F,N,o,1); for(int i=0;i<N;i++)F[i]=1ll*H[i]*F[i]%MOD; ntt(F,N,oi,-1); printf("%lld",1ll*F[n-1]*fac[n-2]%MOD); return 0; }
这道题有个O(n^2)的递推做法,设$h_n$表示n个点无向不连通图个数,那么$h_n=g_n-f_n$,枚举1所在连通块大小,有:
$$h_n=\sum_{i=1}^{n-1}\binom{n-1}{i-1}*f_i*h_{n-i}$$
这个不包含n自己,就可以递推了。
在这里记个无关紧要的笔记……要拆分$2^{k(n-k)}$,乘法不好拆分,必须转化为加法。
其实就是要把nk转化为有关n或k或n-k的加减法的形式,由初中常用套路$nk=\frac{(n-k)^2-n^2-k^2}{2}$就可以了。