自动搬运

查看原文

来自洛谷，原作者为

	Re_Star 博客：https://www.cnblogs.com/Re-Star

搬运于2025-08-24 22:59:29，当前版本为作者最后更新于2025-03-26 11:19:37，作者可能在搬运后再次修改，您可在原文处查看最新版

自动搬运只会搬运当前题目点赞数最高的题解，您可前往洛谷题解查看更多

以下是正文

---

题面

更好的阅读体验

其实是参考了另一篇题解，但感觉那篇题解很多细节讲的不是很清楚，所以自己写了一篇题解加深印象。

我们假设有 $k$ 种合法的 $n$ 个点的图，现在要对每个图染一种颜色，求方案数。所以答案是 $m^k$。

现在我们要求 $k$ 可以枚举每个子图的大小 $d$,那么将 $n$ 个点每 $d$ 个为一组进行分组的方案数为 ${n\choose d,d,\dots}$，但是这样我们是钦定了顺序的，而每个图内染一种颜色显然是不关心顺序的，所以每种方案都算了 $(\frac nd)!$，所以 $k=\sum\limits_{d|n}\frac{n!}{(d!)^{\frac nd}(\frac nd)!}$。

最终答案为 $m^{\sum\limits_{d|n}\frac{n!}{(d!)^{\frac nd}(\frac nd)!}} \bmod 10^9-401$。

跟据费马小定理，我们实际要求 $m^{\sum\limits_{d|n}\frac{n!}{(d!)^{\frac nd}(\frac nd)!}\bmod 10^9-402}\bmod 10^9-401$，关键是指数部分怎么求。首先对 $10^9-401$ 分解质因数，$10^9-402=2\times13\times5281\times7283$，那么我们可以先对每一个质因数算出取模后的答案再用 CRT 合并。

于是我们转化为解决子问题：$\sum\limits_{d|n}\frac{n!}{(d!)^{\frac nd}(\frac nd)!}\bmod p$，这里问题在于 $d!$，$(\frac nd)!$ 中可能含有因子 $p$ 无法做逆元操作，于是我们用类似 ExLucas 的思路处理一下。

$$\frac{n!}{(d!)^{\frac nd}(\frac nd)!}\bmod p=\frac{\frac{n!}{p^i}}{(\frac{d!}{p^j})^{\frac nd}\frac{(\frac nd)!}{p^k}}\times p^{i-j\frac nd-k}\bmod p $$

那么我们现在需要求的就是形如 $\frac{n!}{p^x}\bmod p$，我们首先将 $1$~$n$ 中 $p$ 的倍数提一个 $p$ 的因子,剩下的按模 $p$ 的模数分组，最后有可能还会剩下几个无法分为完整的一组，即 $n!=p^{\lfloor\frac np\rfloor}（\lfloor\frac np\rfloor）!\times(p-1)!^{\lfloor\frac np\rfloor}\times\prod\limits_{i=p\lfloor\frac np\rfloor+1}^{i<=n}i$ 到这一步我们就可以递归处理下去了。除去所以因子 $p$ 后我们就可以求逆元了，那么这道题就做完了。

code

#include<bits/stdc++.h>
#define ll long long
using namespace std;
const ll mod=1e9-401,N=1e4+5,p[5]={0,2,13,5281,7283},M=1e9-402;
ll qp(ll x,ll y,ll P)
{
	ll res=1;
	while(y)
	{
		if(y&1) (res*=x)%=P;
		(x*=x)%=P,y>>=1;
	}
	return res;
}
void exgcd(ll a,ll b,ll &x,ll &y){
	ll t;
	if(!b){x=1,y=0;return;}
	exgcd(b,a%b,x,y);
	t=x,x=y,y=t-a/b*y;
}
ll inv(ll a,ll P)
{
	ll x,y;
	exgcd(a,P,x,y);
	return (x+P)%P;
}
ll g(ll i,ll P){return (i<P)?0:(g(i/P,P)+i/P);}
ll n,m,ans[5],fac[N],res;
ll S(ll now, ll P)
{
	if(now<P) return fac[now]%P;
	ll res=1;
	for(ll i=(now/P)*P+1;i<=now;i++) (res*=i)%=P;
	return qp(fac[P-1]%P,now/P,P)*S(now/P,P)%P*res%P;
}
ll solve(ll now,ll d,ll P)
{
	if(g(now,P)-(now/d)*g(d,P)-g(now/d,P)) return 0;
	return S(now,P)*inv(qp(S(d,P),now/d,P),P)%P*inv(S(now/d,P),P)%P;
}
inline ll rd()
{
	char c;ll f=1;
	while(!isdigit(c=getchar())) if(c=='-')f=-1;
	ll x=c-'0';
	while(isdigit(c=getchar())) x=x*10+(c^48);
	return x*f;
}
int main()
{
	fac[0]=1;
	for(int i=1;i<p[4];i++) fac[i]=fac[i-1]*i%M;
	for(int T=rd();T--;)
	{
		n=rd(),m=rd(),res=0;
		for(int i=1;i<=4;i++)
		{
			ans[i]=0;
			for(ll d=1;d*d<=n;d++) if(n%d==0)
			{
				(ans[i]+=solve(n,d,p[i]))%=p[i];
				if(d*d!=n) (ans[i]+=solve(n,n/d,p[i]))%=p[i];	
			}
		}
		for(int i=1;i<=4;i++)
		{
			ll Mi=M/p[i];
			(res+=ans[i]*Mi%M*inv(Mi,p[i]))%=M;
		}
		cout<<qp(m,res,mod)<<endl;
	}
    return 0; 
}