#include "stdio.h"
#include "math.h"
#include "conio.h"

#define N 3
static double aa[N][N] = {{1,2,-1} , {1,-1,5} , {4,1,-2}};
static double bb[N+1] = {3,0,2};

void main()
{
    int i,j;
    double a[N+1][N+1] , b[N+1] ;
 double x[N+1] ;
 double det;
    double gaussl(double a[][N+1] , double b[], double x[]);
    void disp(double a[][N+1] , double b[]);
   
    for( i=1; i <= N; i++ )
    {
        for( j=1; j <= N; j++)
          a[i][j] = aa[i-1][j-1] ; 
         
        b[i] = bb[i-1]; 
    }

 //gaussl(a,b,x);
   
   det = gaussl(a,b,x);
   
   if(det != 0)
    {
        printf("\n方程组的解为:");
       
        for( i=1; i <= N; i++ )
          printf(" x[%d]=%f",i,x[i]);

        printf("\n\n系数矩阵的行列式值=%f",det);
    }
    else printf("\n\n系数矩阵奇异阵，高斯方法失败!:");

    getch();
}

double gaussl(double a[][N+1] , double b[], double x[])
{
  
 double max1 , max2;  /* max 主元 */
 double p,q;
 double m1,m2;
 int i , j , k=1 ;
 int t;

    printf("消元前增广矩阵:\n");
    disp( a ,  b);

   for( i=1; i <= N; i++ )
 {
        max1 = a[1][1];
  if( a[i][1] > max1 )
  {
   max1 = a[i][1];
   t = i;
  }
   }
   for( i=1; i <= N; i++ )
   {
  p = a[1][i] ;
  a[1][i] = a[t][i] ;
        a[t][i] = p;
   }

   q = b[1];
   b[1] = b[t];
   b[t] = q;
 

    printf("交换第1,%d行\n",t);  /* t hang */
    disp( a ,  b);

 m1 = a[t][1]/max1;  /* xiang chu wei yi ,yin wei max1 = a[1][1] */

    for( i=1; i <= N; i++ )
 {
  
  for( j=1; j <= N; j++ )
  {
   if(i != 1)
   {
    a[i][j] = a[i][j] - (a[1][j]*m1);
   }
  }
  if(i != 1)
  {
   b[i] = b[i] - (b[1]*m1);
  }
 }
    printf("第%d次消元:\n",k);
    disp( a ,  b);

    for( i=k+1; i <= N; i++ )
 {
  max2 = a[k+1][k+1];
  if( a[i][k+1] > max2 )
  {
   max2 = a[i][k+1];
   t = i;
  }
 }
 for( i=k+1; i <= N; i++ )
 {
  p = a[k+1][i] ;
  a[k+1][i] = a[t][i] ;
        a[t][i] = p;
 }

  q = b[k+1];
  b[k+1] = b[t];
  b[t] = q;           /* jiang b shu zu yong xun huan jiao huan ,que wei yong dao bian liang i */
 
    printf("交换第%d,%d行\n",k+1,t);
    disp( a ,  b);

    m2 = a[t][k+1]/max2;

    for( i=k+1; i <= N; i++ )
 {
  
  for( j=k+1; j <= N; j++ )
  {
   if(i != k+1)
   {
    a[i][j] = a[i][j] - (a[k+1][j]*m2);
   }
  }
  if(i != k+1)
  {
   b[i] = b[i] - (b[k+1]*m2);
  }
 }
    printf("第%d次消元:\n",k+1);
    disp( a ,  b);

  x[3] = b[3] / a[3][3] ;  /* ci chu ke yong for xun huan shu zu lai gai gai zui hao */
 x[2] = (b[2] - x[3]*a[2][3])/a[2][2];
 x[1] = (b[1] - x[3]*a[1][3] - x[2]*a[1][2])/a[1][1];

 return ( a[1][1]*a[2][2]*a[3][3]);  /* ci chu ke yong for xun huan shu zu lai gai gai zui hao */

}

void disp(double a[][N+1] , double b[]) /*用于显示选主元及消元运算中间增广矩阵结果*/
{ /* 由gaussl()调用  */

 int i,j;

  
   for( i=1; i <= N; i++ )
   {
     for( j=1; j <= N; j++ )
  {
   printf(" ");
   printf("%f",a[i][j]);
  }
     printf(" ");
  printf("%f",b[i]);
  printf("\n");
   }
}