gauss.cpp

本程序实现了列主元高斯消去法

#include<math.h>
#include<stdio.h>
#include<iostream.h>
#define M 10 /*系数矩阵的规模*/
void main()
{
/*******************************以下是全局变量定义区***************************************************/
	double a[M][M]={{0}},
		b[M]={0},x[M]={0.0},
		y[M]={0.0},com,cen;/*com是主元素.为什么无法设为DOUBLE??*/
/*******************************以上是全局变量定义区***************************************************/

/*******************************以下是局部变量定义区***************************************************/
	int i,j,k,
		m=100,n,
		key=0,
		p,q,
		xor[M],
		restart=1,
		min; /*i、m是方程数和行数，j、n是未知数数和列数，
			p,q是主元素的行、列下标，key是重新输入开关和中间变量*/
	for(i=0;i<M;i++) 
		xor[i]=i;
/*********************************以上是局部变量定义区***************************************************/

/*******************************以下是变量输入区*********************************************************/

	while(restart==1)
	{
		key=0;
		while(key==0)
		{
			printf("\nProgram start!\n\n\nEnter the number of x n,take care that 1<=n<=%d:",M);
			scanf("%d",&n);
			printf("How many equations in the question?Enter it:");
			scanf("%d",&m);
			min=(m>n?n:m);
			/*这里规定系数矩阵为方阵*/
			for(i=0;i<m;i++)
			{
				printf("Input %d quotions of x[i] in No.%2d equation:\n",n,i+1);
				for(j=0;j<n;j++)
					scanf("%f",&a[i][j]);
				printf("Input b[%2d]:",i+1);
				scanf("%f",&b[i]);
			}
			printf("Enter complete!\nYour equations' quotion are:\n****************************************************************\n");
			for(j=0;j<n;j++)
			printf(" x%2d ",j+1);
			printf("= b \n");
			for(i=0;i<m;i++)
			{
				for(j=0;j<n;j++)
				printf("%9.5f ",a[i][j]);
				printf(" %9.5f\n",b[i]);
			}
			printf("****************************************************************\nAre they correct?Yes enter 1;No enter 0 to input again:");
			scanf("%d",&key);
			restart=0; 
		}
/*********************************以上是变量输入区***************************************************/


/*********************************以下是高斯消去区***************************************************/
/*先找主元素，再消去*/
		for(k=0;k<min;k++)
		{
			com=a[k][k];p=k;q=k;
			for(i=k;i<m;i++)
				for(j=k;j<n;j++)
					if(a[i][j]>com)
					{
						com=a[i][j];
						p=i;
						q=j;
					} 
			/*主元素是com，它的行下标和列下标是p q*/
			for(j=0;j<n;j++)
			{
				cen=a[k][j];
				a[k][j]=a[p][j];
				a[p][j]=cen;
			}
			cen=b[k];
			b[k]=b[p];
			b[p]=cen; 
			for(i=0;i<m;i++)
			{
				cen=a[i][k];
				a[i][k]=a[i][q];
				a[i][q]=cen;
			}
			key=xor[k];
			xor[k]=xor[q];
			xor[q]=key;
			/*主元选取完毕。XOR已经交换。以下是高斯消去*/
			for(i=k+1;(i<m&&com!=0);i++)
			{
				cen=-a[i][k]/a[k][k];
				for(j=k;j<n;j++)
					a[i][j]=a[i][j]+cen*a[k][j];
				b[i]=b[i]+cen*b[k];
			}
			/*高斯消去完成！*/ 
		}
		printf("Enter complete!\nYour equations' quotion are:\n");
		for(j=0;j<n;j++)
			printf(" x%2d ",j+1);
		printf(" b \n");
		for(i=0;i<m;i++)
		{
			for(j=0;j<n;j++)
				printf("%9.5f ",a[i][j]);
			printf(" %9.5f\n",b[i]);
		}/*检验完成！*/


/*********************************以上是高斯消去区***************************************************/ 
/*********************************以下是结果检验区***************************************************/ 
		key=1;
		for(i=0;i<n;i++)
		{
			cen=0.0;
			for(j=0;j<n;j++)
				cen=cen+fabs(a[i][j]);
			if(fabs(cen)<0.00001&&fabs(b[i])<0.00001)
			{
				key=0;
				break;
			}
			else if(fabs(cen)<0.00001&&fabs(b[i])>0.00001)
			{
				key=-1;
				break; 
			}
		} 
		if(key==0)
		{
			printf("\nThe equations have infinit answer groups because some of the equations you entered are invalid.\nDo you want to enter the equations again?\nEnter 1 to enter again;enter 0 to exit:");
			scanf("%d",&restart);
		}
		else if(key==-1)
		{
			printf("\nThe equations have no answer group because some of the equations you entered are contradictory.\nDo you want to enter the equations again?\nEnter 1 to enter again;enter 0 to exit:");
			scanf("%d",&restart);
		}
	}
	if(key==1)
	{


/*********************************以上是结果检验区***************************************************/ 
/*********************************以下是回代求解区***************************************************/ 
		for(k=n-1;k>=0;k--)
		{
			cen=0;
			for(j=m-1;j>k;j--)
				cen=cen+y[j]*a[k][j];
			y[k]=(b[k]-cen)/a[k][k];
		}
		for(i=0;i<n;i++)
			x[xor[i]]=y[i];


/*********************************以上是回代求解区***************************************************/ 
/*********************************以下是结果输出区***************************************************/ 
		for(i=0;i<n;i++)
			printf("x%2d=%9.5f\n",i+1,x[i]);
		restart=0;key=0;
/*********************************以上是结果输出区***************************************************/ 
	}
}