[BZOJ3456]-城市规划-分治NTT

说在前面

这题有复杂度更优秀的做法，比如：Miskcoo的多项式求逆与wzq_QwQ的多项式求$ln$（os:什么？还可以求ln？）
当然me都不会上面的，于是式子推到一半就直接硬上了hhhhh。然而狂T不止，优化一波常数卡到了18s
（果然还是优秀做法比较友好）

---

题目

BZOJ3456传送门

题目大意

求出n个点带标号的简单无向连通图的个数
答案对 1004535809 取模
$n\leq 130000$

输入输出格式

输入格式：
仅一个数字n，含义如题

输出格式：
输出一个数字，表示答案

---

解法

使用惯用式 图的计数类dp（以下均默认带标号）
假设$f[i]$表示$i$个点的连通图的个数，$g[i]$表示连通图的个数

显然有$g[i] = 2^{C_i^2}$，即每条边都可以选或不选

然后我们来统计$f[i]$，为了不重复计数，需要一点小技巧。我们每次枚举和「编号最小的那个点」连通的点是哪一些，这样总不会算重复

$f[i]=g[i] - \sum_{j=1}^if[j]*C^{j-1}_{i-1}*g[i-j]$（枚举上界为$i-1$也可以，因为当$j$取i时值为0）
补集转化，考虑一定不连通的图。除开「编号最小点」，从$i-1$个点中选出$j-1$个点和「编号最小点」在同一连通块的方案数，记这一部分为$S$，则$|S|=f[j]*C^{j-1}_{i-1}$。剩下的部分，在内部随意连边，记这一部分为$G$，则$|G|=g[i-j]$。
因为$G$和$S$之间没有边，所以统计不会重复；因为枚举了$S$的大小，所以统计不会有遗漏

然后答案就是$f[N]$，我们进一步化简上面的式子，得到$f[i]=g[i]-(i-1)!\sum_{j=1}^{i-1}\frac{f[j]}{(j-1)!} \frac{ 2^{C_{i-j}^2}} {(i-j)!}$
发现它是一个，自己和自己卷积的形式，于是我们使用分治NTT来完成计算

优化常数小技巧

1. 利用循环卷积。记当前处理的区间长度为$l$，我们实际上只需要卷积的$0.5l$ 到 $l$的部分，而其它的部分都是废的…所以我们没有必要把NTT的开到两倍区间长度，而只用开到一倍（$f$的部分是$0.5l$，而$g$的部分是$l$，所以卷积次数最高是$1.5l$。我们把区间开到刚好大于$l$，那么多出去的部分只会和$0$到$0.5l$发生混叠，不会影响到我们需要的部分）
2. 优化取模，在NTT里面主要是加减法运算，而仅仅是加减法是无法超过long long的。所以我们只在有乘法的部分取模。在最后把所有数字都取模即可
3. 小范围暴力，当区间长度很小的时候，没有必要再调用NTT，因为有常数，所以暴力跑得比NTT还快……

---

下面是自带大常数的代码

代码里其实还有很多可以优化的地方
不过懒得搞了…

#include <ctime>
#include <cstdio>
#include <cstring>
#include <algorithm>
using namespace std ;

const int P = 1004535809 , g = 3 ;
int N ;
long long fac[130005] , ifac[130005] , mi2C[130005] ;
long long F[270005] , G[270005] , ans[130005] ;

int s_pow( long long x , int b ){
    long long rt = 1 ;
    while( b ){
        if( b&1 ) rt = rt * x %P ;
        x = x * x %P , b >>= 1 ;
    } return rt ;
}

void getG(){
    int yx[1005] , ycnt = 0 , t = P - 1 ;
    for( int i = 2 ; i * i <= t ; i ++ ){
        if( t % i == 0 ){
            yx[++ycnt] = i ;
            while( !t%i ) t /= i ;
        }
    } if( t != 1 ) yx[++ycnt] = t ;
    for( int i = 2 ; i ; i ++ ){
        bool ok = true ;
        for( int j = 1 ; j <= ycnt ; j ++ )
            if( s_pow( i , ( P - 1 ) / yx[j] ) == 1 ){
                ok = false ; break ;
            }
        if( ok ) return ( void ) printf( "g = %d" , i ) ;
    }
}

void preWork(){
    fac[0] = ifac[0] = 1 ;
    for( int i = 1 ; i <= N ; i ++ ) fac[i] = fac[i-1] * i %P ;
    ifac[N] = s_pow( fac[N] , P - 2 ) ;
    for( int i = N - 1 ; i ; i -- ) ifac[i] = ifac[i+1] * ( i+1 ) %P ;
    for( int i = 1 ; i <= N ; i ++ )
        ans[i] = mi2C[i] = s_pow( 2 , 1LL * i * ( i-1 ) / 2 %( P-1 ) ) ;
}

void NTT( long long *a , int len , int rev ){
    for( int i = 1 , now = 0 , t ; i < len ; i ++ ){
        t = len >> 1 ;
        while( now&t ) now ^= t , t >>= 1 ; now ^= t ;
        if( now > i ) swap( a[now] , a[i] ) ;
    }

    for( int i = 2 ; i <= len ; i <<= 1 ){
        long long g = s_pow( ::g , ( P - 1 ) / i ) , nowg ;
        if( rev == -1 ) g = s_pow( g , P - 2 ) ;
        for( int j = 0 ; j < len ; j += i ){
            nowg = 1 ;
            for( int k = 0 ; k < ( i >> 1 ) ; k ++ ){
                long long t1 = a[j+k+(i>>1)] %P * nowg %P ;
                a[j+k+(i>>1)] = a[j+k] - t1 ;
                a[j+k] += t1 ;
                nowg = nowg * g %P ;
            }
        }
    }for( int i = 0 , t = ( rev == -1 ? s_pow( len , P - 2 ) : 1 ) ; i < len ; i ++ )
        a[i] = a[i] %P * t %P ;
}

void div_and_conquer( int lf , int rg ){
    if( rg - lf + 1 <= 100 ){
        long long tmp ;
        for( int i = lf + 1 ; i <= rg ; i ++ ){
            tmp = 0 ;
            for( int j = lf ; j < i ; j ++ )
                tmp += ifac[j-1] * ans[j] %P * mi2C[i-j] %P * ifac[i-j] %P ;
            ans[i] = ( ans[i] - tmp %P * fac[i-1] %P )%P ;
        } return ;
    } int mid = ( lf + rg ) >> 1 , len , t = rg - lf + 1 ;
    div_and_conquer( lf , mid ) ;

    for( len = 1 ; len <= t ; len <<= 1 ) ;
    memset( F , 0 , sizeof( long long ) * len ) ;
    memset( G , 0 , sizeof( long long ) * len ) ;
    for( int i = 1 , j = lf ; i < t ; i ++ , j ++ ){
        if( j <= mid ) F[i] = ifac[j-1] * ans[j] %P ;
        G[i] = mi2C[i] * ifac[i] %P ;
    } NTT( F , len , 1 ) , NTT( G , len , 1 ) ;
    for( int i = 0 ; i < len ; i ++ )
        F[i] = F[i] * G[i] %P ;
    NTT( F , len , -1 ) ;
    for( int i = mid + 1 ; i <= rg ; i ++ )
        ans[i] = ( ans[i] - fac[i-1] * F[i-lf+1]%P )%P ;

    div_and_conquer( mid+1 , rg ) ;
}

void solve(){
    div_and_conquer( 1 , N ) ;
    printf( "%lld\n" , ( ans[N] + P )%P ) ;
}

int main(){
//  getG() ;
    scanf( "%d" , &N ) ;
//  N = 130000 ;
//  int a = clock() ;
    preWork() ; solve() ;
//  printf( "%d" , clock() - a ) ;
}