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;
}