【题解】P4841 城市规划(指数型母函数+多项式Ln)

【题解】P4841 城市规划

P4841 城市规划

超级弱化版本(DP)：POJ - 1737

两张图不同当且仅当边的分布不一样的时候，带编号最后乘一个阶乘即可，现在最主要的问题就是"联通"这个条件。

我首先考虑的容斥，"随意连不联通"的方案太好算了，\(2^{n(n-1)/2}\)，但是发现不会降低复杂度，因为"联通"和"随意连不联通"不好容斥，他们之间的关系不是简单的关系，所以不行了。

其次考虑DP，仍然发现不行，数据范围太大，之前做过一道POJ - 1737 的\(DP\)是\(O(n^2)\)的，那个式子不会优化，可能可以用\(NTT\)优化一下吧，不会。

但是这些思考给了我们启发，考虑一下之前容斥的时候，"\(f(x)\)联通"和"\(g(x)\)随意连不联通"的关系究竟是什么

一张随意联通图是由很多联通图组成，很多究竟是哪些？究竟要哪些自变量？显然既要枚举有多少联通图，也要确定每张连通图的大小。这就是我们之前不能简单容斥的原因。那么，对于这样的问题，既要枚举多少，也要枚举每个的状态的计数怎么办？

一个很好的手法就是构造母函数，我们知道乘法的本质就是组合(大雾)，因为括号乘以括号就自动帮我们完成了枚举多少和枚举状态的工作。

令\(f(x)\)表示\(x\)个无区别点的连通图方案数，\(g(x)\)表示\(x\)个无区别点的随意联通图方案数

又令

\[G(x)=\sum_{i=0}^{\inf} g(i)x^i \\ F(x)=\sum_{i=0}f(i)x^i \]

则有

\[G(x)=\sum_{i=0}^{\inf} \dfrac {(F(x))^i} {i!} \]

看不懂式子的我只想说取追踪\(G(x)\)的\(x^i\)系数来源就懂了，要除以一个阶乘是因为要消去顺序，\((F(x))^i\)是有序的。

联想到指数型生成函数

知道\(e^x\)的麦克劳林级数

\[e^x=\sum_{i=0}^{\inf} \dfrac {x^i} {i!} \]

观察一下我们\(G(x)\)的式子，解出来了，结果很显然就是

\[G(x)= e^{F(x)} \]

则

\[F(x)=\ln G(x) \]

左转默写板子，其中\(G(x)\)的\([x^i]=2^{i(i-1)/2} \times \dfrac 1 {i!}\)

//@winlere
#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
#include<cmath>

using namespace std;  typedef long long ll;
inline int qr(){
      register int ret=0,f=0;
      register char c=getchar();
      while(c<48||c>57)f|=c==45,c=getchar();
      while(c>=48&&c<=57) ret=ret*10+c-48,c=getchar();
      return f?-ret:ret;
}

const int maxn=1<<18|1;
int jc[maxn],inv[maxn],r[maxn];
const int mod=1004535809;
const int g=3;
inline int ksm(const int&,const ll&);
const int gi=ksm(3,mod-2);
inline int ksm(const int&base,const ll&p){
      register int ret=1;
      for(register int t=p%(mod-1),b=base%mod;t;t>>=1,b=1ll*b*b%mod)
	    if(t&1) ret=1ll*ret*b%mod;
      return ret;
}

inline void getr(const int&len){
      static int savlen=0;
      if(len==savlen)return;
      int cnt=0;
      for(register int t=1;t<len;t<<=1)++cnt;
      for(register int t=0;t<len;++t) r[t]=r[t>>1]>>1|(t&1)<<cnt>>1;
}

inline void NTT(int*a,const int&len,const int&tag){
      getr(len);
      for(register int t=0;t<len;++t)
	    if(r[t]>t)swap(a[t],a[r[t]]);
      int *a0,*a1,s=g;
      if(tag!=1)s=gi;
      for(register int t=1,wn;t<len;t<<=1){
	    wn=ksm(s,(mod-1)/(t<<1));
	    for(register int i=0;i<len;i+=t<<1){
		  a1=(a0=a+i)+t;
		  for(register int k=0,w=1,m;k<t;++k,++a1,++a0,w=1ll*w*wn%mod){
			m=1ll**a1*w%mod;
			*a1=(*a0+mod-m)%mod;
			*a0=(*a0+m)%mod;
		  }
	    }
      }
      if(tag!=1) for(register int t=0,w=ksm(len,mod-2);t<len;++t)
		       a[t]=1ll*a[t]*w%mod;
}

void INV(int*a,int*b,const int&len){
      if(len==1){b[0]=ksm(a[0],mod-2);return;}
      INV(a,b,len>>1);
      static int A[maxn],B[maxn];
      for(register int t=0;t<len<<1;++t) A[t]=B[t]=0;
      for(register int t=0;t<len;++t) A[t]=a[t],B[t]=b[t];
      NTT(A,len<<1,1);NTT(B,len<<1,1);
      for(register int t=0;t<len<<1;++t) A[t]=1ll*A[t]*B[t]%mod*B[t]%mod;
      NTT(A,len<<1,-1);
      for(register int t=0;t<len;++t) b[t]=((b[t]+b[t])%mod-A[t]+mod)%mod;
}

inline void inter(int*a,int*b,const int&len){
      for(register int t=len-1;t;--t)
	    b[t]=1ll*a[t-1]*ksm(t,mod-2)%mod;
      b[0]=0;
}

inline void dev(int*a,int*b,const int&len){
      for(register int t=0;t<len-1;++t)
	    b[t]=1ll*(t+1)*a[t+1]%mod;
      b[len-1]=0;
}

inline void LN(int*a,int*b,const int&len){
      static int A[maxn],B[maxn];
      for(register int t=0;t<len<<1;++t) A[t]=B[t]=0;
      INV(a,A,len);
      dev(a,B,len);
      NTT(A,len<<1,1);    NTT(B,len<<1,1);
      for(register int t=0;t<len<<1;++t) B[t]=1ll*B[t]*A[t]%mod;
      NTT(B,len<<1,-1);
      inter(B,b,len);
}


inline void pre(const int&n){
      jc[0]=inv[0]=1;
      for(register int t=1;t<=n;++t)
	    jc[t]=1ll*jc[t-1]*t%mod;
      inv[n]=ksm(jc[n],mod-2);
      for(register int t=n-1;t;--t)
	    inv[t]=1ll*inv[t+1]*(t+1)%mod;
}


int t1[maxn],t2[maxn],n,k;
int main(){
      n=qr();
      pre(n);
      int k=1;
      while(k<=n+1)k<<=1;
      for(register int t=0;t<k;++t) t1[t]=1ll*ksm(2,1ll*t*(t-1ll)>>1)*inv[t]%mod;
      LN(t1,t2,k);
      printf("%lld\n",1ll*t2[n]*jc[n]%mod);
      return 0;
}