poj1845 逆元 因子和

传送门

主要目的还是记录一下，学习了学长博客http://blog.csdn.net/acdreamers/article/details/8220787，写的比我清楚很多……

题意：求a^b的因子和对9901取余 

思路：一个数的因子和求法：对n素数分解，n=p1^a1*p2^a2*...*pk^ak

因子和=(1+p1+p1^2+...+p1^a1)*(1+p2+p2^2+...+p2^a2)*...(1+pk+...pk^ak)
           =(p1^(a1+1)-1)/(p1-1)+...+(pk^(ak+1)-1)/(pk-1)
（这个公式还是搞了很久才明白咋回事，这个式子拆开会发现正好是n的各个因子，正是由于式子的每一项都是n的素因子相乘得到的）
那么这个题求a^b的因子和，我们没有必要求出a^b具体值，只需要进行素因子分解的时候把b幂次加入就可以了，那么式子就化作了
 a^b因子和=(p1^(a1*b+1)-1)/(p1-1)+...+(pk^(ak*b+1)-1)/(pk-1)
    
这里在求的时候运用了公式a/b mod m=a mod (m*b)/b，这个公式比较适用于a是幂的形式，m*b比较小，因为这样可以利用快速幂来求解，此公式无论m是不是素数都适用。

#include<iostream>#include<algorithm>#include<cstdio>#include<cstdlib>#include<cstring>#include<cmath>using namespace std;typedef long long ll;const int mod=9901;const int maxn=1e6+5;bool prime[maxn];int p[maxn];int cnt;void isprime(){    cnt=0;    memset(prime,true,sizeof(prime));    for(int i=2;i<maxn;i++)    {        if(prime[i])        {            p[cnt++]=i;            for(int j=i+i;j<maxn;j+=i)                prime[j]=false;        }    }}ll multi(ll a,ll b,ll m)///由于a太大 会超longlong精度 所以用二分乘法{    ll ans=0;    while(b)    {        if(b&1)        {            ans=(ans+a)%m;        }        b>>=1;        a=(a+a)%m;    }    return ans;}ll quick_mod(ll a,ll b,ll m){    a%=m;    ll ans=1;    while(b)    {        if(b&1)        {            ans=multi(ans,a,m);            b--;        }        b>>=1;        a=multi(a,a,m);    }    return ans;}ll solve(ll a,ll b){    ll ans=1;    for(int i=0;i<cnt&&p[i]*p[i]<=a;i++)    {        if(a%p[i]==0)        {            int num=0;            while(a%p[i]==0)            {                num++;                a/=p[i];            }            ll m=(p[i]-1)*mod;            ll pow;            pow=num*b+1;            ll temp=(quick_mod(p[i],pow,m)-1+m)/(p[i]-1);            ans=ans*temp%mod;        }    }    if(a>1)    {        ll m=(a-1)*mod;        ll pow;        pow=b+1;        ll temp=(quick_mod(a,pow,m)-1+m)/(a-1);        ans=ans*temp%mod;    }    return ans%mod;}int main(){    ll a,b;    isprime();    while(scanf("%lld%lld",&a,&b)!=-1)    {        printf("%lld\n",solve(a,b));    }    return 0;}