diff --git a/GaussJordan.py b/GaussJordan.py new file mode 100644 index 0000000..a4e369d --- /dev/null +++ b/GaussJordan.py @@ -0,0 +1,43 @@ +import numpy as np +# calculate the inverse of the matrix A using Gauss-Jordan algorithm +def GaussJordan (A): + (n, n) = A.shape + B = np.eye(n) + + # Assume A has an inverse matrix + success = 1 + + if np.linalg.det(A) == 0: + # A cannot be reversed + B = np.empty(n) + success = 0 + return + + # merge the two matrices + Ae = np.concatenate((A, B.T), axis=1) + for i in range (0, n): + if Ae[i, i] == 0: + B = np.empty(n) + success = 0 + return + # make the pivot position to one (as in the eye matrix) + Ae[i, :] = Ae[i, :] / Ae[i, i] + + #form zeros above and under the main diagonal in the + # first half and calculate the inverse in the second half + for j in range(0, n): + if i != j: + Ae[j, :] = Ae[j, :] - Ae[i, :] * Ae[j, i] + + #extract the inverse of matrix A from the augmented matrix + B = Ae[:, n : 2 * n ] + return (B, success) + + +# example for a 4x4 matrix +A = np.random.rand(4,4) +(inv, success) = GaussJordan(A) +print "The inverse using Gauss-Jordan is:" +print inv +print "The inverse using inv from linalg is:" +print np.linalg.inv(A) \ No newline at end of file diff --git a/python/Gershgorin.py b/python/Gershgorin.py new file mode 100644 index 0000000..3eb5bb8 --- /dev/null +++ b/python/Gershgorin.py @@ -0,0 +1,41 @@ +import numpy as np +import matplotlib.pyplot as plt + +def Gershgorin (A): + # A must be a square matrix + (m,n) = A.shape + if m != n: + print "You must introduce a square matrix" + return + + plt.axes() + # For each row: + for i in range (0, n): + # the circle has the center in (h, k) where + # h is the real part of A(i,i) and k is the real part of A(i,i) + h = np.real(A[i,i]) + k = np.imag(A[i,i]) + plt.plot(h, k, marker='x', markersize=5, color="blue") + + # the radius of the circle is the sum of + # norm of the elements in the row where i != j + + r = 0 + for j in range (0,n): + if i != j: + r = r + np.linalg.norm(A[i,j]) + # plot the circle + circle = plt.Circle((h, k), r, fill = False) + plt.gca().add_patch(circle) + + eigenval = np.linalg.eigvals(A) + # plot the eigenvalues of the matrix + for x in eigenval: + plt.plot(np.real(x), np.imag(x), marker='o', markersize=5, color="red") + + plt.axis('scaled') + plt.show() +# example for a 3x3 matrix +A = np.random.rand(3,3) +print A +Gershgorin(A) \ No newline at end of file