POJ1845 Sumdiv【因数之和】【二分求等比数列】

题目链接：

http://poj.org/problem?id=1845

题目大意：

对于两个自然数A和B，令s = A^B所有约数和，计算s%9901的结果是多少

思路：

小優的博客对这道题写的非常好，参考博文：http://blog.csdn.net/lyy289065406/article/details/6648539

对A进行素因子分解，A = p1^k1 * p2^k2 * … * pm^km，

则A^B = p1^(k1*B) * p2^(k2*B) * … *pm^(km*B)。

利用乘性函数的性质：

所偶因子和为 s = (1 + p1 + p1^2 + … + p1^(k1*B) ) * (1 + p2 + p2^2 + … + p2^(k2*B)) * … *

                           (1 + pm + pm^2 + … + pm^(km*B))

根据等比数列求和公式可得：s = ( 1+p1^(k1*B+1) )/(p1-1) * (1+p2^(k2*B+1) )/(p2-1) * … * 

  ( 1+pm^(km*B+1) )/(pm-1)

等比数列用二分法来求(这居然都能二分。。。)，p^n用整数快速幂来做。

AC代码：

#include<iostream>
#include<algorithm>
#include<cstdio>
#include<cstring>
#define LL __int64
using namespace std;

LL Power(LL p,LL n)
{
    LL ret = 1;
    while(n > 0)
    {
        if(n&1)
            ret = ret * p % 9901;
        p = p * p % 9901;
        n >>= 1;
    }
    return ret;
}

LL Sum(LL p,LL n)
{
    if(n == 0)
        return 1;

    if(n&1)
        return (Sum(p,n/2)*(1+Power(p,n/2+1)))%9901;
    else
        return (Sum(p,n/2-1)*(1+Power(p,n/2+1)) + Power(p,n/2))%9901;

}

LL P[10000],M[10000];
int main()
{
    int A,B;
    while(cin >> A >> B)
    {
        int ct = 0;
        memset(M,0,sizeof(M));
/*
        for(int i = 2; i<= A; ++i)
        {
            if(A % i == 0)
            {
                P[ct] = i;
                while(A % i == 0)
                {
                    M[ct]++;
                    A /= i;
                }
                ct++;
            }
        }
*/

        for(int i = 2; i*i <= A; )
        {
            if(A % i == 0)
            {
                P[ct] = i;
                while(A % i == 0)
                {
                    M[ct]++;
                    A /= i;
                }
                ct++;
            }
            if(i == 2)
                i++;
            else
                i += 2;
        }

        if(A != 1)
        {
            P[ct] = A;
            M[ct++] = 1;
        }

        LL ans = 1;
        for(int i = 0; i < ct; ++i)
            ans = ans*(Sum(P[i],M[i]*B)%9901)%9901;
        cout << ans << endl;
    }

    return 0;
}