数论poj1845（A^B的约数和）

题目：点击打开链接原文博客：点击打开链接求A^B的约数和对9901取模后的结果(0 <= A,B <= 50000000)；1.根据唯一分解定理将A进行因式分解可得:A = p1^a1 * p2^a2 * p3^a3 * pn^an.
    A^B=p1^(a1*B)*p2^(a2*B)*...*pn^(an*B);
    A^B的所有约数之和sum=[1+p1+p1^2+...+p1^(a1*B)]*[1+p2+p2^2+...+p2^(a2*B)]*[1+pn+pn^2+...+pn^(an*B)]2.
等比数列1+pi+pi^2+pi^3+...+pi^n可以由二分求得(即将需要求解的因式分成部分来求解)
       若n为奇数,一共有偶数项,设p为3,则(1+p)+(p^2+p^3)=(1+p)+p^2(1+p)=(1+p^2)*(1+p)                  1+p+p^2+p^3+........+p^n=(1+p+p^2+....+p^(n/2))*(1+p^(n/2+1));
       若n为偶数,一共有奇数项,设p为4,则(1+p)+p^2+(p^3+p^4)=(1+p)+p^2+p^3(1+p)=(1+p^3)*(1+p)+P^2                  1+p+p^2+p^3+........+p^n=(1+p+p^2+....+p^(n/2-1))*(1+p^(n/2+1));

#include<iostream>
#include<cstdio>
#include<cstring>
#include<algorithm>
using namespace std;
#define Mod 9901
#define LL long long
const int maxn = 10010;
int p[maxn];
void getp()//素数
{
memset(p,0,sizeof(p));

for(int i=2;i<=maxn;i++)
{
if(!p[i])p[++p[0]]=i;
for(int j=1;j<=p[0]&&p[j]<=maxn/i;j++)
{
p[p[j]*i]=1;
if(i%p[j]==0) break;
}
}
}
LL f[maxn][2];
int cnt;
int getf(LL x)//分解
{
cnt = 0;
for(int i = 1; p[i]*p[i] <= x; i++)
{
f[cnt][1] = 0;
if(x % p[i] == 0)
{
f[cnt][0] = p[i];
while(x % p[i] == 0)
{
f[cnt][1]++;
x /= p[i];
}
cnt++;
}
}
if(x != 1)
{
f[cnt][0] = x;
f[cnt++][1] = 1;
}
return cnt;
}
LL qpow(LL a,LL n)//快速幂
{
LL res = 1;
LL tmp = a % Mod;
while(n)
{
if(n & 1)
{
res *= tmp;
res %= Mod;
}
n >>= 1;
tmp *= tmp;
tmp %= Mod;
}
return res;
}
LL sum(LL p, LL n)//计算1+p+p^2+p^3+...+p^n
{
if(p == 0)
return 0;
if(n == 0)
return 1;
if(n & 1)//奇数
return ((1+qpow(p,n/2+1))%Mod*sum(p,n/2)%Mod)%Mod;
else
return ((1+qpow(p,n/2+1))%Mod*sum(p,n/2-1)+qpow(p,n/2)%Mod)%Mod;
}
int main()
{
std::ios::sync_with_stdio(false);
std::cin.tie(0);
int a,b;
getp();
while(cin>>a>>b)
{
getf(a);
LL ans = 1;
for(int i = 0; i < cnt; i++)
{
ans *= (sum(f[i][0],b*f[i][1])%Mod);
ans %= Mod;
}
cout<<ans<<endl;
}
return 0;
}