【BZOJ3456】城市规划 多项式求逆

【BZOJ3456】城市规划

Description

 刚刚解决完电力网络的问题, 阿狸又被领导的任务给难住了.
 刚才说过, 阿狸的国家有n个城市, 现在国家需要在某些城市对之间建立一些贸易路线, 使得整个国家的任意两个城市都直接或间接的连通. 为了省钱, 每两个城市之间最多只能有一条直接的贸易路径. 对于两个建立路线的方案, 如果存在一个城市对, 在两个方案中是否建立路线不一样, 那么这两个方案就是不同的, 否则就是相同的. 现在你需要求出一共有多少不同的方案.
 好了, 这就是困扰阿狸的问题. 换句话说, 你需要求出n个点的简单(无重边无自环)无向连通图数目.
 由于这个数字可能非常大, 你只需要输出方案数mod 1004535809(479 * 2 ^ 21 + 1)即可.

Input

 仅一行一个整数n(<=130000)

Output

 仅一行一个整数, 为方案数 mod 1004535809.

Sample Input

3

Sample Output

4

HINT

 对于 100%的数据, n <= 130000

题解：一直打算做来着，前几天终于抽空把这题过了~

设f[n]表示n个点的简单无向连通图数目，先列出大家都会的DP方程：

$f[n]=2^{n(n-1)\over 2}-\sum\limits_{i=1}^{n-1}f[i]C_{n-1}^{i-1}2^{(n-i)(n-i-1)\over 2}$

哎呀好多n-1太烦了，还是令f[n]表示n+1个点的简单无向连通图数目吧：

$f[n]=2^{n(n+1)\over 2}-\sum\limits_{i=0}^{n-1}f[i]C_n^i2^{(n-i)(n-i-1)\over 2}$

组合数很不优美是不是？那就拆了好啦。

$f[n]=2^{n(n+1)\over 2}-n!\sum\limits_{i=0}^{n-1}{f[i]2^{(n-i)(n-i-1)\over 2}\over i!(n-i)!}$

下面我们想把只含有含有n的放到一起，i放到一起，n-i放到一起，于是疯狂移项得到：

${f[n]\over n!}={2^{n(n+1)\over 2}\over n!}-\sum\limits_{i=0}^{n-1}{f[i]\over i!}\cdot{2^{(n-i)(n-i-1)\over 2}\over (n-i)!}$

明朗很多了是不是？令$F(x)={f[x]\over x!},G(x)={2^{x(x+1)\over 2}\over x!},H(x)={2^{x(x-1)\over 2}\over x!}$，特别地，H(0)=0。

所以$F(x)=G(x)-F(x)*H(x)$，$F(x)={G(x)\over 1+H(x)}$，多项式求个逆就好了。

#include <cstdio>
#include <cstring>
#include <iostream>
#include <algorithm>

using namespace std;
typedef long long ll;
const ll P=1004535809;
const int maxn=(1<<19)+4;
int n;
ll A[maxn],B[maxn],C[maxn],D[maxn],ine[maxn],jc[maxn],jcc[maxn];
inline ll pm(ll x,ll y)
{
	ll z=1;
	while(y)
	{
		if(y&1)	z=z*x%P;
		x=x*x%P,y>>=1;
	}
	return z;
}
inline void NTT(ll *a,int len,int f)
{
	int i,j,k,h;
	for(i=k=0;i<len;i++)
	{
		if(i>k)	swap(a[i],a[k]);
		for(j=len>>1;(k^=j)<j;j>>=1);
	}
	ll wn,w,t;
	for(h=2;h<=len;h<<=1)
	{
		if(f==1)	wn=pm(3,(P-1)/h);
		else	wn=pm(3,P-1-(P-1)/h);
		for(i=0;i<len;i+=h)
		{
			w=1;
			for(j=i;j<i+h/2;j++)	t=w*a[j+h/2]%P,a[j+h/2]=(a[j]-t)%P,a[j]=(a[j]+t)%P,w=w*wn%P;
		}
	}
	if(f==-1)
	{
		t=pm(len,P-2);
		for(i=0;i<len;i++)	a[i]=a[i]*t%P;
	}
}
void get_inv(ll *F,ll *G,int len)
{
	if(len==1)
	{
		G[0]=pm(F[0],P-2);
		return ;
	}
	int i;
	get_inv(F,G,len>>1);
	memcpy(C,F,sizeof(F[0])*len);
	memset(C+len,0,sizeof(F[0])*len);
	NTT(C,len<<1,1),NTT(G,len<<1,1);
	for(i=0;i<len<<1;i++)	G[i]=(G[i]*2-C[i]*G[i]%P*G[i])%P;
	NTT(G,len<<1,-1);
	memset(G+len,0,sizeof(F[0])*len);
}
int main()
{
	scanf("%d",&n);
	int len=1,i;
	while(len<=n)	len<<=1;
	ine[0]=ine[1]=jc[0]=jc[1]=jcc[0]=jcc[1]=1;
	for(i=2;i<=len;i++)	ine[i]=P-(P/i)*ine[P%i]%P,jc[i]=jc[i-1]*i%P,jcc[i]=jcc[i-1]*ine[i]%P;
	for(i=0;i<len;i++)	D[i]=pm(2,i*(i+1ll)/2)*jcc[i]%P;
	for(i=1;i<len;i++)	A[i]=pm(2,i*(i-1ll)/2)*jcc[i]%P;
	A[0]=1;
	get_inv(A,B,len);
	NTT(D,len<<1,1),NTT(B,len<<1,1);
	for(i=0;i<len<<1;i++)	B[i]=D[i]*B[i]%P;
	NTT(B,len<<1,-1);
	printf("%lld",(B[n-1]*jc[n-1]%P+P)%P);
	return 0;
}