View Insert Cell Kemel WidgetsHelp Not TrustedPytho | | Submit | Exercise 2 Let\
ID: 3882090 • Letter: V
Question
View Insert Cell Kemel WidgetsHelp Not TrustedPytho | | Submit | Exercise 2 Let's try implementing and comparing a few different methods for inverting a matrix. In practice, there are many subtleties to implementing an efficient numerically stable inversion routine and it is usually best to make use of those which are provided by the numerical libraries. Even then, for very large matrices with special properties (such as sparsity, etc) there are many specialized routines with better properties. For more information, see the classic text, Numerical Recipes 1. First, let's implement inversion by Cramer's rule as we derived in class. Fill in the following code in the function invcramers below ]: def minor(A, i, j): Return the i,j'th minor of A. n - A.shape[] # list (range()) produces a list [e, 1, # addition of lists concatenates so rows-[e, 1, ,. rows list(range(1))list(range(i+1, n)) cols list(range(list(range(j+1, n)) # this uses 'fancy indexing.. we will discuss # in a later computational Lecture. subMatrix -A[rows, :, cols] return det (subMatrix) ng in more detaft def invCraners (A): Invert the square matrix A using Crer's rule and built in determinant function Cranen's rule: (A1) (ij (1)i+) M(i) / det(A) where M(ji) is the ji'th minor Returns the inverse of A. # assume that the matrix is n × n and get n from the shape n - A.shape[e] HW2-WordExplanation / Answer
Ans: Here is the code for finding the inversion of matrix using cramer's method.
import copy
from fractions import Fraction
def makeIdentity(n):
result = make2dList(n,n)
for i in range(n):
result[i][i] = Fraction(1, 1)
return result
def testMakeIdentity():
print ("Testing makeIdentity..."),
i3 = [[1,0,0],[0,1,0],[0,0,1]]
assert(i3 == makeIdentity(3))
print ("Passed!")
def make2dList(rows, cols):
a=[]
for row in range(rows): a += [[0]*cols]
return a
def multiplyRowOfSquareMatrix(m, row, k):
n = len(m)
rowOperator = makeIdentity(n)
rowOperator[row][row] = k
return multiplyMatrices(rowOperator, m)
def multiplyMatrices(a, b):
# confirm dimensions
aRows = len(a)
aCols = len(a[0])
bRows = len(b)
bCols = len(b[0])
assert(aCols == bRows) # belongs in a contract
rows = aRows
cols = bCols
# create the result matrix c = a*b
c = make2dList(rows, cols)
# now find each value in turn in the result matrix
for row in range(rows):
for col in range(cols):
dotProduct = Fraction(0, 1)
for i in range(aCols):
dotProduct += a[row][i]*b[i][col]
c[row][col] = dotProduct
return c
def addMultipleOfRowOfSquareMatrix(m, sourceRow, k, targetRow):
# add k * sourceRow to targetRow of matrix m
n = len(m)
rowOperator = makeIdentity(n)
rowOperator[targetRow][sourceRow] = k
return multiplyMatrices(rowOperator, m)
def almostEqualMatrices(m1, m2):
# verifies each element in the two matrices are almostEqual to each other
# (and that the two matrices have the same dimensions).
if (len(m1) != len(m2)): return False
if (len(m1[0]) != len(m2[0])): return False
for row in range(len(m1)):
for col in range(len(m1[0])):
if not almostEqual(m1[row][col], m2[row][col]):
return False
return True
def almostEqual(d1, d2):
epsilon = 0.00001
return abs(d1 - d2) < epsilon
def invertMatrix(m):
n = len(m)
assert(len(m) == len(m[0]))
inverse = makeIdentity(n) # this will BECOME the inverse eventually
for col in range(n):
# 1. make the diagonal contain a 1
diagonalRow = col
assert(m[diagonalRow][col] != 0) # @TODO: actually, we could swap rows
# here, or if no other row has a 0 in
# this column, then we have a singular
# (non-invertible) matrix. Let's not
# worry about that for now. :-)
k = Fraction(1,m[diagonalRow][col])
m = multiplyRowOfSquareMatrix(m, diagonalRow, k)
inverse = multiplyRowOfSquareMatrix(inverse, diagonalRow, k)
# 2. use the 1 on the diagonal to make everything else
# in this column a 0
sourceRow = diagonalRow
for targetRow in range(n):
if (sourceRow != targetRow):
k = -m[targetRow][col]
m = addMultipleOfRowOfSquareMatrix(m, sourceRow, k, targetRow)
inverse = addMultipleOfRowOfSquareMatrix(inverse, sourceRow,
k, targetRow)
# that's it!
print(inverse)
return inverse
def testInvertMatrix():
print ("Testing invertMatrix..."),
a = [ [ 1, 2 ], [ 4, 5 ] ]
aInverse = invertMatrix(a)
identity = makeIdentity(len(a))
assert (almostEqualMatrices(identity, multiplyMatrices(a, aInverse)))
a = [ [ 1, 2, 3], [ 2, 5, 7 ], [3, 4, 8 ] ]
aInverse = invertMatrix(a)
identity = makeIdentity(len(a))
assert (almostEqualMatrices(identity, multiplyMatrices(a, aInverse)))
print ("Passed!")
def main():
testMakeIdentity()
testInvertMatrix()
main()
2.Inversion using numpy function:
Compute bit-wise inversion, or bit-wise NOT, element-wise.
numpy.invert(x, /, out=None, *, where=True, casting='same_kind', order='K', dtype=None, subok=True[, signature, extobj]) = <ufunc 'invert'>
x : array_like
Only integer and boolean types are handled.
out : ndarray, None, or tuple of ndarray and None, optional
A location into which the result is stored. If provided, it must have a shape that the inputs broadcast to. If not provided or None, a freshly-allocated array is returned. A tuple (possible only as a keyword argument) must have length equal to the number of outputs.
where : array_like, optional
Values of True indicate to calculate the ufunc at that position, values of False indicate to leave the value in the output alone.
**kwargs
For other keyword-only arguments
Related Questions
Navigate
Integrity-first tutoring: explanations and feedback only — we do not complete graded work. Learn more.