zoukankan      html  css  js  c++  java
  • P4718[模板]PollardRho算法

    正题

    题目链接:https://www.luogu.com.cn/problem/P4718


    题目大意

    给出一个数\(n\),如果它是质数则输出\(Prime\),否则输出它的最大质因子。


    解题思路

    \(\text{Pollard-Rho}\)算法的前置知识是\(\text{Miller-Rabin}\)。在使用\(\text{Miller-Rabin}\)判掉质数之后,\(\text{Pollard-Rho}\)使用基于随机的思想能够较快的求出一个大数的因子之一。

    朴素的随机算法就是随机一个数判断它是不是因子,我们先使用一个较为优秀的随机方式,\(f(x)=f(x-1)^2+c\)(其中\(c\)为一个常数)。

    然后我们利用在这个函数上“跑”的距离来判断,也就是每次拿某两个\(i,j\),判断\(|f(i)-f(j)|\)是否为它的因数。

    但是如果枚举的话\(f\)函数上会出现一些“环”,我们需要快速的判掉“环”的方法。每次拿\(s,t\),令\(t=2s\),若环长为\(c\),那么有\(f(x)=f(x+c)\),当某一时刻\(f(t)=f(s)\)那么环长一定是\(s\)的整数倍。

    然后判到环就退出,如果没有找到就换一个常数重新做,这样的我们的算法雏形就形成了。

    但是这样还是跑的很慢,发现我们在过程中大量调用了\(gcd(d,p)\)导致时间变慢。考虑优化,我们可以每次先做一堆,然后在把这一堆拿过去一起搞定。首先我们有\(gcd(ac,b)=gcd(a,b)\),然后根据\(gcd\)的原理,我们有\(gcd(a,b)=gcd(a\% b,b)\)那么也就是我们有\(gcd(a,b)=gcd(ac\%b,b)\)

    那么假设我们有若干个间隔\(a_1,a_2,a_3,...\)那么我们把这数乘起来模\(p\),然后把得到的结果\(k\)\(p\)\(gcd\)就等价于拿\(a\)中逐个取与\(p\)\(gcd\)

    所以我们的优化方法就是第\(i\)次拿\(2^i\)个间隔去一起与\(p\)判断,但是因为\(i\)后面会很大导致副作用,所以将\(i\)设一个上界\(8\)即可。
    时间复杂度期望是\(O(n^{2.5})\),但跑的飞快

    回到这题来,我们先对\(n\)\(\text{MR}\)判断一次质数,然后跑\(\text{Pr}\)弄出一个因子\(d\),之后将\(n\)的因子\(d\)都去光后分别把\(n\)\(d\)丢下去递归继续跑。可以记录一个目前最大质因子来剪去一些不优状态。


    code

    #include<cstdio>
    #include<cstring>
    #include<algorithm>
    #include<cstdlib>
    #define ll long long
    using namespace std;
    const ll pri[10]={2,3,5,7,11,13,17,19,23,27};
    ll T,n,ans;
    ll ksc(ll a,ll b,ll p){
        ll c=(long double)a*b/p;
        long double ans=a*b-c*p;
        if(ans<0)ans+=p;
        else if(ans>=p)ans-=p;
        return ans;
    }
    ll power(ll x,ll b,ll p){
        ll ans=1;
        while(b){
            if(b&1)ans=ksc(ans,x,p);
            x=ksc(x,x,p);b>>=1;
        }
        return ans;
    }
    bool Mr(ll p){
        if(p==2)return 1;
        if(p<2||!(p&1))return 0;
        ll t=p-1,s=0;
        while(!(t&1))t>>=1,s++;
        for(ll i=0;i<10&&pri[i]<p;i++){
            ll x=power(pri[i],t,p),k;
            for(ll j=0;j<s;j++){
                k=ksc(x,x,p);
                if(k==1&&x!=1&&x!=p-1)
                    return 0;
                x=k;
            }
            if(x!=1)return 0;
        }
        return 1;
    }
    ll gcd(ll a,ll b)
    {return (!b)?a:gcd(b,a%b);}
    ll Pr(ll p){
        ll s=0,t=0,c=1ll*rand()%(p-1)+1;
        for(ll g=1,val=1,d;;g<<=1,s=t,val=1){
            for(ll j=0;j<g;j++){
                t=(ksc(t,t,p)+c)%p;
                val=ksc(val,abs(t-s),p);
                if(j%127==0&&(d=gcd(p,val))>1)
                    return d;
            }
            d=gcd(p,val);
            if(d>1)return d;
        }
        return p;
    }
    void solve(ll n){
        if(n<ans||n<2)return;
        if(Mr(n)){ans=n;return;}
        ll d=0;
        while((d=Pr(n))>=n);
        while(n%d==0)n/=d;
        solve(n);solve(d);
        return;
    }
    signed main()
    {
        srand(998244353);
        scanf("%lld",&T);
        while(T--){
            scanf("%lld",&n);
            if(Mr(n)){
                printf("Prime\n");
                continue;
            }
            ans=0;solve(n);
            printf("%lld\n",ans);
        }
        return 0;
    }
    
  • 相关阅读:
    搭建php环境时解决jpeg6 make: ./libtool:命令未找到
    configure: error: zlib not installed
    南京呼叫中心防火墙配置(备份)
    检查Linux Bash安全漏洞以及各环境修复解决方法
    Linux 内核升级步骤
    CentOS 7没有ifconfig命令处理
    linux kickstart 自动安装
    red hat Linux 使用CentOS yum源更新
    -bash: ./job.sh: /bin/sh^M: bad interpreter: 没有那个文件或目录
    linux LNMP自动安装脚本
  • 原文地址:https://www.cnblogs.com/QuantAsk/p/14253153.html
Copyright © 2011-2022 走看看