全主元素高斯消去法

 

int  rguass(int n,double a[],double b[]){
 int *js,l,k,i,j,is,p,q;
 double d,t;
 js=(int *)malloc(n*sizeof(int));
 l=1;
 for(k=0;k<=n-2;k++){
  d=0.0;
  for(i=k;i<=n-1;i++)
   for(j=k;j<=n-1;j++){
    t=fabs(a[i*n+j]);
    if(t>d)
    {
     d=t;
     js[k]=j;
     is=i;
    }
   
   }
   //两层for循环结束  ，找到最大元素i行j列，j放到js【k】里

// 例如 [1,2             

//    4,3]

//则4为最大，i=1，j=0

   if(d+1.0==1.0)l=0;
   else{

    if(js[k]!=k)            //交换列
     for(i=0;i<=n-1;i++)
     {
      p=i*n+k;q=i*n+js[k];
      t=a[p];a[p]=a[q];a[q]=t;

     
     }
     if(is!=k){             //交换行

      for(j=k;j<=n-1;j++){
       p=k*n+j;q=is*n+j;
       t=a[p];a[p]=a[q];a[q]=t;
      }
      t=b[k];b[k]=b[is];b[is]=t;
   }
   }
   
  
  if(l==0)
  {
   free(js);
   printf("fail\n");
   return 0;
  }
  d=a[k*n+k];         //行，列交换后，a【k】【k】就是最大
   printf("----d=%f\n",d);
  for(j=k+1;j<=n-1;j++)          //归一化
  {
   p=k*n+j;a[p]=a[p]/d;
  }
  b[k]=b[k]/d;

  printf("a[0]=%f,a[1]=%f  b[0]=%f\n",a[0],a[1],b[0]);
  printf("a[2]=%f,a[3]=%f b[1]=%f\n",a[2],a[3],b[1]);
  for (i=k+1;i<=n-1;i++)           //消元
  {
   for (j=k+1;j<=n-1;j++)
   {
    p=i*n+j;
    a[p]=a[p]-a[i*n+k]*a[k*n+j];
   printf("b[%d]=%f\n",i,b[i]);
   }
   b[i]=b[i]-a[i*n+k]*b[k];
  
  }
  
  
 }          //最外面的k循环结束

 d=a[(n-1)*n+n-1];
 if(fabs(d)+1.0==1.0)
 {
  free(js);
  printf("fail\n");
  return 0;
 }
 b[n-1]=b[n-1]/d;   //解向量的最后一个分量
 
 for(i=n-2;i>=0;i--)//回带
 {
  t=0.0;
  for (j=i+1;j<=n-1;j++)
  {
   t=t+a[i*n+j]*b[j];
   
  }
  b[i]=b[i]-t;
 }

 

 

//恢复解向量
 js[n-1]=n-1;
 for (k=n-1;k>=0;k--)
 {
  if(js[k]!=k)
  {
   
   t=b[k];b[k]=b[js[k]];b[js[k]]=t;
  }
 }
 free(js);
 return 1;
 }

 

 

int _tmain(int argc, TCHAR* argv[], TCHAR* envp[])
{
 int nRetCode = 0;

 // initialize MFC and print and error on failure
 if (!AfxWinInit(::GetModuleHandle(NULL), NULL, ::GetCommandLine(), 0))
 {
  // TODO: change error code to suit your needs
  cerr << _T("Fatal Error: MFC initialization failed") << endl;
  nRetCode = 1;
 }
 else
 {
  // TODO: code your application's behavior here.
  CString strHello;
  strHello.LoadString(IDS_HELLO);

  int i=0;

 

//  [1,2

//    4,3]

//      static double a[2][2]={{1,2},
//      {4,3}
//      };
//     static double b[2]={3,7};
//
  static double a[4][4]={{0.2368,0.2471,0.2568,1.2671},
  {0.1968,0.2071,1.2168,0.2271},
  {0.1581,1.1675,0.1768,0.1871},
  {1.1161,0.1254,0.1397,0.1490}
  };
  static double b[4]={1.8471,1.7471,1.6471,1.5471};
 // rguass(2,*a,b);
  if(rguass(4,*a,b)!=0){
   for (i=0;i<=3;i++)
   {
    cout<<"x["<<i<<"]="<<b[i]<< endl;
   }
  }
  
 }

 return nRetCode;
}