bzoj1478/1815[Shoi2006]color 有色图

题目链接：bzoj1478  bzoj1815

题意：
染色图是无向完全图，且每条边可被染成M种颜色中的一种。两个染色图是同构的，当且仅当可以改变一个图的顶点的编号，使得两个染色图完全相同。问N个顶点，M种颜色，本质不同（两两互不同构）的染色图个数（模质数P）。（1<=N<=53，1<=M<=1000，N<P<=109，时间限制10s)

题解：
置换-ploya

双倍经验喔

关于这题，这文档讲得很清楚：http://wenku.baidu.com/view/fee9e9b9bceb19e8b8f6ba7a.html?from=search###
这题想起来挺难的。
首先它是对点的置换，但是是边染上了颜色，就是说实际上是边的置换。所以我们要看一下点置换和边置换之间的关系。
假定一个点置换，把它表示为循环，比如是(a1,a2,....)(b1,b2...)(c1,c2...)...
1、对于不在一个循环里面的点：
比如a1,b1, 那么会有边循环((a1,b1),(a2,b2)...) 设a循环的循环节是l1,b循环的循环节是l2,那么形成的边循环的循环节显然是LCM(l1,l2)。
一共有l1*l2个点对，每个点对都在一个循环节为LCM(l1,l2)的循环上，所以一共有l1*l2/LCM(l1,l2)=GCD(l1,l2)个循环节，所以C(f)=m^GCD(l1,l2)。（回到burnside引理，C为置换之后仍为本身的数目，就是说要循环节里的每条边都一样的颜色）
2、对于在一个循环里面的点：
比如a1、a2。设这个a循环的循环节为l1。
如果l1是奇数，那么循环长度为l1，一共有C(l1,2)个点对，所以是(l1-1)/2个循环节，所以C(f)=m^((l1-1)/2)。
如果l1是偶数，除了上面这种情况之外，还有一种的循环节是l1/2（就是两个点刚好相隔半个周期，而边是双向的），所以一共有(C(l1,2)-l1/2)/l1+1=l1/2个点对。
整理一下：

题解引用hyc的

也许最好什么都预处理一下。比如预处理gcd什么的，效果：时间减半。让我卡时过了..

#include<cstdio>
#include<cstdlib>
#include<cstring>
#include<iostream>
#include<algorithm>
using namespace std;
typedef long long LL;
#define maxn 60
 
LL fac[maxn],l[maxn],gd[maxn][maxn],ans,n,m,mod;
LL gcd(LL a,LL b) {return (b==0)?a:gcd(b,a%b);}
LL qpow(LL x,LL t)
{
    LL ret=1;
    while (t)
    {
        if (t&1) ret=(ret*x)%mod;
        x=(x*x)%mod;t>>=1;
    }return ret;
}
void get_ans(LL cnt)
{
    LL i,j,C=0;
    for (i=1;i<=cnt;i++) C=C+l[i]/2;
    for (i=1;i<cnt;i++)
     for (j=i+1;j<=cnt;j++)
      C=C+gd[l[i]][l[j]];
    LL now=1,tt=0;
    for (i=1;i<=cnt;i++)
    {
        if (l[i]!=l[i-1])
        {
            now=(now*fac[tt])%mod;
            tt=0;
        }tt++;
    }
    now=(now*fac[tt])%mod;
    for (i=1;i<=cnt;i++) now=(now*l[i])%mod;
    LL S=(fac[n]*qpow(now,mod-2))%mod;
    ans=(ans+(S*qpow(m,C))%mod)%mod;
}
void div(LL cnt,LL x,LL ret)
{
    if (ret==0) {get_ans(cnt);return;}
    if (ret<x) return;
    cnt++;
    while (x<=ret)
    {
        l[cnt]=x;
        div(cnt,x,ret-x);
        x++;
    }
}
int main()
{
    LL i,j;
    scanf("%lld%lld%lld",&n,&m,&mod);
    for (i=1;i<=n;i++)
      for (j=i;j<=n;j++) gd[i][j]=gd[j][i]=gcd(i,j);
    fac[0]=1;for (i=1;i<=n;i++) fac[i]=(fac[i-1]*i)%mod;
    ans=0;div(0,1,n);
    printf("%lld\n",(ans*qpow(fac[n],mod-2))%mod);
    return 0;
}