import random

#Es el algoritmo mas empleado, debido a su facilidad de implementacion. Sea p el numero que queremos saber si es primo

#Se calcula b, siendo b el numero de veces que 2 divide a (p - 1), es decir, 2 es la mayor potencia de 2 que divide a (p - 1).

def calculaB(p):
	nuevo = p - 1
	numero = 0
	while nuevo % 2 == 0:
		nuevo /= 2
		numero += 1
	return numero

#Calculamos entonces m, tal que p = 1 + 2^b * m.

def calculaM(p,b):
	return (p - 1)/ 2**b

def RabinMiller(numero):
	p = numero
	print "p",p

	b = calculaB(p)
	print "b",b
	m = calculaM(p,b)
	print "m",m

	#1. Escoger un numero aleatorio a < p.
	a = random.randrange(1,p)
	print "a",a

	#2. Sea j = 0 y z = a^m (mod p).
	j = 0
	z = (a**m)%p
	print "z",z
	
	
	#3. Si z = 1, o z = p - 1, entonces p pasa el test y puede ser primo.
	if z == 1 or z == p -1:
		return True
	
	while True:
		print "z",z
		#4. Si j > 0 y z = 1, p no es primo.
		if j > 0 and z == 1:
			print 4
			return False
		
		#5. Sea j = j + 1. Si j = b y z = p - 1, p no es primo.
		j = j + 1
		
		if j == b and z == p - 1:
			print 5
			return False
		
		#6. Si j < b y z = p - 1, z = z 2 (mod p). Volver al paso (4).
		if j < b and z == p - 1:
			z = z ** 2 % p
			continue
		
		#7. Si j < b y z = p - 1, entonces p pasa el test y puede ser primo.
		if j < b and z == p - 1:
			print 7
			return True
		
		#8. p no es primo.
		print 8
		return False



#La probabilidad de que un numero compuesto pase este algoritmo para un numero a es del 25 %. Esto quiere decir que necesitaremos menos pasos para llegar al mismo nivel de confianza que el obtenido con el algoritmo de Lehmann.


#Los primeros 1000 primos: http://www.utm.edu/research/primes/lists/small/1000.txt