我正在尝试在Python上实现Miller-Rabin算法。 我已经按照教科书的伪代码所说的编码进行编码，但是由于某种原因，它无法正常工作。 要为sepcific，函数“测试”返回有时会返回“ true”。

Question

这里是一种不同的米勒·拉宾实施，它运行良好。 Millerrabin功能可能是您最感兴趣的：

def llinear_diophantinex(a, b, divmodx=1, x=1, y=0, withstats=False, pow_mod_p2=False): 
  origa, origb = a, b 
  r=a  
  q = a//b 
  prevq=1  
  if withstats == True: 
    print(f"a = {a}, b = {b}, q = {q}, r = {r}")   
  while r != 0:  
       prevr = r  
       a,r,b = b, b, r   
       q,r = divmod(a,b) 
       x, y = y, x - q * y 
       if withstats == True: 
         print(f"a = {a}, b = {b}, q = {q}, r = {r}, x = {x}, y = {y}")  
  y = 1 - origb*x // origa - 1 
  x,y=y,x 
  modx = (-abs(x)*divmodx)%origb 
  if withstats == True: 
    print(f"x = {x}, y = {y}, modx = {modx}") 
  if pow_mod_p2==False:   
    return x, y, modx 
  else: 
    if x < 0: return modx 
    else: return x 

def ltrailing(N):
    return len(str(bin(N))) - len(str(bin(N)).rstrip('0'))

def MillerRabin(N, primetest, iterx, powx, withstats=False): 
  primetest = pow(primetest, powx, N) 
  if withstats == True:
     print("first: ",primetest) 
  if primetest == 1 or primetest == N - 1: 
    return True 
  else: 
    for x in range(0, iterx-1): 
       primetest = pow(primetest, 2, N) 
       if withstats == True:
          print("else: ", primetest) 
       if primetest == N - 1: return True 
       if primetest == 1: return False 
  return False 

def sfactorint_isprime(N, withstats=False):
    if N == 2:
      return True
    if N % 2 == 0:
      return False
    if N < 2:
        return False
    iterx = ltrailing(N - 1)
    """ This k test is an algorithmic test builder instead of using
        random numbers. The offset of k, from -2 to +2 produces pow tests
        that fail or pass instead of having to use random numbers and more
        iterations. All you need are those 5 numbers from k to get a 
        primality answer. This is the same as doing:
        pow(N, (1<<N.bit_length()) - 1, 1<<N.bit_length()) but much faster
        using a linear diophantine formula which gives the same result for 
        powers of 2
    """
    k = llinear_diophantinex(N, 1<<N.bit_length(), pow_mod_p2=True) - 1
    t = N >> iterx
    tests = [k-2, k-1, k, k+1, k+2]
    for primetest in tests:
        if primetest >= N:
            primetest %= N
        if primetest >= 2:
            if MillerRabin(N, primetest, iterx, t, withstats) == False:
                return False
    return True

在Python中的Miller-Rabin算法的实现是：

# Check for primality using the Miller-Rabin algorithm
def is_prime(n, k = 40):
    import random
    if n < 2 or (n != 2 and not n & 1):
        return False
    if n < 6:
        return True
    random_gen = random.SystemRandom()
    for _ in range(k):
        a = random_gen.randrange(2, n - 1)
        exp = n - 1
        while not exp & 1:
            exp >>= 1
        if pow(a, exp, n) == 1:
            continue
        while exp < n - 1:
            if pow(a, exp, n) == n - 1:
                break
            exp <<= 1
        else:
            return False
    return True

Answer 1