poj 1845 Sumdiv

Sumdiv

Time Limit: 1000MS		Memory Limit: 30000K
Total Submissions: 19529		Accepted: 4916

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). 

题目大意：求A^B所有约数之和

思路： 
由于乘性函数的性质 
题目等价于求(1+a1+a1^2+a1^3+a1^4+.....+a1^(b1k)) *(1+a2+a2^2+a2^3+a2^4+......+a2^(b2k))+......
=s1*s2*s3......*sn
根据等比数列求和易知
(1-an^(bnk+1))/1-an
an为约数，bn为该约数的个数 

由于出现除法取模，需要求逆元，记得讨论逆元不存在，逆元存在当且仅当gcd(x,MOD)==1

代码如下：

#include<cmath>
#include<cstdio>
#include<cstring>
#include<iostream>
#include<algorithm>
using namespace std;
typedef long long LL;
const LL maxn=100005;
const LL MOD=9901;
LL a[maxn],b[maxn];
LL quickpow(LL n,LL m)
{
	n%=MOD;
	if(n==0)return MOD-1;
	if(n==1)return 0;
	LL ans=1;
	while(m)
	{
		if(m&1)ans*=n,ans%=MOD;
		m>>=1;
		n*=n,n%=MOD;
	}
	ans=(ans%MOD+MOD)%MOD;
	return ans-1;
}
LL getfac(LL n)
{
	LL index=0;
	for(int i=2;i*i<=n;i++)
	{
		while(n%i==0&&n)b[index]++,n/=i;
		if(b[index])a[index++]=i;
	}
	if(n>1)b[index]++,a[index]=n,index++;
	return index;
}
LL x,y;
long long exgcd(long long a, long long b, long long &x, long long &y)
{
	if (a == 0){
		x = 0;
		y = 1;
		return b;
	}
	else
	{
		long long tx, ty;
		long long d = exgcd(b%a, a, tx, ty);
		x = ty - (b / a)*tx;
		y = tx;
		return d;
	}
}
int main()
{
	LL A,B;
	LL index=0;
	while(cin>>A>>B)
	{
		memset(b,0,sizeof(b));
		LL cnt=getfac(A),ans=1;
		for(int i=0;i<cnt;i++)
		{
			//if(a[i]%MOD==0)continue;
			if(a[i]%MOD==1){ans=ans*(b[i]*B+1),ans%=MOD;continue;}
			LL L=quickpow(a[i],b[i]*B+1);
			exgcd(a[i]-1,MOD,x,y);
			LL R=x;
			L=((L%MOD)+MOD)%MOD;
			R=((R%MOD)+MOD)%MOD;
			ans*=((L*R)%MOD+MOD)%MOD,ans%=MOD;
		}
		ans=((ans%MOD)+MOD)%MOD;
		cout<<ans<<endl;
	}
}