#matrix stuff. 

import math

num_rows = 4
num_cols = 5

M = []
for i in range(num_rows):  #each row
    M.append([])     # a new row
    for j in range(num_cols):
        M[i].append(i*num_cols+j)   #each item of this row

#print as a list of lists
print(M)

M = []
for i in range(num_rows):
    M.append([i*num_cols+j for j in range(num_cols)]) #list of items of this row
print(M)

M = [[i*num_cols+j for j in range(num_cols)] for i in range(num_rows)]
print(M)

     
#print each row as a list
for row in M:
    print(row)

#print each item of each row
for i in range(len(M)):
    for j in range(len(M[i])):
        print(M[i][j], end=" ")
    print()
#also:
for row in M:
    for col in row:
        print(col, end=" ")
    print()
    
col_x = 2
print("column", col_x, ":", end="")
for row in M:
    print(row[col_x], end=" ")
print()
    
#print a column as a list
print([row[col_x] for row in M])

#return a COLUMN as a list
def column(M,col):
    return [row[col] for row in M]

print(column(M,2))


#zero matrix of size rows*cols
def zero_matrix(rows,cols):
    return [[0 for j in range(cols)] for i in range(rows)]

Z5 = zero_matrix(5,5)
print(Z5)


#return IDENTITY matrix of size n
def identity_matrix(n):
    M = [[0 for j in range(n)] for i in range(n)]   #zero matrix first
    for i in range(n):
        M[i][i] = 1
    return M

I5 = identity_matrix(5)
print(I5)



#return TRANSPOSE of a matrix
def transpose(M):
    Mt = []
    for j in range(len(M[0])):    #each column
        #Mt.append([row[j] for row in M])  #or use the column function above
        Mt.append(column(M,j))
    return Mt

M = [[i*num_cols+j for j in range(num_cols)] for i in range(num_rows)]
Mt = transpose(M)
print(M)
print(Mt)
print(transpose(I5))


#return SUM matrix of two same-sized matrices
def sum_matrix(A,B):
    return [[A[i][j]+B[i][j] for j in range(len(A[0]))] for i in range(len(A))]

M2 = [[i*num_cols+2*j for j in range(num_cols)] for i in range(num_rows)]
print(M2)
print(sum_matrix(M,M2))


#return MINOR submatrix less row i and column j
def minorMij(M,i,j):
    #n-1 sqr matrix, 0 not needed, just a placeholder
    minorIJ =  [[0 for j in range(len(M)-1)] for i in range(len(M)-1)] 
    newrow = newcol = 0

    for oldrow in range(len(M)):
        if oldrow != i:
            for oldcol in range(len(M)):
                if oldcol != j:
                    minorIJ[newrow][newcol] = M[oldrow][oldcol]
                    newcol += 1
            newrow += 1
            newcol = 0    	
    #print(minorIJ)
    return minorIJ

# return DETERMINANT of square matrix
def det(M):
    if len(M) == 2:
        return M[0][0]*M[1][1] - M[0][1]*M[1][0]
    elif len(M) == 3:
        return M[0][0]*M[1][1]*M[2][2] + M[0][1]*M[1][2]*M[2][0] + M[0][2]*M[1][0]*M[2][1] -\
               M[0][2]*M[1][1]*M[2][0] - M[0][1]*M[1][0]*M[2][2] - M[0][0]*M[1][2]*M[2][1] 
    else:
        sum = 0
        for j in range(len(M)):
            sum += M[0][j] * (-1)**(2+j) * det(minorMij(M,0,j))
        return sum

M2 = [[1,2],[3,4]]
print(det(M2))    #-2
print(det(transpose(M2)))    #-2
M3 = [[-2,2,-3],[-1,1,3],[2,0,-1]]
print(det(M3))    #18
print(det(transpose(M3)))    #18
print(det(I5))   #1
M = [[i*num_cols+j for j in range(5)] for i in range(5)]
print(M)
print(det(M))
print(det(transpose(M)))

M44 = [[2,3,4,5],[0,-1,2,2],[3,4,2,4],[2,3,4,1]]
print(det(M44))   #-40
#print([[0 for j in range(len(M)-1)] for i in range(len(M)-1)] )


######################   MATRIX    ########################
print("MATRIX")

#return MATRIX MULTIPLY of two matrices nm * mp = np product
def matrix_multiply(A, B):
    if len(A[0]) != len(B):
        print("ERROR. #columns of first matrix != #rows of second matrix")
    else:
        AB = [[0 for j in range(len(B[0]))] for i in range(len(A))]
        for i in range(len(A)):  #rows of A
            for j in range(len(B[0])):   #columns of B
                sum = 0
                for k in range(len(B)):
                    sum += A[i][k] * B[k][j]
                AB[i][j] = sum
                
                #AB[i][j] = dot_product(A[i],column(B,j))            
        return AB

print(matrix_multiply(M,I5))
print(matrix_multiply(I5,M))

def adjoint(M):
    #n sqr matrix, 0 not needed, just a placeholder
    adj =  [[0 for j in range(len(M))] for i in range(len(M))] 
    
    for i in range(len(M)):
        for j in range(len(M)):
            if len(M) > 2: 	
                adj[i][j] = (-1)**(i+j+2) * det(minorMij(M,i,j))
            #else:
              #  adj[i][j] = (-1)**(i+j+2) * det(M)
    
    if len(M) == 2:   # a b c d  -> d -c -b a
        adj[0][0] = M[1][1]
        adj[0][1] = -M[1][0]
        adj[1][0] = -M[0][1]
        adj[1][1] = M[0][0]
 
    return transpose(adj)

#given adjoint matrix and determinant, return inverse
def inverse (adj, det):
    #n sqr matrix, 0 not needed, just a placeholder
    inv =  [[0 for j in range(len(adj))] for i in range(len(adj))] 
    
    for i in range(len(adj)):
        for j in range(len(adj)): 	
            inv[i][j] = 1/det * adj[i][j]
                    
    return inv

print(inverse(adjoint(I5),det(I5)))
if det(M2) != 0:
    M2inv = inverse(adjoint(M2),det(M2))
    print(matrix_multiply(M2,M2inv))
    print(matrix_multiply(M2inv,M2))
if det(M3) != 0:
    M3inv = inverse(adjoint(M3),det(M3))
    print(matrix_multiply(M3,M3inv))
    print(matrix_multiply(M3inv,M3))

M44 = [[2,3,4,5],[0,-1,2,2],[3,4,2,4],[2,3,4,1]]
print(inverse(adjoint(M44),det(M44)))
M44i = [[-1.25, 1, 1, .25],[.65,-.8,-.4,-.05],[.075,.1,-.2,.225],[.25,0,0,-.25]]
print(inverse(adjoint(M44i),det(M44i)))