#### Asub-routineforthein-situinversionofasquarematrix When a programmer deals with mathematics, one of the most common routines added to his/her math library is matrix inversion. The following code implements a very efficient in-situ matrix inversion (meaning that no memory is allocated for another matrix to hold the result; the inverse of the matrix is built gradually in the matrix itself!). The code is based on an algorithm published in Numerical Recipes.

##### Classic VB (vb6) code

The code is written in classic vb (vb6). I am afraid that I do not have the time to convert it to .NET. Have fun:

Public Function TLMatrixInvertGaussJordan(mMatrix1() As Double, gDim As Long) As Long

On Error GoTo Er

'Matrix inversion with Gauss-Jordan elimination and full pivoting.

'mMatrix1 is transformed gradually to its inverse for reasons of memory usage :-)

'if you need to keep the original matrix, copy it to another matrix before inverting!

'Code taken from Numerical Recipes and modified for Visual Basic

'For this version we removed the solution vectors, because we only need the inverse of the matrix

Dim N As Long

Dim nP As Long

Dim i As Long

Dim iCol As Long

Dim iRow As Long

Dim J As Long

Dim K As Long

Dim L As Long

Dim lL As Long

Dim Indxc() As Long

Dim Indxr() As Long

Dim Ipiv() As Long

Dim Big As Double

Dim Dum As Double

Dim piVinv As Double

N = gDim

ReDim Indxc(N) As Long

ReDim Indxr(N) As Long

ReDim Ipiv(N) As Long

For i = 1 To N 'do 22 i=1,n

Big = 0

For J = 1 To N 'do 13 j=1,n

If Ipiv(J) <> 1 Then 'if(ipiv(j).ne.1)then

For K = 1 To N 'do 12 k=1,n

If Ipiv(K) = 0 Then 'if (ipiv(k).eq.0) then

If Abs(mMatrix1(J, K)) >= Big Then ' if (abs(a(j,k)).ge.big)then

Big = Abs(mMatrix1(J, K)) 'Big = Abs(a(J, K))

iRow = J

iCol = K

End If

End If

Next K '12

End If

Next J '13

Ipiv(iCol) = Ipiv(iCol) + 1

If iRow <> iCol Then 'If (iRow.ne.iCol) Then

For L = 1 To N 'do 14 l=1,n

Dum = mMatrix1(iRow, L) '  Dum = a(iRow, L)

mMatrix1(iRow, L) = mMatrix1(iCol, L) 'a(irow,l)=a(icol,l)

mMatrix1(iCol, L) = Dum

Next L '14

'do 15 l=1,m

'Dum = b(iRow, L)

'b(iRow, L) = b(iCol, L)

'b(iCol, L) = Dum

'enddo 15

End If

Indxr(i) = iRow

Indxc(i) = iCol

If Abs(mMatrix1(iRow, iCol)) <= 0.000000000000001 Then 'if (a(icol,icol).eq.0.) pause

Debug.Print "ERROR : TLMatrixInvertGaussJordan : Singular Matrix."

Exit Function

End If

piVinv = 1 / mMatrix1(iCol, iCol) 'pivinv=1./a(icol,icol)

mMatrix1(iCol, iCol) = 1 'a(iCol, iCol) = 1

For L = 1 To N 'do 16 l=1,n

mMatrix1(iCol, L) = mMatrix1(iCol, L) * piVinv 'a(icol,l)=a(icol,l)*pivinv

Next L '16

'do 17 l=1,m

'b(iCol, L) = b(iCol, L) * piVinv

'enddo 17

For lL = 1 To N 'do 21 ll=1,n

If lL <> iCol Then 'if(ll.ne.icol)then

Dum = mMatrix1(lL, iCol) 'dum=a(ll,icol)

mMatrix1(lL, iCol) = 0 'a(ll,icol)=0.

For L = 1 To N 'do 18 l=1,n

mMatrix1(lL, L) = mMatrix1(lL, L) - mMatrix1(iCol, L) * Dum 'a(ll,l)=a(ll,l)-a(icol,l)*dum

Next L ' 18

'do 19 l=1,m

'b(lL, L) = b(lL, L) - b(iCol, L) * Dum

'enddo 19

End If

Next lL '21

Next i '22

For L = N To 1 Step -1 'do 24 l=n,1,-1

If Indxr(L) <> Indxc(L) Then 'if(indxr(l).ne.indxc(l))then

For K = 1 To N 'do 23 k=1,n

Dum = mMatrix1(K, Indxr(L)) 'dum=a(k,indxr(l))

mMatrix1(K, Indxr(L)) = mMatrix1(K, Indxc(L)) 'a(k,indxr(l))=a(k,indxc(l))

mMatrix1(K, Indxc(L)) = Dum 'a(k,indxc(l))=dum

Next K '23

End If

Next L '24

TLMatrixInvertGaussJordan = 1

Exit Function

Er:

Debug.Print "ERROR : TLMatrixInvertGaussJordan : Unspecified error."

End Sub