poj 1845 Sumdiv (等比求和+逆元)

题目链接：http://poj.org/problem?id=1845

题目大意：给出两个自然数a,b,求a^b的所有自然数因子的和模上9901 (0 <= a,b <= 50000000)

解题思路：我们先利用唯一分解定理，将a分解成(p1^q1)*(p2^q2)……(pk^qk)的形式，则a^b=((p1^q1)*(p2^q2)……(pk^qk))^b=(p1^q1b)*(p2^q2b)……(pk^qkb)

a^b的因子和就会等于(1+p1+p1^2+……p1^q1b)*(1+p2+p2^2+……p2^q2b)*……(1+pk+pk^2+……pk^qkb)

然后我们可以用等差求和公式转化为((p1^(q1b+1)-1)/(p1-1))*((p2^(q2b+1)-1)/(p2-1))……((pk^(qkb+1)-1)/(pk-1))

对于求逆元：

(a/b)%mod=(a%(mod*b))/b%mod。对B*mod取余，剩余的值必定是B的倍数，这种方法是用于mod和B小的时候,用在这题就刚好了。

#include<iostream>
using namespace std;
typedef long long ll;
const int MAXN=500005;
const int mod=9901;
ll a,b,prime[MAXN],tot;
void getPrime(int N){  //筛素数 
    for(int i=1;i<=N;i++) prime[i]=1;
    for(int i=2;i<=N;i++){
        if(prime[i])
            prime[tot++]=i;
        for(int j=0;j<tot&&prime[j]*i<=N;j++){
            prime[i*prime[j]]=0;
            if(i%prime[j]==0)break;
        }
    }
}
ll qmul(ll a,ll b,ll m){
    ll res=0;
    while(b){
        if(b&1) res=(res+a)%m;
        b>>=1;
        a=(a+a)%m;
    }
    return res;
}
ll qpow(ll a,ll b,ll m){
    ll res=1;
    while(b){
        if(b&1) res=qmul(res,a,m);  //直接相乘会爆，可以一个一个加 
        a=qmul(a,a,m);
        b>>=1;
    }
    return res;
}
ll solve(ll x,ll y){
    ll ans=1;
    for(int i=0;prime[i]*prime[i]<=x;i++){
        if(x%prime[i]==0){
            int cnt=0;
            while(x%prime[i]==0){
                cnt++;
                x/=prime[i];
            }
            ll M=(prime[i]-1)*mod;
            ans=ans*(qpow(prime[i],cnt*y+1,M)-1+M)%M/(prime[i]-1)%mod;
        }
    }
    if(x>1){
        ll M=(x-1)*mod;
        ans=ans*(qpow(x,y+1,M)-1+M)%M/(x-1)%mod;
    }
    return ans;
}
int main(){
    cin>>a>>b;
    getPrime(500000);
    cout<<solve(a,b)<<endl;
    return 0;
}