zoukankan      html  css  js  c++  java
  • 米勒罗宾大素数测定

    View Code
    /*
    判断素数和分解合数
    Miller-Rabin测试素数,Pollard_Rho分解质因子
    由于算法本身基于概率,所以存在TLE、WA的可能,多次提交即可
    */
    //////////////模板开始//////////////
    #include <ctime>
    #include <iostream>
    #include <algorithm>
    #include <cmath>
    using namespace std;
    #define M 1510
    #define N 100010
    #define inf 200000000
    long long factor[100], fac_top = -1;
    //计算两个数的gcd
    long long gcd(long long a,long long b) {
        if(a==0) return b;
        long long c;
        while(b!=0) {
            c=b;
            b=a%b;
            a=c;
        }
        return a;
    }
    //ret=(a*b)%n (n<2^62)
    long long muti_mod(long long a,long long b,long long n) {
        long long exp=a%n,res=0;
        while(b) {
            if(b&1) {
                res+=exp;
                if(res>n) res-=n;
            }
            exp<<=1;
            if(exp>n)
                exp-=n;
            b>>=1;
        }
        return res;
    }
    //ret=(a^b)%n
    long long mod_exp(long long a,long long p,long long m) {
        long long exp=a%m,res=1; //
        while(p>1) {
            if(p&1) res=muti_mod(res,exp,m);
            exp=muti_mod(exp,exp,m);
            p>>=1;
        }
        return muti_mod(res,exp,m);
    }
    //miller-rabin法测试素数,time 测试次数
    bool miller_rabin(long long n,int times) {
        if(n==2) return 1;
        if(n==1||(n!=2&&!(n%2))||(n!=3&&!(n%3))||(n!=5&&!(n%5))||(n!=7&&!(n%7)))
            return 0;
        long long a,u=n-1,x,y;
        int t=0;
        while(u%2==0) {
            ++t;
            u/=2;
        }
        srand(time(0));
        for(int i=0; i<times; i++) {
            a=rand()%(n-1)+1;
            x=mod_exp(a,u,n);
            for(int j=0; j<t; j++) {
                y=muti_mod(x,x,n);
                if(y==1&&x!=1&&x!=n-1) return false;
                x=y;
            }
            if(y!=1) return false;
        }
        return true;
    }
    
    long long pollard_rho(long long n,int c) {
        //找出一个因子
        long long x,y,d,i=1,k=2;
        srand(time(0));
        x=rand()%(n-1)+1;
        y=x;
        while(true) {
            ++i;
            x=(muti_mod(x,x,n)+c)%n;
            d=gcd(y-x,n);
            if(1<d&&d<n) return d;
            if(y==x) return n;
            if(i==k) {
                y=x;
                k<<=1;
            }
        }
    }
    
    void findFactor(long long n,int k) {
        //二分找出所有质因子,存入factor
        if(n==1) return;
        if(miller_rabin(n,6)) {
            factor[++fac_top]=n;
            return;
        }
        long long p=n;
        while(p>=n)
            p=pollard_rho(p,k--);//k值变化,防止死循环
        findFactor(p,k);
        findFactor(n/p,k);
    }
    //////////////模板完毕//////////////
    
    const int MAXN = 100010;
    int it, prime_num[MAXN], plist[MAXN] = {0};
    void getPrime() {
        for(it = 2; it <= MAXN; it++) {
            if(!plist[it])
                for(int j = it * 2; j <= MAXN; j += it)
                    plist[j] = 1;
        }
        it = 0;
        for(int i = 2; i <= MAXN; i++)
            if(!plist[i])
                prime_num[it++] = i;
    }
    
    int main() {
        int t;
        getPrime();
        scanf("%d", &t);
        while (t--) {
            int num = 0;
            long long n, ans = 0, now, pos;
            scanf("%I64d", &n);
            if (n <= 1000000000) {
                now = n;
                long long limit = sqrt((double)n);
                for (int i = 0; i < it && prime_num[i] <= limit; i++) {
                    if (n % prime_num[i] == 0) {
                        pos = prime_num[i];
                        num++;
                        long long tmp = prime_num[i];
                        n /= prime_num[i];
                        while (n % prime_num[i] == 0) {
                            tmp *= prime_num[i];
                            n /= prime_num[i];
                        }
                        ans += tmp;
                    }
                }
                if (n != 1)
                    num++, ans += n;
                if (num == 1) {
                    printf("%d %I64d\n", num, now/pos);
                    continue;
                }
                printf("%d %I64d\n", num, ans);
                continue;
            }
            fac_top = -1;
            if (miller_rabin(n, 10)) { //随机测试10次
                printf("1 %I64d\n", n);
            } else {
                findFactor(n, 30);
                sort(factor, factor+fac_top+1);
                long long tmp = factor[0];
                for (int i = 1; i <= fac_top; i++) {
                    if (factor[i] != factor[i-1]) {
                        ans += tmp;
                        tmp = factor[i];
                        num++;
                    } else {
                        tmp *= factor[i-1];
                    }
                }
                num++;
                ans += tmp;
                if (num == 1) {
                    printf("%d %I64d\n", num, n/factor[0]);
                    continue;
                }
                printf("%d %I64d\n", num, ans);
            }
        }
        return 0;
    }
  • 相关阅读:
    NLB网路负载均衡管理器详解
    Nginx配置详解
    Nginx代理功能与负载均衡详解
    .Net使用RabbitMQ详解
    说说面向服务的体系架构SOA
    .Net中的RealProxy实现AOP
    搭建自己的Nuget服务器
    VMware虚拟网络连接模式详解(NAT,Bridged,Host-only)
    JsonUtils
    Linux三剑客
  • 原文地址:https://www.cnblogs.com/gray035/p/3045788.html
Copyright © 2011-2022 走看看