import numpy as np
import scipy.linalg

print("Section 1.3: #30")

coeff = (
    [[ 1, 1, 1 ],
     [ 1, 2, 2 ],
     [ 2, 3, -4]]
    )
print(
np.linalg.solve(
    coeff,
    [6, 11, 3])
)

print(
np.linalg.solve(
    coeff,
    [7, 10, 3])
)

print("Section 1.3: #32")

dimension = 3
tot = [0]*3
ct = 10000
for x in range(ct):
    P, L, U = scipy.linalg.lu(np.random.rand(dimension, dimension))
    for i in range(dimension):
        tot[i] += U[i][i]
print([x/ct for x in tot]);

print("Section 1.4: #21")

A = [[.5, .5], [.5, .5]]
B = [[1, 0], [0, -1]]
C = np.matmul(A, B)
print("A^2")
print(np.linalg.matrix_power(A, 2))
print("A^3")
print(np.linalg.matrix_power(A, 3))
print("A^k = A");
print("B^2")
print(np.linalg.matrix_power(B, 2))
print("B^3")
print(np.linalg.matrix_power(B, 3))
print("B^{2k} = B^2. B^{2k+1} = B.")
print("C^2")
print(np.linalg.matrix_power(C, 2))
print("C^3")
print(np.linalg.matrix_power(C, 3))
print("C^k = 0")

print("Section 1.4: #59")

print("A * v = [3, 4, 5]")
print("v' * v = 50")
print(f"v * A = {np.matmul([3,4,5], np.identity(3))}")

print("Section 1.4: #60")

print("Av = [4. 4. 4. 4.]")
print("Bw = [10. 10. 10. 10.]")

print("Section 1.6: #12")

print("Section 1.6: #32")

print("the inverse of the first matrix is")
print(np.linalg.inv(5*np.identity(4) - np.ones((4,4))))
print("(meaning a = .4, b = .2)")
print("ones(4)*ones(4) =")
print(np.matmul(np.ones((4,4)),np.ones((4,4))))
print("so we can solve this like the following equation: (a*eye +\n\
b*ones)(c*eye + d*ones) = ac*eye + (ad+bc+dimension*bd)*ones = eye")
print("c = 1/a. d = -b/a/(a+dimension*b), giving for the next matrix")
print("a = 1/6, d = 6/(1/6+5*-1)") # WRONG

print("Section 1.6: #47")

print("I'm not sure. This is Python. It says (in either case):")
print("numpy.linalg.LinAlgError: Singular matrix")
'''
print(np.linalg.solve(np.ones((4, 4)), np.random.rand(4, 1)))
print(np.linalg.solve(np.ones((4, 4)), np.ones((4, 1))))
'''

''' #68
dimension = 500
northwest = np.random.rand(dimension, dimension)
for x in range(0, dimension):
    for y in range(x+1, dimension):
        northwest[y][x] = 0
print(northwest)
'''
print("Section 1.6: #69")


from time import perf_counter
matrix = np.random.rand(500, 500)
start = perf_counter()
np.linalg.inv(matrix)
first_time = perf_counter() - start
print("500x500 inverse:", str(first_time) + "s")
matrix = np.random.rand(1000, 1000)
start = perf_counter()
np.linalg.inv(matrix)
second_time = perf_counter() - start
print("1000x1000 inverse:", str(second_time) + "s")
print("factor:", second_time / first_time)

print("Section 1.6: #70")

dimension = 1000
I = np.identity(dimension)
A = np.random.rand(dimension, dimension)
U = np.random.rand(dimension, dimension)
for x in range(0, dimension):
    for y in range(x+1, dimension):
        U[y][x] = 0

start = perf_counter()
np.linalg.inv(U)
print("inverse of U:", str(perf_counter() - start) + "s")
start = perf_counter()
np.linalg.solve(U, I)
print("U\I:", str(perf_counter() - start) + "s")
start = perf_counter()
np.linalg.inv(A)
print("inverse of A:", str(perf_counter() - start) + "s")
start = perf_counter()
np.linalg.solve(A, I)
print("A\I:", str(perf_counter() - start) + "s")

#print("Section 1.6: #71")

print("Section 2.2: #33")

print("For the first matrix:")
A = [[1, 3, 3], [2, 6, 9], [-1, -3, 3]]
b = [[1, 5, 5]]
# TODO
print(scipy.linalg.null_space(A))
print("For the second matrix:")
A = [[1, 3, 1, 2], [2, 6, 4, 8], [0, 0, 2, 4]]
b = [[1, 3, 1]]

print("Section 2.2: #35")

print("For the first system:")
A = [[1, 2], [2, 4], [2, 5], [3, 9]]
print(sympy.Matrix.rref(A)) # ew is there no automated way to do this ?

print("Section 2.2: #36")

print("(a)")
A = np.array([[1, 2, 1], [2, 6, 3], [0, 2, 5]])
print("Basis of the column space:")
print(scipy.linalg.orth(A))
print("Multiplying the rows by this gives zero (for this matrix, \
equivalent to [0 0 0]^T)")
print(scipy.linalg.null_space(A.transpose()))

print("(b)")
A = np.array([[1, 1, 1], [1, 2, 4], [2, 4, 8]])
print("Basis of the column space:")
print(scipy.linalg.orth(A))
print("Multiplying the rows by this gives zero (for this matrix, \
equivalent to [0 2 -1]^T)")
print(scipy.linalg.null_space(A.transpose()))

print("Section 2.3: #2")

print("Section 2.3: #5")

print("Section 2.3: #13")

print("Section 2.3: #16")

print("Section 2.3: #18")

print("Section 2.3: #24")