zoukankan      html  css  js  c++  java
  • Pollard_Rho 整数分解法

    如果要对比较大的整数分解,显然之前所学的筛选法和是试除法都将不再适用。所以我们需要学习速度更快的Pollard_Rho算法

    pollard_rho 算法流程

    Pollard_rho算法的大致流程是 先判断当前数是否是素数(Miller_rabin)了,如果是则直接返回。如果不是素数的话,试图找到当前数的一个因子(可以不是质因子)。然后递归对该因子和约去这个因子的另一个因子进行分解。

    那么自然的疑问就是,怎么找到当前数n的一个因子?当然不是一个一个慢慢试验,而是一种神奇的想法。其实这个找因子的过程我理解的不是非常透彻,感觉还是有一点儿试的意味,但不是盲目的枚举,而是一种随机化算法。我们假设要找的因子为p,他是随机取一个x1,由x1构造x2,使得{p可以整除x1-x2 && x1-x2不能整除n}则p=gcd(x1-x2,n),结果可能是1也可能不是1。如果不是1就找寻成功了一个因子,返回因子;如果是1就寻找失败,那么我们就要不断调整x2,具体的办法通常是x2=x2*x2+c(c是自己定的)直到出现x2出现了循环==x1了表示x1选取失败重新选取x1重复上述过程

    通过这个方式,我们只需要知道x1和c就可以生成一系列数值,c可以取任意给定值,一般取c=1。

    若序列出现循环,则退出。

    计算p=gcd(xi-1-xi, n), 若p=1,返回上一步,直到p>1为止;若p=n,则n为素数,否则p为一个约数并递归分解pn/p

      1 #include <iostream>
      2 #include <time.h>
      3 #include <algorithm>
      4 #include <stdio.h>
      5 
      6 typedef long long LL;
      7 
      8 using namespace std;
      9 
     10 const int times = 20;
     11 LL fac[1001];
     12 int cnt;
     13 
     14 LL mul(LL a,LL b,LL mod){
     15     LL ans = 0;
     16     while (b){
     17         if (b & 1){
     18             ans = (ans + a) % mod;
     19         }
     20         a = (a<<1) % mod;
     21         b >>= 1;
     22     }
     23     return ans;
     24 }
     25 
     26 
     27 LL pow(LL a,LL b,LL mod){
     28     LL ans = 1;
     29     while (b){
     30         if (b & 1){
     31             ans = mul(ans,a,mod);
     32         }
     33         b >>= 1;
     34         a = mul(a,a,mod);
     35     }
     36     return ans;
     37 }
     38 
     39 
     40 bool witness(LL a,LL n){
     41     LL temp = n - 1;
     42     int j = 0;
     43     while (temp % 2 == 0){  // 其实就是得到 m
     44         j++;
     45         temp /= 2;
     46     }
     47     LL x = pow(a,temp,n);
     48     if (x == 1 || x == n-1){   // 判断a^m
     49         return true;
     50     }
     51     while (j--){
     52         x = mul(x,x,n);  // 进一步判断 a^(2m)  a^(4m) ...
     53         if (x == n-1)
     54             return true;
     55     }
     56     return false;
     57 }
     58 
     59 bool miller_rabin(LL n){
     60     if (n == 2){ //  如果是2肯定是素数
     61         return true;
     62     }
     63     if (n<2 || n % 2 == 0){ //如果小于2或者是大于2的偶数肯定不是素数
     64         return false;
     65     }
     66     for (int i=0;i<times;i++){ //随机化检验
     67         LL a = rand() % (n-1) + 1;
     68         if (!witness(a,n))
     69             return false;
     70     }
     71     return true;
     72 }
     73 
     74 LL gcd(LL a,LL b){ // 这里的gcd和一般的gcd不一样
     75     if (a == 0){  // pollard_rho的需要
     76         return 1;
     77     }
     78     if (a < 0){ // 可能有负数
     79         return gcd(-a,b);
     80     }
     81     while (b){
     82         LL t = a % b;
     83         a = b;
     84         b = t;
     85     }
     86     return a;
     87 }
     88 
     89 LL pollard_rho(LL n,LL c){ // 找因子
     90     LL i = 1,k = 2; // 用来判断是否成环
     91     LL xx = rand() % n,y = xx;
     92     while (1){
     93         i++;
     94         xx = (mul(xx,xx,n) + c) % n;
     95         LL d = gcd(y-xx,n);
     96         if (1<d && d<n){ // 找到一个因数
     97             return d;
     98         }
     99         if (y == xx){ // 出现循环,那么查找失败
    100             return n;
    101         }
    102         if (i == k){ // 相当一个优化?
    103             y = xx;
    104             k <<= 1;
    105         }
    106     }
    107 }
    108 
    109 void find(LL n){ // 通过因数来找质因子
    110     if (miller_rabin(n)){
    111         fac[cnt++] = n; // 记录质因子
    112         return ;
    113     }
    114     LL p = n;
    115     while (p >= n)
    116         p = pollard_rho(p,rand() % (n-1) + 1); //如果转了一圈还是p那么继续
    117     find(p);
    118     find(n/p);
    119 }
    120 
    121 int main(){
    122     srand(time(NULL));
    123     int t;
    124     scanf("%d",&t);
    125     while (t--){
    126         LL x;
    127         scanf("%lld",&x);
    128         if (miller_rabin(x)){
    129             printf("Prime
    ");
    130             continue;
    131         }
    132         cnt = 0;
    133         find(x);
    134         LL ans = fac[0];
    135         for (int i=1;i<cnt;i++){
    136             if (fac[i]<ans)
    137                 ans = fac[i];
    138         }
    139         printf("%lld
    ",ans);
    140     }
    141     return 0;
    142 }
  • 相关阅读:
    Python解析Yahoo的XML格式的天气预报数据
    如何卸载wineQQ?
    纪念我的第一篇
    hihocoder1062 最近公共祖先·一
    hihocoder1055 刷油漆(树形DP)
    hihocoder1050 树中的最长路径
    hihocoder1049 根据二叉树的先序序列和中序序列求后序序列
    hihocoder1044 状态压缩
    hihocoder1043 完全背包
    hihocoder1038 01背包
  • 原文地址:https://www.cnblogs.com/-Ackerman/p/11354920.html
Copyright © 2011-2022 走看看