I feel old, how the times have changed. Old document I written in 1991(Matrix inverse). I was a student, I could not afford matlab, I thought I’ll be smart and write my own matrix code in C++ just for fun.
julia> a = rand(1000,1000);
julia> @time ia =inv(a);
0.169260 seconds (6 allocations: 8.126 MiB)
old document about c++ object I written myself back in 1991
Document version V1.0
Umatrix (The Universal Matrix object in C++) by Steven Siew.
Table of contents
1) Introduction
2) List of functions
3) Examples
4) Common Errors and FAQ
------------------------------------------------------------------------------
Section (1) Introduction
Umatrix was written out of frustration when trying to implement many
C programs which require matrix operations. Instead of writting lots
of "for" loops.I have decided to write a universal matrix object once and
for all.
Using umatrix it's possible to write the following
#include <stdio.h>
#include <math.h>
#include "umatrix.h"
int main()
{ matrix A,B,C;
A.create(4,3);
B.create(3,6);
A.fill(2.4);
B.fill(4.5);
C= A*B;
C.show();
return 0;
}
As you can see matrix object is a dynamic object which can be created in
runtime. Also there is no predefined limit on the size of the matrix object.
The only practical limits are the amount of memory on the computer. All the
usual operators "+","-","*" have been overloaded.
I hope you will have fun with umatrix, I know I did (Damn those bloody
bugs!) and I know it will make you life easier.
Time taken to find the inverse of a 1000 by 1000 matrix
is 1:15 (1 hour 15 minutes). Of which 74 minutes is system time and
40 seconds is user time. Note: This is the unoptimise version.
Using g++ -O2 option. Time taken to find the inverse of a 1000 by 1000
matrix is 33 minutes. Of which 31 minutes is system time and 23 seconds
is user time.
Using g++ -m486 -O2 option on a Pentium II 233Mhz PC running Linux with
128MB of memory with no other processes running, the time taken to find
the inverse of a 1000 by 1000 matrix is 240 seconds or 4 minutes exactly.
Using egcs -mcpu=686 -O2 option on a Pentium II 333Mhz PC running Red Hat
Linux 6.1 with 128MB of memory with gnome running, the time taken to find
the inverse of a 1000 by 1000 matrix is 188 seconds or 3 minutes 8 seconds.
------------------------------------------------------------------------------
Section (2) List of functions
2.1 matrix(Mrow,Mcol)
This is the matrix constructor. for example
matrix A(3,3);
This creates a matrix object called A and dynamically allocate
3*3*sizeof(double) bytes in memory to store the contents of the
matrix A.
2.2 .create(Mrow,Mcol)
This dynamically allocate memory as storage for the contents of the
matrix object.
matrix A;
A.create(3,3);
2.3 .show()
This display the contents of the matrix on the stdout. Warning! This
is not intended as a universal means of outputing the matrix. It's
included as a QUICK & DIRTY method of displaying the contents of a
matrix. YOU SHOULD WRITE YOUR OWN METHOD OF DISPLAYING A MATRIX.
matrix A(4,5);
A.fill(666.0);
A.show();
2.4 .set()
This obtains the contents of a matrix from the user through the stdin.
Warning! This is not intended as a universal means of obtaining the
contents of a matrix from the user. It's included as a QUICK & DIRTY
method for doing so. YOU SHOULD WRITE YOUR OWN METHOD OF OBTAINING A
MATRIX FROM THE USER.
2.5 .destroy()
This deallocates the memory allocated to a matrix in storing it's
(numeric) contents. Use this the same way you would use the free()
function.
matrix A(7,8); // allocated 7*8*sizeof(double) bytes
A.destory(); // deallocate them
2.6 .fill()
This fills a matrix up with a number.
matrix A(2,5);
A.fill(42.0);
2.7 .ivector() and .dvector()
This are two internal functions provided for POWER USERS only. Do not
use this unless you know what you are doing.
2.8 .mread()
This returns the value of an element at row R and col C.
matrix A(2,3);
A(0,0,1); A(0,1,2); A(0,2,3);
A(1,0,4); A(1,1,5); A(1,2,6);
printf("row 0 col 1 is %f which is 2\n",A.mread(0,1));
2.9 .mwrite()
This puts a value into the element at row R and col C.
matrix A(2,3);
A(0,0,1); A(0,1,2); A(0,2,3);
A(1,0,4); A(1,1,5); A(1,2,6);
A.mwrite(0,1,666);
printf("row 0 col 1 is %f which is 666\n",A.mread(0,1));
2.10 operator +
This is the overloaded + operator.
2.11 operator -
This is the overloaded - operator.
2.12 operator *
This is the overloaded * operator.
2.13 operator /
This is the overloaded + operator.
2.14 ()
This is overloaded operator () used to a quick means of obtaining the
value of an element in the matrix.
matrix A(2,3);
A(0,0,1); A(0,1,2); A(0,2,3);
A(1,0,4); A(1,1,5); A(1,2,6);
printf("row 0 col 1 is %f which is 2\n",A(0,1));
2.15 ()
This is overloaded operator () used to a quick means of setting the
value of an element in the matrix.
matrix A(2,3);
A(0,0,1); A(0,1,2); A(0,2,3);
A(1,0,4); A(1,1,5); A(1,2,6);
A(0,1,666); // setting a new value into row 0 col 1
printf("row 0 col 1 is %f which is 666\n",A.mread(0,1));
2.16 =
This is the overloaded = operator
2.17 <<
This is the overloaded << operator. THIS HAS A SPECIAL MEANING.
2.18 .appfunc()
This applies a (double) func to every element in the matrix. For example
if you want to replace every element in the matrix with its squart root.
matrix A(2,3);
A(0,0,1); A(0,1,2); A(0,2,3);
A(1,0,4); A(1,1,5); A(1,2,6);
A.appfunc(sqrt);
A.show;
2.19 .copy()
This creates a copy of the matrix.
matrix A(2,3),B;
A.fill(8.1);
B=A.copy();
2.20 .transpose()
This creates the transpose of a matrix.
matrix A(2,3),B;
A(0,0,1); A(0,1,2); A(0,2,3);
A(1,0,4); A(1,1,5); A(1,2,6);
B=A.transpose();
Note: Warning! Warning! Do not do this!!! You cannot force an X by Y matrix
to be a Y by X matrix.
A=A.transpose();
If you need A to be its own transpose, do this instead.
Temp=A.transpose();
A.destroy();
A=Temp;
2.21 .minormatrix()
This creates the minormatrix of a matrix. Check your maths textbook for
the definition of a minormatrix as I'm too lazy to do it for you. Be careful
as ".minormatrix(r,c)" requires r be between 0..(row-1) and c to be between
0..(col-1). When in doubt, read the source Luke!
matrix A(3,3),B;
A(0,0,1); A(0,1,2); A(0,2,3);
A(1,0,4); A(1,1,5); A(1,2,6);
A(2,0,7); A(2,1,8); A(2,2,9);
B=A.minormatrix(0,0);
2.22 .cofactor()
This creates the cofactor of a matrix. Check your maths textbook for
the definition of a cofactor as I'm too lazy to do it for you. Be careful
as ".cofactor(r,c)" requires r be between 0..(row-1) and c to be between
0..(col-1). When in doubt, read the source Luke!
matrix A(3,3);
double cofactor_result;
A(0,0,1); A(0,1,2); A(0,2,3);
A(1,0,4); A(1,1,5); A(1,2,6);
A(2,0,7); A(2,1,8); A(2,2,9);
cofactor_result=A.cofactor(0,0);
2.23 .det()
This creates the determinant of a matrix. Surely everyone knows what a
determinant of a matrix is right? Girls?
matrix A(3,3);
double det_result;
A(0,0,1); A(0,1,2); A(0,2,3);
A(1,0,4); A(1,1,5); A(1,2,6);
A(2,0,7); A(2,1,8); A(2,2,9);
det_result=A.det();
2.24 .adjoint()
2.25 .divmatrix()
2.26 .inv()
This creates the inverse of a matrix. Surely everyone knows what a
determinant of a matrix is right? Girls?
matrix A(3,3),B;
A(0,0,1); A(0,1,2); A(0,2,3);
A(1,0,4); A(1,1,5); A(1,2,6);
A(2,0,7); A(2,1,8); A(2,2,9);
B=A.inv();
2.27 .gaussj()
2.28 .ludcmp()
2.29 .ludecomp()
2.30 .lubksb()
2.31 .ludsoln()
2.32 .ludinv()
2.34 .luselfinv()
This is same as .inv() however it is currently the fastest method of
finding the inverse of a square matrix. This is what I use to for producing
benchmark results.