import numpy as np
import numpy.linalg as alg

E = np.array([0,1,0,-1,0,0,0,0,0]).reshape((3,3))
a,b,c = 125,31,57
A = np.array([0,a,-b,-a,0,c,b,-c,0]).reshape((3,3))
I = np.eye(3)

# matrices M_n

def M(n,X) :
    return alg.matrix_power(I+(1/n)*X,n)

for k in range(1,5) :
    print(M(k,E),M(k,A))
    
# matrices N_n
    
def N(n,X) :
    B = M(n,X)
    return (B.T).dot(B)
    
for k in range(6) :
    print(N(10**(k+1),A))