/* The following code was written by M. Devi */
/* matrix.h file */
/* To run this program create a header file matrix.h then copy and paste the header code in matrix.h
then create the source file matrix.cpp and copy the source code into this file. Compile and run, this code was tested*/
/* Pure Abstract class MATRIX and derived classes: MARIX_SUM adds two matrices, MATRIX_PROD multiplies two matrices, MATRIX_K multiplies matrix by a constant , MEXP handles matrix expressions */
#ifndef MATRIX_HPP
#define MATRIX_HPP
struct MATRIX
{
virtual int Rows()=0;
virtual int Cols()=0;
virtual float Element(int row,int col)=0;
/* Memory Functions */
virtual MATRIX *New()=0;
virtual void Delete()=0;
};
class MATRIX_SUM:public MATRIX
{
MATRIX *left,*right;
public:
MATRIX_SUM(MATRIX &x,MATRIX &y)
{ left=&x; right=&y; }
int Rows() { return left->Rows(); }
int Cols() { return left->Cols(); }
float Element(int row,int col)
{ return left->Element(row,col)+right->Element(row,col); }
MATRIX *New();
void Delete();
};
class MATRIX_PROD:public MATRIX
{
MATRIX *left,*right;
public:
MATRIX_PROD(MATRIX &x,MATRIX &y)
{
left=&x; right=&y;
}
int Rows() { return left->Rows(); }
int Cols() { return right->Cols(); }
float Element(int row,int col);
MATRIX *New();
void Delete();
};
class MATRIX_K:public MATRIX
{
MATRIX *mat; float k;
public:
MATRIX_K(float f,MATRIX &a) { mat=&a; k=f; }
int Rows() { return mat->Rows(); }
int Cols() { return mat->Cols(); }
float Element(int row,int col)
{ return k*mat->Element(row,col); }
MATRIX *New();
void Delete();
};
inline MATRIX_SUM operator +(MATRIX &a,MATRIX &b)
{ return MATRIX_SUM(a,b); }
inline MATRIX_PROD operator *(MATRIX &a,MATRIX &b)
{ return MATRIX_PROD(a,b); }
inline MATRIX_K operator *(float f,MATRIX &b)
{ return MATRIX_K(f,b); }
class MEXP:public MATRIX // Matrix Expression
{
MATRIX *expression;
public:
MEXP(); ~MEXP();
int Rows();
int Cols();
float Element(int row,int col);
MATRIX *New();
MATRIX &operator=(MATRIX &);
void Delete() { }
};
void Print(MATRIX &b);
#endif
/* Maps an array into a matrix */
#ifndef EMATRIX_H
#define EMATRIX_H
#include "matrix.h"
class MATRIX_ARRAY:public MATRIX
{
float *array; int row,col;
public:
MATRIX_ARRAY(int r,int c,float *ar)
{ row=r; col=c; array=ar; }
int Rows() { return row; }
int Cols() { return col; }
float Element(int r,int c)
{ return array[col*r+c]; }
MATRIX *New();
void Delete();
};
#endif EMATRIX_H
/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
/*Source file matrix.cpp corresponding to the matrix.h */
/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
#include <iostream.h>
#include "matrix.h"
float MATRIX_PROD::Element(int row,int col)
{
float fl=0; int cols=left->Cols();
for(int i=0;i<cols;++i)
fl += left->Element(row,i) * right->Element(i,col);
return fl;
}
void Print(MATRIX &a)
{
int rows=a.Rows();
int cols=a.Cols();
cout<<"Matrix Printout:\n";
for(int i=0;i<rows;++i)
{
for(int j=0;j<cols;++j)
{
cout<<a.Element(i,j) << " ";
}
cout<<"\n";
}
}
MATRIX *MATRIX_SUM::New()
{
MATRIX *l=left->New();
MATRIX *r=right->New();
return new MATRIX_SUM(*l,*r);
}
MATRIX *MATRIX_PROD::New()
{
MATRIX *l=left->New();
MATRIX *r=right->New();
return new MATRIX_PROD(*l,*r);
}
MATRIX *MATRIX_K::New()
{
MATRIX *m=mat->New();
return new MATRIX_K(k,*m);
}
void MATRIX_SUM::Delete()
{
left->Delete();
right->Delete();
delete left; left=0;
delete right; right=0;
}
void MATRIX_PROD::Delete()
{
left->Delete(); right->Delete();
delete left; left=0;
delete right; right=0;
}
void MATRIX_K::Delete()
{
mat->Delete();
delete mat; mat=0;
}
MEXP::MEXP():expression(0) { }
MEXP::~MEXP()
{
if(expression!=0) { expression->Delete(); delete expression; }
}
int MEXP::Rows()
{
if(expression!=0) return expression->Rows();
else return 0;
}
int MEXP::Cols()
{
if(expression!=0) return expression->Cols();
else return 0;
}
float MEXP::Element(int r,int c)
{
if(expression!=0)
return expression->Element(r,c);
else return 0;
}
MATRIX *MEXP::New()
{
if(expression!=0) return expression->New();
else return 0;
}
MATRIX &MEXP::operator=(MATRIX &a)
{ MATRIX *b;
b=a.New();
if(expression!=0)
{
expression->Delete();
delete expression;
}
expression=b;
return *expression;
}
MATRIX * MATRIX_ARRAY::New()
{ return new MATRIX_ARRAY(row,col,array); }
void MATRIX_ARRAY::Delete() {}
/////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
//EXAMPLE
////////////////////////////////////////////////////////////////////////////////////////////////////////////////////////
int main()
{
float fl[9]={1,2,3,4,5,6,7,8,9}; //initialize array
MATRIX_ARRAY ar(3,3,fl); // create a 3 by 3 matrix object from array fl.
MEXP me; //declare an expression matrix
me=ar;
/* multiply me by itself 13 times.
for(int i=0;i<13;++i)
{
me = me*ar;
}
cout<<"Value: "<<me.Element(0,0)<<"\n"; //show the result of (0,0) value of me raised to 13th power.
cout<<"Value: "<<me.Element(1,1)<<"\n"; //show the result of (1,1) value of me raised to 13th power.
Print(me); // show the entire matrix me raised to the 13th power.
return 0;
}
/* remark: take out Print(me); line then run it again. Notice how much faster it executes. That is because each entry
in the result is obtained without storing any intermediate results. Actually you can compute any entry of expressions of matrices of any size since the intermediate results are never stored. */