row_principal_element_elimination.c

主元素消去法解线性方程组

/************************************************************************ 
# Row principal element elimination solution system of linear equations
# The N value may change.
*************************************************************************/
#define N 3
#include <stdio.h>
#include <stdlib.h>
double absv(double x){
	double z;
	z = (x >= 0 ? x : -1*x);
	return z;
}
void main(){
	int k,i,j,r;
	double a[N][N],b[N],x[N],t,temp,max,s=0;
	printf ("Please input the coefficient matrix a[N][N]:\n");
	for (i=0; i<N; i++)
		for (j=0; j<N; j++)
			scanf ("%lf",&a[i][j]);
	printf ("Please input the constant vector b[N]:\n");
	for (i=0; i<N; i++)
			scanf ("%lf",&b[i]);
	printf ("\n");
	for (k=0; k<N-1; k++){
		max = absv(a[k][k]);
		r = k;
		for (i=k+1; i<N; i++){
			if (max < absv(a[i][k])){
				max = absv(a[i][k]);
				r = i;
			}
		}
		if (k==r && 0==a[k][k]){
			printf ("Since det(a) == 0.\n");
			printf ("The solution is not unique!\n\n");
			exit(1);
		}
		else{
			if (k==r && a[k][k]!=0){
				for (i=k+1; i<N; i++){
					t = a[i][k] / a[k][k];
					for (j=k; j<N; j++)
						a[i][j] -= a[k][j] * t;
					b[i] -= b[k] * t;
				}
			}
			else{
				for (j=k; j<N; j++){
					temp = a[k][j];
					a[k][j] = a[r][j];
					a[r][j] = temp;
				}
				temp = b[k];
				b[k] = b[r];
				b[r] = temp;
				for (i=k+1; i<N; i++){
					t = a[i][k] / a[k][k];
					for (j=k; j<N; j++)
						a[i][j] -= a[k][j] * t;
					b[i] -= b[k] * t;
				}
			}
		}
	}
	for (i=0; i<N; i++){
		for (j=0; j<N; j++)
			printf ("%lf  ",a[i][j]);
		printf ("%lf\n",b[i]);
	}
	for (i=N-1; i>=0; i--){
		for (j=N-1; j>i; j--)
			s += a[i][j] * x[j];
		x[i] = (b[i] - s) / a[i][i];
		s = 0;
	}
	printf ("\nThe solution of linear equations system is:\n");
	for (i=0; i<N; i++)
		printf ("x%d = %g\n",i+1,x[i]);
	printf ("\n");
}