POJ 1845 Sumdiv

算是一道数学题。

首先，由唯一分解定理，可以将任意一个数分解为a=p1^a1*p2^a2*....*pk^ak。那么，此时对应的所有因子的和应该是：(1+p1+p1^2+p1^3+....+p1^a1)*.....*(1+pk+pk^2+...+pk^ak)。这个式子还是很好理解的，对于a的一个因子，它对pi一共有ai+1种选择方式。

之后这题目的关键就是求1+pk+pk^2+...+pk^ak。

主要有两种方法：

1.二分求解，定义cal(pk,ak)=1+pk+pk^2+...+pk^ak，当ak是奇数的时候cal(pk,ak)=cal(pk,ak/2)*(1+pow_mod(pk,ak/2+1));当ak是偶数的时候，

cal(pk,ak)=cal(pk,ak-1)+pow_mod(pk,ak)。

2.利用乘法逆元。注意到该数列是等比数列，那么cal(p,k)＝(p^(ak+1)-1)/(p-1)，那么就只要求p-1关于mod的乘法逆元就可以解决了，但求逆元的时候要注意p-1

是mod的倍数的情况。

#include "iostream"
#include "cstdio"
#include "cmath"
#include "algorithm"
#include "map"
#define LL long long
using namespace std;
#define MOD 9901

//返回ab mod n
LL mul_mod(LL a,LL b,LL c)
{
LL ret=0,tmp=a%c;
while(b)
{
if(b&0x1)if((ret+=tmp)>=c)ret-=c;
if((tmp<<=1)>=c)tmp-=c;
b>>=1;
}
return ret;
}
//快速幂
LL pow_mod(LL a,LL p,LL n){
if(p==0)  return 1%n;
LL ans=pow_mod(a,p/2,n);
ans=mul_mod(ans, ans, n);
if(p%2==1)  ans=mul_mod(ans, a, n);
return ans%n;
}

LL prime[1111];
LL t[1111];

LL solve(LL p,LL t){
if(t==0)  return 1;
if(t==1)  return (1+p)%MOD;
if(t%2==1){
return solve(p,t/2)*(pow_mod(p,t/2+1,MOD)+1)%MOD;
}
else  return (solve(p,t-1)+pow_mod(p,t,MOD))%MOD;
}
int main(){
LL a,b;
cin>>a>>b;
int cnt=0;
if(a==0){
printf("0\n");
return 0;
}
if(b==0){
printf("1\n");
return 0;
}
LL tem=a;
LL high=(LL)sqrt((double)a)+10;
for(int i=2;i<=high;i++){
if(tem==1)  break;
if(tem%i==0){
prime[cnt]=i;
while(tem%i==0){
tem/=i;
t[cnt]++;
}
cnt++;
}
}
if(tem!=1){
prime[cnt]=tem;
t[cnt++]=1;
}
LL ans=1;
//printf("%lld\n",prime[0]);
//printf("%lld\n",solve(2,3));
//for(int i=0;i<cnt;i++){
//printf("%lld %lld\n",prime[i],t[i]);
//}
for(int i=0;i<cnt;i++){
ans=(ans*solve(prime[i]%MOD,t[i]*b))%MOD;
}
printf("%lld\n",(ans+MOD)%MOD);
return 0;
}