【NOIP2017模拟9.2A组】春思

Description

Input

Output

Sample Input
2 3
Sample Output
15

这题有一个很高级的公式，就是计算因子和的公式。
要快速计算x的因子个数和因子和，首先将x分解质因数。
$X=p1^{q1}×p2^{q2}×p3^{q3}×…×pn^{qn}=∏pi^{qi}$
$X的因子和=（1+p1^1+p1^2+…+p1^{q1}）×（1+p2^1+p2^2+…+p2^{q2})×…×（1+pn^1+p1^2+…+pn^{qn}）=∏（1+pi^1+pi^2+⋯+pi^{qi}）$
$（1+pi^1+pi^2+⋯+pi^n）$那我们可以拆出1来，然后算出后面再加上1就好了。
那么我们考虑怎么算后面的；
①如果n是偶数，$Sn=pi^1+pi^2+⋯+pi^{\frac{n}{2}}+pi^{\frac{n}{2}+1}+⋯+pi^n$
我们把前面抽出来$pi^1+pi^2+⋯+pi^{\frac{n}{2}}$，后面的提出$pi^{\frac{n}{2}}$来，就变成了$pi^{\frac{n}{2}}*（pi^1+pi^2+⋯+pi^{\frac{n}{2}}）$，于是原式=$（1+pi^{\frac{n}{2}}）*S\frac{n}{2}$，于是用递归思想就解决了。(当然，$pi^{\frac{n}{2}}$用快速幂算)
②如果n是奇数，$Sn=pi^1+pi^2+⋯+pi^{\frac{n-1}{2}}+pi^{\frac{n-1}{2}+1}+⋯+pi^n$
我们把前面抽出来$pi^1+pi^2+⋯+pi^{\frac{n-1}{2}}$，后面的提出$pi^{\frac{n-1}{2}}$来，就变成了$pi^{\frac{n-1}{2}}*（pi^1+pi^2+⋯+pi^{\frac{n-1}{2}}$，于是原式=$(1+pi^{\frac{n-1}{2}}）*S\frac{n-1}{2}$,但我们发现，$\frac{n-1}{2}+\frac{n-1}{2}=n-1$，所以还要加上一个$pi^n$。
这样我们先给a分解质因数，然后每个的个数都乘上b，就可以套用这个公式去算啦。

#include<cstdio>
#include<iostream>
using namespace std;
long long num[1001][3];
long long a,b,ans,aa;
long long ksm(long long a,long long b)
{
    long long t=1,y=a%9901;
    while (b) 
    {
        if (b&1==1) t=t*y%9901;
        (y*=y)%=9901;
        b=b>>1;
    }
    return t;
}
long long find(long long x,long long y)
{
    if (y==1) return x;
    if (y%2==0) return (find(x,y/2)%9901)*(1+ksm(x,y/2));
    if (y%2==1) return ((find(x,(y-1)/2)%9901)*(1+ksm(x,(y-1)/2))%9901)+ksm(x,y);
}
int main()
{
    scanf("%lld%lld",&a,&b);
    long long i,j;
    j=2;aa=a;
    if (b==0) 
    {
        a=1;
        b=1;
    }
    ans=1;
    for (j=2;j*j<=aa;++j)
    {
        if (a%j==0) 
        {
            ++num[0][0];
            while (a%j==0)
            {
                num[num[0][0]][1]=j;
                ++num[num[0][0]][2];
                a/=j;   
            }
            num[num[0][0]][2]*=b;
            ans=(ans*(find(num[num[0][0]][1],num[num[0][0]][2])+1))%9901;
        }
    }
    if (a>1)
    {   
        num[num[0][0]][1]=a;
        num[num[0][0]][2]=b;
        ans=(ans*(find(num[num[0][0]][1],num[num[0][0]][2])+1))%9901;
    }
    printf("%lld\n",ans);
}