找到素数的概率(使用 miller-rabin 检验)
Probability of finding a prime (using miller-rabin test)
我已经实施了 Miller-Rabin 素数测试,并且每个函数似乎都可以单独正常工作。但是,当我尝试通过生成 70 位随机数来找到素数时,我的程序在找到通过 Miller-Rabin 测试(10 步)的数字之前平均生成超过 100000 个数字。这很奇怪,对于一个小于 70 位的随机奇数,它是质数的概率应该非常高(根据 Hadamard-de la Vallée Poussin 定理,大于 1/50)。我的代码可能有什么问题?随机数生成器是否有可能以非常低的概率抛出质数?我猜不是...非常欢迎任何帮助。
import random
def miller_rabin_rounds(n, t):
'''Runs miller-rabin primallity test t times for n'''
# First find the values r and s such that 2^s * r = n - 1
r = (n - 1) / 2
s = 1
while r % 2 == 0:
s += 1
r /= 2
# Run the test t times
for i in range(t):
a = random.randint(2, n - 1)
y = power_remainder(a, r, n)
if y != 1 and y != n - 1:
# check there is no j for which (a^r)^(2^j) = -1 (mod n)
j = 0
while j < s - 1 and y != n - 1:
y = (y * y) % n
if y == 1:
return False
j += 1
if y != n - 1:
return False
return True
def power_remainder(a, k, n):
'''Computes (a^k) mod n efficiently by decomposing k into binary'''
r = 1
while k > 0:
if k % 2 != 0:
r = (r * a) % n
a = (a * a) % n
k //= 2
return r
def random_odd(n):
'''Generates a random odd number of max n bits'''
a = random.getrandbits(n)
if a % 2 == 0:
a -= 1
return a
if __name__ == '__main__':
t = 10 # Number of Miller-Rabin tests per number
bits = 70 # Number of bits of the random number
a = random_odd(bits)
count = 0
while not miller_rabin_rounds(a, t):
count += 1
if count % 10000 == 0:
print(count)
a = random_odd(bits)
print(a)
这在 python 2 而不是 python 3 中起作用的原因是两者处理整数除法的方式不同。在 python 2 中,3/2 = 1
,而在 python 3 中,3/2=1.5
.
看起来你应该在 python 3 中强制进行整数除法(而不是浮点数除法)。如果您更改代码以强制整数除法 (//
),如下所示:
# First find the values r and s such that 2^s * r = n - 1
r = (n - 1) // 2
s = 1
while r % 2 == 0:
s += 1
r //= 2
无论您使用什么 python 版本,您都应该看到正确的行为。
我已经实施了 Miller-Rabin 素数测试,并且每个函数似乎都可以单独正常工作。但是,当我尝试通过生成 70 位随机数来找到素数时,我的程序在找到通过 Miller-Rabin 测试(10 步)的数字之前平均生成超过 100000 个数字。这很奇怪,对于一个小于 70 位的随机奇数,它是质数的概率应该非常高(根据 Hadamard-de la Vallée Poussin 定理,大于 1/50)。我的代码可能有什么问题?随机数生成器是否有可能以非常低的概率抛出质数?我猜不是...非常欢迎任何帮助。
import random
def miller_rabin_rounds(n, t):
'''Runs miller-rabin primallity test t times for n'''
# First find the values r and s such that 2^s * r = n - 1
r = (n - 1) / 2
s = 1
while r % 2 == 0:
s += 1
r /= 2
# Run the test t times
for i in range(t):
a = random.randint(2, n - 1)
y = power_remainder(a, r, n)
if y != 1 and y != n - 1:
# check there is no j for which (a^r)^(2^j) = -1 (mod n)
j = 0
while j < s - 1 and y != n - 1:
y = (y * y) % n
if y == 1:
return False
j += 1
if y != n - 1:
return False
return True
def power_remainder(a, k, n):
'''Computes (a^k) mod n efficiently by decomposing k into binary'''
r = 1
while k > 0:
if k % 2 != 0:
r = (r * a) % n
a = (a * a) % n
k //= 2
return r
def random_odd(n):
'''Generates a random odd number of max n bits'''
a = random.getrandbits(n)
if a % 2 == 0:
a -= 1
return a
if __name__ == '__main__':
t = 10 # Number of Miller-Rabin tests per number
bits = 70 # Number of bits of the random number
a = random_odd(bits)
count = 0
while not miller_rabin_rounds(a, t):
count += 1
if count % 10000 == 0:
print(count)
a = random_odd(bits)
print(a)
这在 python 2 而不是 python 3 中起作用的原因是两者处理整数除法的方式不同。在 python 2 中,3/2 = 1
,而在 python 3 中,3/2=1.5
.
看起来你应该在 python 3 中强制进行整数除法(而不是浮点数除法)。如果您更改代码以强制整数除法 (//
),如下所示:
# First find the values r and s such that 2^s * r = n - 1
r = (n - 1) // 2
s = 1
while r % 2 == 0:
s += 1
r //= 2
无论您使用什么 python 版本,您都应该看到正确的行为。