题目问$A^B$的所有因数和。

根据唯一分解定理将A进行因式分解可得:A = p1^a1 * p2^a2 * p3^a3 * pn^an.

A^B=p1^(a1*B)*p2^(a2*B)*...*pn^(an*B);

A^B的所有约数之和sum=[1+p1+p1^2+...+p1^(a1*B)]*[1+p2+p2^2+...+p2^(a2*B)]*[1+pn+pn^2+...+pn^(an*B)]

知道这个，问题就变成求出A的所有质因数pi以及个数n，然后$\prod(1+p_i+p_i^2+\cdots+p_i^{n-1}+p_i^n)$就行了。可以构造矩阵来求：

记$S_n=p_i+p_i^2+\cdots+p_i^{n-1}+p_i^n$

$$ \begin{bmatrix} p_i & 1 \\ 0 & 1 \end{bmatrix} \times \begin{bmatrix} S_n \\ p_i \end{bmatrix} = \begin{bmatrix} S_{n+1} \\ p_i \end{bmatrix} $$

 

$$ \begin{bmatrix} S_n \\ p_i \end{bmatrix} = \begin{bmatrix} p_i & 1 \\ 0 & 1 \end{bmatrix} ^n \times \begin{bmatrix} S_0 \\ p_i \end{bmatrix} $$

 

A忘了$\pmod {9901}$，爆intWA到头疼= =

1 #include
     
       2 #include
      
        3 using namespace std; 4 struct Mat{ 5     int m[2][2]; 6 }; 7 Mat operator*(const Mat &m1,const Mat &m2){ 8     Mat m={
   0}; 9     for(int i=0; i<2; ++i){10         for(int j=0; j<2; ++j){11             for(int k=0; k<2; ++k){12                 m.m[i][j]+=m1.m[i][k]*m2.m[k][j];13                 m.m[i][j]%=9901;14             }15         }16     }17     return m;18 }19 int calu(int a,int n){20     a%=9901;21     Mat e={
   1,0,0,1},x={a,1,0,1};22     while(n){23         if(n&1) e=e*x;24         x=x*x;25         n>>=1;26     }27     return (e.m[0][1]*a+1)%9901;28 }29 bool isPrime(int n){30     if(n<2) return 0;31     for(int i=2; i*i<=n; ++i){32         if(n%i==0) return 0;33     }34     return 1;35 }36 int main(){37     int a,b;38     scanf("%d%d",&a,&b);39     if(isPrime(a)){40         printf("%d",calu(a,b));41         return 0;42     }43     int res=1;44     for(int i=2; i*i<=a; ++i){45         if(a%i) continue;46         if(isPrime(i)){47             int cnt=0,tmp=a;48             while(tmp%i==0){49                 ++cnt;50                 tmp/=i;51             }52             res*=calu(i,cnt*b);53             res%=9901;54         }55         if(i!=a/i && isPrime(a/i)){56             int cnt=0,tmp=a;57             while(tmp%i==0){58                 ++cnt;59                 tmp/=i;60             }61             res*=calu(a/i,cnt*b);62             res%=9901;63         }64     }65     printf("%d",res);66     return 0;67 }