BZOJ[1815]男人一上午写一道题，机房人不解；看到内容，众人沉默，网友：长见识了

传送门ber~

详细题解

知道点的顺序，边的顺序也能确定下来，所以对于点置换可以算出相应的边置换计算
可我们总不能$n!$枚举全排列，我们可以进一步发现，结构相同的点置换对应的边置换肯定是一样的(废话)，可以搜出所有不同结构的点置换，发现$n=53$时点置换个数不到$30w$

当前搜到点置换$L_1\le L_2\le L_3\le ...\le L_k$，那么需要计算边置换循环个数
连接不同点循环循环i,j的边：一个循环覆盖$Lcm(L_i,L_j)$条边，所以$c=\frac{L_i\ast L_j}{Lcm(L_i,L_j)}=gcd(L_i,L_j)$个循环(画画图更容易理解一点)
点循环$i$内部的边：循环数为$\frac{L_i-1}{2}$(因为一半会重)
所以总循环数$C=\sum \lfloor \frac{L_i}{2}\rfloor +\sum \sum gcd(L_i,L_j)$
$L_1\le L_2\le L_3\le ...\le L_k$循环数量容易推得为$\frac{n!}{L_1{L}_2...L_n{B}_1!B_2!...B_n!}$，其中$B_i$表示循环长度为$i$的循环个数
剩下的套Polya就行了

代码如下：

#include<bits/stdc++.h>
#define int long long
#define N 120
using namespace std;
inline int read(){
    int x=0,f=1;char c;
    do c=getchar(),f=c=='-'?-1:f; while(!isdigit(c));
    do x=(x<<3)+(x<<1)+c-'0',c=getchar(); while(isdigit(c));
    return x*f;
}
typedef long long LL;
int n,m,p,s[N],a[N],jc[N],gd[N][N];
LL ans;
LL gcd(LL a,LL b){
    return b?gcd(b,a%b):a;
}
inline LL qpow(LL x,int k){
    LL sum=1;
    while(k){
        if(k&1) sum=sum*x%p;
        x=x*x%p;
        k>>=1;
    }
    return sum;
}
void dfs(int x,int cnt,int pre){
    if(cnt==n){
        LL t=0,pp=1;
        for(int i=1;i<x;i++)
            t+=s[i]/2;
        for(int i=1;i<x;i++)
            for(int j=i+1;j<x;j++)
                t+=gd[s[i]][s[j]];
        for(int i=1;i<=n;i++) (pp*=jc[a[i]])%=p;
        for(int i=1;i<x;i++) (pp*=s[i])%=p;
       // cout<<pp<<" "<<t<<endl;
        (ans=ans+jc[n]*qpow(pp,p-2)%p*qpow(m,t)%p)%=p;
        return;
    }
    for(int i=pre;i<=n-cnt;i++){
        s[x]=i;a[i]++;
        dfs(x+1,cnt+i,i);
        a[i]--;
    }
}
main(){
    n=read();m=read();p=read();
    for(int i=jc[0]=1;i<=n;i++) jc[i]=jc[i-1]*i%p;
    for(int i=1;i<=n;i++)
        for(int j=i;j<=n;j++) gd[i][j]=gd[j][i]=gcd(i,j);
    dfs(1,0,1);
    printf("%lld\n",ans*qpow(jc[n],p-2)%p);
    return 0;
}