Gauss Seidel - Finite Element Method

390 Views Asked by At

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];

}