POJ 1845 Sumdiv 【推公式 | 逆元】

题意： 

          求 a^b 所有因子和。

题解：

        a={p1^x1*p2^x2.....pn^xn}

        sum(a^b)= (p1^0+p1^1+...+p1^x1).......;

        就是算多个等比数列和。

       等比数列求和公式 (a^n+1-1)/(a-1)

       除法需要逆元，

      逆元如果余数c和除数b不大 可以直接

       a/b%modc =a%(b*c)/b%c;

       除了逆元还有种二分求解的方法，每次整个能从一半推出，分奇偶。

#include<cstdio>
#include<iostream>
#include<cmath>
#include<algorithm>
#include<cstring>
#define ll long long
using namespace std;
const int maxn=1e5+10;
const int maxs=1e7+10;
const int inf=0xfffff;
const ll mod=9901;
ll a,b,c;
int tot_p=0;
int p[maxn];
bool v[maxn];
void init(){
    v[1]=true;
    for(int i=2;i<maxn;++i){
        if(!v[i]) p[++tot_p]=i;
        for(int j=1;j<=tot_p&&i*p[j]<maxn;++j){
            v[i*p[j]]=true;
            if(i%p[j]==0) break;
        }
    }
}
ll mul(ll a,ll b,ll p)//大数乘法，直接相乘会爆int64，需要逐位相乘
{
    ll s=0;
    while(b)
    {
        if(b&1)
            s=(s+a)%p;
        a=(a*2)%p;
        b/=2;
    }
    return s;
}
ll kru(ll x,ll y,ll p){
    ll ret=1;
    while(y){
        if(y&1){
            ret=mul(ret,x,p);
        }
        y/=2;
        x=mul(x,x,p);
    }
    return ret;
}
ll sum(ll a,ll b)
{
    if(b==1)return 1;
    if(b&1)return (sum(a,b/2)*(1+kru(a,b/2+1,9901))+kru(a,b/2,9901))%9901;
    else return sum(a,b/2)*(1+kru(a,b/2,9901))%9901;
}
int main(){
    init();
    int pos=1;
    scanf("%lld %lld",&a,&b);
    ll s,t,ans=1,f,nu;
    while(a&&pos<tot_p){
        if(a==1) break;
        if(a%p[pos]!=0) {pos++;continue;}
        s=0;
        t=p[pos];
        while(a%p[pos]==0){
            a/=p[pos];
            s++;
        }
        if((t-1)%mod==0) ans=ans*(kru(t,b*s+1,mod*(t-1))/(t-1))%mod;
        else ans=ans*(kru(t,b*s+1,mod)-1)*kru(t-1,mod-2,mod)%mod;
    }
    if(a>1){
        if((a-1)%mod==0) ans=ans*(kru(a,b+1,mod*(a-1))/(a-1))%mod;
        else ans=ans*(kru(a,b+1,mod)-1)*kru(a-1,mod-2,mod)%mod;
    }
    printf("%lld\n",(ans+mod)%mod);
    return 0;
}