poj

题意：求\(a^b\)的所有约数之和\(mod9901\)的值.(\(1<=a,b<=5*10^7\))

分析：根据算数基本定理(任何一个大于1的正整数都能唯一分解为有限个质数的乘积)，把a分解质因数，表示为\(p_1^{c_1}*p_2^{c_2}*...*p_n^{c_n}\),再由算数基本定理的"约数和"推论可得,\(a^b\)的所有约数之和为\((1+p_1+p_1^2+...+p1^{b*c_1})*...*(1+p_n+p_n^2+...+pn^{b*c_n})\)

又因为上面每一项都是等比数列,由等比数列求和公式得\(1+p_1+p_1^2+...+p1^{b*c_1}=(p_1^{b*c_1+1}-1)/(p_1-1)\)

所以我们可以直接快速幂计算\(p_1^{b*c_1+1}-1\).

因为9901是一个质数,如果\(p_1-1\)不是9901的倍数,就可以直接计算\(p_1-1\)的乘法逆元inv,用乘inv代替除以\(p_1-1\);而如果是倍数的话,一定满足\(p1\)mod\(9901=1\),所以在mod9901意义下,\(1+p_1+p_1^2+...+p1^{b*c_1}=1+1+1^2+..+1^{b*c_1}=b*c_1+1\),也是可以直接计算的.

const int mod=9901;
int a,b,m,ans=1;
int p[10005],c[10005];
void divide_prime(int n){
    for(int i=2;i*i<=n;i++){
        if(n%i==0){
            p[++m]=i;
            while(n%i==0)n/=i,c[m]++;
        }
    }
    if(n>1)p[++m]=n,c[m]=1;
}
int quickpow(int a,long long b){
    int cnt=1;
    while(b){
        if(b&1)cnt=((long long)cnt*a)%mod;
        a=((long long)a*a)%mod;
        b>>=1;
    }
    return cnt%mod;
}
int main(){
    a=read(),b=read();
    divide_prime(a);
    for(int i=1;i<=m;i++){
        if((p[i]-1)%mod==0){
            ans=(ans*((long long)b*c[i]+1)%mod)%mod;
            continue;
        }
        int x=quickpow(p[i],(long long)b*c[i]+1);
        x=(x-1+mod)%mod;
        int y=quickpow(p[i]-1,mod-2);
        ans=((long long)ans*x*y)%mod;
    }
    printf("%d\n",ans%mod);
    return 0;
}