bzoj1488 [HNOI2009]图的同构（群论+组合数学+polya）

感觉我的排列组合需要回炉重造一下了qaq

求不同构的n个点的简单图的种类数。

因为每条边有存在或不存在两种选择，所以可以转化为完全图上对边的2染色计数。但是我们只能处理点的置换（就是个n次对称群），而对于边的置换我们就很gg了。我们考虑对于一个点的置换，能否直接得出这种情况下的边的置换。

每条边无非两种情况：
1、两点处在同一循环节中
2、两点处在不同循环节中
我们考虑一个长度为x的循环节内部，边的循环节个数为n/2.
（个人理解是$\lceil n*(n-1)/2/n\rceil$，即边的个数/n上取整，可以画画图找规律）
我们考虑两个长度分别为x,y的循环节之间的边，循环节个数为gcd(x,y)
（个人理解是$xy/lcm(x,y)$，即边数/边循环节长度）

因此我们可以直接由一个点置换得到边置换的循环节个数m，那么不动点个数就是$2^m$，因此根据polya原理，最后的答案就是$\frac{\sum2^m}{n!}$

我们考虑枚举点置换，是O(n!)的，gg

我们发现边置换的循环节个数其实只与点置换的循环节长度有关，因此我们考虑枚举点置换的循环节长度，即枚举n的划分，我们钦定它有序，从小到大枚举，假设我们现在枚举了n的一个划分为l1,l2,…lk,长度为i的循环节个数有si个，则这种划分所代表的点置换的个数为

$$\frac{n!}{l_1*l_2*...*l_k*s_1!*s_2!...s_q!}$$
除以$l_i$是因为每一个循环节的重复圆排列，除以$s_i$是因为相同大小的循环节的重复排列。
然后此题就解决了orz，复杂度大概是n的有序整数拆分数，极限情况f(60)=966467,大概还行。

#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
using namespace std;
#define inf 0x3f3f3f3f
#define ll long long
#define mod 997
#define N 70
inline char gc(){
    static char buf[1<<16],*S,*T;
    if(S==T){T=(S=buf)+fread(buf,1,1<<16,stdin);if(T==S) return EOF;}
    return *S++;
}
inline int read(){
    int x=0,f=1;char ch=gc();
    while(ch<'0'||ch>'9'){if(ch=='-')f=-1;ch=gc();}
    while(ch>='0'&&ch<='9') x=x*10+ch-'0',ch=gc();
    return x*f;
}
int n,inv[N],fac[N],ifac[N],len[N],tot=0,num[N],ans=0;
inline void init(){
    inv[1]=1;fac[1]=1;ifac[1]=1;
    for(int i=2;i<=n;++i) inv[i]=inv[mod%i]*(mod-mod/i)%mod;
    for(int i=2;i<=n;++i) fac[i]=fac[i-1]*i%mod,ifac[i]=ifac[i-1]*inv[i]%mod;
}
inline int ksm(int x,int k){
    int res=1;for(;k;k>>=1,x=x*x%mod) if(k&1) res=res*x%mod;return res;
}
inline int gcd(int x,int y){return y?gcd(y,x%y):x;}
inline void dfs(int val,int left){//枚举n的划分n1+n2+...=n,钦定n1<=n2<=...,该填val了,还剩left可以划分
    if(left==0){
        int m=0,cnt=fac[n];
        for(int i=1;i<=tot;++i){
            m+=len[i]/2*num[i]+num[i]*(num[i]-1)/2*len[i];
            for(int j=i+1;j<=tot;++j)
                m+=num[i]*num[j]*gcd(len[i],len[j]);
            (cnt*=ksm(inv[len[i]],num[i])*ifac[num[i]])%=mod;
        }ans+=ksm(2,m)*cnt;ans%=mod;return;
    }if(val>left) return;dfs(val+1,left);
    for(int i=1;i*val<=left;++i){
        len[++tot]=val;num[tot]=i;dfs(val+1,left-i*val);--tot;
    }
}
int main(){
//  freopen("a.in","r",stdin);
    n=read();init();dfs(1,n);
    printf("%d\n",ans*ifac[n]%mod);
    return 0;
}