matrix.cpp

模式识别的作业代码

// Matrix.cpp: implementation of the CMatrix class.
//
//////////////////////////////////////////////////////////////////////
#include "Stdafx.h"
#include "matrix.h"
#include <math.h>

//////////////////////////////////////////////////////////////////////
// Construction/Destruction
//////////////////////////////////////////////////////////////////////

CMatrix::CMatrix()
{
	m_nNumColumns=1;     //矩阵列数
	m_nNumRows=1;        //矩阵行数
	m_pData=NULL;
    Init(m_nNumRows,m_nNumColumns);
}


CMatrix::CMatrix(int nRows,int nCols)
{
	m_nNumColumns=nCols;     //矩阵列数
	m_nNumRows=nRows;        //矩阵行数
	m_pData=NULL;
    Init(m_nNumRows,m_nNumColumns);
}


CMatrix::CMatrix(int nRows,int nCols,double value[])
{
    m_nNumColumns=nCols;     //矩阵列数
	m_nNumRows=nRows;        //矩阵行数
	m_pData=NULL;
    Init(m_nNumRows,m_nNumColumns);
	SetData(value);
}

CMatrix::CMatrix(int nSize)
{
    m_nNumColumns=nSize;     //矩阵列数
	m_nNumRows=nSize;        //矩阵行数
	m_pData=NULL;
    Init(nSize,nSize);
}

CMatrix::CMatrix(int nSize,double value[])
{
	m_nNumColumns=nSize;     //矩阵列数
	m_nNumRows=nSize;        //矩阵行数
	m_pData=NULL;
    Init(m_nNumRows,m_nNumColumns);
	SetData(value);
}

CMatrix::CMatrix(const CMatrix& other)
{
	m_nNumColumns=other.GetNumColumns();
    m_nNumRows=other.GetNumRow();
	m_pData=NULL;
	Init(m_nNumRows,m_nNumColumns);
	memcpy(m_pData,other.m_pData,sizeof(double)*m_nNumRows*m_nNumColumns);
}

CMatrix::~CMatrix()
{
	if(m_pData)
	{
		delete[] m_pData;
		m_pData=NULL;
	}

}


BOOL CMatrix::Init(int nRows,int nCols)
{
	if(m_pData)
	{
		delete[] m_pData;
		m_pData=NULL;
	}
	m_nNumRows=nRows;
	m_nNumColumns=nCols;
	int nSize=nRows*nCols;
	if(nSize<0)
		return FALSE;

	//动态堆分配内存
	m_pData=new double[nSize];
	if(m_pData==NULL)
		return FALSE;
	//将各元素置0
	memset(m_pData,0,sizeof(double)*nSize);
	return TRUE;
}


BOOL CMatrix::MakeUnitMatrix(int nSize)
{
	if(!Init(nSize,nSize))
		return FALSE;
	for(int i=0;i<nSize;i++)
		for(int j=0;j<nSize;j++)
			if(i==j)
				SetElement(i,j,1);
	return TRUE;
}


BOOL CMatrix::SetElement(int nRows,int nCols,double value)
{
	if(nCols<0||nCols>=m_nNumColumns||nRows<0||nRows>=m_nNumRows)
		return FALSE;
	if(m_pData==NULL)
		return FALSE;
    m_pData[nCols+nRows*m_nNumColumns]=value;
	return TRUE;

}

double CMatrix::GetElement(int nRows,int nCols) const
{
	ASSERT(nCols>=0||nCols<m_nNumColumns||nRows>=0||nRows<m_nNumRows);
    ASSERT(m_pData);
	return m_pData[nCols+nRows*m_nNumColumns];
}

int CMatrix::GetNumColumns() const
{
	return m_nNumColumns;
}


int CMatrix::GetNumRow() const
{
	return m_nNumRows;
}


double* CMatrix::GetData() const
{
	return m_pData;
}


void CMatrix::SetData(double value[])
{
	memset(m_pData,0,sizeof(double)*m_nNumRows*m_nNumColumns);
	memcpy(m_pData,value,sizeof(double)*m_nNumRows*m_nNumColumns);
}


//-------数学运算符重载--------------------------------------
CMatrix& CMatrix::operator=(const CMatrix& other)
{
	if(&other !=this)
	{
		m_nNumColumns=other.m_nNumColumns;
		m_nNumRows=other.m_nNumRows;
		Init(m_nNumRows,m_nNumColumns);
		memcpy(m_pData,other.m_pData,sizeof(double)*m_nNumRows*m_nNumColumns);
	}
	return *this;
}


BOOL CMatrix::operator==(const CMatrix& other) const
{
	//首先检查行列数是否相等
	if(m_nNumColumns!=other.GetNumColumns()||m_nNumRows!=other.GetNumRow())
		return FALSE;
	for(int i=0;i<m_nNumRows;i++)
		for(int j=0;j<m_nNumColumns;j++)
		{
			if(GetElement(i,j)!=other.GetElement(i,j))
				return FALSE;
		}
	return TRUE;
}


BOOL CMatrix::operator!=(const CMatrix& other) const
{
	return !(*this==other);
}

CMatrix CMatrix::operator+(const CMatrix& other)const
{
	ASSERT(m_nNumColumns==other.GetNumColumns()&&m_nNumRows==other.GetNumRow());
	//构造结果矩阵
	CMatrix result(*this);
	for(int i=0;i<m_nNumRows;i++)
		for(int j=0;j<m_nNumColumns;j++)
			result.SetElement(i,j,result.GetElement(i,j)+other.GetElement(i,j));
    return result;
}


CMatrix CMatrix::operator-(const CMatrix& other)const
{
	ASSERT(m_nNumColumns==other.GetNumColumns()&&m_nNumRows==other.GetNumRow());
	//构造结果矩阵
	CMatrix result(*this);
	for(int i=0;i<m_nNumRows;i++)
		for(int j=0;j<m_nNumColumns;j++)
			result.SetElement(i,j,result.GetElement(i,j)-other.GetElement(i,j));
    return result;
}


CMatrix CMatrix::operator*(double value)const
{
	CMatrix result(*this);
	for(int i=0;i<m_nNumRows;i++)
		for(int j=0;j<m_nNumColumns;j++)
			result.SetElement(i,j,this->GetElement(i,j)*value);
    return result;
}

CMatrix CMatrix::operator*(const CMatrix& other)const
{
	ASSERT(m_nNumColumns==other.GetNumRow());
	CMatrix result(m_nNumRows,other.GetNumColumns());
	double value;
	for(int i=0;i<result.GetNumRow();i++)
	{
		for(int j=0;j<other.GetNumColumns();j++)
		{
			value=0.0;
			for(int k=0;k<m_nNumColumns;k++)
			{
				value+=GetElement(i,k)*other.GetElement(k,j);
			}
			result.SetElement(i,j,value);
		}
	}
	return result;
}

double CMatrix::operator[](int index) const
{
	if( index <0 && index > m_nNumColumns*m_nNumRows)
	{
		AfxMessageBox("数组越界!");
		return 1.0e+100;
	}
	return m_pData[index];
}

BOOL CMatrix::InvertGaussJordan()
{
	int *pnRow,*pnCol,i,j,k,l,u,v;
	double d=0,p=0;
	pnRow=new int[m_nNumColumns];
	pnCol=new int[m_nNumColumns];
	if(pnRow==NULL||pnCol==NULL)
		return FALSE;
	//消元
	for(k=0;k<=m_nNumColumns-1;k++)
	{
		d=0.0;
		for(i=k;i<=m_nNumColumns-1;i++)
		{
			for(j=k;j<=m_nNumColumns-1;j++)
			{
				l=i*m_nNumColumns+j;
				p=fabs(m_pData[l]);
				if(p>d)
				{
					d=p;
					pnRow[k]=i;
					pnCol[k]=j;
				}
			}
		}
		//失败
		if(d==0.0)
		{
			delete[] pnRow;
			delete[] pnCol;
			return FALSE;
		}
		if(pnRow[k]!=k)
		{		
			for(j=0;j<=m_nNumColumns-1;j++)		
			{			
				u=k*m_nNumColumns+j;			
				v=pnRow[k]*m_nNumColumns+j;			
				p=m_pData[u];            
				m_pData[u]=m_pData[v];			
				m_pData[v]=p;		
			}	
		}

		if(pnCol[k]!=k)
		{		
			for(i=0;i<=m_nNumColumns-1;i++)		
			{			
				u=i*m_nNumColumns+k;			
				v=i*m_nNumColumns+pnCol[k];			
				p=m_pData[u];           
				m_pData[u]=m_pData[v];			
				m_pData[v]=p;		
			}	
		}
		l=k*m_nNumColumns+k;
		m_pData[l]=(1.0/m_pData[l]);

		for(j=0;j<=m_nNumColumns-1;j++)	
		{	
			if(j!=k)	
			{		
				u=k*m_nNumColumns+j;		
				m_pData[u]=m_pData[u]*m_pData[l];
			}
		}
		for(i=0;i<=m_nNumColumns-1;i++)
		{
			if(i!=k)			
			{				
				for(j=0;j<=m_nNumColumns-1;j++)				
				{
					if(j!=k)					
					{						
						u=i*m_nNumColumns+j;						
						m_pData[u]=m_pData[u]-m_pData[i*m_nNumColumns+k]*m_pData[k*m_nNumColumns+j];					
					}				
				}		
			}
		}
		for(i=0;i<=m_nNumColumns-1;i++)
		{
			if(i!=k)
			{
				u=i*m_nNumColumns+k;
				m_pData[u]=-m_pData[u]*m_pData[l];
			}
		}	
	}
	//调整恢复行列次序
	for(k=m_nNumColumns-1;k>=0;k--)
	{
		if(pnCol[k]!=k)
		{
			for(j=0;j<=m_nNumColumns-1;j++)
			{
				u=k*m_nNumColumns+j;
				v=pnCol[k]*m_nNumColumns+j;
				p=m_pData[u];
                m_pData[u]=m_pData[v];
			    m_pData[v]=p;
			}
		}
        if(pnRow[k]!=k)
		{
			for(i=0;i<=m_nNumColumns-1;i++)
			{
				u=i*m_nNumColumns+k;
				v=i*m_nNumColumns+pnRow[k];
				p=m_pData[u];
                m_pData[u]=m_pData[v];
			    m_pData[v]=p;
			}
		}
	}
	//清理内存
	delete[] pnRow;
	delete[] pnCol;
	return TRUE;
}

//////////////////////////////////////////////////////////////////////
// 矩阵的转置
//
// 参数：无
//
// 返回值：CMatrix型，指定矩阵转置矩阵
//////////////////////////////////////////////////////////////////////
CMatrix CMatrix::Transpose() const
{
	// 构造目标矩阵
	CMatrix	Trans(m_nNumColumns, m_nNumRows);

	// 转置各元素
	for (int i = 0 ; i < m_nNumRows ; ++i)
	{
		for (int j = 0 ; j < m_nNumColumns ; ++j)
			Trans.SetElement(j, i, GetElement(i, j)) ;
	}

	return Trans;
}

////////////////////////////
//GAUSS －－扫描矩阵  (袁伟添加)
// 参数：	IP ， 对角线的元素位置索引
//			A :   待扫描的矩阵（为方阵）
// 参数返回 ： A (扫描后的矩阵)
void CMatrix::WGAU(int IP, CMatrix &A) 
{
	double CQ;
	int KP;
	CMatrix B(A);

	int i,j;
	double tmp;

	KP = A.GetNumColumns(); //对于方阵而言，行列相等
	if(KP != A.GetNumRow() )
	{
		TRACE("不是方阵");
		return;
	}
	
	CQ=1.0/(A.GetElement(IP,IP)+1.e-15); //计算对角线上的元素的倒数（防止为0，加了一个很小的数）

	for(i=0;i<KP;i++)
	for(j=0;j<KP;j++)//消元
	{
			tmp = A.GetElement(i,j) - A.GetElement(IP,j)*CQ*A.GetElement(i,IP);
			B.SetElement(i,j,tmp);
	}