q-analog of KWMF by exhaustive search.ipynb

Kernel: SageMath (stable)

from sage.all import *
from Weight import Weight
from PartitionTree import PartitionTree
import __future__
import sys
import resource
import time


count = []

# Constants Declaration
bChangeA4 = matrix([[1, 0, 0, 0, 0],
                    [-1, 1, 0, 0, 0],
                    [0, -1, 1, 0, 0],
                    [0, 0, -1, 1, 0],
                    [0, 0, 0, -1, 1]]).inverse();

bChangeG2 = matrix([[2, 1, 0],
[1, 0, 0],
[1, 1, 1]])

bChangeF4 = matrix([[1, 1, 0, 0],
[2, 1, 1, 0],
[3, 1, 1, 1],
[2, 0, 0, 0]])

bChangeE6 = matrix([[0,   0,    0,    0,    0,   -2,    0,    0],
[ 1/2,  1/2,  1/2,  1/2,  1/2, -3/2,    0,    0],
[-1/2,  1/2,  1/2,  1/2,  1/2, -5/2,    0,    0],
[   0,    0,    1,    1,    1,   -3,    0,    0],
[   0,    0,    0,    1,    1,   -2,    0,    0],
[   0,    0,    0,    0,    1,   -1,    0,    0],
[   0,    0,    0,    0,    0,   -1,    1,    0],
[   0,    0,    0,    0,    0,    1,    0,    1]])

bChangeE7 = matrix([[   0,    0,    0,    0,    0,    0,   -2,    0],
[ 1/2,  1/2,  1/2,  1/2,  1/2,  1/2,   -2,    0],
[-1/2,  1/2,  1/2,  1/2,  1/2,  1/2,   -3,    0],
[   0,    0,    1,    1,    1,    1,   -4,    0],
[   0,    0,    0,    1,    1,    1,   -3,    0],
[   0,    0,    0,    0,    1,    1,   -2,    0],
[   0,    0,    0,    0,    0,    1,   -1,    0],
[   0,    0,    0,    0,    0,    0,    1,    1]])

bChangeE8 = matrix([[   0,    0,    0,    0,    0,    0,    0,    2],
[ 1/2,  1/2,  1/2,  1/2,  1/2,  1/2,  1/2,  5/2],
[-1/2,  1/2,  1/2,  1/2,  1/2,  1/2,  1/2,  7/2],
[   0,    0,    1,    1,    1,    1,    1,    5],
[   0,    0,    0,    1,    1,    1,    1,    4],
[   0,    0,    0,    0,    1,    1,    1,    3],
[   0,    0,    0,    0,    0,    1,    1,    2],
[   0,    0,    0,    0,    0,    0,    1,    1]])

basisChangeMatrices = {"G2": bChangeG2, "F4": bChangeF4, "E6": bChangeE6, "E7": bChangeE7, "E8": bChangeE8, "A4":bChangeA4}
numOfRoots = {"G2": 6, "F4": 4, "E6": 36, "E7": 63, "E8": 120}
vectorDimension = {"G2": 3, "F4": 4, "E6": 8, "E7": 8, "E8": 8}
num = 0

def getPartition(name, root):
    # used to change the basis from the standard basis of R^n to simple roots
    changeBasis = basisChangeMatrices[name]
    
    # this collect all of the positive roots for use in partitioning
    poset = [vector(list(Eval(x))) for x in RootSystem(name).ambient_space().positive_roots()]
    alphaPoset = [Weight(list(changeBasis * x)) for x in poset]
    print("Positive roots in terms of simple roots: ")
    for x in alphaPoset:
        print(x)
        
    print("\n")
    
    tree = PartitionTree(Weight(root), alphaPoset, 0, 0);
    tree.cleanTree()
    
    parts = []
    tree.getPartitions(parts);
    print("The partitions are:")
    for x in parts:
        print(x)
        
    print(len(parts))

def Eval(string):
    return eval(compile(str(string), '<string>', 'eval', __future__.division.compiler_flag))

def computeMq(name, lamb, mu):
    num = 0
    success = 0
    # initialize constants and vector space for the lie algebra
    lie_algebra = RootSystem(name).ambient_space()
    weyl_group = WeylGroup(name, prefix = "s")
    
    # equation is a vector that will hold the coefficients of Mq
    equation = vector([0 for i in range(0, numOfRoots[name] + 25)])
    
    # used to change the basis from the standard basis of R^n to simple roots
    changeBasis = basisChangeMatrices[name]
    
    # this collect all of the positive roots for use in partitioning
    poset = [vector(list(Eval(x))) for x in RootSystem(name).ambient_space().positive_roots()]
    alphaPoset = [Weight(list(changeBasis * x)) for x in poset]
    
    #sys.exit()
    
    # if lambda is not specified, the highest root is used
    if lamb == "none":
        lamb = lie_algebra.highest_root()
    else:
        lamb = changeBasis.inverse() * vector(lamb)
        lamb = weyl_group.domain()(list(eval(str(lamb))))
    
    # if mu is not specified, 0 vector is used
    if mu == "none":
        mu = weyl_group.domain()([0 for i in range(0, vectorDimension[name])])
    else:
        mu = changeBasis.inverse() * vector(mu)
        mu = weyl_group.domain()(list(eval(str(mu))))
    
    # rho is 1/2 the sum of the positive roots
    rho = lie_algebra.rho()
    
    for elm in weyl_group:
        
        # expression in partition function
        res = elm.action(lamb + rho) - (mu + rho)
        
        #change basis from standard basis to simple roots
        res = vector(list(Eval(res)))
        res = changeBasis * res
        
        # make a root object out of the result, check if it can be partitioned
        r = Weight(list(eval(str(res))))
        
        if not (r.isNegative() or r.hasFraction()):
            # root can be partitioned, find its Pq, add it to the sum
            tree = PartitionTree(r, alphaPoset, 0, 0)
            tree.cleanTree()
            
            Pq = [0 for x in range (0, numOfRoots[name] + 25)]
            tree.generatePq(Pq)
            
            # save results to a table
            count.append([str(elm), str(elm.length()), str(r), str(Pq)])
            
            Pq = ((-1)**elm.length()) * vector(Pq)
            
            equation = equation + Pq
            num += ((-1)**elm.length()) * tree.countPartitions()
            
            #print("Success: " + str(res))
            
            del(tree)
            del(Pq)
        #else:
            #print("Failure: " + str(res))
            
        del(res)
        del(r)
        
    return equation

start = time.time()
eq = computeMq("E7", "none", "none")
end = time.time()
print(str(end - start))

from Weight import Weight
import __future__
import sys
from sage.misc.parser import Parser
import time

# Here we define the variables that will appear in the expansion, since we are working with the exceptional Lie algebras we will need at most 8 A's and q
a, b, c, d, e, f, g, h, q = var('a, b, c, d, e, f, g, h, q')
# We build terms based on the simple roots, i.e. alpha1->A1, so we use this dictionary to easily make the correlation
variables = {'A1':a, 'A2':b, 'A3':c, 'A4':d, 'A5':e, 'A6':f, 'A7':g, 'A8':h, 'q':q}
P = Parser(make_var=var)

count = []

# Constants Declaration
bChangeA4 = matrix([[1, 0, 0, 0, 0],
                    [-1, 1, 0, 0, 0],
                    [0, -1, 1, 0, 0],
                    [0, 0, -1, 1, 0],
                    [0, 0, 0, -1, 1]]).inverse();

bChangeG2 = matrix([[2, 1, 0],
[1, 0, 0],
[1, 1, 1]])

bChangeF4 = matrix([[1, 1, 0, 0],
[2, 1, 1, 0],
[3, 1, 1, 1],
[2, 0, 0, 0]])

bChangeE6 = matrix([[0,   0,    0,    0,    0,   -2,    0,    0],
[ 1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2, -3.0/2,    0,    0],
[-1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2, -5.0/2,    0,    0],
[   0,    0,    1,    1,    1,   -3,    0,    0],
[   0,    0,    0,    1,    1,   -2,    0,    0],
[   0,    0,    0,    0,    1,   -1,    0,    0],
[   0,    0,    0,    0,    0,   -1,    1,    0],
[   0,    0,    0,    0,    0,    1,    0,    1]])

bChangeE7 = matrix([[   0,    0,    0,    0,    0,    0,   -2,    0],
[ 1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,   -2,    0],
[-1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,   -3,    0],
[   0,    0,    1,    1,    1,    1,   -4,    0],
[   0,    0,    0,    1,    1,    1,   -3,    0],
[   0,    0,    0,    0,    1,    1,   -2,    0],
[   0,    0,    0,    0,    0,    1,   -1,    0],
[   0,    0,    0,    0,    0,    0,    1,    1]])

bChangeE8 = matrix([[   0,    0,    0,    0,    0,    0,    0,    2],
[ 1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,  5.0/2],
[-1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,  7.0/2],
[   0,    0,    1,    1,    1,    1,    1,    5],
[   0,    0,    0,    1,    1,    1,    1,    4],
[   0,    0,    0,    0,    1,    1,    1,    3],
[   0,    0,    0,    0,    0,    1,    1,    2],
[   0,    0,    0,    0,    0,    0,    1,    1]])

basisChangeMatrices = {"G2": bChangeG2, "F4": bChangeF4, "E6": bChangeE6, "E7": bChangeE7, "E8": bChangeE8, "A4":bChangeA4}
numOfRoots = {"G2": 6, "F4": 4, "E6": 36, "E7": 63, "E8": 120}
vectorDimension = {"G2": 3, "F4": 4, "E6": 8, "E7": 8, "E8": 8}
num = 0

def Eval(string):
    return eval(compile(str(string), '<string>', 'eval', __future__.division.compiler_flag))

def convertToAs(weight):
    x = 1
    for i in range(0, len(weight)):
        for j in range(0, weight[i]):
            x = x * variables["A" + str(i+1)]
    return 1/(1 - q*x)

def convertToDs(weight):
    x = 1
    for i in range(0, len(weight)):
        for j in range(0, weight[i]):
            x = x * variables["A" + str(i+1)]
    return x


def computeMq(name, lamb, mu):
    # initialize constants and vector space for the lie algebra
    lie_algebra = RootSystem(name).ambient_space()
    weyl_group = WeylGroup(name, prefix = "s")

    # start mq = 0 then build it term by term
    equation = 0

    # used to change the basis from the standard basis of R^n to simple roots
    changeBasis = basisChangeMatrices[name]

    # this collect all of the poitive roots and builds the general terms for the summations
    poset = [vector(list(Eval(x))) for x in RootSystem(name).ambient_space().positive_roots()]
    alphaPoset = [changeBasis * x for x in poset]
    termsForSums = [convertToAs(list(x)) for x in alphaPoset]
    mainEq = 1
    for x in termsForSums:
        mainEq *= x

    # if lambda is not specified, the highest root is used
    if lamb == "none":
        lamb = lie_algebra.highest_root()
    else:
        lamb = changeBasis.inverse() * vector(lamb)
        lamb = weyl_group.domain()(list(eval(str(lamb))))

    # if mu is not specified, 0 vector is used
    if mu == "none":
        mu = weyl_group.domain()([0 for i in range(0, vectorDimension[name])])
    else:
        mu = changeBasis.inverse() * vector(mu)
        mu = weyl_group.domain()(list(eval(str(mu))))

    # rho is 1/2 the sum of the positive roots
    rho = lie_algebra.rho()

    for elm in weyl_group:
        # expression in partition function
        res = elm.action(lamb + rho) - (mu + rho)

        #change basis from standard basis to simple roots
        res = vector(list(Eval(res)))
        res = changeBasis * res
        r = Weight(list(res))

        if not (r.isNegative() or r.hasFraction()):
            #print(str(res) + " - " + str(elm))

            answer = mainEq

            for i in range(0, len(res)):
                    answer = answer.series(variables["A" + str(i+1)], res[i] + 1).truncate()
                    answer = answer.coefficient(variables["A" + str(i+1)], res[i])
                    #print(answer)

            answer = answer.expand()
            #print("FINAL: " + str(answer))

            # get the term that we are looking for
            #desiredterm = convertToDs(res)

            """
            NOTICE: The following code is just searching through the string representation of the polynomial to find the term we are interested in, and then stripping the A's so that we have the result of Pq, which we then use to compute the contribution to m_q
            """
            # terms in the expression are separated by +'s, split the string into an array of terms by separating at the +'s
            #terms = str(answer).split("+")
            #ar = []

            #for x in terms:
                # get the potential terms by checking to see if our desired term is present. Note that our term may be in a larger term that we do not want.
                #if str(desiredterm) in x:
                    #ar.append(x)

            # remove the A's so that we just have the result of Pq for the weight
            #for i in range(0, len(ar)):
                #ar[i] = ar[i].replace(str(desiredterm), '').replace("*", '').strip()

            # if an A is still present, then we know that our desired term was present in an unwanted larger term, so we get rid of it
            # notice that sage orders polynomials so that the term with the highest degree is at the beginning, so we do this operation in reversed order so that we get the terms with the lowest degree at the beginning, which is consistent with our other tables
            #f_ar = []
            #for i in range(len(ar)-1, -1, -1):
                #if not any(x in ar[i] for x in 'abcdefgh'):
                    #f_ar.append(ar[i])

            # saving to later print to a file
            count.append([str(elm), str(elm.length()), str(r), str(answer)])

            # add contribution to our running sum
            equation += answer * ((-1)**elm.length())
            #del(f_ar)
        del(res)
        del(r)

    return equation


start = time.time()
eq = computeMq("E6", "none", "none")
end = time.time()
print(str(end - start))

72

from Weight import Weight
import __future__
import sys
from sage.misc.parser import Parser
import time

# Constants Declaration
bChangeA4 = matrix([[1, 0, 0, 0, 0],
                    [-1, 1, 0, 0, 0],
                    [0, -1, 1, 0, 0],
                    [0, 0, -1, 1, 0],
                    [0, 0, 0, -1, 1]]).inverse();

bChangeG2 = matrix([[2, 1, 0],
[1, 0, 0],
[1, 1, 1]])

bChangeF4 = matrix([[1, 1, 0, 0],
[2, 1, 1, 0],
[3, 1, 1, 1],
[2, 0, 0, 0]])

bChangeE6 = matrix([[0,   0,    0,    0,    0,   -2,    0,    0],
[ 1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2, -3.0/2,    0,    0],
[-1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2, -5.0/2,    0,    0],
[   0,    0,    1,    1,    1,   -3,    0,    0],
[   0,    0,    0,    1,    1,   -2,    0,    0],
[   0,    0,    0,    0,    1,   -1,    0,    0],
[   0,    0,    0,    0,    0,   -1,    1,    0],
[   0,    0,    0,    0,    0,    1,    0,    1]])

bChangeE7 = matrix([[   0,    0,    0,    0,    0,    0,   -2,    0],
[ 1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,   -2,    0],
[-1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,   -3,    0],
[   0,    0,    1,    1,    1,    1,   -4,    0],
[   0,    0,    0,    1,    1,    1,   -3,    0],
[   0,    0,    0,    0,    1,    1,   -2,    0],
[   0,    0,    0,    0,    0,    1,   -1,    0],
[   0,    0,    0,    0,    0,    0,    1,    1]])

bChangeE8 = matrix([[   0,    0,    0,    0,    0,    0,    0,    2],
[ 1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,  5.0/2],
[-1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,  7.0/2],
[   0,    0,    1,    1,    1,    1,    1,    5],
[   0,    0,    0,    1,    1,    1,    1,    4],
[   0,    0,    0,    0,    1,    1,    1,    3],
[   0,    0,    0,    0,    0,    1,    1,    2],
[   0,    0,    0,    0,    0,    0,    1,    1]])

basisChangeMatrices = {"G2": bChangeG2, "F4": bChangeF4, "E6": bChangeE6, "E7": bChangeE7, "E8": bChangeE8, "A4":bChangeA4}
numOfRoots = {"G2": 6, "F4": 4, "E6": 36, "E7": 63, "E8": 120}
vectorDimension = {"G2": 3, "F4": 4, "E6": 8, "E7": 8, "E8": 8}
num = 0


def findAltSet(name, lamb, mu):
    # initialize constants and vector space for the lie algebra
    lie_algebra = RootSystem(name).ambient_space()
    weyl_group = WeylGroup(name, prefix = "s")
    s = weyl_group.gens()
        
    altset = [weyl_group.one()]

    # used to change the basis from the standard basis of R^n to simple roots
    changeBasis = basisChangeMatrices[name]

    # if lambda is not specified, the highest root is used
    if lamb == "none":
        lamb = lie_algebra.highest_root()
    else:
        lamb = changeBasis.inverse() * vector(lamb)
        lamb = weyl_group.domain()(list(eval(str(lamb))))

    # if mu is not specified, 0 vector is used
    if mu == "none":
        mu = weyl_group.domain()([0 for i in range(0, vectorDimension[name])])
    else:
        mu = changeBasis.inverse() * vector(mu)
        mu = weyl_group.domain()(list(eval(str(mu))))

    # rho is 1/2 the sum of the positive roots
    rho = lie_algebra.rho()
    length = len(altset)
    i=0
    while i < length:
        #print "i is:", i
        #print "length is:", length
        for simple in s:
            if ((altset[i] == simple)or (altset[i] == altset[i] * simple)):
                continue
            res = (altset[i]*simple).action(lamb + rho) - (rho + mu)
            res = changeBasis * vector(list(Eval(res)))
            res = Weight(res)
            #print(altset[i]*simple)
            #print(res)
            if not (res.isNegative() or res.hasFraction()):
                if not (altset[i]*simple in altset):
                    altset.append(altset[i]*simple)
                    length += 1
        #print("FINISHED: " + str(i))
        i+=1

    print "LENGTH OF ALT SET:", len(altset)
    return altset



start = time.time()
eq = findAltSet("E8", "none", "none")
end = time.time()
print(str(end - start))
fl = open("E8 altset.txt", "w+")
for x in eq:
    fl.write(str(x) + "\n")
fl.close()

LENGTH OF ALT SET: 2318
71.4812610149

from Weight import Weight
import __future__
import sys
from sage.misc.parser import Parser
import time

# Here we define the variables that will appear in the expansion, since we are working with the exceptional Lie algebras we will need at most 8 A's and q
a, b, c, d, e, f, g, h, q = var('a, b, c, d, e, f, g, h, q')
# We build terms based on the simple roots, i.e. alpha1->A1, so we use this dictionary to easily make the correlation
variables = {'A1':a, 'A2':b, 'A3':c, 'A4':d, 'A5':e, 'A6':f, 'A7':g, 'A8':h, 'q':q}
P = Parser(make_var=var)

count = []
fl = open("E8 stuff.txt", "w")
# Constants Declaration
bChangeA4 = matrix([[1, 0, 0, 0, 0],
                    [-1, 1, 0, 0, 0],
                    [0, -1, 1, 0, 0],
                    [0, 0, -1, 1, 0],
                    [0, 0, 0, -1, 1]]).inverse();

bChangeG2 = matrix([[2, 1, 0],
[1, 0, 0],
[1, 1, 1]])

bChangeF4 = matrix([[1, 1, 0, 0],
[2, 1, 1, 0],
[3, 1, 1, 1],
[2, 0, 0, 0]])

bChangeE6 = matrix([[0,   0,    0,    0,    0,   -2,    0,    0],
[ 1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2, -3.0/2,    0,    0],
[-1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2, -5.0/2,    0,    0],
[   0,    0,    1,    1,    1,   -3,    0,    0],
[   0,    0,    0,    1,    1,   -2,    0,    0],
[   0,    0,    0,    0,    1,   -1,    0,    0],
[   0,    0,    0,    0,    0,   -1,    1,    0],
[   0,    0,    0,    0,    0,    1,    0,    1]])

bChangeE7 = matrix([[   0,    0,    0,    0,    0,    0,   -2,    0],
[ 1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,   -2,    0],
[-1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,   -3,    0],
[   0,    0,    1,    1,    1,    1,   -4,    0],
[   0,    0,    0,    1,    1,    1,   -3,    0],
[   0,    0,    0,    0,    1,    1,   -2,    0],
[   0,    0,    0,    0,    0,    1,   -1,    0],
[   0,    0,    0,    0,    0,    0,    1,    1]])

bChangeE8 = matrix([[   0,    0,    0,    0,    0,    0,    0,    2],
[ 1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,  5.0/2],
[-1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,  1.0/2,  7.0/2],
[   0,    0,    1,    1,    1,    1,    1,    5],
[   0,    0,    0,    1,    1,    1,    1,    4],
[   0,    0,    0,    0,    1,    1,    1,    3],
[   0,    0,    0,    0,    0,    1,    1,    2],
[   0,    0,    0,    0,    0,    0,    1,    1]])

basisChangeMatrices = {"G2": bChangeG2, "F4": bChangeF4, "E6": bChangeE6, "E7": bChangeE7, "E8": bChangeE8, "A4":bChangeA4}
numOfRoots = {"G2": 6, "F4": 4, "E6": 36, "E7": 63, "E8": 120}
vectorDimension = {"G2": 3, "F4": 4, "E6": 8, "E7": 8, "E8": 8}
num = 0

def Eval(string):
    return eval(compile(str(string), '<string>', 'eval', __future__.division.compiler_flag))

def convertToAs(weight):
    x = 1
    for i in range(0, len(weight)):
        for j in range(0, weight[i]):
            x = x * variables["A" + str(i+1)]
    return 1/(1 - q*x)

def convertToDs(weight):
    x = 1
    for i in range(0, len(weight)):
        for j in range(0, weight[i]):
            x = x * variables["A" + str(i+1)]
    return x


def computeMq(name, lamb, mu):
    # initialize constants and vector space for the lie algebra
    lie_algebra = RootSystem(name).ambient_space()
    weyl_group = WeylGroup(name, prefix = "s")
    s = weyl_group.gens()
    
    refs = {"1": weyl_group.one()}
    for i in range(0, len(s)):
        refs["s" + str(i+1)] = s[i]

    # start mq = 0 then build it term by term
    equation = 0

    # used to change the basis from the standard basis of R^n to simple roots
    changeBasis = basisChangeMatrices[name]

    # this collect all of the poitive roots and builds the general terms for the summations
    poset = [vector(list(Eval(x))) for x in RootSystem(name).ambient_space().positive_roots()]
    alphaPoset = [changeBasis * x for x in poset]
    termsForSums = [convertToAs(list(x)) for x in alphaPoset]
    mainEq = 1
    for x in termsForSums:
        mainEq *= x

    # if lambda is not specified, the highest root is used
    if lamb == "none":
        lamb = lie_algebra.highest_root()
    else:
        lamb = changeBasis.inverse() * vector(lamb)
        lamb = weyl_group.domain()(list(eval(str(lamb))))

    # if mu is not specified, 0 vector is used
    if mu == "none":
        mu = weyl_group.domain()([0 for i in range(0, vectorDimension[name])])
    else:
        mu = changeBasis.inverse() * vector(mu)
        mu = weyl_group.domain()(list(eval(str(mu))))

    # rho is 1/2 the sum of the positive roots
    rho = lie_algebra.rho()
    inp = open("E8 altset.txt", "r")
    for line in inp:
        ed = line.strip().split("*")
        elm = weyl_group.one()
        for x in ed:
            elm *= refs[x]
        
        # expression in partition function
        res = elm.action(lamb + rho) - (mu + rho)

        #change basis from standard basis to simple roots
        res = vector(list(Eval(res)))
        res = changeBasis * res
        r = Weight(list(res))

        if not (r.isNegative() or r.hasFraction()):
            #print(str(res) + " - " + str(elm))

            answer = mainEq

            for i in range(0, len(res)):
                    answer = answer.series(variables["A" + str(i+1)], res[i] + 1).truncate()
                    answer = answer.coefficient(variables["A" + str(i+1)], res[i])
                    #print(answer)

            answer = answer.expand()
            #print("FINAL: " + str(answer))

            # get the term that we are looking for
            #desiredterm = convertToDs(res)

            """
            NOTICE: The following code is just searching through the string representation of the polynomial to find the term we are interested in, and then stripping the A's so that we have the result of Pq, which we then use to compute the contribution to m_q
            """
            # terms in the expression are separated by +'s, split the string into an array of terms by separating at the +'s
            #terms = str(answer).split("+")
            #ar = []

            #for x in terms:
                # get the potential terms by checking to see if our desired term is present. Note that our term may be in a larger term that we do not want.
                #if str(desiredterm) in x:
                    #ar.append(x)

            # remove the A's so that we just have the result of Pq for the weight
            #for i in range(0, len(ar)):
                #ar[i] = ar[i].replace(str(desiredterm), '').replace("*", '').strip()

            # if an A is still present, then we know that our desired term was present in an unwanted larger term, so we get rid of it
            # notice that sage orders polynomials so that the term with the highest degree is at the beginning, so we do this operation in reversed order so that we get the terms with the lowest degree at the beginning, which is consistent with our other tables
            #f_ar = []
            #for i in range(len(ar)-1, -1, -1):
                #if not any(x in ar[i] for x in 'abcdefgh'):
                    #f_ar.append(ar[i])

            # saving to later print to a file
            rans = str(answer).split("+")
            ans = []
            for i in range(len(rans) - 1, -1, -1):
                step = rans[i].split("^")
                if len(step) > 1:
                    step[1] = "{" + step[1] + "}"
                step = "^".join(step)
                ans.append(step)
                
            fl.write(str([str(elm), str(elm.length()), str([int(i) for i in r.coefficients]), "+".join(ans)]) + "\n")

            # add contribution to our running sum
            equation += answer * ((-1)**elm.length())
            #del(f_ar)
        del(res)
        del(r)

    return equation


start = time.time()
eq = computeMq("E8", "none", "none")
end = time.time()
print(str(end - start))

fl.close()
print(eq)

import __future__
def Eval(string):
    return eval(compile(str(string), '<string>', 'eval', __future__.division.compiler_flag))

def getBasisChange(name):
    simples = RootSystem(name).ambient_space().simple_roots()
    bChange = ([Eval((str(x).replace("(", "[").replace(")", "]"))) for x in simples])
    for i in range(len(bChange), len(bChange[0])):
        bChange.append([1 if j==i else 0 for j in range(0, len(bChange[0]))])
    return matrix(bChange).transpose().inverse()

def geometricSumForPartition(positive_root, translations, q_analog):
    x = 1
    for i in range(0, len(positive_root)):
        for j in range(0, positive_root[i]):
            x = x * translations["A" + str(i+1)]
    return 1/(1 - x) if not q_analog else 1/(1 -q*x)

def calculatePartition(name, weight, positive_roots = [], translations = {}, q_analog = False):
    if positive_roots == []:
        bChange = getBasisChange(name)
        positive_roots = [vector(list(Eval(x))) for x in RootSystem(name).ambient_space().positive_roots()]
        positive_roots = [bChange * x for x in positive_roots]
        
    if translations == {}:
        s = ''
        if q_analog:
            s = 'q, '
        s += 'A1'
        for i in range(1, len(weight)):
            s += ', A' + str(i + 1)
        variables = var(s)
        for i in range(0, len(weight)):
            translations["A" + str(i+1)] = eval("A" + str(i+1))
    
    termsForSum = [geometricSumForPartition(list(x), translations, q_analog) for x in positive_roots]
    answer = 1
    for x in termsForSum:
        answer *= x
        
    for i in range(0, len(weight)):
        answer = answer.series(translations["A" + str(i+1)], weight[i] + 1).truncate()
        answer = answer.coefficient(translations["A" + str(i+1)], weight[i])

    answer = answer.expand()
    return answer
    
def findAltSet(name, lamb = None, mu = None):
    # initialize constants and vector space for the lie algebra
    lie_algebra = RootSystem(name).ambient_space()
    weyl_group = WeylGroup(name, prefix = "s")
    simples = weyl_group.gens()
        
    altset = [weyl_group.one()]

    # used to change the basis from the standard basis of R^n to simple roots
    changeBasis = getBasisChange(name)

    # if lambda is not specified, the highest root is used
    if lamb == None:
        lamb = lie_algebra.highest_root()
    else:
        lamb = changeBasis.inverse() * vector(lamb)
        lamb = weyl_group.domain()(list(eval(str(lamb))))

    # if mu is not specified, 0 vector is used
    if mu == None:
        mu = weyl_group.domain()([0 for i in range(0, len(lie_algebra.simple_roots()))])
    else:
        mu = changeBasis.inverse() * vector(mu)
        mu = weyl_group.domain()(list(eval(str(mu))))

    rho = lie_algebra.rho()
    length = len(altset)
    i=0
    while i < length:
        for simple in simples:
            if ((altset[i] == simple)or (altset[i] == altset[i] * simple)):
                continue
            res = (altset[i]*simple).action(lamb + rho) - (rho + mu)
            res = changeBasis * vector(list(Eval(res)))
            res = Weight(res)
            
            if not (res.isNegative() or res.hasFraction()):
                if not (altset[i]*simple in altset):
                    altset.append(altset[i]*simple)
                    length += 1
        i+=1
        
    return altset
    

def calculateMultiplicity(name, lamb = None, mu = None, q_analog = False):
    mult = 0
    lie_algebra = RootSystem(name).ambient_space()
    weyl_group = WeylGroup(name, prefix = "s")
    
    # used to change the basis from the standard basis of R^n to simple roots
    changeBasis = getBasisChange(name)
    
    positive_roots = [vector(list(Eval(x))) for x in RootSystem(name).ambient_space().positive_roots()]
    positive_roots = [getBasisChange(name) * x for x in positive_roots]

    # if lambda is not specified, the highest root is used
    if lamb == None:
        lamb = lie_algebra.highest_root()
    else:
        lamb = changeBasis.inverse() * vector(lamb)
        lamb = weyl_group.domain()(list(eval(str(lamb))))

    # if mu is not specified, 0 vector is used
    if mu == None:
        mu = weyl_group.domain()([0 for i in range(0, len(lie_algebra.simple_roots()))])
    else:
        mu = changeBasis.inverse() * vector(mu)
        mu = weyl_group.domain()(list(eval(str(mu))))

    rho = lie_algebra.rho()
    altset = findAltSet(name)
    #print(changeBasis * vector(list(Eval(lamb))))
    for elm in altset:
        # expression in partition function
        res = elm.action(lamb + rho) - (mu + rho)

        #change basis from standard basis to simple roots
        res = vector(list(Eval(res)))
        res = changeBasis * res

        term = calculatePartition(name, list(res), positive_roots, q_analog=q_analog)
        
        term *= (-1)**elm.length()
        mult += term
        
    return mult

#calculateMultiplicity("A3")
calculatePartition("A4", [2,3,3,2],[0,] True)
printPartitions("Lie Algebra Name",)
#findAltSet("D4")

---------------------------------------------------------------------------
TypeError                                 Traceback (most recent call last)
<ipython-input-1-aa27e3862b75> in <module>()
    135 
    136 #calculateMultiplicity("A3")
--> 137 calculatePartition("A4", [Integer(2),Integer(3),Integer(3),Integer(2)], True)
    138 printPartitions("Lie Algebra Name",)
    139 #findAltSet("D4")
<ipython-input-1-aa27e3862b75> in calculatePartition(name, weight, positive_roots, translations, q_analog)
     34             translations["A" + str(i+Integer(1))] = eval("A" + str(i+Integer(1)))
     35 
---> 36     termsForSum = [geometricSumForPartition(list(x), translations, q_analog) for x in positive_roots]
     37     answer = Integer(1)
     38     for x in termsForSum:
TypeError: 'bool' object is not iterable

RootSystem("D4").highest_root()

192