sl(n,C).py

import numpy as np

# Confined to square matrices
n = 5
width = n
height = n
dim = n*n - 1

# Lie Bracket
def Lie_Bracket(A,B):
    return np.dot(A,B)-np.dot(B,A)

def E(i,j):
    Matrix = [[0 for x in range(width)] for y in range(height)]
    Matrix[i-1][j-1] = 1
    return Matrix
# For instance, E(1,2) = [ 0 , 1 , 0]
#                        [ 0 , 0 , 0]   (for n = 3)
#                        [ 0 , 0 , 0]

def H(k): # Part of the standard basis for sl(n,C)
    Matrix = np.add( E(k,k) , np.negative( E(k+1,k+1) ) )
    return Matrix
# For instance, H(2) = [ 0 , 0 ,  0]
#                      [ 0 , 1 ,  0]   (for n = 3)
#                      [ 0 , 0 , -1]

def up_trig_order(l): # order : { 1, .., n(n-1)/2 } --> {(i,j) | n>=i>j} bijectively
    # check if 0 < l < n(n-1)/2
    if  ( (0 >= l) | (l > n*(n-1)) ):
        print "Error from the function 'order."
        return
    stage = 1
    while l > (n - stage):
        l = l - (n - stage)
        stage = stage + 1
    j = l + stage
    return(stage,j)

def basis_order(m):
    if m <= n*(n-1):
        if (m%2) == 1:
            dummy = (m+1)/2
            return up_trig_order(dummy)
        else:
            dummy = m/2
            return (up_trig_order(dummy)[1], up_trig_order(dummy)[0])
    elif m > n*(n-1):
        dummy = m - n*(n-1)
        return ('H', dummy)

#  The standard basis -Beta- for sl(n,C): { E(i,j), H(k) | i!=j, 0<k<n}
# is ordered by the function "order" as follows..
#  Print & Store (in Basis) the standard ordered basis:
Basis = ['I would like to start from Basis[1]']
for m in range(1,n*n): # 0 < m < n*n
    if basis_order(m)[0] == 'H':
        print "The %d-th vector in the standard ordered basis is " %m
        print "Basis[%d] = H(%d) = " % (m,basis_order(m)[1])
        Basis.append( H(basis_order(m)[1]) )
        print( np.array(Basis[m]) )
        print "\n"
    else:
        print "The %d-th vector in the standard ordered basis is " % m
        print "Basis[%d] = E(%d,%d) = " %(m, basis_order(m)[0] ,basis_order(m)[1])
        Basis.append( E(basis_order(m)[0] ,basis_order(m)[1]) )
        print( np.array(Basis[m]) )
        print "\n"

# Print the Lie algebra structure wrt the ordered basis
# Lie_Algebra_Structure = [] # Store the Lie algebra structure wrt the ordered basis...
# for i in range(1,n*n):
#     for j in range(1,n*n):
#         # Store the Lie algebra structure wrt the ordered basis ...
#         # ... Lie_Algebra_Structure[i][j] = np.array( Lie_Bracket( Basis[i] , Basis[j] ) )
#         print "[ Basis[%d] , Basis[%d] ] = " %(i,j)
#         print ( np.array( Lie_Bracket( Basis[i] , Basis[j] ) ) )

# Demonstrate the matrix of ad(X) w.r.t. "Basis" for each X in -Beta-.
# This only works for sl(n,C)
def matrix_ad(X): # Only for sl(n,C) wrt the basis above
    # Test if X is an n by n matrix with zero trace.
    # ...
    result = []
    column_temp = []
    for column in range(1,n*n):
        temp = Lie_Bracket( X , Basis[column] ) # we shall expand "temp" wrt the ordered basis
        for row in range(1,n*n):
            if row <= n*(n-1): # Working against the basis elements E(i,j)
                i = basis_order(row)[0]
                j = basis_order(row)[1]
                insert = temp[i-1][j-1] # Only works for sl(n,C)
                column_temp.append(insert)
            elif row > n*(n-1): # Working against the basis elements H(k)
                i = row - n*(n-1)
                insert = temp[i-1][i-1]
                column_temp.append(insert)
                temp[i-1][i-1] = temp[i-1][i-1] - insert # Actually we don't need this step
                temp[i][i] = temp[i][i] + insert
        result.append(column_temp)
        column_temp = []
    result = np.array(result)
    result = np.transpose(result)
    return result

def Killing(X,Y):
    return np.trace( np.dot(matrix_ad(X),matrix_ad(Y)) )

# Print the adjoint representation wrt to "Basis"
for i in range(1,n*n):
    print "matrix_ad(Basis[%d]) = " %i
    print( np.array(matrix_ad(Basis[i])) )

Killing_Matrix = np.array( [[0]*dim]*dim )
for i in range(1,dim+1):
    for j in range(1,dim+1):
        Killing_Matrix[i-1][j-1] = Killing(Basis[i],Basis[j])

# ---

def int_det(X):
#  This function naively and slowly computes the determinant of any
# given integer valued square matrix X.
#  Here, "int()" is crucial because while integer operations have
# unlimited accuracy in 'pure' python, numpy integer operations
# will overflow! This bug has taken me a day to find..
#  You may check.. >>> type(X[0][0]) ... <type 'numpy.int64'>
    result = int(0)
    if len(X) == 1:
        return int(X[0][0])
    else:
        for i in range(len(X)):
            if X[i][0] != 0:
                result = result + ((-1)**i) * int(X[i][0]) * int_det(cofactor(X,i,0))
                print"%ld" %result
    return result

def cofactor(X,i,j): # X: integer valued matrix
# Output: The (i,j)-th cofactor of X
    result = np.delete(X,i,0) # axis-0: delete rows
    result = np.delete(result,j,1) # axis-1: delete columns
    return result

# # ---
#
# def is_prime(integer):
#     if integer < 2:
#         return False
#     for i in range(2 , int(np.sqrt(integer))+1):
#         if integer % i == 0:
#             return False
#     return True
#
# def prime(k): # Output: the k-th prime
#     i = 1
#     found = 0
#     while found < k:
#         if is_prime(i) == True:
#             found = found + 1
#         i = i + 1
#     return i - 1
#
# # ---

### Input a prime list ###
f = open("prime-list.txt","r")
data = f.read()
f.close
print "len(data) = %d" %len(data)

k = 1
primelist = ['null']
block = []
for i in range(len(data)):
    if data[i] != ",":
        block.append(int(data[i]))
    else:
        primelist.append(0)
        length = len(block)
        for j in range(length):
            primelist[k] = primelist[k] + ( (10**j) * block[(length-1)-j] )
        print "primelist[%d] = %d" %(k,primelist[k])
        k = k+1
        block = []
##########################

def factor(integer,primelist):
# e.g. factor(60) = [2,1,1] ,which means 60 = 2^2 * 3^1 * 5^1
# e.g. factor(100) = [2,0,2] ,which means 100 = 2^2 * 3^0 * 5^2
    result = ['null'] # The 0-th result is 'null' by convention
    i = 1
    while (integer**2 != 1) & (integer != 0): # Then factorization starts
        times = 0
        while integer % primelist[i] == 0:
            integer = integer / primelist[i]
            times = times + 1
        result.append(times)
        if result[i] != 0: ##
            print "The %d-th prime appears %d (= result[%d]) times." %(i,i,result[i]) ##
        i = i + 1
    return result

# ---

print "\n The matrix of Killing Form is "
# print(np.array(Killing_Matrix))
for i in range(len(Killing_Matrix)):
    for j in range(len(Killing_Matrix)):
        print"%d   " %Killing_Matrix[i][j],
    print "\r"

print "and its determinant is %d = " %int_det(Killing_Matrix)
# print str(factor(int_det(Killing_Matrix),primelist))

This file was intended to compute the determinant of the Killing form matrix of the Lie algebra sl(3,C). It has been generalized to sl(n,C). However, determinant overflows when n > 4. [ver. 3]