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