tag：polya，组合计数

---

很显然要用到polya。

此题中的变换可以理解为，枚举一个排列 \(p\)，把 \(i\) 变成 \(p_i\)。

那么根据polya，有

\[ans=\sum_{p}m^\text{等价类个数} \]

那么此题中的等价类是什么意思呢。

注意到我们要求的是对边进行染色的方案，所以这里的等价关系指的就是 \((u,v)\) 和 \((p_u,p_v)\)。

---

观察数据范围，根据常用套路，我们考虑将置换拆成若干个循环，那么一对长度为 \(a,b\) 的循环能够贡献多少个等价类呢？

以 \(2,3\) 为例子

\[\begin{matrix}1&2&3&1&2&3\\1&2&1&2&1&2\end{matrix} \]

可以发现一个等价类的大小为 \(lcm(a,b)\)，而一共有 \(ab\) 条边，所以等价类就有 \(\frac{ab}{lcm(ab)}=\gcd(a,b)\) 个。

然后考虑一个循环内部能贡献几个等价类。

以 \(5\) 为例子

\[\begin{matrix}1&2&3&4&5\\2&3&4&5&1\end{matrix} \]

\[\begin{matrix}1&2&3&4&5\\3&4&5&1&2\end{matrix} \]

注意到每个等价类会把 \((1,i+1)\) 和 \((1,n-i+1)\) 归为一类，所以实际上就有 \(\lceil\frac{n-1}2\rceil\) 个等价类。

---

实现的时候 \(53\) 的拆分数是 \(329931\)，所以暴力枚举拆分就好了。

然后分 \(3\) 种情况：

• 长度不同的组：\(cnt_i\cdot cnt_j\cdot\gcd(a_i,a_j)\)
• 长度相同的组：\(\binom{cnt_i}2a_i\)
• 组内贡献：\(cnt_i\lceil\frac{a_i-1}2\rceil\)

然后还要乘上这种拆分对应的方案数：

先选出每组的数，然后乘上组内方案，同长度的组要去重。

\[\frac{n!}{\prod(len_i!)}\cdot\prod(len_i-1)\cdot\prod\frac1{cnt_i!}=\frac{n!}{\prod len_i\prod cnt_i!} \]

\(cnt_i\) 表示长度为 \(i\) 的循环有 \(cnt_i\) 个。

由于最后要除以 \(n!\)，所以这里的 \(n!\) 可以扔掉。

---

#include<bits/stdc++.h>
using namespace std;
   
template<typename T>
inline void Read(T &n){
    char ch; bool flag=false;
    while(!isdigit(ch=getchar()))if(ch=='-')flag=true;
    for(n=ch^48;isdigit(ch=getchar());n=(n<<1)+(n<<3)+(ch^48));
    if(flag)n=-n;
}

int n, m, p;

inline int ksm(int base, int k=p-2){
    int res=1;
    while(k){
        if(k&1)
            res = 1ll*res*base%p;
        base = 1ll*base*base%p;
        k >>= 1;
    }
    return res;
}

int a[54], cnt[54], top, jc[54], invjc[54], inv[54];
int gcd[54][54], ans;

inline void solve(){
    int res=0;
    for(int i=1; i<=top; i++) for(int j=i+1; j<=top; j++) res += cnt[i]*cnt[j]*gcd[a[i]][a[j]];
    for(int i=1; i<=top; i++) res += cnt[i]*(cnt[i]-1)/2*a[i];
    for(int i=1; i<=top; i++) res += cnt[i]*(a[i]/2);
    res = ksm(m,res);
    int tmp=1;
    for(int i=1; i<=top; i++) tmp = 1ll*tmp*invjc[cnt[i]]%p*ksm(inv[a[i]],cnt[i])%p;
    ans = (ans+1ll*res*tmp)%p;
}

void dfs(int rem, int prv){
    if(!rem) return solve();
    if(prv and (rem-prv>=prv or rem==prv)) cnt[top]++, dfs(rem-prv,prv), cnt[top]--;
    top++;
    for(int i=prv+1; i<=rem-i; i++){
        a[top] = i; cnt[top] = 1;
        dfs(rem-i,i);
    }
    if(rem>prv) a[top] = rem, cnt[top] = 1, dfs(0,rem);
    a[top] = 0; cnt[top] = 0; top--;
}

int Gcd(int a, int b){return b?Gcd(b,a%b):a;}

int main(){
    Read(n); Read(m); Read(p);
    for(int i=1; i<=n; i++) for(int j=1; j<=n; j++) gcd[i][j] = Gcd(i,j);
    for(int i=1; i<=n; i++) inv[i] = ksm(i);
    jc[0] = 1; for(int i=1; i<=n; i++) jc[i] = 1ll*jc[i-1]*i%p;
    invjc[n] = ksm(jc[n]); for(int i=n; i; i--) invjc[i-1] = 1ll*invjc[i]*i%p;
    dfs(n,0);
    cout<<ans<<'\n';
    return 0;
}