C语言描述 伽罗华域矩阵求逆

/*伽罗华域求矩阵的逆的方法：完全参照全选主元高斯约当法，只要改变运算规则即可：
1.乘法操作照常进行，注意乘法的结果无所谓正负，这是因为伽罗华域的加和减是一样的效果，
所以将原算法中的-号改为+号；
2.加法和减法均改成异或操作
*/


#include "stdio.h"
#include


int brinv(int a[],int n);//矩阵求逆
void brmul(int a[],int b[],int m,int n,int k,int c[]);//验证结果，A*A(-1)是否是单位阵
int brank(int a[],int m,int n);
main()
{ int i,j,rnk;
static int a[4][4]={ {1,1,0,0},
{0,0,0,1},
{1,0,1,0},
{0,1,0,1} };
static int b[4][4],c[4][4];
for (i=0; i<=3; i++)
for (j=0; j<=3; j++)
b[i][j]=a[i][j];

rnk=brank(a,4,4);
if (rnk!=4)
{
printf("奇异矩阵，不可求逆！");
exit(1);
}
printf("矩阵是满值矩阵，rank= %d,可以对其求逆\n",rnk);

for (i=0; i<=3; i++)
for (j=0; j<=3; j++)
a[i][j]=b[i][j];

i=brinv(a,4);
if (i!=0)
{ printf("MAT A IS:\n");
for (i=0; i<=3; i++)
{ for (j=0; j<=3; j++)
printf("%d ",b[i][j]);
printf("\n");
}
printf("\n");
printf("MAT A- IS:\n");
for (i=0; i<=3; i++)
{ for (j=0; j<=3; j++)
printf("%d ",a[i][j]);
printf("\n");
}
printf("\n");
printf("MAT AA- IS:\n");
brmul(b,a,4,4,4,c);
for (i=0; i<=3; i++)
{ for (j=0; j<=3; j++)
printf("%d ",c[i][j]);
printf("\n");
}
}
}

int brinv(int a[],int n)
{
int *is,*js,i,j,k,l,u,v;
int d,p;
is=malloc(n*sizeof(int));
js=malloc(n*sizeof(int));
for (k=0; k<=n-1; k++)
{ d=0;
for (i=k; i<=n-1; i++)
for (j=k; j<=n-1; j++)
{ l=i*n+j; p=a[l];
if (p>d) { d=p; is[k]=i; js[k]=j;}
}
if (d+1==1)
{ free(is); free(js); printf("err**not inv\n");
return(0);
}
if (is[k]!=k)
for (j=0; j<=n-1; j++)
{ u=k*n+j; v=is[k]*n+j;
p=a[u]; a[u]=a[v]; a[v]=p;
}
if (js[k]!=k)
for (i=0; i<=n-1; i++)
{ u=i*n+k; v=i*n+js[k];
p=a[u]; a[u]=a[v]; a[v]=p;
}
l=k*n+k;
a[l]=a[l];
for (j=0; j<=n-1; j++)
if (j!=k)
{ u=k*n+j; a[u]=a[u]*a[l];}
for (i=0; i<=n-1; i++)
if (i!=k)
for (j=0; j<=n-1; j++)
if (j!=k)
{ u=i*n+j;
a[u]=a[u]^(a[i*n+k]*a[k*n+j]);
}
for (i=0; i<=n-1; i++)
if (i!=k)
{ u=i*n+k; a[u]=a[u]*a[l];}
}
for (k=n-1; k>=0; k--)
{ if (js[k]!=k)
for (j=0; j<=n-1; j++)
{ u=k*n+j; v=js[k]*n+j;
p=a[u]; a[u]=a[v]; a[v]=p;
}
if (is[k]!=k)
for (i=0; i<=n-1; i++)
{ u=i*n+k; v=i*n+is[k];
p=a[u]; a[u]=a[v]; a[v]=p;
}
}
free(is); free(js);
return(1);
}

void brmul(int a[],int b[],int m,int n,in

t k,int c[])
{ int i,j,l,u;
for (i=0; i<=m-1; i++)
for (j=0; j<=k-1; j++)
{ u=i*k+j; c[u]=0;
for (l=0; l<=n-1; l++)
c[u]=c[u]^(a[i*n+l]*b[l*k+j]);
}
return;
}
int brank(int a[],int m,int n)
{
int i,j,k,nn,is,js,l,ll,u,v;
int q,d;
nn=m;
if (m>=n) nn=n;
k=0;
for (l=0; l<=nn-1; l++)
{ q=0;
for (i=l; i<=m-1; i++)
for (j=l; j<=n-1; j++)
{ ll=i*n+j; d=abs(a[ll]);
if (d>q) { q=d; is=i; js=j;}
}
if (q+1==1) return(k);
k=k+1;
if (is!=l)
{ for (j=l; j<=n-1; j++)
{ u=l*n+j; v=is*n+j;
d=a[u]; a[u]=a[v]; a[v]=d;
}
}
if (js!=l)
{ for (i=l; i<=m-1; i++)
{ u=i*n+js; v=i*n+l;
d=a[u]; a[u]=a[v]; a[v]=d;
}
}
ll=l*n+l;
for (i=l+1; i<=n-1; i++)
{ d=a[i*n+l]/a[ll];
for (j=l+1; j<=n-1; j++)
{ u=i*n+j;
a[u]=a[u]^(d*a[l*n+j]);
}
}
}
return(k);
}