How Can We Automate Matrix Inversion for Matrices Near the Identity?

Click For Summary

Discussion Overview

The discussion revolves around algorithms for automating matrix inversion, particularly for matrices that are near the identity matrix. Participants explore various methods, their efficiencies, and practical applications, including numerical routines and programming languages.

Discussion Character

  • Exploratory
  • Technical explanation
  • Debate/contested
  • Mathematical reasoning

Main Points Raised

  • Some participants inquire about algorithms for matrix inversion suitable for automation, expressing dissatisfaction with traditional methods described in linear algebra texts.
  • One participant suggests Gaussian elimination as a method that works for all matrices but notes it is the slowest option available.
  • Another participant mentions the importance of matrix properties (e.g., sparsity, symmetry) and suggests that certain methods, like the conjugate gradient, are better suited for specific types of matrices, even if they do not compute the inverse explicitly.
  • A participant shares a Fortran code snippet for matrix inversion, emphasizing the need for translation to modern programming languages and noting that calculating the inverse is often unnecessary.
  • Some participants express mixed feelings about the use of Fortran, with one acknowledging its efficiency despite its age, while another suggests looking into Numerical Recipes for additional resources.
  • One participant discusses a practical application of matrix inversion in glaze manipulation for pottery, illustrating how matrix equations can be used to derive ingredient proportions based on desired glaze properties.
  • A later reply introduces the idea of using geometric series for matrices of the form I-f, where f is small, suggesting an alternative approach for matrices near the identity.

Areas of Agreement / Disagreement

Participants express differing opinions on the best methods for matrix inversion, with no consensus reached on a single preferred approach. Various methods and their applicability are debated, indicating multiple competing views.

Contextual Notes

Some limitations are noted regarding the assumptions about matrix properties and the context in which different algorithms are applied. The discussion does not resolve the mathematical steps involved in the proposed methods.

Who May Find This Useful

This discussion may be of interest to those involved in numerical analysis, programming, and applications of linear algebra in practical scenarios, such as engineering or materials science.

radagast
Messages
483
Reaction score
1
Anybody know of a link to a page that describes an algorithm for matrix inversion. My old linear algebra book describes a 'by hand' method, but it's unsuitable for automating.
 
Physics news on Phys.org
Gaussian elimination isn't that bad, is it? :smile:


Important concerns include the size of your matrix, any qualitative properties it may have (such as sparse, symmetric, or banded), and what you want to do with the inverse.

For example, http://www.icos.ethz.ch/teaching/ml02/painless-conjugate-gradient.pdf is well suited for solving the equation Ax=b when A is large sparse matrix, but it won't explicitly compute A inverse.
 
Last edited by a moderator:
Yeah, Gaussian elimination is the only method that's gauranteed to work on all matrices, but it is the slowest.

- Warren
 
Rather then do a web search, I searched my book shelf. These routines are in Fortran (thats how old my book shelf is!) Not sure if this is the best format, but you will have to translate to a modern language anyway. I believe that the syntax is simple enough that you should (with fundamental programing knowledge) be able to figure out what is happening. I think I have edited out all of the OCR glitches, my software wanted to read an = as a z, so if there is a weird z that makes no sense floating around replace it with an =.

These were lifted from Elementary Numerical Analysis by Conte and de Boor.
I will do my best to answer any questions.

BTW
In the text preceding the inversion routine the authors say

..., we hasten to point out that there is usually no good reason for ever calculating the inverse.

SUBROUTINE SUBST ( W, IPIVOT, B, N, X )
INTEGER IPIVOT(N) , I,IP,J
REAL B(N) ,W(N,N) ,X(N), SUM
C****** I N PUT ******
C W, IPIVOT, N ARE AS ON OUTPUT FROM F ACT 0 R , APPLIED TO THE
C MATRIX A OF ORDER N.
C B IS AN N-VECTOR, GIVING THE RIGHT SIDE OF THE SYSTEM TO BE SOLVED.
C****** 0 U T PUT ******
C X IS THE N-VECTOR SATISFYING A*X . B.
C****** MET HOD ******
C ALGORITHM 4.4 IS USED, I.E., THE FACTORIZATION OF A CONTAINED IN
C W AND IPIVOT (AS GENERATED IN FACTOR) IS USED TO SOLVE A*X = B
C FOR X BY SOLVING TWO TRIANGULAR SYSTEMS.
C
IF (N .LE. 1) THEN
X(1) = B(1)/W(1,1)
RETURN
END IF
IP = IPIVOT(1)
X(1) = B(IP)
DO 15 I=2,N
SUM = 0.
DO 14 J=I,I-l
14 SUM = W(I,J)*X(J) + SUM
IP z IPIVOT(I)
15 X(I) = B(IP) - SUM
C
X(N) = X(N)/W(N,N)
DO 20 I=N-l,I,-1
SUM = 9.
DO 19 J=I+l,N
19 SUM = W(I,J)*X(J) + SUM
20 X(I) = (X(I) - SUM)/W(I,I)
RETURN
END

SUBROUTINE FACTOR ( W, N, D, IPIVOT, IFLAG )
INTEGER IFLAG,IPIVOT(N) , I,ISTAR,J,K
REAL D(N) ,W(N,N), AWIKOD,COLMAX,RATIO,ROWMAX,TEMP
**** INPUT ******
C W ARRAY OF SIZE (N,N) CONTAINING THE MATRIX A OF ORDER N TO BE
C FACTORED.
C N THE ORDER OF THE MATRIX
C***** WORK AREA ******
C D A REAL VECTOR OF LENGTH N, TO HOLD ROW SIZES
C***** 0UTPUT ******
C W ARRAY OF SIZE (N,N) CONTAINING THE LU FACTORIZATION OF P*A FOR
C SOME PERMUTATION MATRIX P SPECIFIED BY IPIVOT.
C IPIVOT INTEGER VECTOR OF LENGTH N INDICATING THAT ROW IPIVOT(K)
C WAS USED TO ELIMINATE X(K) , K-l,...,N .
C IFLAG AN INTEGER,
C = 1, IF AN EVEN NUMBER OF INTERCHANGES WAS CARRIED OUT,
C = -1,IF AN ODD NUMBER OF INTERCHANGES WAS CARRIED OUT,
C = 0, IF THE UPPER TRIANGULAR FACTOR HAS ONE OR MORE ZERO DIA-
C GONAL ENTRIES.
C THUS, DETERMINANT (A) - IFLAG*W(l,l)*...*W(N,N) .
C IF IFLAG .NE. 0, THEN THE LINEAR SYSTEM A*X 8 B CAN BE SOLVED FOR
C X BY A
C CALL SUBST (W, IPIVOT, B, N, X )
C**** METHOD ******
C THE PROGRAM FOLLOWS ALGORITHM 4.2, USING SCALED PARTIAL PIVOTING.
C
IFLAG = 1
C INITIALIZE IPIVOT, D
DO 9 I-l,N
IPIVOT(I) = I
ROWMAX = 0.
DO 5 J=l,N
5 ROWMAX = AMAXl(ROWMAX,ABS(W(I,J)))
IF (ROWMAX .Eq. 0.) THEN
IFLAG = 0
ROWMAX = 1.
END IF
9 D(I) = ROWMAX
IF (N .LE. 1) RETURN
C FACTORIZATION
DO 20 K=l,N-l
C DETERMINE PIVOT ROW, THE ROW ISTAR.
COLMAX = ABS(W(K,K)/D(K)
ISTAR = K
DO 13 I=K+l,N
AWIKOD = ABS(W(I,K))/D(I)
IF (AWIKOD .GT. COLMAX) THEN
COLMAX = AWIKOD
ISTAR = I
END IF
13 CONTINUE
IF (COLMAX .EO. 0.) THEN
IFLAG = 0
ELSE
IF (ISTAR .GT. K) THEN
C MAKE K THE PIVOT ROW BY INTERCHANGING IT WITH
C THE CHOSEN ROW ISTAR.
IFLAG = -IFLAG
I = IPIVOT(ISTAR)
IPIVOT(ISTAR) = IPIVOT(K)
IPIVOT(K) = I
TEMP = D(ISTAR)
D(ISTAR) = D(K)
D(K) = TEMP
DO 15 J=l,N
TEMP = W(ISTAR,J)
W(ISTAR,J) = W(K,J)
15 W(K,J) = TEMP
END IF
C ELIMINATE X(K) FROM ROWS K+l,...,N.
l6 DO 19 I=K+l,N
W(I,K) = W(I,K)/W(K,K)
RATIO = W(I,K)
DO 19 J=K=1,N
W(I,J) = W(I,J) - RATIO*W(K,J)
19 CONTINUE
END IF
20 CONTINUE
IF (W(N,N) .EQ. 0.) IFLAG -0
RETURN
END

C PROGRAM FOR CALCULATING THE INVERSE OF A GIVEN MATRIX
C CALLS FACT0R, SUBST.
PARAMETER NMAX=30,NMAXSQ=NMAX*NMAX
INTEGER I,IBEG,IFLAG,IPIVOT(NMAX) ,J,N,NSQ
REAL A(NMAXSQ) ,AINV(NMAXSQ) ,B(NMAX)
1 READ 501, N
501 FORMAT(I2)
IF (N .LT. 1 .OR. N .GT. NMAX) STOP
C READ IN MATRIX ROW BY ROW
NSQ = N*N
DO 10 I=I,N
18 READ 510, (A(J) ,J=I,NSQ,N)
510 FORMAT(5EI5.7)
C
CALL FACTOR ( A, N, B, IPIVOT, IFLAG )
IF (IFLAG .EQ. 0) THEN
PRINT 611
611 FORMAT('IMATRIX IS SINGULAR')
GO TO 1
END IF
DO 21 I=l,N
21 B(I) = 0.
IBEG = 1
DO 30 J=I,N
B(J) = 1.
CALL SUBST ( A, IPIVOT, B, N, AINV(IBEG) )
B(J) = 0.
38 IBEG = IBEG + N
PRINT 630
630 FORMAT('ITHE COMPUTED INVERSE IS '//)
DO 31 I=l,N
31 PRINT 631, I, (AINV(J) ,J-I,NSQ,N)
631 FORMAT('0ROW ',I2,8EI5.7/(7X,8EI5.7))
GO TO 1
END
 
Last edited:
Well, Fortran is pretty lame overall.

Numerical Recipes also has a lot of good material on the topic, in both fortran and C++. And, as quite a nice gift to the scientific community, the books are available in their entirety online:

http://ww.nr.com

- Warren
 
Last edited by a moderator:
Warren,
The routines deserve a look inspite of the language, the above text generally presented very nice efficient routines, even if the language is archaic. :)

BTW: you may want to edit your link www usually works better then ww!
 
Thanks to all of you for the info. Do worry about it being in FORTRAN, unfortunately I'm old enough to be able to read it. :smile:
I was after the algorithm, I can always code it up into Java.

I throw and fire pottery for a hobby, and I've always had it in the back of my mind to write a glaze manipulation program. All the important properties of the ingredients are linearly translated into the glaze, so by using AX = B where A is a vector of ingredient percentages, X is the component/properties matrix, and B is the glaze component/properties vector, I can make changes in glaze properties (B), derive X-1, then come up with an new set of proportions corrosponding to the glaze I want.

A mundane application, but it's a real pain when the thermal coefficient of expansion for the glaze and clay body don't match.
 
Originally posted by chroot
Numerical Recipes also has a lot of good material on the topic, in both fortran and C++. And, as quite a nice gift to the scientific community, the books are available in their entirety online:

http://www.nr.com


great site...
I'we wanted a good book on this ever since I studied this course in college (the profeesor was such a looser though I couldn't understand much from him).
My eternal gratitude...:smile:
PS: edited the above link...
 
for matrices near the identity matrix, i.e. of form I-f for small f, what about the geometric series?

aha! i got the center of the earth, or the beginning of time. the goal of all paleontologists!
 

Similar threads

Undergrad How can I convince myself that I can find the inverse of this matrix?
  • · Replies 34 ·
2
Replies
34
Views
4K
Graduate Particle swarm optimization for matrix inversion
  • · Replies 1 ·
Replies
1
Views
2K
Undergrad Getting eigenvalues of an arbitrary matrix with programming
  • · Replies 1 ·
Replies
1
Views
2K
Undergrad What Does the F Matrix Look Like for a Linear Bezier Curve?
  • · Replies 1 ·
Replies
1
Views
2K
High School What unique contributions to math did linear algebra make?
  • · Replies 19 ·
Replies
19
Views
4K
Undergrad Can we retrieve the inverse of matrix A in this example?
  • · Replies 14 ·
Replies
14
Views
2K
Undergrad Want to understand how to express the derivative as a matrix
  • · Replies 8 ·
Replies
8
Views
3K
Undergrad What is the connection between cofactors, determinants, and matrix inverses?
  • · Replies 8 ·
Replies
8
Views
2K
Undergrad Change of Basis Matrix vs Transformation matrix in the same basis....
  • · Replies 12 ·
Replies
12
Views
5K
Graduate Solve Badly Scaled Problem with Newton-Chord Technique
  • · Replies 12 ·
Replies
12
Views
3K