• 数学#素数判定Miller_Rabin+大数因数分解Pollard_rho算法 POJ 1811&2429


    素数判定Miller_Rabin算法详解:

    http://blog.csdn.net/maxichu/article/details/45458569

    大数因数分解Pollard_rho算法详解:

    http://blog.csdn.net/maxichu/article/details/45459533

    然后是参考了kuangbin的模板:

    http://www.cnblogs.com/kuangbin/archive/2012/08/19/2646396.html

    模板如下:

    //快速乘 (a*b)%mod
    //二进制竖式乘法:
    //10101*1011=
    //10101*1+10101*2^1*1+10101*2^2*0+10101*2^3*1
    //位上是1的就加,是0的就不加
    ll mult_mod(ll a,ll b,ll mod)
    {
        a%=mod;
        b%=mod;
        ll res=0;
        while(b)
        {
            if(b&1)
            {
                res+=a;
                res%=mod;
            }
            a<<=1;
            if(a>=mod) a%=mod;
            b>>=1;
        }
        return res;
    }
    
    //快速幂 (x^n)%mod
    ll pow_mod(ll x,ll n,ll mod)
    {
        if(n==1) return x%mod;
        x%=mod;
        ll t=x;
        ll res=1;
        while(n)
        {
            if(n&1) res=mult_mod(res,t,mod);
            t=mult_mod(t,t,mod);//t平方
            n>>=1;//变成n/2
        }
        return res;
    }
    
    //若为合数返回true,不一定返回false
    //a^(n-1)=1(mod n) -> a^((2^t)*x)=-1(mod n) == a^x=1(mod n)
    bool test(ll a,ll n,ll x,ll t)//miller_rabin算法的核心
    {
        ll res=pow_mod(a,x,n);//a^x mod n
        ll last=res;
        for(int i=1;i<=t;i++)
        {
            res=mult_mod(res,res,n);//res=res*res mod n
            if(res==1&&last!=1&&last!=n-1)
                return true;
            last=res;
        }
        if(res!=1) return true;
        return false;
    }
    
    //若为素数(或伪素数)返回true,合数返回false
    bool miller_rabin(ll n)
    {
        if(n<2) return false;
        if(n==2) return true;
        if((n&1)==0) return false; //偶数
        ll x=n-1,t=0;
        while((x&1)==0)//n-1=(2^t)*x;
        {
            x>>=1;
            t++;
        }
        for(int i=0;i<times;i++)//进行随机判定
        {
            ll a=rand()%(n-1)+1;//随机找0~n-1的整数
            if(test(a,n,x,t))//
                return false;
        }
        return true;
    }
    
    ll factor[100];//保存质因数分解结果
    int tot;//记录质因数个数,下标从0开始
    
    ll gcd(ll a,ll b)
    {
        if(a==0) return 1;
        if(a<0) return gcd(-a,b);
        while(b)
        {
            ll c=a%b;
            a=b;
            b=c;
        }
        return a;
    }
    
    ll pollard_rho(ll x,ll c)
    {
        ll i=1,k=2;
        ll x0=rand()%x;
        ll y=x0;
        while(1)
        {
            i++;
            x0=(mult_mod(x0,x0,x)+c)%x;
            ll d=gcd(y-x0,x);
            if(d!=1&&d!=x) return d;
            if(y==x0) return x;
            if(i==k)
            {
                y=x0;
                k+=k;
            }
        }
    }
    
    //对n进行素因子分解
    void find_factor(ll n)
    {
        if(miller_rabin(n))//若为素数
        {
            factor[tot++]=n;
            return ;
        }
        ll p=n;
        while(p>=n)
            p=pollard_rho(p,rand()%(n-1)+1);
        find_factor(p);
        find_factor(n/p);
    }


    POJ 1811  完完全全的模板题

    #include<iostream>
    #include<cstdio>
    #include<cstdlib>
    #include<ctime>
    #include<algorithm>
    #include<cstring>
    using namespace std;
    typedef long long ll;
    
    const int times=20;//随机算法判定次数
    
    //快速乘 (a*b)%mod
    //二进制竖式乘法:
    //10101*1011=
    //10101*1+10101*2^1*1+10101*2^2*0+10101*2^3*1
    //位上是1的就加,是0的就不加
    ll mult_mod(ll a,ll b,ll mod)
    {
        a%=mod;
        b%=mod;
        ll res=0;
        while(b)
        {
            if(b&1)
            {
                res+=a;
                res%=mod;
            }
            a<<=1;
            if(a>=mod) a%=mod;
            b>>=1;
        }
        return res;
    }
    
    //快速幂 (x^n)%mod
    ll pow_mod(ll x,ll n,ll mod)
    {
        if(n==1) return x%mod;
        x%=mod;
        ll t=x;
        ll res=1;
        while(n)
        {
            if(n&1) res=mult_mod(res,t,mod);
            t=mult_mod(t,t,mod);//t平方
            n>>=1;//变成n/2
        }
        return res;
    }
    
    //若为合数返回true,不一定返回false
    //a^(n-1)=1(mod n) -> a^((2^t)*x)=-1(mod n) == a^x=1(mod n)
    bool test(ll a,ll n,ll x,ll t)//miller_rabin算法的核心
    {
        ll res=pow_mod(a,x,n);//a^x mod n
        ll last=res;
        for(int i=1;i<=t;i++)
        {
            res=mult_mod(res,res,n);//res=res*res mod n
            if(res==1&&last!=1&&last!=n-1)
                return true;
            last=res;
        }
        if(res!=1) return true;
        return false;
    }
    
    //若为素数(或伪素数)返回true,合数返回false
    bool miller_rabin(ll n)
    {
        if(n<2) return false;
        if(n==2) return true;
        if((n&1)==0) return false; //偶数
        ll x=n-1,t=0;
        while((x&1)==0)//n-1=(2^t)*x;
        {
            x>>=1;
            t++;
        }
        for(int i=0;i<times;i++)//进行随机判定
        {
            ll a=rand()%(n-1)+1;//随机找0~n-1的整数
            if(test(a,n,x,t))//
                return false;
        }
        return true;
    }
    
    ll factor[100];//保存质因数分解结果
    int tot;//记录质因数个数,下标从0开始
    
    ll gcd(ll a,ll b)
    {
        if(a==0) return 1;
        if(a<0) return gcd(-a,b);
        while(b)
        {
            ll c=a%b;
            a=b;
            b=c;
        }
        return a;
    }
    
    ll pollard_rho(ll x,ll c)
    {
        ll i=1,k=2;
        ll x0=rand()%x;
        ll y=x0;
        while(1)
        {
            i++;
            x0=(mult_mod(x0,x0,x)+c)%x;
            ll d=gcd(y-x0,x);
            if(d!=1&&d!=x) return d;
            if(y==x0) return x;
            if(i==k)
            {
                y=x0;
                k+=k;
            }
        }
    }
    
    //对n进行素因子分解
    void find_factor(ll n)
    {
        if(miller_rabin(n))//若为素数
        {
            factor[tot++]=n;
            return ;
        }
        ll p=n;
        while(p>=n)
            p=pollard_rho(p,rand()%(n-1)+1);
        find_factor(p);
        find_factor(n/p);
    }
    
    int main()
    {
        srand(time(0)); //需要ctime文件,poj的g++无法实现
        int t;
        scanf("%d",&t);
        ll n;
        while(t--)
        {
            scanf("%I64d",&n);
            if(miller_rabin(n))
            {
                printf("Prime
    ");
                continue;
            }
            tot=0;
            find_factor(n);
            sort(factor,factor+tot);
            printf("%I64d
    ",factor[0]);
        }
    }


    POJ 2689 

    这道题的函数是自己重新打了一遍的,因为理解不透彻,打算边打边想,结果还是没想清楚,尤其是后面的pollard_rho算法,而且还因为把long long型变量声明成了int型,导致一直TLE。。。

    AC代码如下:

    #include<iostream>
    #include<cstdio>
    #include<cstdlib>
    #include<ctime>
    #include<algorithm>
    #include<cstring>
    #include<cmath>
    using namespace std;
    typedef long long ll;
    
    const int times=20;//随机算法判定次数
    ll factor[1000],tot,ans,mina,k;
    
    ll gcd(ll a,ll b)
    {
        if(b==0)
            return a;
        return gcd(b,a%b);
    }
    
    ll mult_mod(ll a,ll b,ll mod)
    {
        a%=mod; b%=mod;
        ll res=0;
        while(b)
        {
            if(b&1)
            {
                res+=a;
                res%=mod;
            }
            a<<=1;
            if(a>=mod) a%=mod;
            b>>=1;
        }
        return res;
    }
    
    ll pow_mod(ll x,ll n,ll mod)
    {
        if(n==1) return x%mod;
        x%=mod;
        ll t=x,res=1;
        while(n)
        {
            if(n&1) res=mult_mod(res,t,mod);
            t=mult_mod(t,t,mod);
            n>>=1;
        }
        return res;
    }
    
    bool test(ll a,ll n)
    {
        ll x=n-1,t=0,res,last;
        while((x&1)==0)
        {
            x>>=1;
            t++;
        }
        last=pow_mod(a,x,n);
    
        for(int i=0;i<t;i++)
        {
            res=pow_mod(last,2,n);
            if(res==1&&last!=1&&last!=n-1)
                return true;
            last=res;
        }
        if(res!=1) return true;
        return false;
    }
    
    bool millier_rabin(ll n)
    {
        if(n==2) return true;
        if(n==1||(n&1)==0) return false;
        for(int i=0;i<times;i++)
        {
            ll a=rand()%(n-1)+1;
            if(test(a,n))
                return false;
        }
        return true;
    }
    
    ll pollard_rho(ll x,ll c)
    {
        ll x0,y,i=1,k=2;
        x0=rand()%x;
        y=x0;
        while(1)
        {
            i++;
            x0=(mult_mod(x0,x0,x)+c)%x;
            ll d=gcd(y-x0,x);
            if(d>1&&d<x) return d;
            if(y==x0) break;
            if(i==k)
            {
                y=x0;
                k+=k;
            }
        }
        return x;
    }
    
    void find_fac(ll n,int c)
    {
        if(n==1) return;
        if(millier_rabin(n))
        {
            factor[tot++]=n;
            return;
        }
        ll p=n;
        while(p>=n)
            p=pollard_rho(p,c--);
        find_fac(p,c);
        find_fac(n/p,c);
    }
    
    void dfs(ll i,ll x)
    {
        if(i>=tot)
        {
            if(x>mina&&x<=k)
                mina=x;
            return;
        }
        dfs(i+1,x);
        dfs(i+1,x*factor[i]);
    }
    
    int main()
    {
        ll n,m;
        while(~scanf_s("%I64d%I64d",&n,&m))
        {
            srand(time(0)); //需要ctime文件,poj的g++无法实现
            if(n==m)
            {
                printf("%I64d %I64d
    ",n,m);
                continue;
            }
            tot=0;
            ans=m/n;
            find_fac(ans,107);
            sort(factor,factor+tot);
            int i,j=0;
            for(i=1;i<tot;i++)
            {
                while(factor[i-1]==factor[i]&&i<tot)
                    factor[j]*=factor[i++];
                if(i<tot)
                    factor[++j]=factor[i];
            }
            tot=j+1; mina=1;
            k=(ll)sqrt(ans*1.0);
            dfs(0,1);
            printf("%I64d %I64d
    ",mina*n,ans/mina*n);
        }
        return 0;
    }
  • 相关阅读:
    python3字典删除元素和添加元素的几种方法
    查看ef core 3.1/3.0/2.1.2生成的sql语句
    C#Qrcode生成二维码支持带标题
    .net RabbitMQ 介绍、安装、运行
    Golang之初探
    MySql5.7多实例配置教程
    Centos7 安装MySql 5.7
    sqlserver智能提示插件-sql prompt(9.4.6)的安装及注册流程
    ABP框架(asp.net core 2.X+Vue)模板项目学习之路(二)--切换MySql数据库
    ABP框架(asp.net core 2.X+Vue)模板项目学习之路(一)
  • 原文地址:https://www.cnblogs.com/atmacmer/p/5264093.html
Copyright © 2020-2023  润新知