seven.cpp

这是我数值分析课程的实习题目的源代码(并附有文档),文档包括对问题进行了描述,并对结果进行分析,希望对大家有所帮助.

//选列主元高斯方法解线性方程

#include<iostream.h>
#include <math.h>
#include<iomanip.h>
int Funtion1();
int Funtion(double A[5][5]);
double x1,x2,x3;
main()
{
	Funtion1();  //单精度求根
	// 未扰动系数
	double A[5][5]={0,0,0,0,0,
		0,3.01,6.03,1.99,1,
		0,1.27,4.16,-1.23,1,
		0,0.987,-4.81,9.34,1,
	0,0,0,0,0};
	cout<<"未扰动系数,";
	Funtion(A);//双精度计算
	//扰动系数
     double A1[5][5]={0,0,0,0,0,
		0,3.00,6.03,1.99,1,
		0,1.27,4.16,-1.23,1,
		0,0.99,-4.81,9.34,1,
	0,0,0,0,0};
	 cout<<"扰动系数,";
	 Funtion(A1);
	return 0;
}
//单精度求根
int Funtion1()
{
	float A[5][5]={0,0,0,0,0,
		0,3.01,6.03,1.99,1,
		0,1.27,4.16,-1.23,1,
		0,0.987,-4.81,9.34,1,
	0,0,0,0,0};
	int n=3;
	int R[4];
	float m[5][5];
	for(int i=1;i<=3;i++)
		R[i]=i;  //方程指标
	float temp;
	int r;
	for(int k=1;k<=(n-1);k++)
	{
	    temp=fabs(A[R[k]][k]);
		r=R[k];
		for(int h=k+1;h<=n;h++)
			if(fabs(A[R[h]][k])>temp)
			{
				temp=fabs(A[R[h]][k]);
				r=h;
			}
		if(temp==0)
		{
			cout<<"方程无唯一解!";
			return 0;
		}

		 float  temp1=R[k];
	     R[k]=R[r];
		 R[r]=temp1;
		 for(i=k+1;i<=n;i++)
		 {
			 m[R[i]][k]=A[R[i]][k]/A[R[k]][k];
			 for(int j=k+1;j<=(n+1);j++)
				 A[R[i]][j]=A[R[i]][j]-m[R[i]][k]*A[R[k]][j];

		 }

	}
	if(A[R[n]][n]==0)
	{
			cout<<"方程无唯一解!";
			return 0;
	}
	float x[5];
	x[3]=A[R[n]][n+1]/A[R[n]][n];
	for( i=n-1;i>=1;i--)
	{
		float t=0;
		for(int j=i+1;j<=n;j++)
		{
			
			t=t+A[R[i]][j]*x[j];
		}
		x[i]=(A[R[i]][n+1]-t)/A[R[i]][i];

	}
	cout<<"单精度："<<endl<<setprecision(10)<<"x1="<<x[1]
		<<"\t"<<"x2="<<x[2]<<"\t"<<"x3="<<x[3]<<endl;
	return 0;
	
}
//双精度计算
int  Funtion(double A[5][5])
{
	
	int n=3;
	int R[4];
	double m[5][5];
	for(int i=1;i<=3;i++)
		R[i]=i;  //方程指标
	double temp;
	int r;
	for(int k=1;k<=(n-1);k++)
	{
	    temp=fabs(A[R[k]][k]);
		r=R[k];
		for(int h=k+1;h<=n;h++)
			if(fabs(A[R[h]][k])>temp)
			{
				temp=fabs(A[R[h]][k]);
				r=h;
			}
		if(temp==0)
		{
			cout<<"方程无唯一解!";
			return 0;
		}

		 int  temp1=R[k];
	     R[k]=R[r];
		 R[r]=temp1;
		 for(i=k+1;i<=n;i++)
		 {
			 m[R[i]][k]=A[R[i]][k]/A[R[k]][k];
			 for(int j=k+1;j<=(n+1);j++)
				 A[R[i]][j]=A[R[i]][j]-m[R[i]][k]*A[R[k]][j];

		 }

	}
	if(A[R[n]][n]==0)
	{
			cout<<"方程无唯一解!";
			return 0;
	}
	double x[5];
	x[3]=A[R[n]][n+1]/A[R[n]][n];
	for( i=n-1;i>=1;i--)
	{
		double t=0;
		for(int j=i+1;j<=n;j++)
		{
			
			t=t+A[R[i]][j]*x[j];
		}
		x[i]=(A[R[i]][n+1]-t)/A[R[i]][i];

	}
	cout<<"双精度:"<<endl<<setprecision(10)<<"x1="<<x[1]<<"\t"<<"x2="<<x[2]<<"\t"
		<<"x3="<<x[3]<<endl;
	return 0;
		
}