from random import randrange, SystemRandom from sys import maxsize # For Python 2, maxsize must be replaced by maxint , # As is the case once more below... def sift_primes(primes_list, pc): for dp in primes_list: if pc % dp == 0: return False primes_list.append(pc) return True def smallprimes(max = 65535): """ Compute a list of known primes. """ primes_list = [2, 3] i = 5 while i <= max: sift_primes(primes_list, i) i += 2 return primes_list def is_prime(n, k=192, def_list=[], min_s=0, totient_coprimes=[65537]): """ Test if a number is prime Args: n -- int -- the number to test k -- int -- the number of combined tests def_list -- int list -- definite known primes, if any, gapless and in ascending order Must start from 2 or 3 min_s -- int -- See below... totient_coprimes -- int list -- What the totient must not be divisible by, which is odd return True if n is prime, and n-1 is divisible by 2^min_s And as long as the probability of non-primeness is extremely small Hence, False could result when n is actually prime This probability can be reduced by setting (min_s = 2), thereby also saving time Setting (min_s > 2) or (min_s = 1) has no known benefits Setting (min_s = -1) will reduce the security, and guarantee that False is only returned, when n is not prime (Not recommended for cryptography) """ # Test if n is not even. # Prevent some error-loops. if min_s < 1 and n == 2: return True if min_s < 2 and n == 3: return True if n == 3: return False if n <= 1 or n % 2 == 0: return False """ First, make sure that the tentative prime - 1 will be coprime with a Public-Key Exponent. """ for dp in totient_coprimes: if dp <= 1 or dp % 2 == 0: continue if (n - 1) % dp == 0: return False # find r and s s = 0 r = n - 1 while r & 1 == 0: s += 1 r //= 2 """ A radical safeguard: Reject if s < 2 Even though n could be prime, Just because the available test would be no stronger than the Fermat Test. """ if min_s > 0 and s < min_s: return False """ Next, check against known primes. """ for dp in def_list: if dp >= n: return True if n % dp == 0: return False if len(def_list) > 0: highest_known_prime = def_list[len(def_list) - 1] highest_known_prime *= highest_known_prime if n < highest_known_prime: return True """ Additional safeguard: If a^r mod n == 1 ... The probability turned out to be 1 / (2^s) . If this probability is exceeded greatly, reject n. """ odd_ones = 0 # max_oo = (3 * k) // (1 << (s + 1)) # do k tests n1 = (1 << 32) - 1 for _ in range(k): # Fermat Test will no longer be performed... # Proceed with (modified) Miller-Rabin Test. a = randrange(2, n - 1) # Because the Mersennes Twister RNG is weak, mix it up # between even and odd... if a > n1 and a < n - n1: a += randrange(0, n1) x = pow(a, r, n) if min_s != -1 and x == 1: odd_ones += 1 if odd_ones == k: return False if x != 1 and x != n - 1: j = 1 while j < s and x != n - 1: x = pow(x, 2, n) if x == 1: return False j += 1 if x != n - 1: return False return True def generate_prime_candidate(length, rng, min_s=0): """ Generate an odd integer randomly Args: length -- int -- the length of the number to generate, in bits Minimum: min_s + 1 rng -- object -- System Random Number Generator min_s -- int -- Power of 2, which (prime - 1) must not be divisible by return a integer """ # Find chunk size of SystemRandom # For Python 2, maxsize must be replaced by maxint . maxint_v = 2 maxint_s = 1 while maxint_v < maxsize: maxint_s += 1 maxint_v = maxint_v << 1 maxint_v -= 1 # Create min_s related mask min_s_mask = 0 i = 1 while i < min_s: min_s_mask |= 1 << i i += 1 # generate random bits p = 0 remainder = length % maxint_s for _ in range((length + maxint_s - 1) // maxint_s): p = p << maxint_s p += rng.randint(0, maxint_v) if remainder > 0: p = p >> (maxint_s - remainder) """ If the primeness test was merely to fail, because (p - 1) was divisible by 2^min_s Then, as many requests for the system random number generator would be wasted. Therefore, the cause of such failures can be eliminated, before testing for primeness. """ # Apply min_s_mask if min_s > 1: p |= 1 << max(length - 1, min_s) p -= 1 << min_s p |= min_s_mask p += 2 # apply a mask to set MSB and LSB to 1 p |= (1 << length - 1) | 1 return p def generate_prime_number(length=1024): """ Generate a prime Args: length -- int -- length of the prime to generate, in bits return a prime """ rng = SystemRandom() # Invoke strictness min_s = 2 p = 4 def_list = smallprimes(65535) # keep generating while the primality test fail while not is_prime(p, 256, def_list, min_s): p = generate_prime_candidate(length, rng, min_s) return p print(generate_prime_number())