Sumdiv
Time Limit: 1000MS Memory Limit: 30000K
Description
Consider two natural numbers A and B. Let S be the sum of all natural divisors of A^B. Determine S modulo 9901 (the rest of the division of S by 9901).
Input
The only line contains the two natural numbers A and B, (0 <= A,B <= 50000000)separated by blanks.
Output
The only line of the output will contain S modulo 9901.
Sample Input
2 3
Sample Output
15
Hint
2^3 = 8.
The natural divisors of 8 are: 1,2,4,8. Their sum is 15.
15 modulo 9901 is 15 (that should be output).
Source
Romania OI 2002
解题思路
题意:求 AB A B 的所有约数和模9901
首先我们需要知道约数和公式:
由唯一分解定理,对 ∀A∈Z ∀ A ∈ Z ,可以将 A A 写作:,其中 pi p i 是质数,那么 A A 的约数和为
这道题中,可以先把
A
A
分解质因数,假设,那么
AB=pBk11pBk22⋯pBknn
A
B
=
p
1
B
k
1
p
2
B
k
2
⋯
p
n
B
k
n
,所以答案就是
∏∑j=0Bkipjimod9901
∏
∑
j
=
0
B
k
i
p
i
j
mod
9901
所以现在的问题变成了:
1. 如何分解质因数
将
A
A
不断除2直至除尽,然后不断除3直至除尽……由于是从小到大不断相
除,可以保证不会除到一个合数(因为
A
A
一定已经被这个合数的最小质因子除尽了),除的同时记录的因子及其出现次数即可
int fac = 2;
while(fac * fac <= A){
while(A % fac == 0){
if(p[p[0]] != fac) p[++p[0]] = fac;
k[p[0]]++;
A /= fac;
}
if(fac == 2) fac++;
else fac += 2;
}
if(A != 1) p[++p[0]] = A, k[p[0]] = 1;
2.如何快速地求出 ∑j=0npjmod9901 ∑ j = 0 n p j mod 9901 ,即 (p0+p1+⋯+pn)mod9901 ( p 0 + p 1 + ⋯ + p n ) mod 9901
- 法一:
- 如果
n
n
是奇数,
- 如果
n
n
是偶数,
- 求和部分递归即可,
pn
p
n
和
pn+1
p
n
+
1
快速幂即可
这样就可以 O(log2N) O ( l o g 2 N ) 地计算了
- 如果
n
n
是奇数,
- 法二:
当然,也可以直接用等比数列求和公式, (p0+p1+⋯+pn)=1−pn+11−p ( p 0 + p 1 + ⋯ + p n ) = 1 − p n + 1 1 − p ,但是由于要取模,所以注意乘逆元,另外,又因为 9901 9901 是一个小质数,用快速幂时还需要特判
Code
#include<cstdio>
#include<cstring>
using namespace std;
typedef long long LL;
const int MOD = 9901;
int A, B, p[10001], k[10001], ans;
LL fPow(LL bs, LL idx){
LL res = 1;
bs %= MOD;
while(idx){
if(idx & 1) (res *= bs) %= MOD;
(bs *= bs) %= MOD;
idx >>= 1;
}
return res;
}
LL cal(int p, LL k){
if(k == 0) return 1ll;
if(k & 1) return cal(p, k/2) * (1ll + fPow(p, k/2+1)) % MOD;
else return (fPow(p, k/2) + cal(p, k/2-1) * (1ll + fPow(p, k/2+1)) % MOD) % MOD;
}
int main(){
while(scanf("%d%d", &A, &B) != EOF){
ans = 1;
p[0] = 0;
int fac = 2;
while(fac * fac <= A){
while(A % fac == 0){
if(p[p[0]] != fac) p[++p[0]] = fac;
k[p[0]]++;
A /= fac;
}
if(fac == 2) fac++;
else fac += 2;
}
if(A != 1) p[++p[0]] = A, k[p[0]] = 1;
for(int i = 1; i <= p[0]; i++)
(ans *= cal(p[i], 1ll * k[i] * B)) %= MOD;
printf("%d\n", ans);
}
return 0;
}