jzoj3303. 【集训队互测2013】城市规划

Description

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

Input

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

Output

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

Sample Input

输入1：
3

输入2：
4

输入3：
100000

Sample Output

输出1：
4

输出2：
38

输出3：
829847355

Data Constraint

对于20%的数据, n <= 10
对于40%的数据, n <= 1000
对于60%的数据, n <= 30000
对于80%的数据, n <= 60000
对于100%的数据, n <= 130000

题解

题目疯狂暗示NTT。
然鹅比赛时不会，基本属于弃疗状态。
比赛完后疯狂学习FFT与NTT。
总算有个了解了。
不懂NTT的点这个：https://blog.csdn.net/HiChocolate/article/details/95023398

现在我们了解NTT是什么东东了。
回到这一题：
10%
直接暴力
30%
首先我们考虑DP。
直接DP求答案不是那么好搞，考虑容斥
设$f_{[i]}$表示当前图大小为i的方案数。
方程长这样：
$f_{[n]}=2^{C_n^2}-sum_{i=1}^{n-1}f_{[i]}*C_{n-1}^{i-1}*2^{C_{n-i}^2}$
为什么，请诸君一一听我细细道来：
前面的$2^{C_n^2}$表示当前无向图中所有的边，连接或不连接的方案数。
在sigma里面表示：我们钦定一个点在一个合法的联通图内，然后其他的点随意选择。
十分显然明了。
时间复杂度：$O(n^2)$
话说某xzb强行把$O(n^2)$的方法卡过了，%%%orz。
100%
我们看到上面的柿子。
考虑把那个丑陋的$C_{n-1}^{i-1}$给拆掉
得到：
$f_{[n]}=2^{C_n^2}-sum_{i=1}^{n-1}f_{[i]}*frac{(n-1)!}{(i-1)!(n-i)!}*2^{C_{n-i}^2}$
移项得：
$f_{[n]}=2^{C_n^2}-(n-1)!*sum_{i=1}^{n-1}frac{f_{[i]}}{(i-1)!}*frac{2^{C_{n-i}^2}}{(n-i)!}$
我们设A=$frac{f_{[i]}}{(i-1)!}$B=$frac{2^{C_{n-i}^2}}{(n-i)!}$
那么上面这个柿子可以表示成一个卷积的形式。
NTT打法好！
但是怎么套呢？
两个方法：
cdq分治
我们发现，当前的n都是有前面的i得到的。
那么考虑分治。
每次分治一个区间[l,mid]把这个区间的代价加到[mid+1,r]里面去。
怎么加？
将[l,mid]的A提出来，与一个长度为r-l的B跑一边NTT，得到一个新的C，加入[mid+1,r]里面。
对于一个单点，直接求解答案即可。
思路清晰，实现恶心。
多项式求逆
还是画柿子。
稍微变一下：
$frac{2^{C_n^2}}{(n-1)!}=sum_{i=1}^{n}frac{f_{[i]}}{(i-1)!}*frac{2^{C_{n-i}^2}}{(n-i)!}$
这个是个更加明显的卷积形式。
设$A=frac{2^{C_n^2}}{(n-1)!}B=frac{f_{[i]}}{(i-1)!}C=frac{2^{C_{n-i}^2}}{(n-i)!}$
只不过里面要求的B在卷积里面。
怎么办？
我们换一下即可知道，我们求出C的逆序列即可再次利用NTT。
这叫做多项式求逆。
这里有个blog，自己学。（我不想打）
https://www.cnblogs.com/yoyoball/p/8724115.html

代码

cdq分治的（跑得贼慢）

{$inline on}
uses math;
type
        arr=array[-10..262144] of int64;
var
        i,j,k,l,n,m,lim:longint;
        rev:int64;
        mo:int64=1004535809;
        f,g,w,a,b,jc,ny:arr;
function qsm(a,b:int64):int64;
var
        t,y:int64;
begin
        t:=1;
        y:=a;
        while b>0 do
        begin
                if b and 1=1 then t:=t*y mod mo;
                y:=y*y mod mo;
                b:=b div 2;
        end;
        exit(t);
end;
procedure dft(var a:arr;n:longint);inline;
var
        i,j,k,now:longint;
        x,y:int64;
        r:arr;
begin
        k:=n-1;j:=0;
        while k>0 do
        begin
                k:=k div 2;
                inc(j);
        end;
        r[0]:=0;
        for i:=0 to n-1 do
        begin
                r[i]:=(r[i div 2] div 2) or ((i and 1) shl (j-1));
        end;
        for i:=0 to n-1 do
        begin
                if i<r[i] then
                begin
                        j:=a[i];a[i]:=a[r[i]];a[r[i]]:=j;
                end;
        end;
        now:=2;
        while now<=n do
        begin
                i:=now div 2;
                for j:=0 to i-1 do
                begin
                        k:=j;
                        while k<n do
                        begin
                                x:=a[k];y:=a[k+i]*w[j*(n div now)] mod mo;
                                a[k]:=(x+y) mod mo;
                                a[k+i]:=(x-y+mo) mod mo;
                                k:=k+now;
                        end;
                end;
                now:=now*2;
        end;
end;
procedure fft(var a,b:arr;n:longint);
var
        i:longint;
        j,k:int64;
begin
        rev:=qsm(3,(mo-1) div n);
        w[0]:=1;
        for i:=1 to n do w[i]:=w[i-1]*rev mod mo;
        dft(a,n);
        dft(b,n);
        for i:=0 to n-1 do a[i]:=a[i]*b[i] mod mo;
        for i:=0 to n div 2 do
        begin
                j:=w[i];w[i]:=w[n-i];w[n-i]:=j;
        end;
        dft(a,n);
        k:=qsm(n,mo-2);
        for i:=0 to n-1 do
        begin
                a[i]:=a[i]*k mod mo;
        end;
end;
procedure solve(l,r:longint);
var
        i,j,k,mid:longint;

begin
        if l=r then
        begin
                if r=44824 then
                mid:=mid;
                f[l]:=(g[l]-f[l]*jc[l-1] mod mo+mo) mod mo;
        end
        else
        begin
                mid:=(l+r) div 2;
                solve(l,mid);
                if mid-l+1>r-l then j:=mid-l+1
                else j:=r-l;
                k:=1;
                while k<=j do k:=k*2;
                for i:=0 to 2*k do
                begin
                        a[i]:=0;
                        b[i]:=0;
                end;
                for i:=l to mid do a[i-l+1]:=f[i]*ny[i-1] mod mo;
                for i:=1 to r-l do b[i]:=g[i]*ny[i] mod mo;
                j:=max((mid-l+1)*2,r-l+1);
                fft(a,b,j*2);
                for i:=mid+1 to r do
                begin
                        f[i]:=(f[i]+a[i-l+1]) mod mo;
                end;
                solve(mid+1,r);
        end;
end;
begin
        //assign(input,'4data.in');reset(input);
        g[0]:=1;
        for i:=1 to 262144 do g[i]:=qsm(2,((int64(i-1)*i div 2) mod (mo-1)) mod mo);
        jc[0]:=1;
        for i:=1 to 262144 do jc[i]:=jc[i-1]*i mod mo;
        ny[262144]:=qsm(jc[262144],mo-2);
        for i:=262143 downto 0 do ny[i]:=ny[i+1]*(i+1) mod mo;
        readln(n);
        lim:=1;
        while lim<n+m do lim:=lim*2;
        solve(1,lim);
        writeln(f[n]);
end.
end.