sor.cpp

Jacobi和SOR算法实现

#include<iostream.h>
#include<stdio.h>
#include<math.h>
#define MAXSIZE 20
#define MAX 100
#define epsilon 0.00000001//求解精度
int main()
{
	int n;
	int i,j,k;
	double err,w;
	static double a[MAXSIZE][MAXSIZE],b[MAXSIZE][MAXSIZE],c[MAXSIZE],g[MAXSIZE];
	static double x[MAXSIZE],y[MAXSIZE];
	cout<<"Please input the dim of Ax=c"<<endl;//输入方程组维数
	cin>>n;
	if(n<=0)
	{cout<<"Please input a number between 1 and "<<n;
	return 1;
	}
	cout<<"Now input the matrix a(i,j):"<<endl;//输入A矩阵
	for(i=0;i<n;i++)
		for(j=0;j<n;j++)
			cin>>a[i][j];
		cout<<"Now input the matrix c(i):"<<endl;//输入c矩阵
		for(i=0;i<n;i++)
			cin>>c[i];
		cout<<"Now input the value of w:"<<endl;
			cin>>w;
		if(w<1||w>=2)
		{
			cout<<"w must between 1 and 2"<<endl;
	        return 1;	
		}
		for(i=0;i<n;i++)//变换
			for(j=0;j<n;j++)
			{
				b[i][j]=-a[i][j]/a[i][i];
				g[i]=c[i]/a[i][i];
			}
			for(int m=0;m<MAX;m++)//迭代
			{
				for(j=0;j<n;j++)
					y[j]=g[j];
				for(j=0;j<n;j++)
				{
					for(k=0;k<j;k++)
						y[j]+=b[j][k]*y[k];
					for(k=j+1;k<n;k++)
						y[j]+=b[j][k]*x[k];
					y[j]=(1-w)*x[j]+w*y[j];
				}
				
				err=0;//误差
				for(j=0;j<n;j++)
					if(err<fabs(y[j]-x[j]))err=fabs(y[j]-x[j]);
					for(j=0;j<n;j++)
						x[j]=y[j];
					if(err<epsilon)//输出
					{cout<<"the result is:"<<endl;
					for(i=0;i<n;i++)cout<<x[i]<<" ";
					cout<<endl<<"the iterations: "<<m<<endl;//迭代次数
					return 0;
					}
			}

			cout<<"no result"<<endl;
			return 1;
}