aks.sagews

# Step 1: Finding if n is a perfect power
# for b = 2 to floor(log(n,2))
# check if n^(1/b) is an integer
def is_perf_pow(n):
    for a in xrange(2, int((sqrt(n))) + 1):#fix this. sqrt gives a float
        for b in xrange(2, log(n, a) + 1):
            res = a^b
            if (res > n):
                break
            if (n == res):
                return "composite"
    return "not perfect power"

# Fixes "Python int too large to convert to C long" error for range() and xrange() in Python 2
def is_perf_pow2(n):
    a = 2
    a_lim = floor(sqrt(n)) + 1
    while a < a_lim:
        b = 2
        b_lim = log(n, a) + 1
        while b < b_lim:
            res = a^b
            if res > n:
                break
            if n == res:
                return "composite"
            b = b + 1
        a = a + 1
    return "not perfect power"

# r-1 is the worst case because that is the set of all possible remainders before the order is found
# r < log(n)^5 in general
def ord_mod(n, r):#calculates order for Step 2
    """
    Given the parameters (n, r), returns the order of n modulo r.
    """
    if r < 2:
        return 'mod must be greater than 2'
    if gcd(n, r) != 1:
        return "n and r must be coprime!"

    res = 0 # priming the loop
    i = 1
    while res != 1:
        res = n^i % r
        i = i + 1

    return i - 1


R.<x> = ZZ[]#Let R be the ring
# R.<x> = PolynomialRing(ZZ,1)

# Start where n^k > r. Example: Start at 2^5 when r = 21 (k > log(r, n))
# r - 1 > ord(n, r) > log(n, 2)^2 ... r > log(n, 2)^2 + 1
def smallest_r(n):
    """
    Finding the smallest r such that ord_mod(n, r) > log(n, 2)^2
    """
    ord = 0
    r = floor(log(n, 2)^2)
    #ord_array = []
    # incrementing r doesnt necessarily increase the order
    while ord <= log(n, 2)^2:
        r = r + 1
        if gcd(r,n) == 1:
            ord = ord_mod(n,r)
            #ord_array.append((ord,r))
    return r

# R.<x> = PolynomialRing(Integers(n)) vs power algorithm
# Step 5: Polynomial division/mod
def polynom_mod(r, n):
    R.<x> = Integers(n)[] # does modding at every step
    for a in xrange(1, ZZ(floor((sqrt(euler_phi(r)) * log(n, 2))))):#estimate squareroom in terms of log and figure out casting
        left_pol = (x + a)^n
        #right_pol = x^n + a
        # x^(n%r) + a is the remainder of (x^n + a) / (x^r - 1)
        if left_pol.quo_rem(x^r - 1)[1] != (x^(n%r) + a):
            return "composite"
            #return "composite"
    return "prime"


def aks(n):
    pow = is_perf_pow2(n)
    if pow == "composite":
        return pow
    r = smallest_r(n)
    for a in range(2, min(r, (n/2) + 1)):
        if gcd(a, n) > 1:
            return "composite"
    if n <= r:
        return "prime"
    else:
        p = polynom_mod(r, n)
        return p


#from sage.misc.html import HtmlFragment
#HtmlFragment('<iframe>test</iframe>')

#is_prime(216091)

@interact
def _(n = input_box(2)): print(aks(n))