poj1845 Sumdiv(算数基本定理,质因数分解)

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
",ans%mod);
    return 0;
}