Decomposicao de matriz em LU com pivoteamento

.py

def matrixPrinter(A):
    print('\n'.join([''.join(['{:4}'.format(item) for item in row])
      for row in A]))


def LUPivo(matrix):

    n = len(matrix)

    M = [[0 for x in range(n)]
                for y in range(n)];
    L = [[0 for x in range(n)]
                for y in range(n)];
    U = [[0 for x in range(n)]
                for y in range(n)];
    matrixPivo = [[0 for x in range(n)]
                for y in range(n)];

    # pivoteamento
    for i in range(n - 1):

        pivo = matrix[i][i]
        pivoLine = i

        for j in range(n):
            if (matrix[i][j] > pivo):
                pivo = matrix[j][i]
                pivoLine = j

        matrix[i], matrix[pivoLine] = matrix[pivoLine], matrix[i]

        for j in range(i, n-1):
            M[j+1][i] = matrix[j+1][i]/pivo
            for k in range(n):
                U[j+1][k] = matrix[j+1][k] - (M[j+1][i]*matrix[i][k])
                matrix[j+1][k] = U[j+1][k]
        L = M

    for i in range(n):
        for j in range(n):
            if i == j:
                L[i][j] = 1
            if j > i:
                L[i][j] = 0

    U[0] = matrixPivo[0]

    return L, U


matriz = [[1, 2, 3],
          [2, 3, 1],
          [4, 0, 2]]

A = LUPivo(matriz)



n = len(matriz)

print("Lower           \t\tUpper");

for i in range(n):

    # Lower
    for j in range(n):
        print(round(A[0][i][j], 2), end = "\t")
    print("", end = "\t")

    # Upper
    for j in range(n):
        print(round(A[1][i][j], 2), end = "\t")
    print("")