main.cpp

计算机数值方法中gauss列主元素消去法的演示算法

#include <stdio.h>
#include <math.h>
#include <iostream>
using namespace std;

/**************************************************

函数名：	findMax
功能：		寻找列主元
参数说明：	
			float* M	数组头地址
			int n		需要寻找的列
			int max		数组行数

**************************************************/

int findMax(float* M,int n,int max)
{
	float tem=0;
	int no;//标记主元所在行
	int i;
	for(i=n;i<max;i++)
	{
		if(fabs(M[i*(max+1)+n])>fabs(tem))
		{
			no=i;
			tem=M[i*(max+1)+n];
		}
	}
	if(fabs(tem)<0.000001)
		return 0;//若没有绝对值大于零的主元，返回0并退出

	for(i=0;i<max+1;i++)
	{
		tem=M[no*(max+1)+i];
		M[no*(max+1)+i]=M[n*(max+1)+i];
		M[n*(max+1)+i]=tem;
	}//交换主元至最上 
	return 1;
}

/**************************************************

函数名：	clearElmt
功能：		对特定的列消元
参数说明：	
			float* M		数组头地址
			int row			需要消元的列
			int maxLine		数组行数

**************************************************/

void clearElmt(float* M,int row ,int maxLine)
{
	float tem;//当前行的乘数
	int i,j;
	for( i=row+1;i<maxLine;i++)
	{
		tem=M[i*(maxLine+1)+row]/M[row*(maxLine+1)+row];
		for(j=row;j<maxLine+1;j++)
		{
			M[i*(maxLine+1)+j]-=M[row*(maxLine+1)+j]*tem;//消元
		}
	}
}

/**************************************************

函数名：	result
功能：		对上三角数组求解
参数说明：	
			float* M		数组头地址
			int maxLine		数组行数

**************************************************/

int result(float* M,int maxLine)
{
	int i;
	for(i=0;i<maxLine;i++)
	{
		if(fabs(M[i*(maxLine+1)+i])<0.000001)
			return 0;//若M[n][n]为零，返回零并退出
	}
	M[(maxLine-1)*(maxLine+1)+maxLine]/=M[(maxLine-1)*(maxLine+1)+maxLine-1];//最后一位的解
	for(i=maxLine-1;i>0;i--)
	{
		for(int j=i;j<maxLine;j++)
		{
			M[(i-1)*(maxLine+1)+maxLine]-=M[(j)*(maxLine+1)+maxLine]*M[(i-1)*(maxLine+1)+j];
			//将求出的数回代并累加到b[n]中
				
		}
		M[(i-1)*(maxLine+1)+maxLine]/=M[(i-1)*(maxLine+1)+(i-1)];//b[n]/a[n][n]
	}
	return 1;
}

int main()
{
	int line;
	int i,j;
	float *M;
	
	char ch;
	while(1)
	{
		while(1)
		{
			printf("请输入运行方式（y使用默认数组n输入数组e退出程序）：");
			cin>>ch;
			if(ch=='y')
			{
				line=5;
				float N[30]=
				{
					5,10,15,20,25,15,
						7,21,42,63,84,49,
						6,16,32,56,66,24,
						6,17,37,76,96,27,
						7,19,39,73,93,33
				};
				M=new float[30];
				for(i=0;i<30;i++)
				{
					M[i]=N[i];
				}
				break;
			}
			else if(ch=='n')
			{
				printf("请输入要求解矩阵的维数：");
				cin>>line;
				M=new float[line*(line+1)];//为数组分配地址空间
				
				for(i=0;i<line;i++)//为数组赋值
				{
					for(j=0;j<line+1;j++)
					{
						printf("请输入【%d，%d】",i,j);
						cin>>M[i*(line+1)+j];
					}
				}
				break;
			}
			else if(ch=='e')
			{
				return 1;
			}
			else
			{
				printf("输入有误，请重新输入！！\n");
			}
		}
		printf("输入的数组：\n");
		for(i=0;i<line;i++)//打印要求解的数组
		{
			for(j=0;j<line+1;j++)
			{
				printf("%12f",M[i*(line+1)+j]);
			}
			printf("\n");
		}
		printf("\n");
		for(i=0;i<line-1;i++)
		{
			if(findMax(M,i,line)==0)//找主元
			{
				printf("消元出现问题，题目无解，程序退出。\n");
				break;
			}
			clearElmt(M,i,line);
		}
		if(result(M,line)==0)
		{
			printf("回代出现问题，题目无解，程序退出。\n");
			break;
		}
		printf("运算结果为：\n");
		for(i=0;i<line;i++)//打印运算结果
		{
			printf("x[%d]=%12f\n",i,M[i*(line+1)+line]);
		}
		
	}
	return 1;
}