realmatrix.cpp

这是一个实数矩阵计算的源代码

		while(pData1->nRow<n)
			pData1=pData1->pNext;
	}
	else if(pData1->nRow>n||pData1->nRow==0||pData1==pTail)
	{
		while(pData1->nRow>n||pData1->nRow==0||pData1==pTail)
			pData1=pData1->pPrev;
	}

	pData2=pData1;

	if(pData2->nRow<n1)
	{
		while(pData2->nRow<n1)
			pData2=pData2->pNext;
	}

	//换行
	pData1=pData1->pRight;
	while(pData1!=NULL)
	{
		NewMatrix.Store(n1,pData1->nCol,pData1->data);
		pData1=pData1->pRight;
	}

	pData2=pData2->pRight;
	while(pData2!=NULL)
	{
		NewMatrix.Store(n,pData2->nCol,pData2->data);
		pData2=pData2->pRight;
	}
	pCurrent=pHead->pNext;
	for(int i=1;i<=pHead->nRow;++i)
	{
		if(i==n||i==n1)
		{
			pCurrent=pCurrent->pNext;
			continue;
		}
		pData=pCurrent;
		pData=pData->pRight;
		while(pData!=NULL)
		{
			NewMatrix.Store(i,pData->nCol,pData->data);
			pData=pData->pRight;
		}
		pCurrent=pCurrent->pNext;
	}
	return NewMatrix;
}

CRealMatrix CRealMatrix::ColExchange(int col1,int col2)
{
	CRealMatrix NewMatrix(pHead->nRow,pHead->nCol);
	if(pHead==NULL||col1<0||col1>pHead->nCol
		||col2<0||col2>pHead->nCol)

		return NewMatrix;

	int n=col1<col2?col1:col2;
	int n1=col1>col2?col1:col2;

	CNode* pData1=new CNode;
	CNode* pData2=new CNode;

	pData1=pCurrent;
	if(pData1==pHead)
	{
		for(int i=0;i<n;++i)
			pData1=pData1->pNext;
	}
	if(pData1->nCol<n)
	{
		while(pData1->nCol<n)
			pData1=pData1->pNext;
	}
	if(pData1->nCol>n||pData1->nCol==0||pData1==pTail)
	{
		while(pData1->nCol>n||pData1->nCol==0||pData1==pTail)
			pData1=pData1->pPrev;
	}

	pData2=pData1;

	if(pData2->nCol<n1)
	{
		while(pData2->nCol<n1)
			pData2=pData2->pNext;
	}

	//换行
	pData1=pData1->pDown;
	while(pData1!=NULL)
	{
		NewMatrix.Store(pData1->nRow,n1,pData1->data);
		pData1=pData1->pDown;
	}

	pData2=pData2->pDown;
	while(pData2!=NULL)
	{
		NewMatrix.Store(pData2->nRow,n,pData2->data);
		pData2=pData2->pDown;
	}
	pCurrent=pHead->pNext;
	for(int i=1;i<=pHead->nCol;++i)
	{
		if(i==n||i==n1)
		{
			pCurrent=pCurrent->pNext;
			continue;
		}
		pData=pCurrent;
		pData=pData->pDown;
		while(pData!=NULL)
		{
			NewMatrix.Store(pData->nRow,i,pData->data);
			pData=pData->pDown;
		}
		pCurrent=pCurrent->pNext;
	}
	return NewMatrix;
}

double CRealMatrix::Get_det()
{
	int lp;//标志
	double m_det=0;
	if(pHead!=NULL&&pHead->nRow==pHead->nCol)
	{
		int n=pHead->nRow;
		CRealMatrix src;
		src.pHead=pHead;
		src.pTail=pTail;
		src.pCurrent=pCurrent;
		src.pData=pData;
		if(n==1)
		{
			GetData(1,1,m_det);
			return m_det;
		}

		//递归算法
		double value;
		for(int i=1;i<=n;++i)
		{
			CRealMatrix matrix(src);
			if((i+1)%2)lp=-1;
			else lp=1;
			GetData(1,i,value);
			if(value==0) continue;
			matrix.DeleteRow(1);
			matrix.DeleteCol(i);
			double det;
			if(matrix.pHead->nRow==1)matrix.GetData(1,1,det);
			else det=matrix.Get_det();
			m_det+=lp*value*det;
		}
		src.pTail=NULL;
		src.pData=NULL;
		src.pCurrent=NULL;
		src.pHead=NULL;
	}
	return m_det;
}

BOOL CRealMatrix::Decompose0()
{
	int n=pHead->nRow;
	if(n==pHead->nCol&&n!=0)
	{
		CNode* psrcCurrent=pHead->pNext;
		CNode* psrcData=psrcCurrent;
		while(psrcCurrent!=pTail)
		{
			//计算Lij
			for(INT j=1;j<=psrcCurrent->nRow-1;++j)
			{
				CNode* pNode=psrcCurrent;
				CNode* pNew=pHead->pNext;
				double valueU=GetDownEle(psrcData->nRow,j);
				if(valueU==0)	continue;
				while(pNew->nRow<j)
					pNew=pNew->pNext;
				while(pNew->pRight!=NULL&&pNew->pRight->nCol<=j)
					pNew=pNew->pRight;
				while(pNode->pRight!=NULL&&pNode->pRight->nCol<=j)
					pNode=pNode->pRight;
				pNode->data/=pNew->data;
				
			}
			pData=psrcCurrent;
			while(pData->pRight!=NULL&&pData->pRight->nCol<=psrcCurrent->nRow)
				pData=pData->pRight;
			//计算Dii
			psrcData=pData;
			double value=GetDiaEle(pData->nRow,pData->nRow);

			if(value==0)	return	1;//对角元为0

			//计算Uij
			for(j=psrcCurrent->nRow+1;j<=n;++j)
				double valueU=GetUpEle(psrcCurrent->nRow,j);
			CNode* pNode=psrcData;
			if(pNode->nCol>psrcCurrent->nRow)
				pNode->data/=value;
			while(pNode->pRight!=NULL)
			{
				pNode=pNode->pRight;
				pNode->data/=value;
			}
			psrcCurrent=psrcCurrent->pNext;
			psrcData=psrcCurrent;
		}//#end while
		return 2;
	}
	return 0; //该矩阵为空或不为方阵
}
double CRealMatrix::GetDownEle(INT row,INT col)
{
	double sum=0;
	double valueD,valueL,valueU;
	CNode* pRow=pHead->pNext;
	CNode* pCol=pHead;
	CNode* pTemp=pHead;

	while(pRow->nRow<row)
		pRow=pRow->pNext;
	for(INT k=1;k<=col-1;++k)
	{
		if(pRow->pRight!=NULL&&pRow->pRight->nCol<row)
		{
			pRow=pRow->pRight;
			valueL=pRow->data;
			k=pRow->nCol;
			if(k>col-1) break;
		}
		else
			break;
		pTemp=pTemp->pNext;
		while(pTemp->nRow<k)
			pTemp=pTemp->pNext;
		pCol=pTemp;
		valueD=GetDiaEle(k,pTemp);
		valueU=GetUpEle(col,pCol);
		sum+=valueL*valueD*valueU;
	}
	CNode* pNode=pHead->pNext;
	while(pNode->nRow<row)
		pNode=pNode->pNext;
	while(pNode->pRight!=NULL&&pNode->pRight->nCol<=col)
		pNode=pNode->pRight;
	if(pNode->nCol==col)
	{
		pNode->data-=sum;
		if(pNode->data==0)
		{
			DeleteNode(row,col);
			return 0;
		}
		return pNode->data;
	}
	Add(row,col,-sum);
	return -sum;
}
double CRealMatrix::GetUpEle(INT row,INT col)
{
	double sum=0;
	double valueD,valueL,valueU;
	CNode* pRow=pHead->pNext;
	CNode* pCol=pHead;
	CNode* pTemp=pHead;

	while(pRow->nRow<row)
		pRow=pRow->pNext;
	for(INT k=1;k<=row-1;++k)
	{
		if(pRow->pRight!=NULL&&pRow->pRight->nCol<row)
		{
			pRow=pRow->pRight;
			valueL=pRow->data;
			k=pRow->nCol;
			if(k>row-1)	break;
		}
		else
			break;
		pTemp=pTemp->pNext;
		while(pTemp->nRow<k)
			pTemp=pTemp->pNext;
		pCol=pTemp;
		valueD=GetDiaEle(k,pTemp);
		valueU=GetUpEle(col,pCol);
		sum+=valueL*valueD*valueU;
	}
	CNode* pNode=pHead->pNext;
	while(pNode->nRow<row)
		pNode=pNode->pNext;
	while(pNode->pRight!=NULL&&pNode->pRight->nCol<=col)
		pNode=pNode->pRight;
	if(pNode->nCol==col)
	{
		pNode->data-=sum;
		if(pNode->data==0)
		{
			DeleteNode(row,col);
			return 0;
		}
		return pNode->data;
	}
	Add(row,col,-sum);
	return -sum;
}
double CRealMatrix::GetDiaEle(INT row,INT col)
{
	double sum=0;
	double valueD,valueL,valueU;
	CNode* pRow=pHead->pNext;
	CNode* pCol=pHead;
	CNode* pTemp=pHead;
	while(pRow->nRow<row)
		pRow=pRow->pNext;

	for(INT k=1;k<=row-1;++k)
	{
		if(pRow->pRight!=NULL&&pRow->pRight->nCol<row)
		{
			pRow=pRow->pRight;
			valueL=pRow->data;
			k=pRow->nCol;
			if(k>row-1)	break;
		}
		else
			break;
		pTemp=pTemp->pNext;
		while(pTemp->nRow<k)
			pTemp=pTemp->pNext;
		pCol=pTemp;
		valueD=GetDiaEle(k,pTemp);
		valueU=GetUpEle(row,pCol);
		sum+=valueL*valueD*valueU;
	}//#end for
	if(pData->nCol==row)
	{
		pData->data-=sum;
		return pData->data;
	}
	Add(row,row,-sum);
	return -sum;
}

double CRealMatrix::GetUpEle(INT k,CNode* pCol)
{
	while(pCol->pRight!=NULL&&pCol->pRight->nCol<=k)
		pCol=pCol->pRight;
	if(pCol->nCol==k)
		return pCol->data;
	return 0;
}

double CRealMatrix::GetDiaEle(INT n,CNode* pTemp)
{
	while(pTemp->pRight!=NULL&&pTemp->pRight->nCol<=n)
		pTemp=pTemp->pRight;
	if(pTemp->nRow==pTemp->nCol)
		return pTemp->data;
	return 0;
}

int CRealMatrix::rank()
{
	int m_rank=0;
	int n=pHead->nRow;
	int m=pHead->nCol;
	if(pHead!=NULL&&pHead->nRow!=0)
	{
		double* Temp=new double[n*m];	//矩阵数组
		
		//矩阵数组赋初值
		//赋零
		for(int i=0;i<n*m;++i)
			Temp[i]=0;


		//倒出矩阵元素
		pCurrent=pHead->pNext;
		while(pCurrent!=pTail)
		{
			pData=pCurrent;
			if(pData->pRight==NULL)
			{
				pCurrent=pCurrent->pNext;
				continue;
			}
			pData=pData->pRight;
			while(pData!=NULL)
			{
				Temp[(pData->nRow-1)*n+pData->nCol-1]=pData->data;
				pData=pData->pRight;
			}
			pCurrent=pCurrent->pNext;
		}

		// 选主元并消元
		for(i=0;i<n;++i)
		{
			int p=i;
			double max=fabs(Temp[i*n+i]);
			for(int k=i+1;k<n;++k)
			{
				if(max<fabs(Temp[k*n+i]))
				{
					max=fabs(Temp[k*n+i]);
					p=k;
				}
			}
			if(max==0) return m_rank;
			//交换两行
			++m_rank;
			if(p!=i)
			{
				for(int j=i;j<m;++j)
				{
					max=Temp[i*n+j];
					Temp[i*n+j]=Temp[p*n+j];
					Temp[p*n+j]=max;
				}
			}

			//消元
			for(int j=i+1;j<m;++j)
			{
				for(k=i+1;k<n;++k)
				{
					Temp[k*n+j]-=Temp[k*n+i]*Temp[i*n+j]/Temp[i*n+i];
				}
			}
		}
	delete[]Temp;
	}
	return m_rank;
}

CRealMatrix CRealMatrix::AddRow(int row,CRealMatrix& src)
{
	CRealMatrix newMatrix;
	if(pHead==NULL&&src.pHead->nCol==0&&src.pHead->nRow!=0)
	{
		newMatrix.Create(1,1);
		return newMatrix;
	}
	if(pHead==NULL&&src.pHead->nCol!=0&&src.pHead->nRow!=0)
	{
		newMatrix.Create(1,src.pHead->nCol);

		//copy data
		CNode* srcpData=new CNode;
		srcpData=src.pHead->pNext;
		newMatrix.pData=newMatrix.pHead->pNext;

		while(srcpData->pRight!=NULL)
		{
			CNode* newpData1=new CNode;
			srcpData=srcpData->pRight;
			newpData1->nRow=srcpData->nRow;
			newpData1->nCol=srcpData->nCol;
			newpData1->data=srcpData->data;
			newpData1->pDown=NULL;
			newpData1->pRight=NULL;
			newMatrix.pData->pRight=newpData1;
			newMatrix.pData=newMatrix.pData->pRight;
		}
		return newMatrix;
	}
	if(pHead!=NULL)
	{
		newMatrix.Create(pHead->nRow+1,pHead->nCol);
		if(row>pHead->nRow) row=pHead->nRow+1;
		if(row<1)row=1;

		//copy matrix data to newMatrix
		pData=pHead->pNext;
		pCurrent=pData;

		CNode* newpData=new CNode;
		CNode* newpCurrent=new CNode;
		newpData=newMatrix.pHead->pNext;
		newpCurrent=newpData;
		
		while(pCurrent->nRow<row)
		{
			while(pData->pRight!=NULL)
			{
				CNode* newpData1=new CNode;
				pData=pData->pRight;
				newpData1->nRow=pData->nRow;
				newpData1->nCol=pData->nCol;
				newpData1->data=pData->data;
				newpData1->pDown=NULL;
				newpData1->pRight=NULL;
				newpData->pRight=newpData1;
				newpData=newpData->pRight;


				CNode *pTemp=newMatrix.pHead;
				for(int i=0;i<newpData->nCol;++i)
					pTemp=pTemp->pNext;
				while(pTemp->pDown!=NULL)
					pTemp=pTemp->pDown;
				pTemp->pDown=newpData;

			}
			pCurrent=pCurrent->pNext;
			pData=pCurrent;
			newpCurrent=newpCurrent->pNext;
			newpData=newpCurrent;
		}

		if(src.pHead!=NULL)
		{
			src.pData=src.pHead->pNext;
			while(src.pData->pRight!=NULL)
			{
				CNode* newpData1=new CNode;
				src.pData=src.pData->pRight;
				newpData1->nRow=newpData->nRow;
				newpData1->nCol=src.pData->nCol;
				newpData1->data=src.pData->data;
				newpData1->pDown=NULL;
				newpData1->pRight=NULL;
				newpData->pRight=newpData1;
				newpData=newpData->pRight;

				CNode *pTemp=newMatrix.pHead;
				for(int i=0;i<newpData->nCol;++i)
					pTemp=pTemp->pNext;
				while(pTemp->pDown!=NULL)
					pTemp=pTemp->pDown;
				pTemp->pDown=newpData;
			}
		}
		newpCurrent=newpCurrent->pNext;
		newpData=newpCurrent;

		while(pCurrent!=pTail)
		{
			while(pData->pRight!=NULL)
			{
				CNode* newpData1=new CNode;
				pData=pData->pRight;
				newpData1->nRow=pData->nRow+1;
				newpData1->nCol=pData->nCol;
				newpData1->data=pData->data;
				newpData1->pDown=NULL;
				newpData1->pRight=NULL;
				newpData->pRight=newpData1;
				newpData=newpData->pRight;

				CNode *pTemp=newMatrix.pHead;
				for(int i=0;i<newpData->nCol;++i)
					pTemp=pTemp->pNext;
				while(pTemp->pDown!=NULL)
					pTemp=pTemp->pDown;
				pTemp->pDown=newpData;
			}
			pCurrent=pCurrent->pNext;
			pData=pCurrent;
			newpCurrent=newpCurrent->pNext;
			newpData=newpCurrent;
		}
	}
	return newMatrix;
}

CRealMatrix CRealMatrix::AddCol(int col,CRealMatrix& src)
{
	CRealMatrix newMatrix;
	if(pHead==NULL&&src.pHead->nCol==0&&src.pHead->nRow!=0)
	{
		newMatrix.Create(1,1);
		return newMatrix;
	}

	if(pHead==NULL&&src.pHead->nCol!=0&&src.pHead->nRow!=0)
	{
		newMatrix.Create(src.pHead->nRow,1);

		//copy data to newMatrix 
		CNode* newpData=new CNode;
		newpData=newMatrix.pHead->pNext;
		src.pData=src.pHead->pNext;
		while(src.pData!=pTail)
		{
			if(src.pData->pRight!=NULL&&src.pData->pRight->nCol==1)
			{
				CNode* pTemp=new CNode;
				pTemp->nRow=src.pData->pRight->nRow;
				pTemp->nCol=src.pData->pRight->nCol;
				pTemp->data=src.pData->pRight->data;
				pTemp->pRight=NULL;
				pTemp->pDown=NULL;
				newpData->pRight=pTemp;

				pTemp=newMatrix.pHead;
				for(int i=0;i<newpData->nCol;++i)
					pTemp=pTemp->pNext;
				while(pTemp->pDown!=NULL)
					pTemp=pTemp->pDown;
				pTemp->pDown=newpData;
			}
			src.pData=src.pData->pNext;
			newpData=newpData->pNext;