日渐成熟的Pollard_Rho板子

namespace PR{
    map<ll,int >book;
    ll mul(ll a,ll b,ll MOD)
    {return (a*b-(ll)((long double)a/MOD*b)*MOD+MOD)%MOD;}
    ll pp[]={2,3,5,7,11,13,17};
    ll ksm(ll a,ll b,ll MOD)
    {
        ll ans=1;
        for(;b;b>>=1,a=mul(a,a,MOD))
            if(b&1LL)ans=mul(ans,a,MOD);
        return ans;
    }
    bool ip(ll x)
    {
        if(x==1)return 0;
        for(int i=0;i<7;i++)
            if(x==pp[i])return 1;
            else if(x%pp[i]==0)return 0;
        ll tmp=x-1;
        int cnt=0;
        while((tmp&1)^1)++cnt,tmp>>=1;
        for(int i=0;i<7;i++){
            ll val=ksm(pp[i],tmp,x),nxt;
            for(int op=1;op<=cnt;op++){
                nxt=mul(val,val,x);
                if(nxt==1&&val!=1&&val!=x-1)return 0;
                val=nxt;
            }
            if(val!=1)return 0;
        }
        for(int i=0;i<3;i++){
            ll val=ksm(rand()%(x-1)+1,tmp,x),nxt;
            for(int op=1;op<=cnt;op++){
                nxt=mul(val,val,x);
                if(nxt==1&&val!=1&&val!=x-1)return 0;
                val=nxt;
            }
            
            if(val!=1)return 0;
        }
        return 1;
    }
    ll Abs(ll x)
    {return x<0?-x:x;}
    ll Found(ll x,ll c)
    {
        ll a=rand()%(x-1)+1,b=a,res;
        for(int i=2,j=2;;i++){
            a=(mul(a,a,x)+c)%x;
            if(a==b)return -1;
            if((res=__gcd(x,Abs(a-b)))>1&&res<x)return res;
            if(i==j){
                b=a;
                j<<=1;
            }
        }
    }
    void Div(ll x)
    {
        
        if(x==1||ip(x)){
            book[x]++;
            return;
        }
        ll t,c=107;
        while((t=Found(x,c--))<0);
        Div(t),Div(x/t);
    }
}
posted @ 2019-05-12 19:46 蒟蒻小果冻 阅读(...) 评论(...) 编辑 收藏