poj1845 Sumdiv

题目大意

求$A^B$的约数和。

题目分析

求法是：假如求A的约数和，那么首先把A质因数分解，$A=p_1^{k_1}*p_2^{k_2}*p_3^{k_3}*...*p_n^{k_n}$，然后答案就是

$$ans=(1+p_1+p_1^2+p_1^3+...+p_1^{k_1})*(1+p_2+p_2^2+p_2^3+...+p_2^{k_2})*...*(1+p_n+p_n^2+p_n^3+...+p_n^{k_n})$$

解法1

注意：以下\号为整除。
二分解法，递归求解，对于$1+p_1+p_1^2+p_1^3+...+p_1^{k_1}$：
假如$k_1$为奇数，即有偶数项，那么可以化为$(1+p_1+p_1^2+...+p_1^{k_1/2})*(1+p_1^{k_1/2+1})$
假如$k_1$为偶数，那么可以化为$(1+p_1+p_1^2+...+p1^{k_1/2-1})*(1+p_1^{k_1/2+1})+p_1^{k_1/2}$
这样就可以啦！

#include<iostream>
#include<cstdio>
#include<climits>
#include<algorithm>
#include<cstring>
#include<cmath>
using namespace std;
#define ll long long
int mod=9901,a,b,tot,sum;
ll ksm(ll x,ll y){//快速幂
    ll re=1;
    while(y){
        if(y&1)re=(re*x)%mod;
        y>>=1;x=x*x%mod;
    }
    return re;
}
ll find(ll x,ll y){
    if(y==0)return 1;
    if(y%2)return (find(x,y/2)*(1+ksm(x,y/2+1)))%mod;//见上解析
    else return (find(x,y/2-1)*(1+ksm(x,y/2+1))+ksm(x,y/2))%mod;
}
void bob(){
    int i;ll num=0,ans=1;
    for(i=2;i*i<=a;){//寻找质因数，直接找就行
        if(a%i==0){
            num=0;while(!(a%i)){num++;a/=i;}
            ans*=find(i,num*b);ans%=mod;
        }
        if(i==2)i++;else i+=2;
    }
    if(a!=1){ans*=find(a,b);ans%=mod;}//特判最后一个质因数
    printf("%lld\n",ans);
}
int main()
{
    int i,j;
    while(scanf("%d%d",&a,&b)==2)bob();
    return 0;
}

解法2

咱们都知道根据方差，可以求出：

$$1+p_1+p_1^2+p_1^3+...+p_1^{k_1}= \frac{p_1^{k_1+1}-1}{p_1-1}$$
但是除法是不满足同余定理的，怎么办呢？
关于逆元有这么一条结论， 然而我并不会证明：
$$\frac{a}{b}mod_m=\frac{a_{mod}(m*b)}{b}$$
但是直接写快速幂乘的时候会炸…

#include<iostream>
#include<cstdio>
#include<climits>
#include<algorithm>
#include<cstring>
#include<cmath>
using namespace std;
#define ll long long
int mod=9901;ll a,b,tot,sum;
ll mul(ll x,ll y,ll m){
    ll re=0;
    while(y){
        if(y&1)re=(re+x)%m;
        y>>=1;x=(x+x)%m;
    }
    return re;
}
ll find(ll x,ll y,ll m){
    ll re=1;x%=m;
    while(y){
        if(y&1)re=mul(re,x,m);//直接相乘会溢出qwq
        y>>=1;x=mul(x,x,m);
    }
    return re;
}
void bob(){
    ll i;ll num=0,ans=1,m;
    for(i=2;i*i<=a;){
        if(a%i==0){
            num=0;while(!(a%i)){num++;a/=i;}
            m=(ll)(i-1)*mod;
            ans*=(find(i,num*b+1,m)+m-1)/(i-1);//神奇的公式
            ans%=mod;
        }
        if(i==2)i++;else i+=2;
    }
    if(a!=1){
        m=(ll)(a-1)*mod;
        ans*=(find(a,b+1,m)+m-1)/(a-1);
        ans%=mod;}
    printf("%lld\n",ans);
}
int main()
{
    int i,j;
    while(scanf("%lld%lld",&a,&b)==2)bob();
    return 0;
}