其他分享
首页 > 其他分享> > bzoj3456 城市规划

bzoj3456 城市规划

作者:互联网

传送门
我们令\(g(n)\)表示有n个点的无向图的个数,\(f(n)\)表示有n个点的无向连通图个数。
\(g(n)\)比较好求,因为一共有\(C_n^2\)条边,每条边可以选或者不选,所以自然有\(g(n) = 2^{C_n^2}\).

之后我们换一种方法,在节点1所在的联通块中有多少个节点,剩下的只需要随便连就可以了。那么就有\(g(n) = \sum_{i=1}^nC_{n-1}^{i-1}f(i)g(n-i)\)

把\(g(n) = 2^{C_n^2}\)带进式子,两边再同时除以\((n-1)!\),得到:
\[\frac{2^{C_n^2}}{(n-1)!} = \sum_{i=1}^n\frac{f(i)}{(i-1)!}\frac{2^{C_{n-i}^2}}{(n-i)!}\]

这玩意可以看出来是卷积的形式,那么我们设:
\[F(x) = \sum_{i=1}^\infty\frac{f(i)}{(i-1)}x^i\]
\[G(x) = \sum_{i=0}^\infty\frac{2^{C_i^2}}{i!}x^i\]
\[H(x) = \sum_{i=1}^\infty\frac{2^{C_i^2}}{(i-1)!}\]

那么\(H(x) = F(x)G(x)\),在\((mod\ x^{n+1})\)意义下求出\(G(x)^{-1}\),那么\(F(x) \equiv G(x)^{-1}H(x) (mod \ x^{n+1})\)即为答案。注意最后要乘上\(fac[n-1]\).

注意不要像我一样犯智障错误……我对指数取模导致第一次只有65pts,但是可以对mod-1取模,因为模数是质数,而\(a^{p-1} \equiv 1(mod\ p)\),对答案没有影响。

#include<bits/stdc++.h>
#define rep(i,a,n) for(int i = a;i <= n;i++)
#define per(i,n,a) for(int i = n;i >= a;i--)
#define enter putchar('\n')
#define fr friend inline
#define y1 poj

using namespace std;
typedef long long ll;
const int M = 800005;
const int mod = 1004535809;
const int G = 3;
const int invG = 334845270;
const double eps = 1e-7;

int read()
{
   int ans = 0,op = 1;char ch = getchar();
   while(ch < '0' || ch > '9') {if(ch == '-') op = -1;ch = getchar();}
   while(ch >= '0' && ch <= '9') ans = ans * 10 + ch - '0',ch = getchar();
   return ans * op;
}

int n,f[M],g[M],c[M],invg[M],d[M],rev[M],fac[M],inv[M],ans;

int inc(int a,int b){return (a+b) % mod;}
int mul(int a,int b){return 1ll * a * b % mod;}
int qpow(int a,ll b)
{
   int p = 1;
   while(b)
   {
      if(b & 1) p = mul(p,a);
      a = mul(a,a),b >>= 1;
   }
   return p;
}

void NTT(int *a,int l,int f)
{
   rep(i,0,l-1) if(i < rev[i]) swap(a[i],a[rev[i]]);
   for(int i = 1;i < l;i <<= 1)
   {
      int w1 = qpow(f ? G : invG,(mod-1) / (i<<1));
      for(int j = 0;j < l;j += (i<<1))
      {
     int w = 1;
     rep(k,0,i-1)
     {
        int kx = a[k+j],ky = mul(a[k+j+i],w);
        a[k+j] = inc(kx,ky),a[k+j+i] = inc(kx,mod-ky),w = mul(w,w1);
     }
      }
   }
   if(!f)
   {
      int inv = qpow(l,mod-2);
      rep(i,0,l-1) a[i] = mul(a[i],inv);
   }
}

void init()
{
   fac[0] = inv[0] = inv[1] = 1;
   rep(i,1,n) fac[i] = mul(fac[i-1],i);
   inv[n] = qpow(fac[n],mod-2);
   per(i,n-1,1) inv[i] = mul(inv[i+1],i+1);
}

int C(int n,int m){return mul(mul(fac[n],inv[m]),inv[n-m]);}

void getinv(int *a,int *b,int len)
{
   if(len == 1) {b[0] = qpow(a[0],mod-2);return;}
   getinv(a,b,(len+1)>>1);
   int l = 1,L = 0;
   while(l < (len<<1)) l <<= 1,L++;
   rep(i,0,l-1) rev[i] = (rev[i>>1] >> 1) | ((i&1) << (L-1));
   rep(i,0,len-1) d[i] = a[i];
   rep(i,len,l-1) d[i] = 0;
   NTT(d,l,1),NTT(b,l,1);
   rep(i,0,l-1) b[i] = mul(inc(2,mod-mul(d[i],b[i])),b[i]);
   NTT(b,l,0);
   rep(i,len,l-1) b[i] = 0;
}

int main()
{
   n = read(),init();
   //rep(i,1,n) c[i] = mul(qpow(2,C(i,2)),inv[i-1]);错误写法
   //rep(i,0,n) g[i] = mul(qpow(2,C(i,2)),inv[i]);
   rep(i,1,n) c[i] = mul(qpow(2,(1ll * i * (i-1) / 2)),inv[i-1]);
   rep(i,0,n) g[i] = mul(qpow(2,(1ll * i * (i-1) / 2)),inv[i]);
   getinv(g,invg,n+1);
   int l = 1,L = 0;
   while(l < n<<1) l <<= 1,L++;
   rep(i,0,l-1) rev[i] = (rev[i>>1] >> 1) | ((i&1) << (L-1));
   NTT(invg,l,1),NTT(c,l,1);
   rep(i,0,l-1) c[i] = mul(c[i],invg[i]);
   NTT(c,l,0);
   ans = mul(c[n],fac[n-1]);
   printf("%d\n",ans);
   return 0;
}

标签:城市规划,ch,const,int,sum,bzoj3456,frac,mod
来源: https://www.cnblogs.com/captain1/p/10460076.html