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;
}