mod*_*tos 2 python algorithm primes
我的解决方案是否有问题找到最大的素数少于n什么时候n可以达到~10 ^ 230?对于更好的方法有什么建议吗?
这是我的尝试,在Python中使用以下版本的Miller-Rabin素性测试:
from random import randrange
small_primes = [
2, 3, 5, 7, 11, 13, 17, 19, 23, 29,
31, 37, 41, 43, 47, 53, 59, 61, 67, 71,
73, 79, 83, 89, 97,101,103,107,109,113,
127,131,137,139,149,151,157,163,167,173,
179,181,191,193,197,199,211,223,227,229,
233,239,241,251,257,263,269,271,277,281,
283,293,307,311,313,317,331,337,347,349,
353,359,367,373,379,383,389,397,401,409,
419,421,431,433,439,443,449,457,461,463,
467,479,487,491,499,503,509,521,523,541,
547,557,563,569,571,577,587,593,599,601,
607,613,617,619,631,641,643,647,653,659,
661,673,677,683,691,701,709,719,727,733,
739,743,751,757,761,769,773,787,797,809,
811,821,823,827,829,839,853,857,859,863,
877,881,883,887,907,911,919,929,937,941,
947,953,967,971,977,983,991,997
]
def probably_prime(n, k):
"""Return True if n passes k rounds of the Miller-Rabin primality
test (and is probably prime). Return False if n is proved to be
composite.
"""
if n < 2: return False
for p in small_primes:
if n < p * p: return True
if n % p == 0: return False
r, s = 0, n - 1
while s % 2 == 0:
r += 1
s //= 2
for _ in range(k):
a = randrange(2, n - 1)
x = pow(a, s, n)
if x == 1 or x == n - 1:
continue
for _ in range(r - 1):
x = pow(x, 2, n)
if x == n - 1:
break
else:
return False
return True
Run Code Online (Sandbox Code Playgroud)
我开始测试probably_prime(n)我减少的位置并测试每个值,n直到我得到一个"可能是素数"的数字.当我在n= ~10 ^ 230的值上测试时,我发现质数大约相差20-30个数字.在阅读了关于素数差距的更多信息后,我的结果似乎不太可能,因为我不应该如此频繁地找到素数.我测试k了高达50,000的值,我得到了相同的答案.我做错了什么,有什么建议可以找到更好的解决方案吗?
你是对的,你的代码一旦超出small_primes表格似乎就会遇到困难.仔细观察,这里有一个错误:
for _ in range(r - 1):
x = pow(x, 2, n)
if x == n - 1:
break
else:
return False
Run Code Online (Sandbox Code Playgroud)
你想要返回False(即复合),如果你从未找到x == n-1(或者你可以短路并返回,False如果x == 1,我认为:见这里).这可以通过更改缩进来完成:
for _ in range(r - 1):
x = pow(x, 2, n)
if x == n - 1:
break
else:
return False
Run Code Online (Sandbox Code Playgroud)
(for/else真的是组合for/if-not-break.)
完成此更改后,我得到:
>>> sum(orig(p, 20) for p in range(10**6, 2*10**6))
54745
>>> sum(fixed(p, 20) for p in range(10**6, 2*10**6))
70435
>>> sum(orig(p, 20) for p in range(10**230, 10**230+10**3))
40
>>> sum(fixed(p, 20) for p in range(10**230, 10**230+10**3))
2
Run Code Online (Sandbox Code Playgroud)
哪个是对的.