高斯消去法

#include "stdio.h"
#include "math.h"
#include "stdlib.h"
typedef double arrtype;
//高斯消去法
void GuassScale(arrtype coff[],int n,arrtype *x = NULL,int m = 1)
{
int col = n + m;
arrtype *s = new arrtype[n];
int *nrow = new int[n];
arrtype *base = coff;
for(int i = 0; i<n;i++)
{
  s[i] = 0;
  for(int j = 0;j<n;j++)
  {
   arrtype t = fabs(coff[i*col + j]);
   if(s[i] < t)
    s[i] = t;
  }
  nrow[i] = i;
}
for(i=0;i<n-1;i++)
{
  int p = i;
  arrtype t = fabs(coff[nrow[p]*col+i])/s[nrow[p]];
  for(int j= i+1;j<n;j++)
  {
   arrtype temp= fabs(coff[nrow[j]*col+i])/s[nrow[j]];
   if(t<temp)
   {
    t=temp;
    p=j;
   }    
  }
  if(coff[nrow[p]*col+i]==0)
   return;
  if(nrow[p]!=nrow[i])
  {
   int ncopy = nrow[i];
   nrow[i] = nrow[p];
   nrow[p] =ncopy;
  }
  for(j= i+1;j<n;j++)
  {
   arrtype m = coff[nrow[j]*col+i]/coff[nrow[i]*col+i];
   for(int k=i+1;k<col;k++)
   {
    coff[nrow[j]*col+k]-=coff[nrow[i]*col+k]*m;
   }
  }
}
if(coff[nrow[n-1]*col+n-1]==0)
  return;
for(i=n;i<col;i++)
{
  x[(n-1)*m+(i-n)]=coff[nrow[n-1]*col+i]/coff[nrow[n-1]*col+n-1];
}
for(i=n-2;i>=0;i--)
{
  for(int j = n;j<col;j++)
  {
   for(int k=n-1;k>i;k--)
    coff[nrow[i]*col+j]-=coff[nrow[i]*col+k]*x[k*m+(j-n)];
   x[i*m+(j-n)]=coff[nrow[i]*col+j]/coff[nrow[i]*col+i];
  }
}
}
void main(int argc,char* argv[])
{
FILE *fpi,*fpo;
if(argc > 1)
  fpi=fopen(argv[1],"r");
else
  fpi=stdin;
if(!fpi)return;
if(argc>2)
  fpo = fopen(argv[2],"w");
else
  fpo = stdout;
if(!fpo)return;
int m,n;
printf("请输入ｎ和ｍ的值：\n");
fscanf(fpi,"%d",&n);
fscanf(fpi,"%d",&m);
arrtype *co;
co = new arrtype[n*(n+m)];
for(int i = 0;i<n*(n+m);i++)
{
  float t;
  fscanf(fpi,"%f",&t);
  co[i]=t;
}
fclose(fpi);
arrtype *x = new arrtype[m*n];
GuassScale(co,n,x,m);
delete[] co;
for(i=0;i<n;i++)
{
  fprintf(fpo,"x[%d]=",i);
  for(int j=0;j<m;j++)
   fprintf(fpo,"%f",x[i*m+j]);
  fprintf(fpo,"\n");
}
delete[] x;
fclose(fpo);
}