求初等矩阵的逆矩阵阵,详细步骤

下次自动登录
现在的位置:
& 综合 & 正文
如何在C#去求矩阵的逆矩阵
,例如先求矩阵行列式的值,然后再求每一项的代数余子式,然后按照矩阵逆的公式去计算。但是等他向我求代码的时候,发现做法并不是那么简单,甚至用上面的思路,很难做出来。于是我参考网上求行列式值的算法,在上面的基础上完成了矩阵逆的算法。
[,] ReverseMatrix( double[,] dMatrix, int Level )
double dMatrixValue = MatrixValue( dMatrix, Level );
if( dMatrixValue == 0 ) return null;
double[,] dReverseMatrix = new double[Level,2*Level];
// Init Reverse matrix
for( int i = 0; i & L i++ )
for( int j = 0; j & 2 * L j++ )
if( j & Level )
dReverseMatrix[i,j] = dMatrix[i,j];
dReverseMatrix[i,j] = 0;
dReverseMatrix[i,Level + i ] = 1;
for( int i = 0, j = 0; i & Level && j & L i++, j++ )
if( dReverseMatrix[i,j] == 0 )
for( ; dMatrix[m,j] == 0; m++ );
if( m == Level )
return null;
// Add i-row with m-row
for( int n = n & 2 * L n++ )
dReverseMatrix[i,n] += dReverseMatrix[m,n];
// Format the i-row with "1" start
x = dReverseMatrix[i,j];
if( x != 1 )
for( int n = n & 2 * L n++ )
if( dReverseMatrix[i,n] != 0 )
dReverseMatrix[i,n] /=
// Set 0 to the current column in the rows after current row
for( int s = Level - 1; s &s-- )
x = dReverseMatrix[s,j];
for( int t = t & 2 * L t++ )
dReverseMatrix[s,t] -= ( dReverseMatrix[i,t]* x );
// Format the first matrix into unit-matrix
for( int i = Level - 2; i &= 0; i-- )
for( int j = i + 1 ; j & L j++ )
if( dReverseMatrix[i,j] != 0 )
c = dReverseMatrix[i,j];
for( int n = n & 2*L n++ )
dReverseMatrix[i,n] -= ( c * dReverseMatrix[j,n] );
double[,] dReturn = new double[Level, Level];
for( int i = 0; i & L i++ )
for( int j = 0; j & L j++ )
dReturn[i,j] = dReverseMatrix[i,j+Level];
MatrixValue( double[,] MatrixList, int Level )
double[,] dMatrix = new double[Level, Level];
for( int i = 0; i & L i++ )
for( int j = 0; j & L j++ )
dMatrix[i,j] = MatrixList[i,j];
int k = 1;
for( int i = 0, j = 0; i & Level && j & L i++, j++ )
if( dMatrix[i,j] == 0 )
for( ; dMatrix[m,j] == 0; m++ );
if( m == Level )
// Row change between i-row and m-row
for( int n = n & L n++ )
c = dMatrix[i,n];
dMatrix[i,n] = dMatrix[m,n];
dMatrix[m,n] =
// Change value pre-value
k *= (-1);
// Set 0 to the current column in the rows after current row
for( int s = Level - 1; s &s-- )
x = dMatrix[s,j];
for( int t = t & L t++ )
dMatrix[s,t] -= dMatrix[i,t]* ( x/dMatrix[i,j] );
double sn = 1;
for( int i = 0; i & L i++ )
if( dMatrix[i,i] != 0 )
sn *= dMatrix[i,i];
double[,] dMatrix = new double[3,3]{{0,1,2},{1,0,1},{4,2,1}};
double[,] dReturn = ReverseMatrix( dMatrix, 3 );
if( dReturn != null )
for( int i=0; i & 3; i++ )
Debug.WriteLine( string.Format( "{0} {1} {2}",
dReturn[i,0], dReturn[i,1],dReturn[i,2] ) );
【上篇】【下篇】码农谷 & 版权所有 (C)
& 湘ICP备号-1扫二维码下载作业帮
3亿+用户的选择
下载作业帮安装包
扫二维码下载作业帮
3亿+用户的选择
求3*3矩阵的逆矩阵A= 1
9最好有详细点的解题步骤。
作业帮用户
扫二维码下载作业帮
3亿+用户的选择
为您推荐:
扫描下载二维码求逆矩阵的几种常用方法总结_百度文库
您的浏览器Javascript被禁用,需开启后体验完整功能,
享专业文档下载特权
&赠共享文档下载特权
&100W篇文档免费专享
&每天抽奖多种福利
两大类热门资源免费畅读
续费一年阅读会员,立省24元!
求逆矩阵的几种常用方法总结
阅读已结束,下载本文需要
定制HR最喜欢的简历
你可能喜欢矩阵求逆运算有多种算法:
伴随矩阵的思想,分别算出其伴随矩阵和行列式,再算出逆矩阵;
LU分解法(若选主元即为LUP分解法:
,每步重新选主元),它有两种不同的实现;
A-1=(LU)-1=U-1L-1,将A分解为LU后,对L和U分别求逆,再相乘;
通过解线程方程组Ax=b的方式求逆矩阵。b分别取单位阵的各个列向量,所得到的解向量x就是逆矩阵的各个列向量,拼成逆矩阵即可。
下面是这两种方法的c++代码实现,所有代码均利用常规数据集验证过。
文内程序旨在实现求逆运算核心思想,某些异常检测的功能就未实现(如矩阵维数检测、矩阵奇异等)。
注意:文中A阵均为方阵。
伴随矩阵法C++程序:
1 #include &iostream&
2 #include &ctime&
//用于产生随机数据的种子
4 #define N 3
//测试矩阵维数定义
6 //按第一行展开计算|A|
7 double getA(double arcs[N][N],int n)
return arcs[0][0];
double ans = 0;
double temp[N][N]={0.0};
int i,j,k;
for(i=0;i&n;i++)
for(j=0;j&n-1;j++)
for(k=0;k&n-1;k++)
temp[j][k] = arcs[j+1][(k&=i)?k+1:k];
double t = getA(temp,n-1);
if(i%2==0)
ans += arcs[0][i]*t;
arcs[0][i]*t;
39 //计算每一行每一列的每个元素所对应的余子式,组成A*
getAStart(double arcs[N][N],int n,double ans[N][N])
ans[0][0] = 1;
int i,j,k,t;
double temp[N][N];
for(i=0;i&n;i++)
for(j=0;j&n;j++)
for(k=0;k&n-1;k++)
for(t=0;t&n-1;t++)
temp[k][t] = arcs[k&=i?k+1:k][t&=j?t+1:t];
getA(temp,n-1);
//此处顺便进行了转置
if((i+j)%2 == 1)
ans[j][i] = - ans[j][i];
71 //得到给定矩阵src的逆矩阵保存到des中。
72 bool GetMatrixInverse(double src[N][N],int n,double des[N][N])
double flag=getA(src,n);
double t[N][N];
if(0==flag)
cout&& "原矩阵行列式为0,无法求逆。请重新运行" &&
return false;//如果算出矩阵的行列式为0,则不往下进行
getAStart(src,n,t);
for(int i=0;i&n;i++)
for(int j=0;j&n;j++)
des[i][j]=t[i][j]/
return true;
97 int main()
bool//标志位,如果行列式为0,则结束程序
int row =N;
int col=N;
double matrix_before[N][N]{};//{1,2,3,4,5,6,7,8,9};
//随机数据,可替换
srand((unsigned)time(0));
for(int i=0; i&N ;i++)
for(int j=0; j&N;j++)
matrix_before[i][j]=rand()%100 *0.01;
cout&&"原矩阵:"&&
for(int i=0; i&N ;i++)
for(int j=0; j&N;j++)
//cout && matrix_before[i][j] &&" ";
cout && *(*(matrix_before+i)+j)&&" ";
double matrix_after[N][N]{};
flag=GetMatrixInverse(matrix_before,N,matrix_after);
if(false==flag)
cout&&"逆矩阵:"&&
for(int i=0; i&i++)
for(int j=0; j&j++)
cout &&matrix_after[i][j] &&" ";
//cout && *(*(matrix_after+i)+j)&&" ";
GetMatrixInverse(matrix_after,N,matrix_before);
cout&&"反算的原矩阵:"&&//为了验证程序的精度
for(int i=0; i&N ;i++)
for(int j=0; j&N;j++)
//cout && matrix_before[i][j] &&" ";
cout && *(*(matrix_before+i)+j)&&" ";
LU分解法C++程序:
1 #include &iostream&
2 #include &cmath&
3 #include &ctime&
5 #define N 300
7 //矩阵乘法
8 double * mul(double A[N*N],double B[N*N])
double *C=new double[N*N]{};
for(int i=0;i&N;i++)
for(int j=0;j&N;j++)
for(int k=0;k&N;k++)
C[i*N+j] += A[i*N+k]*B[k*N+j];
//若绝对值小于10^-10,则置为0(这是我自己定的)
for(int i=0;i&N*N;i++)
if(abs(C[i])&pow(10,-10))
34 //LUP分解
35 void LUP_Descomposition(double A[N*N],double L[N*N],double U[N*N],int P[N])
int row=0;
for(int i=0;i&N;i++)
for(int i=0;i&N-1;i++)
double p=0.0d;
for(int j=i;j&N;j++)
if(abs(A[j*N+i])&p)
p=abs(A[j*N+i]);
cout&& "矩阵奇异,无法计算逆" &&
//交换P[i]和P[row]
int tmp=P[i];
P[i]=P[row];
double tmp2=0.0d;
for(int j=0;j&N;j++)
//交换A[i][j]和 A[row][j]
tmp2=A[i*N+j];
A[i*N+j]=A[row*N+j];
A[row*N+j]=tmp2;
//以下同LU分解
double u=A[i*N+i],l=0.0d;
for(int j=i+1;j&N;j++)
l=A[j*N+i]/u;
A[j*N+i]=l;
for(int k=i+1;k&N;k++)
A[j*N+k]=A[j*N+k]-A[i*N+k]*l;
//构造L和U
for(int i=0;i&N;i++)
for(int j=0;j&=i;j++)
L[i*N+j]=A[i*N+j];
L[i*N+j]=1;
for(int k=i;k&N;k++)
U[i*N+k]=A[i*N+k];
109 //LUP求解方程
110 double * LUP_Solve(double L[N*N],double U[N*N],int P[N],double b[N])
double *x=new double[N]();
double *y=new double[N]();
//正向替换
for(int i = 0;i & N;i++)
y[i] = b[P[i]];
for(int j = 0;j &j++)
y[i] = y[i] - L[i*N+j]*y[j];
//反向替换
for(int i = N-1;i &= 0; i--)
x[i]=y[i];
for(int j = N-1;j &j--)
x[i] = x[i] - U[i*N+j]*x[j];
x[i] /= U[i*N+i];
137 /*****************矩阵原地转置BEGIN********************/
139 /* 后继 */
140 int getNext(int i, int m, int n)
return (i%n)*m + i/n;
145 /* 前驱 */
146 int getPre(int i, int m, int n)
return (i%m)*n + i/m;
151 /* 处理以下标i为起点的环 */
152 void movedata(double *mtx, int i, int m, int n)
double temp = mtx[i]; // 暂存
// 当前下标
int pre = getPre(cur, m, n);
while(pre != i)
mtx[cur] = mtx[pre];
pre = getPre(cur, m, n);
mtx[cur] =
166 /* 转置,即循环处理所有环 */
167 void transpose(double *mtx, int m, int n)
for(int i=0; i&m*n; ++i)
int next = getNext(i, m, n);
while(next & i) // 若存在后继小于i说明重复,就不进行下去了(只有不重复时进入while循环)
next = getNext(next, m, n);
if(next == i)
// 处理当前环
movedata(mtx, i, m, n);
178 /*****************矩阵原地转置END********************/
180 //LUP求逆(将每列b求出的各列x进行组装)
181 double * LUP_solve_inverse(double A[N*N])
//创建矩阵A的副本,注意不能直接用A计算,因为LUP分解算法已将其改变
double *A_mirror = new double[N*N]();
double *inv_A=new double[N*N]();//最终的逆矩阵(还需要转置)
double *inv_A_each=new double[N]();//矩阵逆的各列
//double *B
=new double[N*N]();
=new double[N]();//b阵为B阵的列矩阵分量
for(int i=0;i&N;i++)
double *L=new double[N*N]();
double *U=new double[N*N]();
int *P=new int[N]();
//构造单位阵的每一列
for(int i=0;i&N;i++)
//每次都需要重新将A复制一份
for(int i=0;i&N*N;i++)
A_mirror[i]=A[i];
LUP_Descomposition(A_mirror,L,U,P);
inv_A_each=LUP_Solve (L,U,P,b);
memcpy(inv_A+i*N,inv_A_each,N*sizeof(double));//将各列拼接起来
transpose(inv_A,N,N);//由于现在根据每列b算出的x按行存储,因此需转置
return inv_A;
219 int main()
double *A = new double[N*N]();
srand((unsigned)time(0));
for(int i=0; i&N ;i++)
for(int j=0; j&N;j++)
A[i*N+j]=rand()%100 *0.01;
double *E_test = new double[N*N]();
double *invOfA = new double[N*N]();
invOfA=LUP_solve_inverse(A);
E_test=mul(A,invOfA);
//验证精确度
cout&& "矩阵A:" &&
for(int i=0;i&N;i++)
for(int j=0;j&N;j++)
cout&& A[i*N+j]&& " " ;
cout&& "inv_A:" &&
for(int i=0;i&N;i++)
for(int j=0;j&N;j++)
cout&& invOfA[i*N+j]&& " " ;
cout&& "E_test:" &&
for(int i=0;i&N;i++)
for(int j=0;j&N;j++)
cout&& E_test[i*N+j]&& " " ;
两种方法运行时间测试样例(运行环境不同可能会有不同结果,我的主频是2.6GHz,内存8g。时间单位:毫秒ms)
个人认为LU分解法的两个1ms其实是不准确的(实际应远小于1ms,有兴趣可以试试看)。
三种方法的复杂度分析:
伴随矩阵法:此法的时间复杂度主要来源于计算行列式,由于计算行列式的函数为递归形式,其复杂度为O(n2)[参见],而整体算法需要计算每个元素的代数余子式,时间复杂度直接扩大n2倍,变为O(n4)。而递归算法本身是需要占用栈空间的,因此需要注意:当矩阵的维数较大时,随着递归深度的加大,临时占用的空间将会越来越多,甚至可能会出现栈不够用的情况(当然本次实现没有遇到,因为此时的时间开销实在令人难以忍受)!
LU分解法:此法主要是分解过程耗时,求解三角矩阵的时间复杂度是O(n2),分解过程是O(n3),总体来说和高斯消元法差不多,但是避免了高斯消元法的主元素为0的过程。 为了节省空间,A=LU分解的元素存放在A的矩阵中(因为当用过了a[i][j]元素后,便不再用了,所以可以占用原矩阵A的空间)。但是有利就有弊,考虑如果是上千个元素的矩阵,引用传参,这样就改变原矩阵了,因此程序中使用A_mintor作为副本进行使用。另外,可以看出,当矩阵维数超过某值时,内存空间便不够用了(具体是多少没有试验)。还需注意的一点是:程序中未对矩阵是否奇异进行检查,如果矩阵奇异,就不应再进行下去了。
LU分解法中,还可以先分别求出U和L的逆,再相乘,此法其实与常规LU分解法差不多。
文章中用到了矩阵的原地转置算法,具体请参考第4篇文献,这种方法降低了空间复杂度。
需要注意的问题:
本文介绍的方法new了一些指针,未释放,会出现内存泄漏,使用前请释放掉。
本文参考了以下几篇文章:
阅读(...) 评论()}

我要回帖

更多关于 矩阵的逆矩阵怎么求 的文章

更多推荐

版权声明:文章内容来源于网络,版权归原作者所有,如有侵权请点击这里与我们联系,我们将及时删除。

点击添加站长微信