4.2 Pollard p-1算法

• 算法
• Python实现

  1974年，John Pollard发明了p-1算法，其灵感来自于费马小定理。p-1算法一个重要的概念就是B-powersmooth，这个不太好翻译，power，幂，smooth，光滑，就叫幂滑吧。其定义如下，对于质数p，将p-1进行质因数分解，最大的因子 d k d^k dk就是B，而这个质数就是B-powersmooth。举个例子，质数p=1303, p − 1 = 1302 = 2 ⋅ 3 ⋅ 7 ⋅ 31 p-1=1302=2cdot3cdot7cdot31 p−1=1302=2⋅3⋅7⋅31，那么质数1303就是31-powersmooth。再举个例子，质数3697， p − 1 = 3696 = 2 4 ⋅ 3 ⋅ 7 ⋅ 11 p-1=3696=2^4cdot3cdot7cdot11 p−1=3696=24⋅3⋅7⋅11，那么3697就是16-powersmooth。
  Pollard p-1算法基于费马小定理。费马小定理我们已经讲过无数遍了，实际上很多算法都是基于费马小定理的，比如米勒罗宾素性检验。我们从费马小定理开始展开讲讲，其公式如下：
a p − 1 ≡ 1    m o d    p a^{p-1}equiv1;mod;p ap−1≡1modp
  假使我们要对整数n进行质因数分解，他的其中一个质因数是p，数学里的整除符号是 p    ∣    n p;|;n p∣n。对于p，有：
d = a p − 1 − 1 = 0    m o d    p d=a^{p-1}-1=0;mod;p d=ap−1−1=0modp
  所以d是p的倍数，也就是 p    ∣    d p;|;d p∣d，可以得到 p = g c d ( d , n ) p=gcd(d,n) p=gcd(d,n)。讲得这么复杂，其实就是：
g c d ( a p − 1 − 1 , n ) = p gcd(a^{p-1}-1,n)=p gcd(ap−1−1,n)=p
  当然不可能一个个地尝试，那么和原始质因数分解没什么区别。Pollard p-1算法采用了一种随机数的策略。如果随机不出结果，那么就从小于B（这个B代表B-powersmooth）的质数列表里去连成，组成上述公式中的d，直到质数列表被用尽。

# _*_ encoding:utf-8 _*_

import random


def gcd(m, n):

    if m < n:
        m, n = n, m
    while m % n != 0:
        m, n = n, m % n
    return n


def sieve(n):
    mpfs = [0] * (n + 1)
    pn = []
    for i in range(2, n + 1):
        mpf = mpfs[i]
        if mpf == 0:
            pn.append(i)
            mpf = i
        for p in pn:
            if p > mpf:
                break
            x = p * i
            if x > n:
                break
            mpfs[x] = p
    return pn


def pollard(n):
    b = 10_0000
    while b <= 100_0000:
        a = random.randint(2, n - 1)
        g = gcd(a, n)
        if g > 1:
            return g, n // g
        for p in sieve(b):
            if p >= b:
                continue

            p_power = 1
            while p_power * p <= b:
                p_power *= p
            g = gcd(p_power - 1, n)
            if g > 1 and g < n:
                return g, n // g
        b <<= 1
    return 1, n

def factor(n):
    factors = []
    stack = [n]
    while len(stack) > 0:
        x = stack.pop()
        if x == 2:
            factors.insert(0, x)
            continue
        p, q = pollard(x) if x & 1 == 1 else (2, x >> 1)
        if p == 1:
            factors.insert(0, q)
        else:
            stack.append(p)
            stack.append(q)
    return factors


if __name__ == '__main__':
    print(factor(200))

  测试结果：

[2, 2, 2, 5, 5]