miller_rabin_primality_test.sage

import random

def miller_rabin(n):
    nMin1 = n - 1
    a = random.randint(2,nMin1) # generate random integer
    k = 1 
    
    if n % 2 == 0:
        return print(n, "is probably prime")
    
    finding = True
    
    while finding:
        q = int(n/(2^k))
        
        if (2^k*q == nMin1) and (q % 2 != 0): finding = False 
        else: k = k + 1
    
    print("a:",a,"k:", k, "q:", q)
 
    if((a^q) % n == 1):
        print(n, "is probably prime")
        return
     
    for j in range(0, k):
        pow = 2 ^ j
        pow = q * pow
        
        if(a^(pow) % n == nMin1): 
            print(n, "is probably prime")
            return

    return print(n, "is composite")