POJ 2154 Color Polya计数法+欧拉函数

学过Polya计数法的知道答案为∑n^gcd(i, n)，但是n太大了，设gcd(i，n)的值为x，x显然是n的约数，所以可以从1到sqrt(n)枚举x，设i=y*x，i必须要有x这个因子，才能满足gcd(i, n)=x，并且gcd(i/x，n/x) == 1，所以满足的个数为n/x的欧拉函数值。

#include <iostream>
#include <stdlib.h>
#include <algorithm>
#include <stdio.h>
#include <vector>
using namespace std;
const int MAXN = 1000000000;
int isprime[50001];
int prime[8001];
int num,n,p;
void getprime()
{
	num = 0;
	for(int i = 2; i<= 50000; i++)
    {
        if(!isprime[i])
        {
            prime[num++]=i;
            for(int j=1; j*i <= 50000; j++)
            {
                isprime[i*j] = 1;
            }
		}
	}
}
int euler(int x)
{
	int res = x;
	for(int i = 0; i < num && prime[i]*prime[i] <= x; i++)
	{
		if(x%prime[i] == 0)
		{
			res = res/prime[i]*(prime[i]-1);
			while(x%prime[i] == 0)
			{
				x /= prime[i];
			}
		}
	}
	if(x > 1)
        res = res/x*(x-1);
	return res;
}
int expmod(int a,int b,int mod)
{
    int ret = 1;
    a = a%mod;
    while(b > 0)
    {
        if(b&1)
            ret = (ret*a)%mod;
        a = (a*a)%mod;
        b >>= 1;
    }
    return ret;
}

int main()
{
	int Case;
	getprime();
	scanf("%d", &Case);
	while(Case--)
	{
		scanf("%d%d",&n, &p);
		int ans = 0, i;
		for(i = 1; i*i < n; i++)
        {
            if(n%i == 0)
            {
                ans = (ans+euler(i)%p*expmod(n, n/i-1, p)+euler(n/i)%p*expmod(n, i-1, p))%p;
            }
        }
		if(i*i == n)
            ans = (ans+euler(i)*expmod(n, i-1, p))%p;
		cout<<ans<<endl;
	}
    return 0;
}