BZOJ 3456 城市规划 快速傅里叶变换

题目大意：求n个点的无向简单连通图个数，n≤1.3?105

递推式：fi=2C2i?∑i?1j=1fj?Cj?1i?1?2C2i?j

推导戳这里

然后两侧同除(i?1)!得到：

fi(i?1)!=2C2i(i?1)!?∑i?1j=1fj?2C2i?j(j?1)!?(i?j)!

∑ij=1fj?2C2i?j(j?1)!?(i?j)!=2C2i(i?1)!

∑ij=1fj(j?1)!?2C2i?j(i?j)!=2C2i(i?1)!

这显然是一个卷积的形式

然后我们令：

A=∑ni=1fi(i?1)!xi

B=∑ni=02C2ii!xi

C=∑ni=12C2i(i?1)!xi

那么有A?B=C

然后有A≡C?B?1(mod xn+1)

多项式求逆搞一搞就好了

没想到还真有人的FFT比我还慢2333

#include <cstdio>
#include <cstring>
#include <iostream>
#include <algorithm>
#define M 263000
#define MOD 1004535809
#define G 3
using namespace std;
int n,m;
long long Quick_Power(long long x,long long y)
{
    long long re=1;
    while(y)
    {
        if(y&1) (re*=x)%=MOD;
        (x*=x)%=MOD; y>>=1;
    }
    return re;
}
void FFT(int a[],int n,int type)
{
    static int temp[M];
    int i;
    if(n==1) return ;
    for(i=0;i<n;i+=2)
        temp[i>>1]=a[i],temp[i+n>>1]=a[i+1];
    memcpy(a,temp,sizeof(a[0])*n);
    int *l=a,*r=a+(n>>1);
    FFT(l,n>>1,type);FFT(r,n>>1,type);
    long long w=Quick_Power(G,(long long)(MOD-1)/n*type%(MOD-1)),wn=1;
    for(i=0;i<n>>1;i++,(wn*=w)%=MOD)
        temp[i]=(l[i]+wn*r[i])%MOD,temp[i+(n>>1)]=(l[i]-wn*r[i]%MOD+MOD)%MOD;
    memcpy(a,temp,sizeof(a[0])*n);
}
void Get_Inv(int a[],int b[],int n)
{
    static int temp[M];
    int i;
    if(n==1)
    {
        b[0]=Quick_Power(a[0],MOD-2);
        return ;
    }
    Get_Inv(a,b,n>>1);
    memcpy(temp,a,sizeof(a[0])*n);
    memset(temp+n,0,sizeof(a[0])*n);
    FFT(temp,n<<1,1);
    FFT(b,n<<1,1);
    for(i=0;i<n<<1;i++)
        temp[i]=(long long)b[i]*(2-(long long)temp[i]*b[i]%MOD+MOD)%MOD;
    FFT(temp,n<<1,MOD-2);
    long long inv=Quick_Power(n<<1,MOD-2);
    for(i=0;i<n;i++)
        b[i]=temp[i]*inv%MOD;
    memset(b+n,0,sizeof(a[0])*n);
}
int main()
{
    static long long fac[M];
    static int A[M],B[M],C[M],inv_B[M];
    int i;
    cin>>n;
    for(m=1;m<=n;m<<=1);
    for(fac[0]=1,i=1;i<=n;i++)
        fac[i]=fac[i-1]*i%MOD;
    for(i=0;i<=n;i++)
        B[i]=Quick_Power(2,((long long)i*(i-1)>>1)%(MOD-1))*Quick_Power(fac[i],MOD-2)%MOD;
    for(i=1;i<=n;i++)
        C[i]=Quick_Power(2,((long long)i*(i-1)>>1)%(MOD-1))*Quick_Power(fac[i-1],MOD-2)%MOD;
    Get_Inv(B,inv_B,m);
    FFT(inv_B,m<<1,1);
    FFT(C,m<<1,1);
    for(i=0;i<m<<1;i++)
        A[i]=(long long)inv_B[i]*C[i]%MOD;
    FFT(A,m<<1,MOD-2);
    long long inv=Quick_Power(m<<1,MOD-2);
    for(i=1;i<=n;i++)
        A[i]=A[i]*inv%MOD;
    cout<<A[n]*fac[n-1]%MOD<<endl;
    return 0;
}