POJ-1845 Sumdiv(费马小定理+快速幂)

　　原题连接：http://poj.org/problem?id=1845

　　题意：输入两个数A和B（0<=A,B<=50000000）。然后要求我们求出A^B的所有因子的和，并且对S = 9901取模。

　　大致思路：我们可以先将A使用唯一分解定理进行因子分解，即A = p1^a1 * p2^a2 * p3^a3 * p4^a4 .....pk^ak;其中p1,p2,p3,...pk都为素数。

A^B =  p1^(a1*B)  *  p2^(a2*B)  *  p3^(a3*B)  *  p4^(a4*B) .....pk^(ak*B)。

这个式子的所有约数和为(p1^0 + p1^1+p1^2 + ...+ p1^(a1*b)) * (p2^0 + p2^1 + p2^2 + ... + p2^(a2*B)) * .... * (pk^0 + pk^1 + pk^2 + ... + pk(ak*B));

每一项可以看成一个等比数列，每一项的结果可以写为 (pi^(ai*B)-1) / (pi-1) ;然后运用费马小定理可以写为(pi^(ai*B)-1) * (pi-1) mod (S-2)    (注意当p-1 = S时 无法使用 要单独讨论)。

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

typedef long long ll;

const int mod = 9901;

int p[1050],pnum[1050];            //记录因子，因子个数
int vid;        //记录下标 


void init(int s)            //唯一分解 
{
    memset(pnum,0,sizeof(pnum));
    vid = 0;
    for(int i = 2 ; i*i<=s;i++)
    {
        if(s%i==0)
        {
            p[++vid]=i;
            while(s%i==0)
            {
                pnum[vid]++;
                s=s/i;    
            }
        }
    }
    if(s>1)
    {
        p[++vid] = s;
        pnum[vid]=1;
    }
 } 
 
 
 int pows(int x,int y)    //求x^y并且对mod取模 
 {
     x = x%mod;
     if(x==0)
         return mod;
     int ans = 1;
     while(y)
     {
         if(y&1)
             ans = ans *x%mod;
         x = x*x%mod;
         y=y>>1;
    }
     return ans;
 }
 
ll fm(ll a,ll b)        //费马小定理 
 {
     return pows(a,b-2);
 }
 
 void solve(int A,int B)
 {
     ll ans = 1;
     int x,y;
     for(int i = 1 ; i <=vid;i++)
     {
         if(p[i]%mod==1)                    //单独讨论当 p-1 = mod时无法使用费马小定理的情况 
             ans = ans *(pnum[i]*B+1)%mod;            //此时 p[i]^(pnum[i]*B)的所有因子对mod取模都是1 故直接乘以因子的个数然后取模  
        else{
            ans = ans*(pows(p[i],pnum[i]*B+1)-1)%mod;        //利用费马小定理对等比数列求和  (p^n-1)/(p-1)%mod  = (p^n-1)*(p-1)%(mod-2)
            x = fm(p[i]-1,mod);
            ans = ans * x%mod;
        }
     }
     cout<<ans<<endl;
 }
 
 int main()
 {
     int A,B;
     while(scanf("%d%d",&A,&B)!=EOF)
     {
         if(A==0)
         {
             cout<<0<<endl;
             continue;
         }
        if(A==1||B==0)
        {
            cout<<1<<endl;
            continue;
        }
        init(A);
        solve(A,B);
     }
     return 0;
 }

 

这题还有一种用分治的方法对等比数列求和

pi^0+pi^1+pi^2+pi^3+…+pi^ai*b 可以从中部分成两段看成 
pi^0+pi^1+pi^2+pi^3+…+pi^(ai*b/2) 
pi^(ai*b/2+1)+pi^(ai*b/2+2)+pi^(ai*b/2+3)+…+pi^(ai*b/2)

合并后变成(1+pi^(ai*b/2+1))*？

这里的问号要分奇偶讨论