# coding=utf-8
"""
Input: n > 2, an odd integer to be tested for primality;
       k, a parameter that determines the accuracy of the test
Output: composite if n is composite, otherwise probably prime
write n − 1 as 2s·d with d odd by factoring powers of 2 from n − 1
LOOP: repeat k times:
   pick a randomly in the range [2, n − 1]
   x ← ad mod n
   if x = 1 or x = n − 1 then do next LOOP
   for r = 1 .. s − 1
      x ← x2 mod n
      if x = 1 then return composite
      if x = n − 1 then do next LOOP
   return composite
return probably prime
"""


def _try_composite(a, d, n, s):
	if pow(a, d, n) == 1:
		return False
	for i in range(s):
		if pow(a, 2 ** i * d, n) == n - 1:
			return False
	return True  # n  is definitely composite


def is_prime(n, _precision_for_huge_n=16):
	if n in _known_primes or n in (0, 1):
		return True
	if any((n % p) == 0 for p in _known_primes):
		return False
	d, s = n - 1, 0
	while not d % 2:
		d, s = d >> 1, s + 1
	# Returns exact according to http://primes.utm.edu/prove/prove2_3.html
	if n < 1373653:
		return not any(_try_composite(a, d, n, s) for a in (2, 3))
	if n < 25326001:
		return not any(_try_composite(a, d, n, s) for a in (2, 3, 5))
	if n < 118670087467:
		if n == 3215031751:
			return False
		return not any(_try_composite(a, d, n, s) for a in (2, 3, 5, 7))
	if n < 2152302898747:
		return not any(_try_composite(a, d, n, s) for a in (2, 3, 5, 7, 11))
	if n < 3474749660383:
		return not any(_try_composite(a, d, n, s) for a in (2, 3, 5, 7, 11, 13))
	if n < 341550071728321:
		return not any(_try_composite(a, d, n, s) for a in (2, 3, 5, 7, 11, 13, 17))
	# otherwise
	return not any(_try_composite(a, d, n, s)
				   for a in _known_primes[:_precision_for_huge_n])


_known_primes = [2, 3]
_known_primes += [x for x in range(5, 1000, 2) if is_prime(x)]