POJ 2154 Color（polya）

Description
用n种颜色给一个有n个珠子的项链染色，旋转和翻转相同看作同一方案，问一共有多少种不同的方案，结果模p
Input
第一行为一整数T表示用例组数，每组用例占一行为两个整数n和p
(T<=3500，1<=n<=10^9, 1<=p<=30000)
Output
对于每组用例，输出方案数%p
Sample Input
5
1 30000
2 30000
3 30000
4 30000
5 30000
Sample Output
1
3
11
70
629
Solution
polya，先考虑旋转，有n种置换，旋转k个元素后循环节即轮换数为gcd(k,n)，故方案数为，然后是翻转，对称轴显然为直径，当n为偶数时有两种对称轴，一种过两点，另一种过两个中点，第一种对称轴有n/2个，翻转180度的轮换数为n/2+1，第二种对称轴也有n/2个，翻转180度的轮换数为n/2，方案数为n/2*n^(n/2)+n/2*n^(n/2+1)；当n为奇数时对称轴只有一种，即过一点以及一个中点，共n个，翻转180度的轮换数为n/2，方案数为n*n^(n/2)，所以翻转也有n种置换，总置换数为2*n。
Code

#include<cstdio>
#include<iostream>
#include<algorithm>
#include<cstring>
using namespace std;
#define maxn 55555
int prime[maxn],is_prime[maxn],res;
void get_prime(int n)
{
    res=0;
    memset(is_prime,0,sizeof(is_prime));
    for(int i=2;i<n;i++)
        if(!is_prime[i])
        {
            prime[res++]=i;
            for(int j=i;j<n;j+=i)
                is_prime[j]=1;
        }
}
int mod_pow(int a,int b,int p)
{
    int ans=1;
    a%=p;
    while(b)
    {
        if(b&1) ans=(ans*a)%p;
        a=(a*a)%p;
        b>>=1;
    } 
    return ans;
}
int get_euler(int n)
{
    int ans=n;
    for(int i=0;i<res&&prime[i]*prime[i]<=n;i++)
        if(n%prime[i]==0)
        {
            ans=ans/prime[i]*(prime[i]-1);
            while(n%prime[i]==0)
                n/=prime[i];
        }
    if(n>1)ans=ans/n*(n-1);
    return ans;
} 
int main()  
{  
    get_prime(maxn);
    int T,n,p;   
    scanf("%d",&T);  
    while(T--)  
    {  
        scanf("%d%d",&n,&p);
        int ans=0;
        for(int i=1;i*i<=n;i++)
            if(n%i==0)
            {
                if(i*i!=n)
                    ans=(ans+get_euler(i)%p*mod_pow(n,n/i-1,p)%p+get_euler(n/i)%p*mod_pow(n,i-1,p)%p)%p;
                else 
                    ans=(ans+get_euler(i)%p*mod_pow(n,n/i-1,p)%p)%p;
            }
        printf("%d\n",ans);
    }  
    return 0;  
}