poj sumdiv（数学）

将a进行质因数分解，则

        a=(p1^k1)*(p2^k2)*...*(pn^kn);

   故  a^b= p1^(a1*B) * p2^(a2*B) *...* pn^(an*B)；

根据公式，所有因子之和为

sum=(1+p1+p1^2+...p1^k1)*(1+p2+p2^2+...p2^k2)*...*(1+pn+pn^2+...+pn^kn)

计算1+p+p^2+...p^n可以利用二分进行加速

 

当n为奇数时,例如n=5

则1+p+p^2+p^3+p^4+p^5=(1+p+p^2)+p^3*(1+p+p^2)=(1+p+p^2)*(1+p^3)

可以发现1+p+p^2+...p^n=(1+p+p^2+...+p^(k/2))*(1+p^(k/2+1))

 

当n为偶数时，例如n=4

则1+p+p^2+p^3+p^4=(1+p)+p^3*(1+p)+p^2=(1+p)*(1+p^3)+p^2

可以发现1+p+p^2+...p^n=(1+p+p^2+...+p^(k/2-1))*(1+p^(k/2+1))+p^(k/2)

#include<string.h>
#include<cstdio>
const int N=7100;
#define M 9901
int e[N],factor[N];
int cnt;
__int64 a,b;
__int64 pow(__int64 p,__int64 n)        //p^n
{
	__int64 ret=1,s=p;
	while(n)
	{
		if(n&1)  
			ret=(ret*s)%M;
		s=(s*s)%M;
		n>>=1;
	}
	return ret;
}
__int64 get(__int64 p,__int64 k)       //1+p+p^2+...+p^k 
{
	if(k==0)
		return 1;
	if(k&1)        //如果k是奇数     
		return ((1+pow(p,k/2+1))%M*get(p,k/2)%M)%M;
	else
		return ((1+pow(p,k/2+1))%M*get(p,k/2-1)%M+pow(p,k/2)%M)%M;
}
int main()
{
	int i,x=0;
	__int64 ans;
	memset(e,0,sizeof(e));
	scanf("%I64d%I64d",&a,&b);
	for(i=2;i*i<=a;i++)             //求各个因子
		if(a%i==0)
		{
			factor[x]=i;
			while(a%i==0)
			{
				a/=i;
				e[x]++;
			}
			x++;
		}
	if(a>1)
	{
		factor[x]=a;
		e[x++]=1;
	}
	ans=1;
	for(i=0;i<x;i++)
	{
		__int64 tmp=get(factor[i],e[i]*b);
		ans=(ans*tmp)%M;
	}
	printf("%I64d\n",ans);
	return 0;
}

转载于:https://www.cnblogs.com/You-Change/p/3527932.html