逆矩阵法求解矩阵方程AX=B

C程序

/*************************************************************************
* 逆矩阵法求解矩阵方程AX=B，其中A是N*N的矩阵,B是N*N矩阵
* 输入: n----方阵A的行数
*       a----矩阵A
*       m----矩阵B的列数
*       b----矩阵B
* 输出: det----矩阵A的行列式值
*       a----A的逆矩阵
*       b----矩阵方程的解X
**************************************************************************/
double gaussian_jodan_solve(int n,double a[N][N],int m,double b[N][M])
{
double det,f[N];
int i,j,k;

det = gaussian_jodan(n,a);
if(det==0) return (0);
for(k=0;k {
for(i=0;i {
f=0.0;
for(j=0;j f=f+a[j]*b[j][k];
}
for(i=0;i b[k]=f;
}
return(det);
}

调用到的求逆矩阵的子函数

/*************************************************************************
* 高斯--约当列主元素法求矩阵方程A的逆矩阵，其中A是N*N的矩阵
* 输入: n----方阵A的行数
*       a----矩阵A
* 输出: det--A的行列式的值
*       a----A的逆矩阵
**************************************************************************/
double gaussian_jodan(int n,double a[N][N])
{
int i,j,k,mk;
int p[N];  /*记录主行元素在原矩阵中的位置*/
double det,m,f;

det = 1.0;
for(k=0;k {
m=a[k][k];  /*选第K列主元素*/
mk=k;
for(i=k+1;i if(fabs(m)[k]))
{
m=a[k];
mk=i;
}
if(fabs(m) if(mk!=k)
{
for(j=0;j {
f=a[k][j];
a[k][j]=a[mk][j];
a[mk][j]=f;
}
p[k]=mk;
det = -det;
}
else
p[k]=k;
det=det*m;
for(j=0;j if(j!=k)
a[k][j]=a[k][j]/a[k][k];
a[k][k]=1.0/a[k][k];
for(i=0;i {
if(i!=k)
{
for(j=0;j if(j!=k)
a[j]=a[j]-a[k]*a[k][j];
a[k]=-a[k]*a[k][k];
}
}
}
for(k=n-2;k>=0;k--) /*按主行在原矩阵中的位置交换列*/
{
if(p[k]!=k)
for(i=0;i {
f=a[k];
a[k]=a[p[k]];
a[p[k]]=f;
}
}
return(det);
}