poj 2409 Let it Bead (polya)

#include<iostream>#include<math.h>#include<stdio.h>#include<algorithm>#include<string.h>#include<vector>#include<queue>#include<map>#include<set>#define B(x) (1<<(x))using namespace std;typedef long long ll;void cmax(int& a,int b){ if(b>a)a=b; }void cmin(int& a,int b){ if(b<a)a=b; }void cmax(ll& a,ll b){ if(b>a)a=b; }void cmin(ll& a,ll b){ if(b<a)a=b; }void add(int& a,int b,int mod){ a=(a+b)%mod; }void add(ll& a,ll b,ll mod){ a=(a+b)%mod; }const int oo=0x3f3f3f3f;const ll MOD=1000000007;ll quick_pow(ll a,ll k){    ll ans=1;    while(k){        if(k&1) ans=(ans*a);        a=(a*a);        k>>=1;    }    return ans;}int gcd(int a,int b){    return b ? gcd(b,a%b) : a;}int main(){    int c,n;    while(scanf("%d %d",&c,&n)!=EOF){        if(n==0&&c==0)break;        if(n==0){            printf("0\n");            continue;        }        ll ans=0;        if(n&1){            ans=quick_pow(c,(n-1)/2+1)*n;            for(int i=1;i<=n;i++)                ans+=quick_pow(c,gcd(n,i));        }else{            ans=quick_pow(c,(n-2)/2+2)*n/2;            ans+=quick_pow(c,n/2)*n/2;            for(int i=1;i<=n;i++)                ans+=quick_pow(c,gcd(n,i));        }        cout<<ans/n/2<<endl;    }    return 0;}