matrixlyq.cpp

来自「GPS 定位授时源码」· C++ 代码 · 共 541 行

CPP
541
字号
// MATRIX1.cpp: implementation of the MATRIX class.
//
//////////////////////////////////////////////////////////////////////

#include "stdafx.h"
#include "MatrixLYQ.h"
#include "math.h"

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


MATRIX::MATRIX()
{
	this->nRow		= 1;
	this->nCol      = 1;
	this->pMATRIX   =  new double*[1];
	this->pMATRIX[0] = new double[1];
}

MATRIX::~MATRIX()
{
	for(int i=0;i<this->nRow;i++)
	{
		delete[this->nCol] this->pMATRIX[i];
	}
	delete[this->nRow] this->pMATRIX;
	this->pMATRIX = NULL;
}

MATRIX::MATRIX(int nRow, int nCol, double val)
{
	this->nRow = nRow;
	this->nCol = nCol;

	this->pMATRIX = new double*[this->nRow];
	for(int i=0;i<this->nRow;i ++)
	{
		this->pMATRIX[i] = new double[this->nCol];
	}

	for(i = 0;i<this->nRow;i++)
	{
		for(int n =0;n<this->nCol;n++)
		{
			this->pMATRIX[i][n] = val;
		}
	}
}
MATRIX MATRIX::operator = (const MATRIX&  _m)
{

	if (this != & _m)
	{
		for(int n = 0 ;n < this->nRow;n++)
		{
			delete [this->nCol]this->pMATRIX[n];
		}
		delete [nRow]this->pMATRIX;

		this->nRow = _m.nRow;
		this->nCol = _m.nCol;

		this->pMATRIX = new double*[this->nRow];
		for(int l=0;l < this->nRow; l ++)
		{
			this->pMATRIX[l] = new double[this->nCol];
		}
		for(int i =0;i <_m.nRow;i++)
		{
			for(int n =0;n < _m.nCol;n ++)
			{
				this->pMATRIX[i][n] = _m.pMATRIX[i][n];
			}
		}

	}

	return *this;
}
/*

MATRIX  MATRIX::operator +(const MATRIX &_m)
{
	MATRIX N(_m.nRow,_m.nCol);

	if (this->nRow == _m.nRow && this->nCol == _m.nCol)
	{

		
		
		for(int i=0;i < this->nRow;i ++)
		{
			for(int n = 0; n < this->nCol; n ++)
			{
				N.pMATRIX[i][n] = this->pMATRIX[i][n] + _m.pMATRIX[i][n];
			}
		}
	}
	return  N;
}
*/


MATRIX::MATRIX(MATRIX &_m)
{
	this->nRow = _m.nRow;
	this->nCol = _m.nCol;
	this->pMATRIX = new double*[this->nRow];
	for(int i=0;i < this->nRow;i ++)
	{
		this->pMATRIX[i] = new double[this->nCol];
	}
	for(i = 0 ;i < this->nRow;i++)
	{
		for(int n=0;n < this->nCol;n++)
		{
			this->pMATRIX[i][n] = _m.pMATRIX[i][n];
		}
	}
}

bool MATRIX::operator ==(const MATRIX &_m)
{
	if (nRow != _m.nRow || nCol != _m.nCol)
	{
		return false;
	}
	else
	{
		for(int i = 0;i < this->nRow;i++)
		{
			for(int n=0;n < this->nCol; n++)
			{
				if (pMATRIX[i][n] != _m.pMATRIX[i][n])
				{
					return false;
				}
			}
		}
	}
	return true;
}
MATRIX  operator +(const MATRIX &_m,const MATRIX &_n)
{
	MATRIX N(_m.nRow,_m.nCol);

	if (_n.nRow == _m.nRow && _n.nCol == _m.nCol)
	{
		
		
		
		for(int i=0;i < _n.nRow;i ++)
		{
			for(int n = 0; n < _n.nCol; n ++)
			{
				N.pMATRIX[i][n] = _n.pMATRIX[i][n] + _m.pMATRIX[i][n];
			}
		}
	}

	return  N;

}

MATRIX MATRIX::operator !()
{
	MATRIX temp(this->nCol,this->nRow);
	for(int i = 0;i < this->nRow;i ++)
	{
		for(int n=0;n< this->nCol;n++)
		{
			temp.pMATRIX[n][i] = this->pMATRIX[i][n];
		}
	}
	return temp;
}

MATRIX MATRIX::operator ~()
{
	MATRIX  m(*this);
	MATRIX  temp(m.nRow,m.nCol);

	
	int i,j;
	for(i = 0;i < m.nRow;i ++)
	{
		for(j = 0 ;j < m.nCol; j ++)
		{
			if (fabs(m.pMATRIX[i][j]) < 1e-10)
			{
				m.pMATRIX[i][j] = 0;
			}
		}
	}
	if (!IsSqare(m))
	{
//		AfxMessageBox("不是方阵不能求逆");
		return *this;
	}
	
	for(i = 0 ;i < temp.nRow ; i ++)
	{
		for(j = 0 ; j < temp.nCol; j ++)
		{
			if (i == j)
			{
				temp.pMATRIX[i][j] = 1.0;
			}
			else
			{
				temp.pMATRIX[i][j] = 0.0;
			}
			
		}
	}

	int max;
	int n = 0;
	double nMid;
	
	for(i = 0;i < temp.nCol; i ++)//Row
	{
		//逐列选出列的主元
		max		= i;
		for(j =i + 1; j < temp.nRow ;j ++)//可能有问题n
		{
			if(fabs(m.pMATRIX[j][i]) > fabs(m.pMATRIX[max][i]))
			{
				max		= j;
			}

		}
		//max是第i列中最大元素对应的行号
		j = max;
		if (fabs(m.pMATRIX[j][i]) < 1e-20)
		{
//			AfxMessageBox("矩阵是奇异阵,不能求逆");
			for(int u=0;u < this->nRow; u ++)
			{
				for(int v=0;v < this->nCol; v ++)
				{
					this->pMATRIX[u][v] = 0;
				}
			}
			return *this;
		}
		if (j != i)
		{//交换两行元素的位置
			for( n = 0;n < m.nCol;n ++)
			{
				nMid	= m.pMATRIX[i][n];
				m.pMATRIX[i][n]		= m.pMATRIX[j][n];
				m.pMATRIX[j][n]		= nMid;

				nMid	= temp.pMATRIX[i][n];
				temp.pMATRIX[i][n]		= temp.pMATRIX[j][n];
				temp.pMATRIX[j][n]		= nMid;
			}
		}
		
		for(n = i + 1; n < m.nRow; n ++)//
		{
			if (m.pMATRIX[n][i] == 0)
			{
				continue;
			}
			nMid		= (-1.0) * m.pMATRIX[n][i] / m.pMATRIX[i][i];//j;
			for(j = 0 ;j < m.nCol ; j ++)
			{
				m.pMATRIX[n][j]		+= nMid * m.pMATRIX[i][j];
				temp.pMATRIX[n][j]		+= nMid * temp.pMATRIX[i][j];
				
			}
		}
		
	}
	
	
	for(i = m.nCol -1 ; i > 0 ; i --)
	{
		for(j = i -1 ; j >= 0; j --)
		{
			if (m.pMATRIX[j][i] == 0)
			{
				continue;
			}
			nMid		= (-1) * m.pMATRIX[j][i] /m.pMATRIX[i][i];
			for(n = 0 ; n < m.nCol ; n ++)
			{			
				m.pMATRIX[j][n]		 += nMid * m.pMATRIX[i][n];
				temp.pMATRIX[j][n]       += nMid * temp.pMATRIX[i][n];
			}
		}
	}

	for(i = 0; i < m.nRow; i ++)
	{
		for(j = 0;j < m.nCol ; j ++)
		{
			temp.pMATRIX[i][j]		/= m.pMATRIX[i][i];
		}
		
	}
	return temp;
}

bool MATRIX::IsSqare(const MATRIX&  _m)
{
	if (_m.nRow == _m.nCol)
	{
		return true;
	}
	else
	{
		return false;
	}
}

MATRIX operator * (const MATRIX &_m, const MATRIX & _n)
{
	MATRIX _L(_m.nRow,_n.nCol);

	if (_m.nCol != _n.nRow)
	{
//		AfxMessageBox("两个矩阵不能相乘");
		return _L;
	}
	int i,j,k;
	for(i = 0;i < _m.nRow;i ++)
	{
		for(j = 0 ;j < _n.nCol;j ++)
		{
			_L.pMATRIX[i][j] = 0.0;
			for(k = 0; k < _m.nCol; k ++)
			{
				_L.pMATRIX[i][j]  += _m.pMATRIX[i][k] * _n.pMATRIX[k][j];
/*
								if (fabs(_L.pMATRIX[i][j]) < 1e-10)
								{
									_L.pMATRIX[i][j] = 0;
								}*/
				
			}
		}
	}
	return _L;
}


MATRIX operator - (const MATRIX & _m,const MATRIX & _n)
{
	MATRIX  _L(_m.nRow,_m.nCol);
	if (_m.nRow != _n.nRow || _m.nCol != _n.nCol)
	{
//		AfxMessageBox("两个矩阵非同型不能相减");
		return _L;
	}
	int i,j;
	for(i = 0; i < _m.nRow ; i ++)
	{
		for(j = 0 ; j < _m.nCol ; j ++)
		{
			_L.pMATRIX[i][j]	= _m.pMATRIX[i][j] - _n.pMATRIX[i][j];
		}
	}
	return _L;
}

MATRIX operator - (const MATRIX & _m,double dx)
{
	MATRIX  _L(_m.nRow,_m.nCol);
	int i,j;
	for(i = 0; i < _m.nRow ; i ++)
	{
		for(j = 0 ; j < _m.nCol ; j ++)
		{
			_L.pMATRIX[i][j]	= _m.pMATRIX[i][j] - dx;
		}
	}
	return _L;
}
MATRIX operator - (double dx,const MATRIX & _m)
{
	MATRIX  _L(_m.nRow,_m.nCol);
	int i,j;
	for(i = 0; i < _m.nRow ; i ++)
	{
		for(j = 0 ; j < _m.nCol ; j ++)
		{
			_L.pMATRIX[i][j]	=  dx - _m.pMATRIX[i][j];
		}
	}
	return _L;
}
MATRIX operator * (double dx,const MATRIX & _m)
{
	MATRIX  _L(_m.nRow,_m.nCol);
	int i,j;
	for(i = 0; i < _m.nRow ; i ++)
	{
		for(j = 0 ; j < _m.nCol ; j ++)
		{
			_L.pMATRIX[i][j]	=  dx * _m.pMATRIX[i][j];
		}
	}
	return _L;
}

MATRIX operator * (const MATRIX & _m,double dx)
{
	MATRIX  _L(_m.nRow,_m.nCol);
	int i,j;
	for(i = 0; i < _m.nRow ; i ++)
	{
		for(j = 0 ; j < _m.nCol ; j ++)
		{
			_L.pMATRIX[i][j]	= _m.pMATRIX[i][j]  * dx;
		}
	}
	return _L;
}
MATRIX operator ^(const MATRIX &_m,int nPow)
{
	MATRIX _L ;
	_L =  _m;
	for(int i = 0;i < nPow - 1 ; i ++)
	{
		_L = _L * _L;
	}
	return _L;
}
/*

MATRIX & MATRIX::operator =(double val)
{
	for(int n = 0 ;n < this->nRow;n++)
	{
		for(int i =0;i < this->nCol;i++)
		{
			this->pMATRIX[n][i]		= val;
		}
	}
	return *this;
}
*/
MATRIX MATRIX::MallocSpace(int nRow, int nCol)
{
	for(int i=0;i<this->nRow;i++)
	{
		delete[this->nCol] this->pMATRIX[i];
	}
	delete[this->nRow] this->pMATRIX;

	this->nRow = nRow;
	this->nCol = nCol;
	
	this->pMATRIX = new double*[this->nRow];
	
	for( i=0;i<this->nRow;i ++)
	{
		this->pMATRIX[i] = new double[this->nCol];
	}
	
	for(i = 0;i<this->nRow;i++)
	{
		for(int n =0;n<this->nCol;n++)
		{
			this->pMATRIX[i][n] = 0;
		}
	}
	return *this;
}
MATRIX MATRIX::ExpandSpace(int nRow,int nCol) 
{
	int    nRcsRow;
	int	   nRcsCol;
	MATRIX Rut;
	if (nRow == 0 || nCol == 0)
	{
		for(int i = 0;i < this->nRow;i ++)
		{
			delete(this->pMATRIX[i]);
		}
		delete(this->pMATRIX);
		this->pMATRIX = NULL;
		this->nRow = 0;
		this->nCol = 0;
	}
	else
	{
		if (this != NULL)
		{
			nRcsCol = this->nCol;
			nRcsRow = this->nRow;
		}
		else
		{
			nRcsCol		= 0;
			nRcsRow		= 0;
		}
		Rut.MallocSpace(nRow,nCol);
		for(int i=0;i < nRow && i < nRcsRow; i ++)
		{
			for(int j =0;j < nCol && j < nRcsCol;j ++)
			{
				Rut.pMATRIX[i][j]	= this->pMATRIX[i][j];
			}
		}
		*this  = Rut;
	}
	
	
	return *this;
}



//***************************************************************************




















⌨️ 快捷键说明

复制代码Ctrl + C
搜索代码Ctrl + F
全屏模式F11
增大字号Ctrl + =
减小字号Ctrl + -
显示快捷键?