<更新提示>
<第一次更新>
<正文>
sumdiv(POJ 1845)
Description
给定两个自然数A和B,S为A^B的所有正整数约数和,编程输出S mod 9901的结果。
Input Format
只有一行,两个用空格隔开的自然数A和B(0<=A,B<= 50000000)。
Output Format
只有一行,即S mod 9901的结果。
Sample Input
2 3
Sample Output
15
解析
这是一道数学推导+分治的简单运用,大体思路如下。
由算数基本定理可得:
[A=p_1^{a_1}*p_2^{a_2}*...*p_k^{a_k}=prod_{i=1}^kp_i^{a_i}
]
那么
[A=p_1^{a_1B}*p_2^{a_2B}*...*p_k^{a_kB}=prod_{i=1}^kp_i^{a_iB}
]
易知(A^B)的约数之和就是:
[ans=(1+p_1+p_1^2+...+p_1^{a_1B})*...*(1+p_k+p_k^2+...+p_k^{a_kB})\=prod_{i=1}^{k}(1+p_i+p_i^2+...+p_i^{a_iB})\=prod_{i=1}^{k}sum_{j=0}^{a_iB}p_i^j
]
分解质因数我们是可以简单做到的,由于涉及到取模,在不用逆元的情况下,我们不能直接用等比数列求和公式。
所以我们的问题就转化成了等比数列求和,这是可以利用分治在(log)时间内实现的(不涉及除法)。
对于求解(sum(p,c)=1+p+...+p^c),可以分解一下:
- 对于c为奇数
[sum(p,c)=1+p+...+p^c
\=(1+p+...+p^{frac{c-1}{2}})+p^{frac{c+1}{2}}*(1+p+...+p^{frac{c-1}{2}})
\=(1+p^{frac{c+1}{2}})*sum(p,frac{c-1}{2})
]
- 对于c为偶数
[sum(p,c)=1+p+...+p^c
\=(1+p+...+p^{frac{c}{2}-1})+p^{frac{c}{2}}*(1+p+...+p^{frac{c}{2}-1})+p^c
\=(1+p^{frac{c}{2}})*sum(p,frac{c}{2}-1)+p^c
]
再代回约数和公式,直接利用递归来做分治即可。
(Code:)
#include<bits/stdc++.h>
using namespace std;
const int Mod=9901;
long long A,B,p[50],k[50],n,ans=1;
inline void input(void)
{
scanf("%lld%lld",&A,&B);
}
inline void decompose(void)
{
long long temp=A;
for(int i=2;i<=temp;i++)
{
if((temp%i)==0)
{
n++;
p[n]=i;k[n]++;
temp/=i;
while((temp%i)==0)
{
k[n]++;
temp/=i;
}
}
}
for(int i=1;i<=n;i++)
k[i]*=B;
}
inline long long power(long long a,long long b)
{
long long res=1;
for(;b;b>>=1)
{
if(1&b)res=res*a%Mod;
a=a*a%Mod;
}
return res;
}
inline long long sum(long long x,long long y)
{
if(y==1)return x+1;
if(y==0)return 1;
if(y%2==1)
return ((1+power(x,(y+1)/2))%Mod*sum(x,(y-1)/2)%Mod)%Mod;
else return ((1+power(x,y/2))%Mod*sum(x,y/2-1)%Mod+power(x,y)%Mod)%Mod;
}
inline void solve(void)
{
for(int i=1;i<=n;i++)
{
ans*=sum(p[i],k[i])%Mod;
ans%=Mod;
}
}
int main(void)
{
input();
decompose();
solve();
printf("%lld
",ans%Mod);
}
<后记>