综合网上的经验实现了一个double型的Matrix C++类。
不想弄成template了。
支持各种基本运算+-*/ 和求逆矩阵
其中求逆矩阵是矩阵运算的难点,也是最耗时的模块。
#ifndef __CMATRIX_H__
#define __CMATRIX_H__

#include <iostream>
using namespace std;

class CMatrix


{
public:
//default constructor
CMatrix(void);

//unit Matrix constructor
CMatrix(long n);

//one dimension CMatrix constructor
CMatrix(double *pSourceValue, long nWidth);

//two dimension CMatrix constructor
CMatrix(double *pSourceValue, long nWidth,long nHeight);

//copy constructor
CMatrix(const CMatrix &);

//default destructor
~CMatrix(void);
public:


/**//***************basic matrix operator*******************************/

/**//*******************************************************************
* Function: isVector
* Describe: judge the matrix is a number
* Param: NULL
* Return Value:
false: is a number
true: is a matrix or a vector
* Author: Saha
* Date Time: 2010-8-18 10:56:45
*******************************************************************/
bool isVector();

/**//*******************************************************************
* Function: GetDeterminant
* Describe: get the value of the matrix's Determinant.
* Param: NULL
* Return Value:
double: the value of the matrix's Determinant
* Author: Saha
* Date Time: 2010-8-18 10:58:13
*******************************************************************/
double GetDeterminant();


/**//*******************************************************************
* Function: isPositive
* Describe: judge the matrix is a Positive matrix
* Param: NULL
* Return Value:
true: is
false: is not
* Author: Saha
* Date Time: 2010-8-18 11:00:18
*******************************************************************/
bool isPositive();


/**//*******************************************************************
* Function: Transpose
* Describe: get the Transposed matrix
* Param: NULL
* Return Value:
CMatrix: the Transposed matrix
* Author: Saha
* Date Time: 2010-8-18 11:01:40
*******************************************************************/
CMatrix Transpose();


/**//*******************************************************************
* Function: GetsubMatrix
* Describe: get the subMatrix.
example: 1 0 0 2 its subMatrix is 1 0 while offset is 2
2 1 3 1 2 1
3 3 3 3
* Param:
long offset:the subMatrix's width and height
* Return Value:
CMatrix: the subMatrix
* Author: Saha
* Date Time: 2010-8-18 11:04:34
*******************************************************************/
CMatrix GetsubMatrix(long offset);

/**//*******************************************************************
* Function: GetInverse
* Describe: Get Inverse Matrix using adjoint matrix
* Param:
CMatrix &m: the Inversed matrix
* Return Value:
true: has Inversed matrix
false:the matrix is a singular matrix and hasn't Inversed matrix
* Author: Saha
* Date Time: 2010-8-18 14:04:37
*******************************************************************/
bool GetInverse(CMatrix &m);
public:

/**//***********Overloaded Part*****************/
//matrix add
CMatrix operator+(CMatrix &);

//matrix sub
CMatrix operator-(CMatrix &);

//matrix multiply with matrix
CMatrix operator*(CMatrix &);

//matrix multiply with real number
CMatrix operator*(double alpha);

//real number multiply with matrix
friend CMatrix operator*(double alpha, CMatrix &m);
//matrix divide with matrix
CMatrix operator/(CMatrix &);

CMatrix operator+=(CMatrix &);

CMatrix operator-=(CMatrix &);

CMatrix operator*=(double alpha);

//matrix evaluate operation
CMatrix &operator=(CMatrix &);

//matrix suffix operation
double *operator[](long heightPos);

//matrix out stream
friend ostream & operator << (ostream &, CMatrix &);

public:
long getHeight();
long getWidth();

void GetValue(double *p);

private:
double GetAlgebraicCofactor(long nX, long nY);

//the pointer to store the matrix data
double * pMatrix;
long m_nWidth;
long m_nHeight;
};

#endif

#include "Matrix.h"
#include <assert.h>

CMatrix::CMatrix(void)


{
pMatrix = new double[1];
m_nWidth = 0;
m_nHeight = 0;
}

CMatrix::CMatrix(long n)


{
m_nHeight = m_nWidth = n;
pMatrix = new double[n*n];

for(long i = 0; i < n; i++)

{
for(long j = 0; j < n; j++)

{
if(i == j)

{
*(pMatrix + n*i + j)=1;
}
else

{
*(pMatrix + n*i + j)=0;
}
}
}
}

CMatrix::CMatrix(double *pSourceValue,long nWidth)


{
long nHeight = 1;
pMatrix = new double[nWidth*nHeight];

if (NULL == pSourceValue)

{
for(long i = 0; i < nHeight; i++)

{
for(long j = 0; j < nWidth; j++)

{
*(pMatrix + nWidth*i + j) = 0;
}
}
}
else

{
for(long i = 0; i < nHeight; i++)

{
for(long j = 0; j < nWidth; j++)

{
*(pMatrix + nWidth*i + j) = *(pSourceValue + nWidth*i + j);
}
}
}

m_nWidth = nWidth;
m_nHeight = nHeight;
}

CMatrix::CMatrix(double *pSourceValue,long nWidth,long nHeight)


{
pMatrix = new double[nWidth*nHeight];
if (NULL == pSourceValue)

{
for(long i = 0; i < nHeight; i++)

{
for(long j = 0; j < nWidth; j++)

{
*(pMatrix + nWidth*i + j) = 0;
}
}
}
else

{
for(long i = 0; i < nHeight; i++)

{
for(long j = 0; j < nWidth; j++)

{
*(pMatrix + nWidth*i + j) = *(pSourceValue + nWidth*i + j);
}
}
}

m_nWidth = nWidth;
m_nHeight = nHeight;
}

//copy constructor
CMatrix::CMatrix(const CMatrix & m)


{
m_nHeight = m.m_nHeight;
m_nWidth = m.m_nWidth;

pMatrix = new double[m_nHeight*m_nWidth];

for(long i = 0; i < m_nHeight; i++)

{
for(long j = 0; j < m_nWidth; j++)

{
*(pMatrix + m_nWidth*i + j) = *(m.pMatrix + m_nWidth*i + j);
}
}
}

CMatrix::~CMatrix(void)


{
delete []pMatrix;
}

long CMatrix::getWidth()


{
return m_nWidth;
}

long CMatrix::getHeight()


{
return m_nHeight;
}

void CMatrix::GetValue(double *p)


{
for(long i = 0; i < m_nHeight; i++)

{
for(long j = 0; j < m_nWidth; j++)

{
*(p + m_nHeight*j + i) = *(pMatrix + m_nWidth*i + j);
}
}
}

bool CMatrix::isVector()


{
return !(m_nWidth == 1 && m_nHeight == 1);
}

CMatrix CMatrix::GetsubMatrix(long offset)


{
assert(m_nHeight == m_nWidth && offset <= m_nWidth && offset >= 0);
double * pTmp = new double[offset*offset];

for(long i = 0; i < offset; i++)

{
for(long j = 0; j < offset; j++)

{
*(pTmp + offset*i + j)=*(pMatrix + m_nWidth*i + j);
}
}
CMatrix m(pTmp, offset, offset);
delete []pTmp;
return m;
}


double CMatrix::GetDeterminant()


{
long i, j, m;
long lop = 0;
double result = 0;
double mid = 1;

if(1 != m_nWidth)

{
lop = (m_nWidth == 2)? 1 : m_nWidth;
for(m = 0; m < lop; m++)

{
mid=1;
//order add
for(i = 0, j = m; i < m_nWidth; i++, j++)

{
mid = mid * (*(pMatrix + i*m_nWidth + j%m_nWidth));
}
result += mid;
}

for(m = 0; m < lop; m++)

{
mid = 1;
//inverse order sub
for(i = 0, j = m_nWidth - 1 - m + m_nWidth; i < m_nWidth; i++, j--)

{
mid = mid * (*(pMatrix + i*m_nWidth + j%m_nWidth));
}
result -= mid;
}
}
else

{
result = *pMatrix;
}
return (result);
}

bool CMatrix::isPositive()


{
assert(m_nWidth == m_nHeight);
bool result = true;
CMatrix m;

for(long i = 1; i <= m_nHeight; i++)

{
m = this->GetsubMatrix(i);
if(m.GetDeterminant() <= 0)

{
result = false;
break;
}
}
return result;
}

CMatrix CMatrix::Transpose()


{
double * pTmp = new double[m_nWidth*m_nHeight];
for(long i = 0; i < m_nHeight; i++)

{
for(long j = 0; j < m_nWidth; j++)

{
*(pTmp + m_nHeight*j + i) = *(pMatrix + m_nWidth*i + j);
}
}
CMatrix m(pTmp, m_nHeight, m_nWidth);
delete []pTmp;
return m;
}

double CMatrix::GetAlgebraicCofactor(long nX, long nY)


{
assert(nX < m_nWidth && nY < m_nHeight);
assert(m_nHeight == m_nWidth);

long tmpHeight = m_nHeight - 1;
long tmpWidth = m_nWidth - 1;
double *ptmp = new double[tmpWidth*tmpHeight];
double *p = ptmp;

for (int i = 0; i < m_nHeight; i++)

{
for (int j = 0; j < m_nWidth; j++)

{
if (i != nX && j != nY)

{
*p = *(pMatrix + i*m_nWidth + j);
p++;
}
}
}

CMatrix m(ptmp, tmpWidth, tmpHeight);
int nGene = (nX + nY) % 2 == 0 ? 1 : -1;
return (double)nGene * m.GetDeterminant();
}

bool CMatrix::GetInverse(CMatrix &m)


{
double detA = GetDeterminant();
if (0 == detA)

{
return false;
}
if (1 == m_nWidth && 1 == m_nHeight)

{
CMatrix mIns(1);
m = 1/detA * mIns;
return true;
}
double *pTmp = new double[m_nWidth*m_nHeight];

for (int i = 0; i < m_nHeight; i++)

{
for (int j = 0; j < m_nWidth; j++)

{
*(pTmp + i*m_nWidth + j) = GetAlgebraicCofactor(j, i);
}
}
CMatrix mIns(pTmp, m_nWidth, m_nHeight);
detA = 1/detA;
mIns = detA * mIns;
m = mIns;
delete []pTmp;
return true;
}

CMatrix CMatrix::operator +(CMatrix &m1)


{
assert(m1.m_nHeight == m_nHeight && m1.m_nWidth == m_nWidth);

long tmpHeight = m1.m_nHeight;
long tmpWidth = m1.m_nWidth;
double *pTmp = new double[tmpWidth*tmpHeight];

for(long i = 0; i < tmpHeight; i++)

{
for(long j = 0; j < tmpWidth; j++)

{
*(pTmp + tmpWidth*i + j) = *(m1.pMatrix + tmpWidth*i + j) + *(pMatrix + tmpWidth*i + j);
}
}
CMatrix m(pTmp, tmpWidth, tmpHeight);
delete []pTmp;
return m;
}

CMatrix CMatrix::operator -(CMatrix &m1)


{
assert(m1.m_nHeight == m_nHeight && m1.m_nWidth == m_nWidth);

long tmpHeight = m1.m_nHeight;
long tmpWidth = m1.m_nWidth;
double * pTmp = new double[tmpWidth*tmpHeight];

for(long i = 0; i < tmpHeight; i++)

{
for(long j = 0; j < tmpWidth; j++)

{
*(pTmp + tmpWidth*i + j) = *(pMatrix + tmpWidth*i + j) - *(m1.pMatrix + tmpWidth*i + j);
}
}
CMatrix m(pTmp, tmpWidth, tmpHeight);
delete []pTmp;
return m;
}

CMatrix CMatrix::operator *(CMatrix &m1)


{
if(!this->isVector() && m1.isVector())//left number ,right matrix

{
CMatrix m;
m= pMatrix[0] * m1;
return m;
}
else if(this->isVector() && !m1.isVector())//left matrix, right number

{
CMatrix m;
m = (*this) * m1[0][0];
return m;
}
else if(!this->isVector() && !m1.isVector())//left and right are numbers

{
double * pTmp = new double[1];
pTmp[0] = pMatrix[0] * m1[0][0];
CMatrix m(pTmp, 1, 1);
delete []pTmp;
return m;
}
else if(this->isVector() && m1.isVector() && m_nWidth == m1.m_nHeight)//left and right are matrix

{
double sum;
double * pTmp=new double[m_nHeight*m1.m_nWidth];

for(long i = 0; i < m_nHeight; i++)

{
for(long j = 0; j < m1.m_nWidth; j++)

{
sum = 0;
for(long k = 0; k < m_nWidth; k++)

{
sum += (*(pMatrix + m_nWidth*i + k))*(m1[k][j]);
}
*(pTmp + m1.m_nWidth*i + j) = sum;
}
}
CMatrix m(pTmp, m1.m_nWidth, m_nHeight);
delete []pTmp;
return m;
}
else//unknown operation

{
assert(0);
return *this;
}
}

CMatrix operator*(double alpha,CMatrix & m1)


{
CMatrix m = m1;
for(long i = 0; i < m.m_nHeight; i++)

{
for(long j = 0; j < m.m_nWidth; j++)

{
m[i][j] = alpha*m1[i][j];
}
}
return m;
}

CMatrix CMatrix::operator/(CMatrix & m1)


{
CMatrix m2;
m1.GetInverse(m2);
return *this * m2;
}

CMatrix CMatrix::operator*(double alpha)


{
return alpha*(*this);
}

CMatrix CMatrix::operator+=(CMatrix &m)


{
return *this + m;
}

CMatrix CMatrix::operator-=(CMatrix &m)


{
return *this - m;
}

CMatrix CMatrix::operator *=(double alpha)


{
return (*this) * alpha;
}


CMatrix & CMatrix::operator =(CMatrix & m)


{
if(&m == this)

{
return *this;
}
m_nHeight = m.m_nHeight;
m_nWidth = m.m_nWidth;
delete []pMatrix;
pMatrix = new double[m_nHeight*m_nWidth];

for(long i = 0; i < m_nHeight; i++)

{
for(long j = 0; j < m_nWidth; j++)

{
*(pMatrix + m_nWidth*i + j) = *(m.pMatrix + m_nWidth*i + j);
}
}
return *this;
}

double * CMatrix::operator [](long heightPos)


{
assert(heightPos >= 0 && heightPos < m_nHeight);
return pMatrix + m_nWidth*heightPos;
}

ostream & operator<<(ostream & os,CMatrix & m)


{
os << "CMatrix:m_nHeight=" << m.m_nHeight << endl;
os << "CMatrix:m_nWidth=" << m.m_nWidth << endl;

for(long i = 0; i < m.m_nHeight; i++)

{
for(long j = 0; j < m.m_nWidth; j++)

{
os << m[i][j] << " ";
}
os << endl;
}
return os;
}


测试:
#include "Matrix.h"


void main()


{

double arr[3][3]=
{
{1, 0, 1},
{2, 1, 0},
{-3, 2, -5}};
CMatrix m1(&arr[0][0], 3, 3), m2(&arr[0][0], 3, 3);
CMatrix m3;
cout << m3 << endl;
m2.GetInverse(m3);
cout << m1 << endl;
cout << 4*m2 << endl;
cout << m3 << endl;
cout << m1 + m3 << endl;
cout << m1 - m3 << endl;
cout << m1 * m3 << endl;
system("pause");
return;
}
posted on 2010-08-18 17:06
saha 阅读(741)
评论(0) 编辑 收藏 引用