# Alfonso Cevallos. May 5 2015.
# Please report any mistakes.
# This week we work with Python numbers, and built-in elementary operations.
# Therefore the algorithms will NOT be as fast as in our analysis.
# The goal here is simply to understand the higher level idea of these algorithms.

import math
import os
import random
import time

def trim(A): #removes any leading zeros
    while A[len(A)-1] == 0: A.pop()
    return

def exgcd(a,b): #computes (g,x,y) such that gcd(a,b) = g = xa + yb
    if a<b:
        g,y,x = exgcd(b,a)
        return g,x,y
    if b==0: return a,1,0

    r = a%b
    q = (a-r)/b
    g,x,y = exgcd(b,r)
    return g,y,x-y*q

def inverse(a,M):    #computes the inverse of a, modulo M
    g,x,y = exgcd(a,M)
    return x%M

def DFT(A,w,M):
    # Computes Y = DFT_w (f), where A=(a_0,...,a_n-1)
    # represents the polynomial f(x) = sum_i a_i x^i.
    # We assume n = len(A) is a power of 2, and w is a primitive
    # n-th root of unity; and we work modulo M.

    n = len(A)
    if n == 1: return A

    AE, AO = [], []     #we split the polynomial into odd and even coefficients
    for i in range(n/2):
        AE.append(A[2*i])
        AO.append(A[2*i+1])
    w2 = w*w
    YE = DFT(AE,w2,M)
    YO = DFT(AO,w2,M)

    h = 1
    Y = [0]*n
    for i in range(n/2):
        Y[i] = (YE[i] + h*YO[i]) % M
        Y[i+n/2] = (YE[i] - h*YO[i]) % M
        h = (h*w) % M
    return Y

def MultMod(A,B,w,M):
    # Computes the coefficient vector C, where h(x) = f(x)g(x) = sum_i c_i*x^i
    # and f(x) = sum_i a_i*x^i, g(x) = sum_i b_i*x^i.
    # We assume n=len(A)=len(B) is a power of 2, w is a primitive
    # n-th root of unity, and we work mod M.

    n=len(A)
    YF, YG, YH = DFT(A,w,M), DFT(B,w,M), []
    for i in range(n): YH.append((YF[i]*YG[i]) % M) #coordinate-wise multipl.

    C = DFT(YH,inverse(w,M),M)  #perform inverse DFT, then multiply by n inverse.
    nn = inverse(n,M)
    for i in range(n): C[i] = (nn * C[i]) % M
    return C

def Mult(A,B):
    # Computes the coefficient vector C, where h(x)=f(x)g(x)=sum_i c_i*x^i
    # and f(x)=sum_i a_i*x^i, g(x)=sum_i b_i*x^i.
    # Here we think of f, g and h as coming from Z[x] (not Z_M [x]).
    # We focus on finding the right values of M, n, and w,
    # and then we simply use the algorithm above.
    # Check the notes on exercise 3 to better understand this algorithm.

    trim(A)    #we remove unnecessary zeros
    trim(B)
    deg = (len(A)-1) + (len(B)-1)   #deg is the max degree of the product h
    n = 1               #n is the lowest power of 2 such that n > deg
    while n <= deg: n = n*2
    A = A + [0]*(n-len(A))
    B = B + [0]*(n-len(B))  #here we correct the lengths of A and B

    BB = 0      # this is the bound on the absolute value of the coefficients of f and g
    for a in A:
        if BB < abs(a): BB = abs(a)
    for b in B:
        if BB<abs(b): BB = abs(b)
    theRange = 2*(deg+1)*BB*BB    #this is the range of possible values of each coefficient in h

    K = n/2
    L = K
    while 2**L < theRange: L += K   #we define K and L as needed, with K dividing L, and M > theRange
    M = 2**L + 1
    w = 2**(L/K)    #this is a primitive n-th root of unity, as we proved in exercise 2

    C = MultMod(A,B,w,M)    #we use the modular multiplication
    for i in range(n):      # finally we map C back to Z[x]:
        if 2*C[i]>M: C[i]=C[i]-M    #large coefficients are mapped to negative values.
    trim(C)                         #and we remove unnecessary zeros
    return C

def Tests(): #we test with the polynomials given in exercise 4
    A = [3,-4,3,5,0,0,0,0]
    B = [-2,7,-5,2,0,0,0,0]
    C = MultMod(A,B,2,17)
    CC = Mult(A,B)

    print "For the polynomials A = ", A
    print "and B = ", B
    print "The product in Z_17[x] is ", C
    print "And the product in Z[x] is ", CC

Tests()