zoukankan      html  css  js  c++  java
  • 大数因数分解Pollard_rho 算法详解

    有一类问题,要求我们将一个正整数x,分解为两个非平凡因子(平凡因子为1与x)的乘积x=ab。

    显然我们需要先检测x是否为素数(如果是素数将无解),可以使用Miller-Rabin算法来进行测试。

    大数分解最简单的思想也是试除法,就是从2到sqrt(n),一个一个的试验,直到除到1或者循环完,最后判断一下是否已经除到1了即可。(当然这是幼稚做法,复杂度是相当高的,不然我就是想打试除法多方便呀

    Pollard Rho原理

    生日悖论

    如果一年只有365天(不计算闰年),且每个人的生日是任意一天的概率均相等,那么只要随机选取23人,就有50%以上的概率有两人同一天生日

    解释:第一个人不会和前面的人重生日(因为前面没有人),第二个人不重生日的概率为364/365,第三个人363/365……以此类推,那么只要到第23个人,就有,说明这时就有50%以上的概率有两人同生日

    更多的,当一年有n天时,只要人数到达Θ(sqrt(n))的数量级,有至少两个人同一天生日的概率就可以达到50%以上

    图象表达:

     

    利用生日悖论

    利用生日悖论来因数分解,重要的思想就是随机。

    已知我们随机地从[2,N-1]中选择出一个数是N的因数的概率是极小的,这也就意味着需要重复随机选择来提高正确率。

    那么如果我们不是只选择一个数,而是选择k个数,保证其中两个数的差值是N的因数。

    而如果其中两个数x,y的差值为N的因数,就一定会有gcd(|x-y|,N)>1

    假设N只有两个因数(除去自己和1)p和q的情况下,那么就意味着此时只有这两个数能整除N

    但是如果我们要求的是有多少个数x保证gcd(x,N)>1,此时答案就很多了,有p+q-2个

    于是我们就得出了一个简单的策略:

    1.在区间[2,N-1]中随机选取k个数,x1~k
    
    2.判断是否存在gcd(xi-xj,N)>1,若存在,则显然gcd(xi-xj,N)是N的一个因数

    同时也出现了一个问题,就是我们需要选取大约(N¼)个数,内存显然是不可能够的,那么又要如何解决这个问题呢?于是Pollard-rho算法就出现了。

    Pollard Rho算法

    因为数字太多内存不够,所以Pollard-rho只把连续的两个数存在内存中。

    也就是说,我们并不会选出k个数,而是一个一个地生成伪随机数,并检查连续的两个数是否符合条件。

    我们会用一个函数来生成伪随机数,就是这个→f(x)=(x2+a)%N

    其中的a可以自己指定或rand()生成,但是这样也伴随了一个问题就是这个函数生成的伪随机数还是有规律的,会无限循环。

    于是就会出现像希腊字母ρ一样的情况,这也是为什么这个算法名字中含有rho

    那么我们又要如何避免这种情况呢?

    首先为了保证没有答案可能被遗漏,那么至少要把这个环完整的扫一遍。

    联想一下一个比较常见的问题,就是小学数学题做过的两个人在环形操场上跑步,在同时起跑的情况下,当速度快的那个人追上速度慢的那个人的时候,一定已经多跑了一圈,也就是说此时两人肯定都至少跑完了一圈,恰好符合我们的要求。

    那么就是说我们要用两个变量来存储,一个用v的速度扫描环,一个用2v的速度,如果当两个变量相等时还没有找到答案,就退出这个环,重新取随机数,再次带入上面提到的函数中。

    这里有一点要说明一下,就是为什么快的速度一定要是慢的速度的两倍而不能更大。

    假设快的速度为kv(k>2),整个环的路程为s,快的追上慢的所需时间为t,那么我们可以列出式子:

    kvt−vt = s ⇒ (k−1)vt = s ⇒ vt = sk−1

    因为 k > 2,所以 k-1 > 1,那么就有 s/(k1) s

    ,也就意味着此时用来判断答案是否符合条件的(也就是速度较慢的那个)还没有扫描完整个环,那么就有可能会有答案被遗漏。

    以上就是完整的Pollard-rho算法过程,接下来上代码

    void Pollard_rho(long long N){
      long long a=rand()%N+1;//随机生成常数a
      long long x1,x2,d;
      x1=x2=rand()%N+1;
      while(1){
        x1=count(x1,cc);x2=count(count(x2,cc),cc);
        if(x1==x2) {x1=x2=rand()%N+1;cc=rand()%N+1;}
        d=gcd(abs(x2-x1),N);
        if(d>1&&d<N){p=d;q=N/p;return;}//p,q记录N的因数
      }
      return;
    }
    

      

  • 相关阅读:
    【luogu P7112】【模板】行列式求值(数学)(线性代数)(高斯消元)
    【luogu P3506】MOT-Monotonicity 2(线段树)
    【luogu AT5160】Numbers on a Circle(贪心)(堆)
    【luogu AT5159】RGB Balls(贪心)
    【luogu P5435】基于值域预处理的快速 GCD(数学)
    Day-20 基础模块1
    F5学习笔记
    Day-19 面向对象06 多继承 MRO C3算法 super()(重点)
    Day-18 面向对象05 约束、异常处理、日志处理
    MPLS(多协议卷标交换)
  • 原文地址:https://www.cnblogs.com/jaszzz/p/12693526.html
Copyright © 2011-2022 走看看