juzhendiedai.cpp

利用c++编写的矩阵迭代算法

#include<math.h>
#include<iostream.h>
#include<stdio.h>
#define EPS 0.0001  //精度控制量
#define max 100 //最大迭代次数
#define N  10  //方程组的最大阶数
double norm_inf(double x[], int n) //求向量无穷范数
		  {
			  double norm;
			  int i;
			 norm=fabs(x[0]);
			 for(i=1;i<n;i++)
			 {
				 if(fabs(x[i])>norm)
					 norm=fabs(x[i]);
				 
			 }
			 return(norm);
		  }
void jacobi(double a[N][N],double g[N], int n)//jacobi 迭代法
{
 //b表示迭代矩阵，x1_x0表示迭代前后向量x0与x1之差向量
 double b[N][N]={0},x0[N],x1[N],x1_x0[N],norm,temp;
 int i, j,k;
   for(i=0;i<n;i++)
   {
      g[i]=g[i]/a[i][i];
	  for(j=0;j<n;j++)
	  {
		 if(i==j)
		 continue;
		 b[i][j]=-a[i][j]/a[i][i];
	  }
   }
	for(i=0;i<n;i++)       //选取初始向量并构成迭代条件
	{
		x0[i]=0;
		x1[i]=1;
		x1_x0[i]=x1[i]-x0[i];
	}
	k=0;
	norm=norm_inf(x1_x0,n);        //求差向量的无穷范数
	while((norm>=EPS)&&(k<max))  //开始迭代
	{
		for(i=0;i<n;i++)
			x0[i]=x1[i];
		for(i=0;i<n;i++)
			{
				temp=0;
				for(j=0;j<n;j++)
						temp=temp+b[i][j]*x0[j];
					x1[i]=temp+g[i];
					x1_x0[i]=x1[i]-x0[i];
			}
		norm=norm_inf(x1_x0,n);            //求差向量的无穷范数
		k++;                            //累加迭代次数
	}
	for(i=0;i<n;i++)
		printf("x[%d]=%lf\n", i,x1[i]); //输出计算结果
		printf("%d  迭代次数\n",k ); //输出实际迭代次数
			
}
void sater(double a[N][N],double g[N], int n)//jacobi 迭代法
{
	//b表示迭代矩阵，x1_x0表示迭代前后向量x0与x1之差向量
 double b[N][N]={0},x0[N],x1[N],x1_x0[N],norm,temp;
 int i, j,k;
 for(i=0; i<n;i++)
 {
	g[i]=g[i]/a[i][i];
	for(j=0;j<n;j++)
	{
	  if(i==j)
		continue;
	  b[i][j]=-a[i][j]/a[i][i];
	}
 }
 for(i=0;i<n;i++)       //选取初始向量并构成迭代条件
 {
	x0[i]=0;
	x1[i]=1;
	x1_x0[i]=x1[i]-x0[i];
 }
 k=0;
 norm=norm_inf(x1_x0,n);        //求差向量的无穷范数
 while((norm>=EPS)&&(k<max))  //开始迭代
 {
	for(i=0;i<n;i++)
		x0[i]=x1[i];
	for(i=0;i<n;i++)
	{
		 temp=0;
		 for(j=0;j<=i-1;j++)
		{ 
			 temp=temp+b[i][j]*x1[j];}
			 for(j=i+1;j<n;j++)
				{
				 temp=temp+b[i][j]*x0[j];}
				 x1[i]=temp+g[i];
				 x1_x0[i]=x1[i]-x0[i];
				 }
		norm=norm_inf(x1_x0,n);            //求差向量的无穷范数
		k++;                            //累加迭代次数
 }
			  for(i=0;i<n;i++)
				  printf("x[%d]=%lf\n", i,x1[i]); //输出计算结果
			   printf("%d  迭代次数\n",k ); //输出实际迭代次数
			
		  }
void main()
{
	double a[N][N]={ {5,2,3,2},{2,4,1,-2},{1,-3,4,3},{3,2,2,8}};
	double	g[N]={-1,5,4,-6};
	double     G[N]={-1,5,4,-6};
	cout<<"用雅可比方法求解:"<<endl;
	jacobi(a,g,4);
	cout<<"用高斯-赛德尔方法求解:"<<endl;
	sater( a,G,4);
}