题目地址:http://acm.hdu.edu.cn/showproblem.php?pid=3923
因为要取模,所以要乘法逆元
#include<iostream>
#include<cstdio>
#include<cmath>
using namespace std;
typedef long long LL;
const int P=1000000007;
int gcd(int a,int b)
{
if(b==0) return a;
return gcd(b,a%b);
}
int euler_phi(int n)
{
int res=1;
for(int i=2;i*i<=n;i++)
if(n%i==0) { //说明i|n
n/=i,res*=i-1;
while(n%i==0) n/=i,res*=i; //说明i^2|n
}
if(n>1) res*=n-1;
return res%P;
}
LL PowMod(LL x,LL n,LL p) //x^n 对Max取模
{
LL result=1;
LL base = x%p;
while(n>0)
{
if(n & 1)
result=result*base%p; //当二进制不为0时,就要乘个
base=base*base%p; //累乘此时二进制的权值
n>>=1;
}
return result;
}
LL polya(int m,int n) //m color ,n number
{
LL tot=0; //方案数
for(int i=1;i*i<=n;i++) //1~sqrt(n)
{
if(n%i) continue; //当i不是n的约数时就进入下一次循环
tot+=euler_phi(i)*PowMod(m,n/i,P); //d=gcd(n,i) d为n的因数,且有euler_phi(n/i)个
if(i*i!=n) tot+=euler_phi(n/i)*PowMod(m,i,P); //当i*i==n时,不必算两次
}
tot%=P;
if(n%2!=0) tot+=PowMod(m,(n+1)/2,P)*n; //odd
else tot+=(PowMod(m,n/2,P)+PowMod(m,n/2+1,P))*n/2;
tot%=P;
return tot*PowMod(2*n,P-2,P)%P;
//等于(tot/n/2)%P ,tot为未取模的大小
}
int main()
{
int T,kase=0;
cin>>T;
while(T--)
{
int m,n;
cin>>m>>n;
cout<<"Case #"<<++kase<<": "
<<polya(m,n)<<endl;
}
return 0;
}