一道burnside，做了几次之后感觉比较经典了。然后发现网上居然没有题解。

枚举正整数拆分，既每个质换的长度。然后枚举gcd算每类置换的贡献。再用错排公式来算每类置换的数量。然后求和。

#include <cstdio>
#include <cstring>
#include <algorithm>

using namespace std;

#define _l (long long int)

const int maxn = 63;

int n, m, mod, ans, g[maxn][maxn];
int sq[maxn], tq, fac[maxn], finv[maxn], vinv[maxn];

int modPow(int a, int x) {
	int s(1);
	for (; x; x >>= 1, a = _l a * a % mod)
		if (x & 1)
			s = _l s * a % mod;
	return s;
}

int gcd(int a, int b) {
	while (a) {
		register int c(b % a);
		b = a;
		a = c;
	}
	return b;
}

void addAns() {
	int t(0);
	for (int i = 1; i <= tq; ++ i) {
		t += sq[i] >> 1;
	//	if (sq[i] > 1)
	//		++ t;
		for (int j = 1; j < i; ++ j)
			t += g[sq[i]][sq[j]];
	}
	t = modPow(m, t);
	for (int i = 1, j; i <= tq; i = j) {
		for (j = i; j <= tq && sq[i] == sq[j]; ++ j)
			t = _l t * vinv[sq[j]] % mod;
		t = _l t * finv[j - i] % mod;
	}
	ans = (ans + t) % mod;
}

void dfs(int v, int l) {
	if (!v)
		addAns();
	else {
		++ tq;
		for (int i = 1; i <= l && i <= v; ++ i) {
			sq[tq] = i;
			dfs(v - i, i);
		}
		-- tq;
	}
}

int main() {
#ifndef ONLINE_JUDGE
	freopen("in.txt", "r", stdin);
#endif

	scanf("%d%d%d", &n, &m, &mod);
	fac[0] = finv[0] = 1;
	for (int i = 1; i <= n; ++ i) {
		fac[i] = _l fac[i - 1] * i % mod;
		finv[i] = modPow(fac[i], mod - 2);
		vinv[i] = modPow(i, mod - 2);
		for (int j = 1; j <= n; ++ j)
			g[i][j] = gcd(i, j);
	}
	tq = 0;
	dfs(n, n);
	printf("%d\n", ans);
}

首先是burnside，对于旋转x下都要当作一个置换，所以总共有n个。然后第i个的不动点数量等于n^gcd(n,i)，然后统计一下个数是phi(gcd(n,i))个。

然后发现要用奇奇怪怪的东西来求phi。我用的方法是分解质因数然后DFS，跑得飞快。

#include <cstdio>
#include <cstring>
#include <algorithm>

using namespace std;

#define _l (long long int)

const int maxn = 50009;

int n, mod;
int tp, pn[maxn], phi[maxn];
bool pr[maxn];

void pre() {
	memset(pr, 0, sizeof(pr));
	tp = 0;
	phi[1] = 1;
	for (int i = 2; i < maxn; ++ i) {
		if (!pr[i]) {
			pn[tp ++] = i;
			phi[i] = i - 1;
		}
		for (int j = 0; j < tp && i * pn[j] < maxn; ++ j) {
			pr[i * pn[j]] = 1;
			if (i % pn[j] == 0) {
				phi[i * pn[j]] = phi[i] * pn[j];
				break;
			}
			else
				phi[i * pn[j]] = phi[i] * phi[pn[j]];
		}
	}
}

int modPow(int a, int x) {
	int s(1);
	for (a %= mod; x; x >>= 1, a = _l a * a % mod)
		if (x & 1)
			s = _l s * a % mod;
	return s;
}

int fv[maxn], fc[maxn], tv, ans;

void dfs(int p, int pf, int pv) {
	if (p == tv)
		ans = (ans + _l modPow(n, pv - 1) * pf) % mod;
	else {
		dfs(p + 1, pf, pv);
		pf *= fv[p] - 1;
		pv /= fv[p];
		dfs(p + 1, pf, pv);
		for (int i = 2; i <= fc[p]; ++ i) {
			pf *= fv[p];
			pv /= fv[p];
			dfs(p + 1, pf, pv);
		}
	}
}

int calc() {
	int x(n);
	tv = 0;
	for (int i = 0; i < tp && x > 1; ++ i)
		if (x % pn[i] == 0) {
			fv[tv] = pn[i];
			fc[tv] = 0;
			while (x % pn[i] == 0) {
				++ fc[tv];
				x /= pn[i];
			}
			++ tv;
		}
	if (x > 1) {
		fv[tv] = x;
		fc[tv] = 1;
		++ tv;
	}
	ans = 0;
	dfs(0, 1, n);
	return ans;
}

int main() {
	//freopen("in.txt", "r", stdin);
	int t;
	pre();
	scanf("%d", &t);
	while (t --) {
		scanf("%d%d", &n, &mod);
		if (mod == 1)
			puts("0");
		else
			printf("%d\n", calc());
	}
}