poj 2154 Color

            Polya + Euler优化。由于n太大，所以需要用Euler优化。设L = n  / gcd(n， i)，t = i / gcd(i, n)，则gcd(L,  t) = 1.又有i < n所以，t  < L，且t和L互质。于是满足gcd(n, i) = n / L的i的个数为Euler(L)。对于题目原本的求和式应该是sum(pow(n, gcd(i, n) - 1))。所以我们枚举n的所有因子L，然后求和式就变成了

          sum(Euler(L) * pow(n, n / L - 1))。

#include<algorithm>
#include<iostream>
#include<cstring>
#include<cstdio>
#include<vector>
#include<cmath>
#include<set>
#include<map>
#define LL long long
#define CLR(a, b) memset(a, b, sizeof(a))
#define REP(i, n) for(int i = 0; i < n; i ++)
using namespace std;
const int N = 100100;
bool isp[N];
vector<int> p;

void get_P()
{
    CLR(isp, true);p.clear();
    for(int i = 2; i < N; i ++)
    {
        if(isp[i])
        {
            p.push_back(i);
            if(i < 1111) for(int j = i * i; j < N; j += i)
            {
                isp[j] = false;
            }
        }
    }
}

int Euler_phi(int n)
{
    int ret = n;
    for(int i = 0; p[i] * p[i] <= n; i ++) if(n % p[i] == 0)
    {
        ret = ret / p[i] * (p[i] - 1);
        while(n % p[i] == 0) n /= p[i];
    }
    if(n > 1) ret = ret / n * (n - 1);
    return ret;
}

int Pow(int a, int b, int mod)
{
    int ret = 1;a %= mod;
    while(b)
    {
        if(b & 1) ret = ret * a % mod;
        a = a * a % mod;
        b >>= 1;
    }
    return ret;
}

int main()
{
    int n, ans, mod, x;
    get_P();
    scanf("%d", &x);
    while(x --)
    {
        scanf("%d%d", &n, &mod);ans = 0;
        for(int i = 1; i * i <= n; i ++)
        {
            if(n % i != 0) continue;
            ans += Euler_phi(i) % mod * Pow(n, n / i - 1, mod) % mod;
            if(n / i != i)
                ans += Euler_phi(n / i) % mod * Pow(n, i - 1, mod) % mod;
            ans %= mod;
        }
        printf("%d\n", ans);
    }
}