import numpy as np

from matrix import Matrix
from vector import Vector

np.random.seed(0)
###############
# Testing the vector class
print("\n\nTesting the vector class -----------------------------------\n\n")
###############

### 1a Initialization
print("Start 1a initialization")
# from list
x1 = [1, 2, 3, 4, 5]
vec_x1 = Vector(x1)
# from numpy array
x2 = np.random.uniform(0, 1, 10000)
vec_x2 = Vector(x2)
print("End 1a\n")

### 1b __str__ function, string representation
print("Start 1b string representation")
print(f"vec_x1 values: {str(vec_x1)}")
print(f"vec_x2 values: {str(vec_x2)}")
print("End 1b\n")

### 1c shape and transpose
print("Start 1c shape and transpose")
# shape property of the vector
print(f"vec_x1 has shape {vec_x1.shape()} | must be (5,1)")  # Hint: numpy.reshape
# transposition (property) of the vector, only "cosmetic" change
print(f"vec_x2.T() has shape {vec_x2.T().shape()} | must be (1,10000)")
print(f"vec_x2.T().T() has shape {vec_x2.T().T().shape()} | must be (10000,1)")
print("End 1c\n")

### 1d addition and subtraction
print("Start 1d addition and subtraction")
# initialization
a = Vector([1, 3, 5, 7, 9])
b = Vector([-2, 5, 1, 0, -3])
# computation
x = a + b
print(f"a + b = {str(x)} | must be {np.array([1, 3, 5, 7, 9]) + np.array([-2, 5, 1, 0, -3])}")
y = a - b
print(f"a - b = {str(y)} | must be {np.array([1, 3, 5, 7, 9]) - np.array([-2, 5, 1, 0, -3])}")
try:
    x_false = a + b.T()
except ValueError as e:
    print("The correct result is this error: " + str(e))
else:
    print("It is required to raise a system error, e. g., ValueError, since dimensions mismatch!")
try:
    x_false = a - b.T()
except ValueError as e:
    print("The correct result is this error: " + str(e))
else:
    print("It is required to raise a system error, e. g., ValueError, since dimensions mismatch!")
print("End 1d\n")

### 1e multiplication
print("Start 1e multiplication")
# initialization
a = Vector([1, 3, 5, 7, 9])
b = Vector([-2, 5, 1, 0, -3])
# computation with vectors
x = a.T() * b  # scalar
print(f"a.T() * b = {str(x)} | must be {np.sum(np.array([1, 3, 5, 7, 9]) * np.array([-2, 5, 1, 0, -3]))}")
y = a * b  # vector
print(f"a * b = {str(y)} | must be {np.array([1, 3, 5, 7, 9]) * np.array([-2, 5, 1, 0, -3])}")
z = b * a.T()  # matrix
print(f"b * a.T() = \n{str(z)} | must be \n{np.outer(np.array([-2, 5, 1, 0, -3]), np.array([1, 3, 5, 7, 9]))}")
# computation with scalars
x = a * 5
print(f"a * 5 = {str(x)} | must be {np.array([1, 3, 5, 7, 9]) * 5}")
y = 0.1 * b.T()
print(f"0.1 * b.T()= {str(y)} | must be {0.1 * np.array([-2, 5, 1, 0, -3])}")
print("End 1e\n")

### 1f division
print("Start 1f division")
# initialization
a = Vector([1, 3, 5, 7, 9])
b = Vector([-2, 5, 1, 7, -3])
# computation with vectors
x = a / b
print(f"a / b = str{x} | must be {np.array([1, 3, 5, 7, 9]) / np.array([-2, 5, 1, 7, -3])}")
y = a / 5
print(f"a / 5 = str{y} | must be {np.array([1, 3, 5, 7, 9]) / 5} ")
print("End 1f\n")

### 1g norm
print("Start 1g norm")
# initialization
a = Vector([1, 3, 5, 7, 9])
a_norm = a.norm()
a_normalized = a.normalize()
print(f"a_norm = {a_norm} | must be {np.linalg.norm([1, 3, 5, 7, 9])}")
print(f"a_normalize = {str(a_normalized)} | must be {np.array([1, 3, 5, 7, 9]) / np.linalg.norm([1, 3, 5, 7, 9])}")
print("End 1g\n")

### 1h negation
print("Start 1h negation")
# initialization
a = Vector([1, 3, 5, 7, 9])
x = -a
print(f"-a = {str(x)} | must be {-np.array([1, 3, 5, 7, 9])}")
print("End 1h\n")

### 1i manipulation
print("Start 1i manipulation")
# initialization
a = Vector([1, 3, 5, 7, 9])
print(
    f"a[{str([1, 2, 4])}] = {str(np.array(a[1, 2, 4]).reshape(3, ))} | must be {np.array([1, 3, 5, 7, 9]).reshape(5, 1)[np.array([1, 2, 4])].reshape(3, )}")
a[1, 2, 4] = [-1, -1, -1]
print(f"a = {str(a)} | must be {np.array([1, -1, -1, 5, 7, -1])}")
print("End 1i\n")

###############
# Testing the matrix class
print("\n\nTesting the matrix class -----------------------------------\n\n")
###############

### 1a Initialization
print("Start 2a initialization")
a_list = np.array([[1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11], 2 * np.array([1, 2, 3, 4, 5, 6, 7, 8, 9, 10, 11])])
A = Matrix(a_list)
B = Matrix([-1, 2, -1], structure="tridiagonal", n=11)
c_list = [[(i + 1) / (index + 1) for index in range(10)] for i in range(10)]
C = Matrix(c_list, shape=(10, 10))
print("End 2a\n")

### 1b __str__ function, string representation
print("Start 2b string representation")
# print(B.__str__(full = True))
print(str(A))
print(str(B))
print(str(C))
print("End 2b\n")

### 1c shape and transpose
print("Start 2c shape and transpose")
# Initialization
A = Matrix(np.array([i for i in range(12)]).reshape(-1, 1), shape=(4, 3))
print(f"A has shape {A.shape()} | must be (4,3)")
print(f"A.T() has shape {A.T().shape()} | must be (3,4)")
print(f"A.T().T() has shape {A.T().T().shape()} | must be (4,3)")
print("End 2c\n")

### 1d addition and subtraction
print("Start 2d addition and subtraction")
# Initialization 
A = Matrix(structure="diagonal", data=[3], offset=0, n=10)
print(str(A))
A21 = Matrix(structure="diagonal", data=[-1], offset=-1, n=10)
print(str(A21))
A12 = Matrix(structure="diagonal", data=[-1], offset=+1, n=10)
print(str(A12))
B = Matrix(structure="diagonal", data=[1], offset=0, n=10)
print(str(B))
# computation
C = A + A21 + A12 - B
print(str(C) + f"must be\n{Matrix(structure='tridiagonal', data=[-1, 2, -1], n=10)}")
print(str(5 + A - 3))
print("End 2d\n")

### 1e multiplication 
print("Start 2e multiplication")
# initialization
a_mat = [[(i0 + 1) / (i1 + 1) for i1 in range(3)] for i0 in range(10)]
b_mat = [[(i0 + 1) / (i1 + 1) for i1 in range(10)] for i0 in range(3)]
A = Matrix(a_mat)
B = Matrix(b_mat)
# computation
# print(f"A * B =\n{str(A*B)}must be\n{str(np.round(np.array(a_mat) @ np.array(b_mat),decimals=3))}")
print(f"Norm of (A*B - np.(A*B)) is {(A * B - Matrix(np.array(a_mat) @ np.array(b_mat))).norm()} | must be < 1e-8")
# print(f"A.T() * A =\n{str(A.T()*A)}must be\n{str(np.round(np.array(a_mat).T@np.array(a_mat),decimals=3))}")
print(
    f"Norm of (A.T()*A - np(A.T()*A)) is {(A.T() * A - Matrix(np.array(a_mat).T @ np.array(a_mat))).norm()} | must be < 1e-8")
# print(f"A * A.T() =\n{str(A*A.T())}must be\n{str(np.round(np.array(a_mat)@np.array(a_mat).T,decimals=3))}")
print(
    f"Norm of (A*A.T() - np(A*A.T())) is {(A * A.T() - Matrix(np.array(a_mat) @ np.array(a_mat).T)).norm()} | must be < 1e-8")
# print(f"B.T() * A.T() =\n{str(B.T()*A.T())}must be\n{str(np.round(np.array(b_mat).T @ np.array(a_mat).T,decimals=3))}")
print(
    f"Norm of (B.T()*A.T() - np.(B.T()*A.T())) is {(B.T() * A.T() - Matrix(np.array(b_mat).T @ np.array(a_mat).T)).norm()} | must be < 1e-8")
print("End 2e\n")

### 1f division 
print("Start 2f division")
# initialization
A = Matrix(a_mat)
# computation
print(f"Norm of (A/5 - np.(A/5)) is {(A / 5 - Matrix(np.array(a_mat) / 5)).norm()} | must be < 1e-8")
print("End 2f\n")

### 1g norm
print("Start 2g norm")
A = Matrix(structure="tridiagonal", n=50, data=[-1, 2, -1])
print(f"Frobenius norm of tridiagonal matrix: {A.norm('frobenius')} | must be 17.263")
print(f"Row sum norm of tridiagonal matrix: {A.norm('row sum')} | must be 4")
print(f"Col sum norm of tridiagonal matrix: {A.norm('col sum')} | must be 4")
print("End 2g\n")

### 1h negation
print("Start 2h negation")
A = Matrix(structure="tridiagonal", n=50, data=[-1, 2, 1])
print(f"Norm of (A + (-A)) is {(A + (-A)).norm('frobenius')} | must be < 1e-8")
print("End 2h\n")

### 1i manipulation
print("Start 2i manipulation")
A = Matrix(structure="tridiagonal", n=10, data=[-1, 2, 1])
A[1, 1] = 4
A[[1, 2, 3], 2] = [-5, -10, 100]
print(str(A))
print("End 2i\n")