【POJ2154】Color-Polya定理+欧拉函数

测试地址：Color
题目大意：有N种颜色的珠子任意多个，要用N(≤109)颗珠子做成环形的项链，绕中心旋转后相同的方案视为本质相同，问本质不同的方案数对P取模的值。
做法：前几天POJ炸了，所以拖到今天才写这一篇……
这一题需要用到Polya定理（真名不是这个o，但是那个字符不会打QAQ）+欧拉函数。
首先看Polya定理的内容：对于一个置换群G，里面有|G|种置换，表示一种染色方案经过该置换后得到的另一方案与其本质相同，若颜色数为m，那么本质不同的方案数为：

1|G|∑i=1|G|mci

其中ci为第i种置换的循环数，如置换(1,2,3,4,5)−>(3,1,2,5,4)，循环就有2个，即(1,2,3)和(4,5)。
那么再看这一题，容易看出|G|=m=N，那么关键就是求ci了。我们不妨将所有珠子顺时针依次编号为1到N，记第i(1≤i≤N)个置换为所有珠子顺时针位移i个位置形成的置换，那么可以看出ci=gcd(N,i)，所以总方案数就是1N∑Ni=1mgcd(N,i)。但是N很大，O(N)的时间复杂度不足以通过这道题，所以我们要考虑优化。
设di=gcd(N,i)，可以看出di一定是N的因子，那么我们可不可以枚举N的所有因子，然后统计有多少di等于这个因子呢？答案是可以的。对于N的一个因子L，我们要求有多少个i(1≤i≤N)满足di=L，也就是满足gcd(N,i)=L，即满足gcd(N/L,i/L)=1。注意到满足这个式子的i的个数就是欧拉函数φ(N/L)，所以要求的式子就可以转化为：

1N∑L|Nφ(NL)NL

考虑到模意义下乘以一个分数比较困难，所以我们直接把1/N乘进去得到：∑L|Nφ(N/L)NL−1，这样就不用乘分数了。按照上面的方法优化，时间复杂度就降了很多了，只要在枚举L的时候实时算出φ(N/L)的值即可。
以下是本人代码：

#include <cstdio>
#include <cstdlib>
#include <cstring>
#include <iostream>
#include <algorithm>
#define ll long long
using namespace std;
ll T,N,P,tot=0,fac[25],num[25],ans;

void find_factor()
{
  ll s=N;
  tot=0;
  for(ll i=2;i*i<=s;i++)
    if (!(s%i))
    {
      fac[++tot]=i;
      num[tot]=0;
      while(!(s%i)) s/=i,num[tot]++;
    }
  if (s>1) fac[++tot]=s,num[tot]=1;
}

ll phi(ll x)
{
  ll s=x,ss=x;
  for(ll i=2;i*i<=s;i++)
    if (!(s%i))
    {
      ss=ss/i*(i-1);
      while(!(s%i)) s/=i;
    }
  if (s>1) ss=ss/s*(s-1);
  return ss;
}

ll power(ll a,ll b)
{
  ll s=a,sum=1;
  while(b)
  {
    if (b&1) sum=(sum*s)%P;
    b>>=1;s=(s*s)%P;
  }
  return sum;
}

void dfs(int i,ll d)
{
  if (i>tot)
  {
    ans=(ans+phi(N/d)*power(N,d-1))%P;
    return;
  }
  for(int j=0;j<=num[i];j++)
  {
    dfs(i+1,d);
    d*=fac[i];
  }
}

int main()
{
  scanf("%lld",&T);
  while(T--)
  {
    scanf("%lld%lld",&N,&P);
    find_factor();
    ans=0;
    dfs(1,1);
    printf("%lld\n",ans);
  }

  return 0;
}