1. 对您的代码的注释
我将在下面提出的一些观点在其他答案中得到了说明，但将它们放在一起似乎很有用。
在该部分
s = 0
r = n - 1
# factor n - 1 as 2^(r)*s
while r % 2 == 0:
s = s + 1
r = r // 2 # floor
您已经交换了r和s的角色：您实际上已将n - 1 分解为 2 s r。如果您想坚持使用 MathWorld 表示法，那么您必须在这部分代码中交换r和：s
# factor n - 1 as 2^(r)*s, where s is odd.
r, s = 0, n - 1
while s % 2 == 0:
r += 1
s //= 2
在行
for i in range(k):
该变量i未使用：通常命名此类变量_。
你选择一个介于 1 和n - 1之间的随机基数：
a = random.randrange(1, n)
这是它在 MathWorld 文章中所说的，但那篇文章是从数学家的角度写的。事实上，选择基数 1 是没有用的，因为 1 s = 1 (mod n ) 你会浪费试验。同样，选择底数n - 1 是没有用的，因为s是奇数，因此 ( n - 1) s = -1 (mod n )。数学家不必担心浪费的试验，但程序员会这样做，所以改为：
a = random.randrange(2, n - 1)
（n需要至少为 4 才能使此优化起作用，但是我们可以通过在nTrue = 3时返回函数顶部来轻松安排它，就像您为n = 2 所做的那样。）
正如其他回复中所述，您误解了 MathWorld 文章。当它说“ n通过测试”时，它意味着“ n通过了基数a的测试”。关于素数的显着事实是它们通过了所有碱基的测试。因此，当您发现a s = 1 (mod n ) 时，您应该做的是绕过循环并选择下一个要测试的基础。
# a^(s) = 1 (mod n)?
x = pow(a, s, n)
if x == 1:
continue
这里有优化的机会。我们刚刚计算的值x是2 0 s (mod n )。所以我们可以立即对其进行测试并为自己节省一次循环迭代：
# a^(s) = ±1 (mod n)?
x = pow(a, s, n)
if x == 1 or x == n - 1:
continue
在计算2 j s (mod n )的部分中，这些数字中的每一个都是前一个数字的平方 (modulo n )。当您可以将先前的值平方时，从头开始计算每个值是浪费的。所以你应该把这个循环写成：
# a^(2^(j) * s) = -1 (mod n)?
for _ in range(r - 1):
x = pow(x, 2, n)
if x == n - 1:
break
else:
return False
在尝试 Miller-Rabin 之前测试小素数的可分性是个好主意。例如，在拉宾 1977 年的论文中，他说：
在实现算法时，我们结合了一些省力的步骤。首先，我们测试任何素数p < N的可分性，例如N = 1000。
2. 修改代码
把所有这些放在一起：
from random import randrange
small_primes = [2, 3, 5, 7, 11, 13, 17, 19, 23, 29, 31] # etc.
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