Question

十万火急！C/C++高手请进！

分别用Jacobi(雅可比）迭代法和Gauss-Seidel迭代法求解100阶线性方程组：Ax=b，其中矩阵A为：a(ii)=5051-i,(i=1,2,...100);a(ij)=j,(j=1,2,...100且i不等于j）,b(i)=338350+5051*i-2*i*i,(i=1,2,...100).要求解的误差不超过10^（-8）(10的-8次幂）注意，我要的是求解该问题的C或者C++程序。（是两种方法哦）好的追加100~150分。有劳高手了，在此不胜感激！

Answer 1

//以下为Gauss-Seidel迭代法，雅可比迭代法正在调试中.....
#include <iostream>
#include <cmath>
using namespace std;

#define N 100

double norm_inf(double x1_x0[N], int n)
{
double ret = 0;
for(int i=0;i<n;i++)
if( fabs( x1_x0[i] ) > ret)
ret = fabs( x1_x0[i] );

return ret;
}

void seidel(double a[N][N],double g[N],int n)
{
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>=1e-8)&&(k<500))
{
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++)
{
cout << "x["<<(i+1)<<"]="<< x1[i] <<endl;
}

}
void main()
{
double a[N][N], b[N], x[N];
for(int i=1; i<=N; ++i)
{
b[i-1] = 338350+5051*i-2*i*i;
a[i-1][i-1] = 5051-i;
}
for(i=1; i<=N; ++i)
{
for(int j=1; j<=N; ++j)
{
if(i == j) continue;
a[i-1][j-1] = j;
}
}

int n = N;
seidel(a, b, n);
}

//雅可比迭代法：
#include <iostream>
#include <cmath>
using namespace std;

#define N 100
int Jacobidiedai(int n, double a[][N], double b[] , double x[])
{
int i,j;
double *x0,m=0,eps;
x0 = new double[N];
for(i=0;i<n;i++)
x0[i]=i+1;
int k = 0;
while(k<500)
{
for(i=0;i<n;i++)
{
for(j=0;j<n;j++)
m+=a[i][j]*x0[j];
m-=a[i][i]*x0[i];
x[i]=(b[i]-m)/a[i][i];
m=0;
}
eps=fabs(x[0]-x0[0]);
for(i=1;i<n;i++)
{
if(eps<fabs(x[i]-x0[i]))
eps=fabs(x[i]-x0[i]);
}
if(eps<1e-8) return 1;
for(i=0;i<n;i++)
x0[i]=x[i];
++k;
}
delete[] x0;
}
void main()
{
double a[N][N], b[N], x[N];
for(int i=1; i<=N; ++i)
{
b[i-1] = 338350+5051*i-2*i*i;
a[i-1][i-1] = 5051-i;
}
for(i=1; i<=N; ++i)
{
for(int j=1; j<=N; ++j)
{
if(i == j) continue;
a[i-1][j-1] = j;
}
}

int n = N;
Jacobidiedai(n, a, b, x);
for(i=0;i<N;i++)
{
cout << "x["<<(i+1)<<"]="<< x[i] <<endl;
}
}

Answer 2

不会啊