POJ1845-Sumdiv大数约数和

题目链接：http://poj.org/problem?id=1845

题目大意：

　　求A^B的所有约数和s。A和B都很大（0<=A,B<=50000000).

题目分析：

　　这道题让我学会了很多东西，很多我之前没有见识过的专属于数学的技巧之类的。看了别人的博客，受益匪浅。现总结如下：

　　先来分析一下这道题目，再来讲技巧吧。

　　这道题真是大数中的大数。刚看到这道题就完全没有头脑。全凭小优的指点（http://blog.csdn.net/lyy289065406）　。

　　（1） 整数唯一分解定理　

　　任何整数都可以且尽可以分解成若干个素数相乘的形式，如下：

　　A=p1^k1 * p2^k2 *……* pn^kn；其中的pi各表示一个素数.

　　（2）约数和公式

　　按照整数的分解定理，我们有一个整数的约数之和公式为：

　　S=（1+p1+p1^2+……p1^k1) * (1+p2+p2^2+……+p2^k2)  *……*  (1+pn+pn^2+……+pn^kn).

　　（3） 同余模公式

     　      (a+b)%mod=(a%mod + b%mod)%mod;

　　　　 (a*b)%mod=(a%mod * b%mod)%mod;

　　*以上的三个数学知识是解这道题的关键。请先消化理解以后再往下面看*

 

　　第二部份呢，我们来利用上面的理论基础分析一下这道题。

　　如果有：A=p1^k1 * p2^k2 *……* pn^kn；

　　那么有：A^B=p1^(k1*B) * p2^(k2*B) *……* pn^(kn*B);

　　因而   ： S=[1+p1+p1^2+……p1^(k1*B)]  *  [1+p2+p2^2+……+p2^(k2*B)]  *……*  [1+pn+pn^2+……+pn^(kn*B)].

　　至此，答案已经很明了了，下面的事情就是如何来用程序实现上面的过程了，注意，最终的S是要mod 9901 的哦！

　　*还有什么不明白的吗？可以在讨论区提问哦*

 

　　第三部份，我们来尝试用程序实现上述公式：

　　（1）质因数分解

　　在这里我用的根号法+递归法（名字是小优那里借鉴的），起初还没想明白这个方法，试图用素数表的方式来解决问题，看到五千万那么大就没有敲了，如果感兴趣，可以自己试一试哦。看一看下面的“根号法+递归法”吧。

for(int i=2;i*I<=A;) //根号法的体现
{
    int k=0;//这是p[ ]和n[ ]的指针。p是存素因子的，n是存其指数
    if (A%i==0)
    {
        p[k]=i;
        n[k]=0;
        while(A%i) //这是递归法的体现，找到一个素因子就用此法计算其个数
        {
            n[k]++;
            A/=i;
        }
        k++;
    }
    if(i==2) //这是所谓的奇偶法，除了2其余素数都是奇数哦
         i++;
    else
        i+=2;
}

if(A!=1)//常规来讲，这时候的A已经被分解剩下1了，除非A本身就是素数
{
    p[k]=A;
    n[k++]=1;
}

　　* 乍看上去，这也太慢了吧，怎么一个个去试，会不会出现重复啊？会不会把非质数当成是质因子啊？这个问题呢就留给大家自己思考啊。*

　　

　　（2）二分法求等比数列的和

　　解决了A的分解问题，自然需要来解决一下S的求解问题，很明显S是一系列以pi为公比的等比数列和 之积。只要能解决等比数列和的问题那这道题就迎刃而解了啊。最最直接的方法是利用求和公式，但是别忘了，我们的S可是还需要对9901取模的，[pi^(ki*B)-1]/(pi-1)这个结果中pi-1未必和9901互素！因而，解决这个问题，就只好用二分法：

　　　对于一个等比数列求和　　S=1+q+q^2+……q^n

　　　如果n为奇数，那么一共就有偶数个项了，

　　   S=[1+q+……+q^(n/2)] * [1+q^(n/2+1)]

　　　如果n为偶数，那么一共就是有奇数个项了，

　　　S=[1+q+……+q^(n/2-1] * [1+q^(n/2+1)]+q^(n/2)

　　*如果记不住上面的公式的话，就举个例子自己算算，结果自然就很清楚了。*

 

　　现在所有的理论问题都解决了，就看代码吧。

 

#include<iostream>
#include<cstdio>
using namespace std;
#define size 10000
#define mod 9901

int p[size],n[size];

long long pow3(long long int a,long long int b )//快速幂
{
    long long int r = 1, base = a;
    while( b != 0 )
    {
        if( b & 1 )
            r =(r * base) % mod;
        base =(base * base) % mod;
        b >>= 1;
    }
    return r;
}

long long sum(long long p,long long n)//二分递归求解等比数列之和
{
    if(n==0)
        return 1;
    if(n%2)
        return ((sum(p,n/2) % mod) * (1+pow3(p,n/2+1))% mod )% mod;
    else
        return (((sum(p,n/2-1) % mod) * (1+pow3(p,n/2+1))% mod ) % mod + pow3(p,n/2) %mod)%mod;
}

int main()
{
    int A,B;
    while(scanf("%d%d",&A,&B)!=EOF)
    {
        int js=0;
        for(int i=2;i*i<=A;)//质因数分解A，根号法和递归法
        {
            if(A%i==0)
            {
                p[js]=i;
                n[js]=0;
                while(A%i==0)
                {
                    n[js]++;
                    A/=i;
                }
                js++;
            }
            if(i==2)
                i++;
            else
                i+=2;
        }
        if(A!=1)//A本身就是素数的
        {
            p[js]=A;
            n[js++]=1;
        }

        long long ans=1;
        for(int i=0;i<js;i++)
        {
            ans=( ans * sum(p[i],n[i]*B)%mod )%mod;
        }
        cout<<ans<<endl;
    }
    return 0;
}