Friday, September 23, 2011

Is there a computer alegebra package(mathematica, etc.) that could take the inverse of a general matrix?

This is, if you have a matrix without explicit dimensions, labeled only as mn, and have rule for the value of any element ij, would it be possible to get a rule for the value of any element in the inverse matrix?|||I'll assume that m = n here so that the matrix is square (for more general matrices, there is some choice involved in defining an inverse).





So suppose A is n-by-n and you want the row i, column j entry of its inverse. Ask your computer to do the following:





(1) Compute the determinant of A, call that d.





(2) Replace the ith column of A with the n-by-1 column vector that has a 1 in the jth entry and zeros everywhere else.





(3) Evaluate the determinant of this new matrix, and divide it by d. That's what you want.





[This works because the jth column of the inverse of A, is the vector v that solves the equation Av = e_j (where e_j denotes the vector with a 1 in the jth entry and zeros everywhere else). Cramer's Rule then tells you that the ith entry of v is given by the above procedure.]





Mathematica can certainly do this, although I am not familiar with its syntax. It should not be too different from how I can do it in Maple:





----------


with(linalg);


minv := proc (n,i,j)


local A,B,d,k;


A:=matrix(n,n);


d:=det(A);


B:=array(identity,1..n,1..n);


for k from 1 to n do A[k,i]:=B[k,j]; od;


det(A)/d;


end;


----------





After I load this into Maple, executing minv(n,i,j); returns the row i, column j entry of the inverse of a symbolic matrix A. Maple writes the entries of A as A[1,1], A[1,2], et cetera, so for example if I execute minv(2,1,2); I get





A[1, 2]


- ---------------------------------


A[1, 1] A[2, 2] - A[1, 2] A[2, 1]





Which a human mind find more familiar if we let A[1,1] = a, A[1,2] = b, A[2,1] = c, and A[2,2] = d; so that A is the matrix with entries a,b,c,d. The program then tells us that the row 1, column 2 entry the inverse is -b/(ad-bc).





For large n these formulas involve so many terms that they have almost no value in solving problems. You can certainly get Maple to split out the formulas for the entries of a general 5 by 5 matrix, but once you see the output you'll understand what I mean.





A numerical analyst could also tell you that the formulas are poorly suited to numerical calculation (one reason is the way roundoff errors can build up if a machine uses these formulas to calculate with finite precision numbers).





So it's almost never the case that one inverts an n-by-n matrix by taking the formula for the "general" n-by-n matrix inverse and plugging numbers into it.





Nevertheless, I think these formulas are a lot of fun. Enjoy them!





Relevant Wikipedia stuff:





Cramer's Rule: http://en.wikipedia.org/wiki/Cramer%27s_鈥?/a>





Determinants: http://en.wikipedia.org/wiki/Determinant





Numerical analysis: http://en.wikipedia.org/wiki/Numerical_a鈥?/a>

No comments:

Post a Comment