From 32f4af5f369fa9f0b2988ecad7797f4bec3661c3 Mon Sep 17 00:00:00 2001 From: Holden Rohrer Date: Tue, 21 Sep 2021 17:12:46 -0400 Subject: notes and homework --- li/matlab_hw.py | 183 ++++++++++++++++++++++++++++++++++++++++++++++++++++++++ 1 file changed, 183 insertions(+) create mode 100644 li/matlab_hw.py (limited to 'li/matlab_hw.py') diff --git a/li/matlab_hw.py b/li/matlab_hw.py new file mode 100644 index 0000000..71b89b6 --- /dev/null +++ b/li/matlab_hw.py @@ -0,0 +1,183 @@ +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") -- cgit