poj1286 Necklace of Beads（polya）

polya计数。

一个置换是长这样的，与列顺序无关：

$$\left(\begin{matrix} 1&2&3&4\\ 2&3&4&1 \end{matrix} \right)$$
群，自己学习下吧，懒得写了(逃

由拉格朗日定理，如果A是S的一个子群，则$|A| \ | \ |S|$，且A的所有陪集大小均相等，为|A|。

由轨道-稳定化子定理 $|orbit(x)|*|stab(x)|=|G|$

由burnside引理：集合M关于置换群G的轨道数，等于G中每个置换下不动点的个数的算术平均数。即

$$\frac{1}{|G|}\sum\limits_{f\in G}|\{x\in M|f(x)==x\}|$$
由polya计数法：对一个置换p，记c(p)表示p中的循环节个数
设G是定义在X上的一个置换群，如果用m种颜色对X中的元素进行染色，本质不同的方案数为
$$\frac{\sum\limits_{\pi\in G}m^{c(\pi)}}{|G|}$$

然后对于此题，用3种颜色给n个点的环染色，旋转同构，翻转同构。

我们搞出置换群：
首先是对于旋转同构：显然有n个置换，对于循环移位k位形成的置换，循环节个数为gcd(n,k)。我们可以$O(nlogn)$的计算这n个置换的贡献，也可以考虑循环节个数为k的有几个置换，也就是求gcd(n,x)==k的x的个数，显然等于$\phi(n/k)$。因此就是$\sum\limits_{k|n}m^k\phi(n/k)$。我们O(n)的先处理出phi函数，然后用$O(\sqrt nlogn )$的复杂度计算贡献。

然后是对于翻转同构：分奇偶讨论
n为奇数：有n条对称轴，循环节个数分别为n/2+1，可以直接O(logn)计算。
n为偶数：有n条对称轴，对称轴在点上时，循环节个数为n/2+1，不在点上时，循环节个数为n/2，也可以直接O(logn)计算。

综上，复杂度为$O(n+T*\sqrt nlogn)$

这题保证了答案最大不超过int，然而中途计算还是会超 囧

#include <cstdio>
#include <cstring>
#include <algorithm>
#include <cmath>
using namespace std;
#define inf 0x3f3f3f3f
#define ll long long
#define N 30
inline char gc(){
    static char buf[1<<16],*S,*T;
    if(S==T){T=(S=buf)+fread(buf,1,1<<16,stdin);if(T==S) return EOF;}
    return *S++;
}
inline int read(){
    int x=0,f=1;char ch=gc();
    while(ch<'0'||ch>'9'){if(ch=='-')f=-1;ch=gc();}
    while(ch>='0'&&ch<='9') x=x*10+ch-'0',ch=gc();
    return x*f;
}
int n,phi[N],prime[N],tot=0;
bool notprime[N];
inline void getprime(){
    notprime[1]=1;phi[1]=1;
    for(int i=2;i<24;++i){
        if(!notprime[i]) prime[++tot]=i,phi[i]=i-1;
        for(int j=1;prime[j]*i<24;++j){
            notprime[prime[j]*i]=1;
            if(i%prime[j]==0){phi[i*prime[j]]=phi[i]*prime[j];break;}
            phi[i*prime[j]]=phi[i]*phi[prime[j]];
        }
    }
}
inline ll ksm(ll x,int k){
    ll res=1;for(;k;k>>=1,x*=x) if(k&1) res*=x;return res;
}
int main(){
//  freopen("a.in","r",stdin);
    getprime();
    while(1){
        n=read();if(n==-1) break;if(!n){puts("0");continue;}
        ll ans=0;int x=1;
        for(;x*x<n;++x) if(n%x==0) ans+=ksm(3,x)*phi[n/x]+ksm(3,n/x)*phi[x];
        if(x*x==n) ans+=ksm(3,x)*phi[n/x];if(n&1) ans+=ksm(3,n+1>>1)*n;
        else ans+=ksm(3,n/2+1)*n/2+ksm(3,n/2)*n/2;printf("%lld\n",ans/(n*2));
    }return 0;
}