Efficient Matrix Transformations for C++ Programming

[DrKareem]

I'm making a matrix class that allows the user to manipulate with matrices, on C++ that is.

I'm having difficulties finding the algorithms of transforming a matrix to it's row echelon form or reduced row echelon form. The last one could help in finding the inverse of a matrix.

My teacher advised me not to develop the algorithms and get them online. So I'm wondering if anybody knows anything about those that could help.

 

[dduardo]

Here is the strategy for reduced row echelon form. Follow it carefully

Say you have a 3x3 integer array

[ 2 , 21 , 3 ]
[ 4 , 15 , 18 ]
[ 8 , 3 , 9 ]

Definition: (Row) Pivot = first non-zero digit in the row

Step 1

Pivot = (0,0) = 2
Multiplier = (1, 0) / Pivot = 4/2 = 2

(1,0) = (1, 0) - Multiplier * (0,0)
(1,1) = (1,1) - Multiplier * (0,1)
(1,2) = (1, 2) - Multiplier * (0,2)

Or simply

Row 1 = Row1 - Multiplier * Row 0

[ 2 , 21 , 3 ]
[ 0 , -27 , 12 ]
[ 8 , 3 , 9 ]

Step 2

Pivot = (0,0) = 2
Multiplier = (2, 0) / Pivot = 8/2 = 4

(2,0) = (2, 0) - Multiplier * (0,0)
(2,1) = (2,1) - Multiplier * (0,1)
(2,2) = (2, 2) - Multiplier * (0,2)

Or simply

Row 2 = Row2 - Multiplier * Row 0

[ 2 , 21 , 3 ]
[ 0 , -27 , 12 ]
[ 0 , -81 , -3 ]

Step 3

Pivot = (1,1) = -27
Multiplier = (2, 1) / Pivot = -81/-27 = 3

(2,1) = (2,1) - Multiplier * (1,1)
(2,2) = (2, 2) - Multiplier * (1,2)

Or simply

Row 2 = Row2 - Multiplier * Row 1

[ 2 , 21 , 3 ]
[ 0 , -27 , 12 ]
[ 0 , 0 , -33 ]

Step 4

Pivot = (2,2) = -33
Multiplier = (1, 2) / Pivot = -(12/33)

(1,2) = (1,2) - Multiplier * (2,2)

Or simply

Row 1 = Row1 - Multiplier * Row 2

[ 2 , 21 , 3 ]
[ 0 , -27 , 0 ]
[ 0 , 0 , -33 ]

Step 5

Pivot = (2,2) = -33
Multiplier = (0, 2) / Pivot = 3/-33 = -(1/11)

(0,2) = (0,2) - Multiplier * (2,2)

Or simply

Row 0 = Row0 - Multiplier * Row 2

[ 2 , 21 , 0 ]
[ 0 , -27 , 0 ]
[ 0 , 0 , -33 ]

Step 6

Pivot = (1,1) = -27
Multiplier = (0, 1) / Pivot = 21/-27 = -(7/9)

(0,1) = (0,1) - Multiplier * (1,1)

Or simply

Row 0 = Row0 - Multiplier * Row 1

[ 2 , 0 , 0 ]
[ 0 , -27 , 0 ]
[ 0 , 0 , -33 ]

Step 7

Pivot = (0,0) = 2

Row 0 = Row 0 / Pivot

[ 1 , 0 , 0 ]
[ 0 , -27 , 0 ]
[ 0 , 0 , -33 ]

Step 8

Pivot = (1,1) = -27

Row 1 = Row 1 / Pivot

[ 1 , 0 , 0 ]
[ 0 , 1 , 0 ]
[ 0 , 0 , -33 ]

Step 9

Pivot = (2,2) = -33

Row 2 = Row 2 / Pivot

[ 1 , 0 , 0 ]
[ 0 , 1 , 0 ]
[ 0 , 0 , 1 ]

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

If your interested in finding the Inverse matrix, just augment the one above with the identiy matrix

[ 2 , 21 , 3 | 1 , 0 , 0 ]
[ 4 , 15 , 18 | 0 , 1 , 0 ]
[ 8 , 3 , 9 | 0 , 0 , 1 ]

Then do the same procedure but when doing the multiplier row operations, make sure to do it to the right side of the augmented matrix as well.

Here is the first step with the the augmented matrix:

[ 2 , 21 , 3 | 1 , 0 , 0 ]
[ 0 , -27 ,12 | -2 , 1, 0 ]
[ 8 , 3 , 9 | 0 , 0 , 1 ]

 

[chroot]

Have a look at Numerical Recipes:

http://www.library.cornell.edu/nr/cbookcpdf.html

cf. Chapter 2.

- Warren

 

[DrKareem]

dduardo, the algorithm that you used is really nice, but in certain cases, when you have for example huge matrices is perhaps inefficient, since it would require a huge numbre of iterations. It's easy to figure out an algorithm for square matrices, but for rectangular ones it becomes harder.
Chroot, i haven't looked on the lib you gave yet, but i will see it later and comment on it.

 

[dduardo]

Are you looking for a fast matrix inversion algorithm or are you looking for a reduced row echelon algorithm that you can then apply to do matrix inversion? As you probable already know, the method I showed you in my previous post was a basic Gaussian elimination alogrithm that is O(n^3). If you just want to invert matrices you have to ask yourself what type of matrices you intend to invert? Are they dense, sparse, symmetric, banded or anything? Gaussian elimination will handle all the cases you can throw at it, but is slower than other methods which target specific cases. You could speed up gaussian elimination a little bit with strassen matrix multiplication, making it an O(n^ln(7)) algorithm. Unfortunetly it is really only good for square matrices.

 

[Hurkyl]

The other thing is that do you really need the inverse? There are algorithms that will solve the equation Ax=b without ever finding the inverse of A. (e.g. conjugate gradients)

 

[DrKareem]

The only purpose of doing this is to practive my programming skills Hurkyl.

And dduardo, my knowledge about Linear Algebra is kinda limited.

My first impression was that if i found the algorithm for reduced row echelon form, i could use the algorithm easily to do the inverse algorithm as well since they are extermely close to each other.

Since i finished my exams and now i got some free time again, i will try to implement your algorithm on C++ code, and tell you guys what i came up with :).



Also, a recursive method of doing the task came to my mind, but i read the recursive functions are ineffecient, especially if the data inserted is huge...