Determinant

/*
行列式を計算するプログラム
*/

#include 
#include 
#define N 4
double mdeterm(double a[][N],int n);
void mcopy(double a[][N],double b[][N],int m,int n);
double mdet_row_trans(double b[][N],int n,int i,double d);

main()
{
 double a[][N]={{1,1,0,2},{3,-1,1,0},{0,2,2,1},{2,0,1,3}},d;
 d=mdeterm(a,N);
 printf("|A|=%10.4f\n",d);
}

double mdeterm(double a[][N],int n)
{
 double b[N][N],d;
 int i;
 d=1;
 mcopy(a,b,n,n);
 for(i=0;i<=n-1;i++){
   d=mdet_row_trans(b,n,i,d);
   if(d==0)break;
   }
 return d;
}

void mcopy(double a[][N],double b[][N],int m,int n)
{
  int i,j;
  for(i=0;i<=m-1;i++)
    for(j=0;j<=n-1;j++)b[i][j]=a[i][j];
}

double mdet_row_trans(double b[][N],int n,int i,double d)
{
 int i1,j,k,change;
 double b1,b2;
 i1=i-1;
 do{
   i1++;
   b1=b[i1][i];
 }
 while(b1==0 && i1＜n-1);
 if(b1==0)d=0;
 else{
   change=1;
   if(i!=i1){
     change=-1;
     for(j=i;j＜=n-1;j++){
       b2=b[i][j];
       b[i][j]=b[i1][j];
       b[i1][j]=b2;}
     }

   b1=b[i][i];
   d=b1*d*change;
   for(j=i;j＜=n-1;j++) b[i][j]=b[i][j]/b1;
   for(j=i+1;j＜=n-1;j++){
     b1=b[j][i];
     for(k=i;k＜=n-1;k++)b[j][k]=b[j][k]-b1*b[i][k];
   }}
 return d;
}

[出力結果]

|A|=  -24.0000

Back to C Language