I am solving an equation using finite element method, and for that I have to use Gauss Seidel to invert a matrix. In Gauss Seidel I am using a "while" which breaks if the absolute error reaches the error that I want but what happens is that the absolute error is always increasing along iterations...what could it be my error?
I think it is the initial guess for Gauss Seidel and so I tried some different initial guess than (0,...,0) such as (1, ..., 1) and (10, ..., 10).
Here is a part of my code:
double u_h[d], XO[d];
//initial guess
for(int i=0; i<d; i++) XO[i]=1.;
bool erro=false;
double erro_abs=100;
double erro_teorico=pow (10.0, -20.);
while (erro==false){
for(int i=0; i<d; i++){
for(int j=0; j<=i-1; j++) soma_j-=A[i][j]*u_h[j];
for(int j=i+1; j<d; j++) soma_j-=A[i][j]*XO[j];
u_h[i]=(soma_j+B[i])/(A[i][i]);
}
//error calculation:
//1st i will calculate the largest absolute error of this iteration:
double erro_abs_max=fabs(u_h[0]-XO[0]);
for(int i=0; i<d; i++) {
if (fabs(u_h[i]-XO[i])>erro_abs_max) erro_abs_max=fabs(u_h[i]-XO[i]);
}
erro_abs=erro_abs_max;
if (erro_abs<erro_teorico) erro=true;
for(int i=0; i<d; i++) XO[i]=u_h[i];
}