InPlace
Matrix Inversion by Modified GaussJordan Algorithm
D. DasGupta, M.Tech., M.ASCE,
P.E., MCP
Former
V.P.Development, LEAP Software, Inc., Tampa, FL, USA
Former
Principal Consultant, McDonnell Douglas Automation Co., St.
Louis, MO, USA
Former
Assistant Director, Central Water and Power Commission, New
Delhi, India
REF:
Applied
Mathematics, 2013, 4, 13921396
Published Online
October 2013
ABSTRACT
The classical GaussJordan method for matrix
inversion involves
augmenting the matrix with a unit matrix and requires a workspace twice
as large as the original matrix as well as
computational operations to be performed on
oth the original and the unit matrix. A modified version of the method
for performing the inversion without explicitly generating the unit
matrix by replicating its functionality
within the original matrix space for more efficient utilization of
computational resources is presented in this
article.
Although the algorithm described here picks the
pivots solely from the
diagonal which, therefore, may not contain a zero, it did not pose any
problem for the author because he used it
to invert structural stiffness matrices which met this requirement.
Techniques such as row/column swapping to handle
offdiagonal pivots are also applicable to this method but are beyond
the scope of this article.
Keywords: Numerical Methods,
GaussJordan, Matrices, Inversion, InPlace, InCore, Structural
Analysis
Excerpt
of the document:
ADDENDUM
MicrosoftÂ® Visual BasicÂ® Routine for InPlace Matrix Inversion
'(A) Sample calling routine
'Sample matrix size Dim
iNmat As Integer iNmat
= 3
'Sample matrix elements ReDim
dblMatrix(iNmat, iNmat) As Double dblMatrix(1, 1) = 1: dblMatrix(1, 2) = 1:
dblMatrix(1, 3) = 1 dblMatrix(2, 1) = 1.2: dblMatrix(2, 2) = 1:
dblMatrix(2, 3) = 1.6 dblMatrix(3, 1) = 0.4: dblMatrix(3, 2) = 0:
dblMatrix(3, 3) = 0.2
'InPlace Inversion InPlace iNmat, dblMatrix()
'The matrix dblMatrix will now contain the
inverse
'(B) InPlace Matrix Inversion Subroutine
Public
Sub InPlace(iNmat As Integer, arrMat() As Double)
'Legend: 'iNmat =
Matrix size 'arrMat =
Matrix array
Dim
iRow As Integer, jCol As Integer, iCycle As Integer Dim
dPivot As Double, dAbsDiag As Double, iBig As Integer ReDim
blnRowProcessed(iNmat) As Boolean
'Preset entire Index vector to Active For
iRow = 1 To iNmat: blnRowProcessed(iRow) = False: Next iRow
'Loop over degrees of freedom For
iCycle = 1 To iNmat
'
Preset Pivot to Zero
dPivot = 0!
' Find
the row number and value of largest diagonal element of '
active rows (Index = 1)
' Loop
over matrix rows For
iRow = 1 To iNmat ' If
Index for this row is YES If
blnRowProcessed(iRow) = False Then '
Store ABS value of its diagonal element
dAbsDiag = Abs(arrMatrix(iRow, iRow)) '
If ABS(Diagonal) > Current Pivot then
If dAbsDiag > dPivot Then '
Store current row number and diagonal value as pivot '
Note: Since all diagonal elements
must be positive, '
this test must succeed at least once
iBig = iRow: dPivot = dAbsDiag
End If
End If
Next
' The
number of active row with highest diagonal has been determined ' and
the value of its diagonal element has been stored as Pivot
'
Reset Index of the row with highest diagonal to Inactive and its '
diagonal element to Unity.
blnRowProcessed(iBig) = True
dPivot = arrMatrix(iBig, iBig)
arrMatrix(iBig, iBig) = 1!
'
Divide the entire Pivotal row with Pivot For
jCol = 1 To iNmat
arrMatrix(iBig, jCol) = arrMatrix(iBig, jCol) / dPivot
Next
'
This will make the Pivotal Row diagonal = unity ' the
reciprocal of its original value?
'
Loop over all rows For
iRow = 1 To iNmat
'
Skip the Pivotal row which has already been processed
If iRow <> iBig Then
'
Retrieve the Pivotal column element of the current row and '
reset it to zero
dPivot = arrMatrix(iRow, iBig) arrMatrix(iRow, iBig) = 0!
'
Subtract from each element of the current row the '
Pivotal column element times the Pivotal column element in '
the row equal to current column
For jCol = 1 To iNmat
arrMatrix(iRow, jCol) =
arrMatrix(iRow, jCol)  _
dPivot * arrMatrix(iBig, jCol)
Next
End If
Next
Next
End Sub
From 'Matrix Inversion' to
home
From 'Matrix
Inversion' to Linear Algebra
