綜合網上的經驗實現了一個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 閱讀(744)
評論(0) 編輯 收藏 引用