gauss2.c

列主元高斯消去：解方程组用

#define N 3
#include <stdio.h>
#include<math.h>
void main()
{
int i,j,k,t=0;
double a[N][N]={{1,0,0},{0,0,0},{0,0,1}},b[N]={1,1,1},m[N][N],s1,s2,s;

/*第一步：按列选主元素*/
for(k=0;k<N-1;k++)
{
	for(i=k;i<N;i++)
		if(fabs(a[t][k])<fabs(a[i][k]))
			t=i;
/*第二步：检验是否奇异，并输出增广矩阵*/
	if(a[t][k]==0)
	{
		printf("A奇异\n");
		  for(i=0;i<N;i++)
	{
	    for(j=0;j<N;j++)
		printf("%4f  ",a[i][j]);
	printf("\n");
	}

    for(i=0;i<N;i++)
	{
        printf("b[%d]=%f",i,b[i]);
        printf("\n");
	}
		return;
	}
/*第三步：换行*/
	else if(t!=k)
	{
		for(j=k;j<N;j++)
		{
			s1=a[k][j];
			a[k][j]=a[t][j];
			a[t][j]=s1;
		}
			s2=b[k];
            b[k]=b[t];
			b[t]=s2;
	}
/*第四步：计算乘子*/
    for(i=k+1;i<N;i++)
		m[i][k]=a[i][k]/a[k][k];
/*第五步：消元计算*/
	for(i=k+1;i<N;i++)
	{
		for(j=k;j<N;j++)
			a[i][j]=a[i][j]-m[i][k]*a[k][j];
		b[i]=b[i]-m[i][k]*b[k];
	}
}
	/*输出增广矩阵*/
    for(i=0;i<N;i++)
	{
	    for(j=0;j<N;j++)
		printf("%4f  ",a[i][j]);
	printf("\n");
	}

    for(i=0;i<N;i++)
	{
        printf("b[%d]=%f",i,b[i]);
        printf("\n");
	}
/*第六步：回代求解*/    
for(i=N-1;i>=0;i--)
{
	if(i==N-1)
	b[i]=b[i]/a[i][i];
	else
	{
		for(j=i+1,s=0;j<N;j++)
			s=s+a[i][j]*b[j];
		b[i]=(b[i]-s)/a[i][i];
	}
}
/*输出结果*/
printf("\n");
for(i=0;i<N;i++)
{
printf("x[%d]=%f",i,b[i]);
printf("\n");
}
}