在Python中使用大型素数

Question

使用Python处理大质数的有效方法是什么？你在这里或谷歌上搜索,你会发现许多不同的方法... sieves,primality test algorithms ...哪些方法适用于较大的素数？

Answer 1

为了确定一个数字是否是素数,有一个筛子和素数测试.

# for large numbers, xrange will throw an error.
# OverflowError: Python int too large to convert to C long
# to get over this:

def mrange(start, stop, step):
    while start < stop:
        yield start
        start += step

# benchmarked on an old single-core system with 2GB RAM.

from math import sqrt

def is_prime(num):
    if num == 2:
        return True
    if (num < 2) or (num % 2 == 0):
        return False
    return all(num % i for i in mrange(3, int(sqrt(num)) + 1, 2))

# benchmark is_prime(100**10-1) using mrange
# 10000 calls, 53191 per second.
# 60006 function calls in 0.190 seconds.

这似乎是最快的.not any你看到还有另一个版本,

def is_prime(num)
    # ...
    return not any(num % i == 0 for i in mrange(3, int(sqrt(num)) + 1, 2))

然而,在基准测试中,我接受70006 function calls in 0.272 seconds.了all 60006 function calls in 0.190 seconds.测试时的使用情况100**10-1.

如果您需要找到下一个最高素数,此方法将不适合您.您需要进行素性测试,我发现Miller-Rabin算法是一个不错的选择.Fermat方法稍微慢一点,但对伪赝假体更准确.使用上述方法在此系统上需要+5分钟.

Miller-Rabin 算法:

from random import randrange
def is_prime(n, k=10):
    if n == 2:
        return True
    if not n & 1:
        return False

    def check(a, s, d, n):
        x = pow(a, d, n)
        if x == 1:
            return True
        for i in xrange(s - 1):
            if x == n - 1:
                return True
            x = pow(x, 2, n)
        return x == n - 1

    s = 0
    d = n - 1

    while d % 2 == 0:
        d >>= 1
        s += 1

    for i in xrange(k):
        a = randrange(2, n - 1)
        if not check(a, s, d, n):
            return False
    return True

Fermat algoithm:

def is_prime(num):
    if num == 2:
        return True
    if not num & 1:
        return False
    return pow(2, num-1, num) == 1

获得下一个最高素数:

def next_prime(num):
    if (not num & 1) and (num != 2):
        num += 1
    if is_prime(num):
        num += 2
    while True:
        if is_prime(num):
            break
        num += 2
    return num

print next_prime(100**10-1) # returns `100000000000000000039`

# benchmark next_prime(100**10-1) using Miller-Rabin algorithm.
1000 calls, 337 per second.
258669 function calls in 2.971 seconds

使用该Fermat测试,我们得到了一个基准测试45006 function calls in 0.885 seconds.,但你的假冒伪劣率更高.

因此,如果只需要检查数字是否为素数,那么第一种方法is_prime就可以了.如果你使用它的话,它是最快的mrange.

理想情况下,您可能希望存储由此生成的素数next_prime并从中读取.

例如,使用next_prime与所述Miller-Rabin算法:

print next_prime(10^301)

# prints in 2.9s on the old single-core system, opposed to fermat's 2.8s
1000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000000000000000000000000000000000000000000000000000000000000000
0000000000000000000000000000000000000000000000000000000000000000000000000000
00000000000000000000000000000000000000000000000000000000000000000000000531

你无法return all(num % i for i in mrange(3, int(sqrt(num)) + 1, 2))及时做到这一点.我甚至不能在这个旧系统上做到这一点.

并且为了确保next_prime(10^301)并Miller-Rabin产生正确的值,还使用Fermat和Solovay-Strassen算法进行了测试.

请参阅:fermat.py,miller_rabin.py和solovay_strassen.py上gist.github.com

编辑:修正了一个错误 next_prime