POJ - 2154 Color - （Ploya定理，欧拉函数，1~n的gcd(n,i)之和）

 

Description

Beads of N colors are connected together into a circular necklace of N beads (N<=1000000000). Your job is to calculate how many different kinds of the necklace can be produced. You should know that the necklace might not use up all the N colors, and the repetitions that are produced by rotation around the center of the circular necklace are all neglected. 

You only need to output the answer module a given number P. 

Input

The first line of the input is an integer X (X <= 3500) representing the number of test cases. The following X lines each contains two numbers N and P (1 <= N <= 1000000000, 1 <= P <= 30000), representing a test case.

Output

For each test case, output one line containing the answer.

Sample Input

5
1 30000
2 30000
3 30000
4 30000
5 30000

Sample Output

1
3
11
70
629

题意：给出一个数n，问由n种颜色组成长为n的项链的方案数，这里认为旋转后相同的为同一种方案,最后结果mod p

明显的Polya定理，注意这里的置换只有旋转置换，那么由公式：

不同染色的方案数L=1/|G|*[m^p(a1)+m^p(a2)+....+m^p(ak)].其中p(ai)是某个置换的循环节数，|G|即是置换的个数，m是颜色数，这里|G|和m都等于n.旋转i个位置的置换循环节数=gcd(i,n);

可由下面代码完成计算：

 

1. for(int i=1;i<=n;i++) //旋转i对应置换的循环个数为gcd（n,i）  
2.         sum+=pow(1.0*n,gcd(n,i));

但是题目的n非常大，这个循环显然超时，而且麻烦在要mod p，所以上面代码不行的

 

但是我们又有方法可以快速求得1~n的gcd（n，i）值之和，原理如下：

由于要求1~n的gcd（n，i）值之和，而且这些gcd的值只能等于n的因子，假设n有t个因子p1~pt，即这些gcd的值只能等于p1~pt中的某一个，我们可以按n的因子的值将其分组，那么第k组的gcd值之和即为(n^pk)*（第k组的个数），那么第k组的个数即为gcd(n,i)=k的i的个数，化为gcd（n/k，i/k）=1的个数，即求euler(n/k),这里euler为欧拉函数，euler（n）等于小于n的数中与n互质的数的个数。通过枚举n的因子i，那么这一组的gcd之和等于euler(n/i)*(n^i)。

而mod p的问题，由于最后要除以n而除法不能直接取模，但是我们在计算euler(n/i)*pow_mod(n,i)时可以直接计算euler(n/i)*pow_mod(n,i-1)即直接少乘一个n，之后就不用除了

Polya定理见：http://blog.csdn.net/sdau20163942/article/details/78981839

代码：

 

//只考虑旋转置换
#include<iostream>
#include<string>
#include<cstdio>
#include<algorithm>
#include<cmath>
#include<math.h>
#include<iomanip>
#include<queue>
#include<cstring>
#include<map>
using namespace std;
typedef long long ll;
#define N 1000000000
#define M 32000

int n,mod;

int euler_phi(int n)
{
    int ans=n;
    for(int i=2;i*i<=n;i++)
        if(n%i==0)
    {
        ans-=ans/i;
        while(n%i==0)
            n/=i;
    }
    if(n>1)
        ans-=ans/n;
    return ans;
}

int pow_mod(int n,int k) //快速幂
{
    if(n>=mod)
        n=n%mod;
    int res=1;
    while(k>0)
    {
        if(k&1)
            res=res*n%mod;
        n=n*n%mod;
        k>>=1;
    }
    return res;
}

int Polya()
{
    int i,res,ans=0;
    double t=sqrt(1.0*n);

    for(i=1;i<=t;i++)
    {
        if(n%i==0)
        {
            res=((euler_phi(n/i)%mod)*pow_mod(n,i-1))%mod; //处理因子i，即gcd（n，x）=i
            ans=(ans+res)%mod;
            if(i*i!=n)
            {
                res=((euler_phi(i)%mod)*pow_mod(n,n/i-1))%mod; //处理因子 n/i
                ans=(ans+res)%mod;
            }
        }
    }
    return ans;
}
int main()
{
    int T;
    scanf("%d",&T);
    while(T--)
    {
        scanf("%d%d",&n,&mod);

        printf("%d\n",Polya());
    }
    return 0;
}