poj1845 数论

 1 //Accepted    204K    16MS
 2 //约数和
 3 //n=p1^e1*p2^e2***pk^ek
 4 //约数和为:(p1^0+p1^1+..+p1^e1)*(p2^0+p2^1+..+p2^e2)*..(pk^0+pk^1+..pk^ek)
 5 //现考虑: S=p1^1+p1^2+..p1^e1
 6 // 令t=e1/2
 7 // if (e1%2==0)
 8 // S=(p1^1+p1^2+..+p1^t)+p1^t*(p1^1+p1^2+..+p1^t)
 9 // if (e1%2==1)
10 // S=(p1^1+p1^2+..+p1^t)+p1^t*(p1^1+p1^2+..+p1^t)+p1^e1
11 //由此可递归求解
12 #include <cstdio>
13 #include <cstring>
14 #include <iostream>
15 using namespace std;
16 const int imax_n = 10000;
17 const int pp = 9901;
18 int pri[imax_n];
19 int cnt;
20 void prime()
21 {
22     cnt=0;
23     for (int i=2;i<imax_n;i++)
24     {
25         for (int j=i*i;j<imax_n;j+=i)
26         pri[j]=1;
27     }
28     for (int i=2;i<imax_n;i++)
29     if (pri[i]==0)
30     pri[cnt++]=i;
31 }
32 int exp_mod(int a,int b)
33 {
34     int res=1;
35     a=a%pp;
36     while (b)
37     {
38         if (b&1) res=res*a%pp;
39         a=a*a%pp;
40         b>>=1;
41     }
42     return res;
43 }
44 int getSum(int a,int k)
45 {
46     if (k==0) return 1;
47     if (k==1) return (1+a)%pp;
48     if (k==2) return (1+a%pp+a%pp*a%pp)%pp;
49     int temp=getSum(a,(k+1)/2-1);
50     if (k%2==1)
51     return (temp+exp_mod(a,k/2+1)*temp%pp)%pp;
52     return (temp+exp_mod(a,k/2)*temp%pp+exp_mod(a,k))%pp;
53 }
54 int split(int n,int m)
55 {
56     int t=0;
57     int ans=1;
58     for (int i=0;i<cnt && (__int64 )pri[i]*pri[i]<=(__int64 )n;i++)
59     {
60         if (n%pri[i]==0)
61         {
62             t=0;
63             while (n%pri[i]==0)
64             {
65                 t++;
66                 n/=pri[i];
67             }
68             ans=getSum(pri[i],t*m)*ans%pp;
69         }
70     }
71     if (n!=1)
72     {
73         ans=getSum(n,m)*ans%pp;
74     }
75     return ans;
76 }
77 int main()
78 {
79     prime();
80     int a,b;
81     scanf("%d%d",&a,&b);
82     printf("%d\n",split(a,b));
83     return 0;
84 }