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}^nfrac{f(i)}{(i-1)!}frac{2^{C_{n-i}^2}}{(n-i)!} ]

这玩意可以看出来是卷积的形式，那么我们设:

[F(x) = sum_{i=1}^inftyfrac{f(i)}{(i-1)}x^i ]

[G(x) = sum_{i=0}^inftyfrac{2^{C_i^2}}{i!}x^i ]

[H(x) = sum_{i=1}^inftyfrac{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('
')
#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
",ans);
   return 0;
}