[转帖]矩阵的各种运算_转置,相乘,行列式,求逆

鼯鼠

最近写一个程序，要用过矩阵的各种运算，于是在网上搜索了一下，发现代码还是很多的，比较了一下，发现“程序世界”的百度空间总结的最全，于是使出拿来主义，根据“程序世界”百度空间上的网友的评论，在求逆函数中增加了奇异处理的代码，稍微完善了一下。

程序代码：

//hbl编写于2008.5.8
//编译环境Dev-C++4.9.9.2 + xpSP2
//转置,相乘,求行列式,求逆函数的绝大部分代码来自"程序世界"的百度空间:
//http://hi.baidu.com/vcmfc/blog/item/f386632ad368f53a5243c169.html,对此表示感谢


#include <iostream>
//#include <iomanip> //setw函数
using namespace std;

//--显示矩阵，形参m为行，n为列
void MatrixDisplay(double *A,int m,int n)
{
for(int i=0;i<m;i++)
{
for(int j=0;j<n;j++)
cout<<A[i*n + j]<<" ";
cout<<endl;
}
cout<<endl;
}

//--求矩阵转置，形参m为行，n为列,A转置后存为B
//-------------------------------------------------
void MatrixInverse(double *A,double *B, int m,int n)
{
for (int i=0;i<n;i++)
{
for (int j=0;j<m;j++)
B[i * m + j] = A[j*n + i];
}
}
//-------------------------------------------------

//--求矩阵相乘,A矩阵为[m,p],B矩阵为[p,n],C为[m,n]
//-------------------------------------------------------------------
void MatrixMultiply(double *A,double *B,double *C ,int m,int p,int n)
{
for (int i = 0; i < m; i++)
{
for (int j = 0; j < n; j++)
{
double sum = 0;
for (int k = 0; k < p; k++)
{
sum = sum + A[i * p + k] * B[k * n + j];
}
C[i*n +j] = sum;
}
}
}
//-------------------------------------------------------------------

//--求矩阵行列式
//-------------------------------------------------
double MatrixDet(double *A,int m)
{
int i = 0, ii = 0, j = 0, jj = 0, k = 0, t = 0, tt = 1;
double det = 1, mk = 0;

double *pA = new double [m*m];
double *pB = new double [m];
for (i = 0; i < m; i++)
{
pB = 0;
for (j = 0; j < m; j++)
{
pA[i*m + j] = A[i*m +j];
}
}

for (k = 0; k < m; k++)
{
for (j = k; j < m; j++)
{
if (pA[k*m + j])
{
for (i = 0; i < m; i++)
{
pB = pA[i*m + k];
pA[i*m + k] = pA[i*m + j];
pA[i*m + j] = pB;
}

if (j-k)
{
tt = tt * (-1);
}

t = t + 1;

break ;
}

}

if ( t )
{
for (ii = k + 1; ii < m; ii++)
{
mk = (-1)*pA[ii*m + k]/pA[k*m + k];
pA[ii*m + k] = 0;

for (jj = k + 1; jj < m ;jj++)
{
pA[ii*m +jj] = pA[ii*m + jj] + mk*pA[k*m + jj];
}
}

det = det * pA[k*m + k];
t = 0;
}
else
{
det = 0;
break ;
}

}

det = det*tt;

delete pA;
delete pB;
pA = NULL;
pB = NULL;

return det;
}
//-------------------------------------------------

//--求矩阵逆,求A的逆矩阵C
//-------------------------------------------------
void MatrixInv(double *A,double *C, int m)
{
int i,j,x,y;
double X = 0;

double *SP=new double [m*m];
double *AB=new double [m*m];
double *B=new double [m*m];

X=MatrixDet(A,m);
//奇异处理
if(X+1.0==1.0)
{
cout<<"Err!The Matrix is singular!"<<endl;
system("pause");
exit(1);
}
X=1/X;
for(i=0;i<m;i++)
{
for(j=0;j<m;j++)
{
for(x = 0;x < m; x++)
{
for (y = 0; y < m; y++)
{
B[x*m + y]=A[x*m + y];
}
}
for(x=0;x<m;x++)
B[x*m+j]=0;
for(y=0;y<m;y++)
B[i*m+y]=0;

B[i*m+j]=1;
SP[i*m+j]=MatrixDet(B,m);
SP[i*m+j] = SP[i*m+j];
AB[i*m+j]=X*SP[i*m+j];
}
}
MatrixInverse(AB,C,m,m);

delete SP;
delete AB;
delete B;
SP = NULL;
AB = NULL;
B = NULL;

}
//-------------------------------------------------

int main()
{
cout<<"以下运算结果是正确的,已通过matlab验证"<<endl;
//验证矩阵转置
double testA[4][3]={1,2,3,4,5,6,7,8,9,10,11,12};
double testB[3][4];
MatrixInverse(&testA[0][0],&testB[0][0],4,3);
cout<<"原矩阵A:"<<endl;
MatrixDisplay(&testA[0][0],4,3);
cout<<"A转置后的矩阵B:"<<endl;
MatrixDisplay(&testB[0][0],3,4);

//验证矩阵相乘
double testC[4][4];
MatrixMultiply(&testA[0][0],&testB[0][0],&testC[0][0],4,3,4);
cout<<"A,B相乘后的矩阵C:"<<endl;
MatrixDisplay(&testC[0][0],4,4);

//验证矩阵求逆
double testD[3][3]={12,52,3,4,5,6,7,8,9};
double testE[3][3];
MatrixInv(&testD[0][0],&testE[0][0],3);
cout<<"原矩阵D:"<<endl;
MatrixDisplay(&testD[0][0],3,3);
cout<<"D求逆后的矩阵E:"<<endl;
MatrixDisplay(&testE[0][0],3,3);

//奇异处理的情况
double testF[4][4];
cout<<"C求逆后的结果:"<<endl;
MatrixInv(&testC[0][0],&testF[0][0],4);
system("pause");
}

运行结果：

---------------

注：如果你觉得我的表述有误或还不够好的话，请通过此邮箱与我联系：2008.hbl@gmail.com，提出您宝贵的意见，我们一起来完善它。您的一言一语，都是对本人的莫大支持。

--------------