The hardware and bandwidth for this mirror is donated by dogado GmbH, the Webhosting and Full Service-Cloud Provider. Check out our Wordpress Tutorial.
If you wish to report a bug, or if you are interested in having us mirror your free-software or open-source project, please feel free to contact us at mirror[@]dogado.de.

Matrix inversion by elementary row operations

Michael Friendly

2024-10-02

The following examples illustrate the steps in finding the inverse of a matrix using elementary row operations (EROs):

These have the properties that they do not change the inverse. The method used here is sometimes called the Gauss-Jordan method, a form of Gaussian elimination. Another term is (row-reduced) echelon form.

Steps:

  1. Adjoin the identity matrix to the right side of A, to give the matrix \([A | I]\)
  2. Apply row operations to this matrix until the left (\(A\)) side is reduced to \(I\)
  3. The inverse matrix appears in the right (\(I\)) side

Why this works: The series of row operations transforms \[ [A | I] \Rightarrow [A^{-1} A | A^{-1} I] = [I | A^{-1}]\]

If the matrix is does not have an inverse (is singular) a row of all zeros will appear in the left (\(A\)) side.

Load the matlib package

library(matlib)

Create a 3 x 3 matrix

   A <- matrix( c(1, 2, 3,
                  2, 3, 0,
                  0, 1,-2), nrow=3, byrow=TRUE)

Join an identity matrix to A

   (AI <-  cbind(A, diag(3)))
##      [,1] [,2] [,3] [,4] [,5] [,6]
## [1,]    1    2    3    1    0    0
## [2,]    2    3    0    0    1    0
## [3,]    0    1   -2    0    0    1

Apply elementary row operations to reduce A to an identity matrix.

The right three cols will then contain inv(A). We will do this three ways:

  1. first, just using R arithmetic on the rows of AI
  2. using the ERO functions in the matlib package
  3. using the echelon() function

1. Using R arithmetic

    (AI[2,] <- AI[2,] - 2*AI[1,])     # row 2 <- row 2 - 2 * row 1
## [1]  0 -1 -6 -2  1  0
    (AI[3,] <- AI[3,] + AI[2,])       # row 3 <- row 3 + row 2
## [1]  0  0 -8 -2  1  1
    (AI[2,] <- -1 * AI[2,])           # row 2 <- -1 * row 2
## [1]  0  1  6  2 -1  0
    (AI[3,] <-  -(1/8) * AI[3,])        # row 3 <- -.25 * row 3
## [1]  0.000  0.000  1.000  0.250 -0.125 -0.125

Now, all elements below the diagonal are zero

    AI
##      [,1] [,2] [,3] [,4]   [,5]   [,6]
## [1,]    1    2    3 1.00  0.000  0.000
## [2,]    0    1    6 2.00 -1.000  0.000
## [3,]    0    0    1 0.25 -0.125 -0.125
      #--continue, making above diagonal == 0
    AI[2,] <- AI[2,] - 6 * AI[3,]     # row 2 <- row 2 - 6 * row 3
    AI[1,] <- AI[1,] - 3 * AI[3,]     # row 1 <- row 1 - 3 * row 3
    AI[1,] <- AI[1,] - 2 * AI[2,]     # row 1 <- row 1 - 2 * row 2

    AI
##      [,1] [,2] [,3]  [,4]   [,5]   [,6]
## [1,]    1    0    0 -0.75  0.875 -1.125
## [2,]    0    1    0  0.50 -0.250  0.750
## [3,]    0    0    1  0.25 -0.125 -0.125
   #-- last three cols are the inverse
  (AInv <- AI[,-(1:3)])
##       [,1]   [,2]   [,3]
## [1,] -0.75  0.875 -1.125
## [2,]  0.50 -0.250  0.750
## [3,]  0.25 -0.125 -0.125
   #-- compare with inv()
  inv(A)
##       [,1]   [,2]   [,3]
## [1,] -0.75  0.875 -1.125
## [2,]  0.50 -0.250  0.750
## [3,]  0.25 -0.125 -0.125

2. Do the same, using matlib functions rowadd(), rowmult() and rowswap()

   AI <-  cbind(A, diag(3))

   AI <- rowadd(AI, 1, 2, -2)        # row 2 <- row 2 - 2 * row 1
   AI <- rowadd(AI, 2, 3, 1)         # row 3 <- row 3 + row 2
   AI <- rowmult(AI, 2, -1)          # row 1 <- -1 * row 2
   AI <- rowmult(AI, 3, -1/8)        # row 3 <- -.25 * row 3

   # show result so far
   AI
##      [,1] [,2] [,3] [,4]   [,5]   [,6]
## [1,]    1    2    3 1.00  0.000  0.000
## [2,]    0    1    6 2.00 -1.000  0.000
## [3,]    0    0    1 0.25 -0.125 -0.125
    #--continue, making above-diagonal == 0
   AI <- rowadd(AI, 3, 2, -6)        # row 2 <- row 2 - 6 * row 3
   AI <- rowadd(AI, 2, 1, -2)        # row 1 <- row 1 - 2 * row 2
   AI <- rowadd(AI, 3, 1, -3)        # row 1 <- row 1 - 3 * row 3
   AI
##      [,1] [,2] [,3]  [,4]   [,5]   [,6]
## [1,]    1    0    0 -0.75  0.875 -1.125
## [2,]    0    1    0  0.50 -0.250  0.750
## [3,]    0    0    1  0.25 -0.125 -0.125

3. Using echelon()

echelon() does all these steps row by row, and returns the result

   echelon( cbind(A, diag(3)))
##      [,1] [,2] [,3]  [,4]   [,5]   [,6]
## [1,]    1    0    0 -0.75  0.875 -1.125
## [2,]    0    1    0  0.50 -0.250  0.750
## [3,]    0    0    1  0.25 -0.125 -0.125

It is more interesting to see the steps, using the argument verbose=TRUE. In many cases, it is informative to see the numbers printed as fractions.

   echelon( cbind(A, diag(3)), verbose=TRUE, fractions=TRUE)
## 
## Initial matrix:
##      [,1] [,2] [,3] [,4] [,5] [,6]
## [1,]  1    2    3    1    0    0  
## [2,]  2    3    0    0    1    0  
## [3,]  0    1   -2    0    0    1  
## 
## row: 1 
## 
##  exchange rows 1 and 2
##      [,1] [,2] [,3] [,4] [,5] [,6]
## [1,]  2    3    0    0    1    0  
## [2,]  1    2    3    1    0    0  
## [3,]  0    1   -2    0    0    1  
## 
##  multiply row 1 by 1/2
##      [,1] [,2] [,3] [,4] [,5] [,6]
## [1,]   1  3/2    0    0  1/2    0 
## [2,]   1    2    3    1    0    0 
## [3,]   0    1   -2    0    0    1 
## 
##  subtract row 1 from row 2
##      [,1] [,2] [,3] [,4] [,5] [,6]
## [1,]    1  3/2    0    0  1/2    0
## [2,]    0  1/2    3    1 -1/2    0
## [3,]    0    1   -2    0    0    1
## 
## row: 2 
## 
##  exchange rows 2 and 3
##      [,1] [,2] [,3] [,4] [,5] [,6]
## [1,]    1  3/2    0    0  1/2    0
## [2,]    0    1   -2    0    0    1
## [3,]    0  1/2    3    1 -1/2    0
## 
##  multiply row 2 by 3/2 and subtract from row 1
##      [,1] [,2] [,3] [,4] [,5] [,6]
## [1,]    1    0    3    0  1/2 -3/2
## [2,]    0    1   -2    0    0    1
## [3,]    0  1/2    3    1 -1/2    0
## 
##  multiply row 2 by 1/2 and subtract from row 3
##      [,1] [,2] [,3] [,4] [,5] [,6]
## [1,]    1    0    3    0  1/2 -3/2
## [2,]    0    1   -2    0    0    1
## [3,]    0    0    4    1 -1/2 -1/2
## 
## row: 3 
## 
##  multiply row 3 by 1/4
##      [,1] [,2] [,3] [,4] [,5] [,6]
## [1,]    1    0    3    0  1/2 -3/2
## [2,]    0    1   -2    0    0    1
## [3,]    0    0    1  1/4 -1/8 -1/8
## 
##  multiply row 3 by 3 and subtract from row 1
##      [,1] [,2] [,3] [,4] [,5] [,6]
## [1,]    1    0    0 -3/4  7/8 -9/8
## [2,]    0    1   -2    0    0    1
## [3,]    0    0    1  1/4 -1/8 -1/8
## 
##  multiply row 3 by 2 and add to row 2
##      [,1] [,2] [,3] [,4] [,5] [,6]
## [1,]    1    0    0 -3/4  7/8 -9/8
## [2,]    0    1    0  1/2 -1/4  3/4
## [3,]    0    0    1  1/4 -1/8 -1/8

These binaries (installable software) and packages are in development.
They may not be fully stable and should be used with caution. We make no claims about them.
Health stats visible at Monitor.