# C++ my Gauss–Seidel method doesn’t work wery well

I want to solve system of equations (3×3) in C++ useing Gauss–Seidel method. I read some of pages about this: Wiki Page2 and Page3. By pattern oneself on those pages I wrote this code:

```#include <iostream>
#include <cmath>
#define epsilon 0.001

int main()
{
//declaration and initialization of variables
long double a, b, x, xo;

a = 3;   a = -6.5; a = 3.33; b = 7;
a = 4.7; a = 8;    a = -17;  b = 5;
a = 6;   a = 2;    a = 8;    b = -5.2;
//By calculator: x=0,5063353908 y=-1,215147148 z=-0,725964756

int n = 3;

//initialization of x array
for (int i = 0; i < n; i++)
{
x[i] = b[i] / a[i][i];
}

//start Gauss-Seidel method
do
{
for (int i = 0; i < n; i++)
{
//xo is useing to chceck accuracy
xo[i] = x[i];
for (int j = 0; j < n; j++)
{
if (i != j) x[i] -= (a[i][j] * x[j]) / a[i][i];
}
}
//chceck the accuracy
} while (fabs(x - xo) < epsilon &&
fabs(x - xo) < epsilon &&
fabs(x - xo) < epsilon);

std::cout << "x = " << x << std::endl;
std::cout << "y = " << x << std::endl;
std::cout << "z = " << x << std::endl;

return 0;
}
```

This code doesn’t work very well – the result of the program isn’t the same like solving by calculator. I really don’t have an idea what is going on.

## Answer

There are several issues in your code:

1. The accuracy checking is not correct
2. In calculation of the solution, a term `b[i]` is missing
3. The matrix is not diagonal dominant, far from it: the algorithm cannot converge

In addition, to avoid getting `nan` when the algorithm doesn’t converge, I have added a condition to limit the number of iterations.

I have also added a final check, i.e. calculation of `ax`, to see if a good result is effectively obtained.

```#include <iostream>
#include <cmath>
#define epsilon 0.001

int main()
{
//declaration and initialization of variables
long double a, b, x, xo;

a = 3;   a = -6.5; a = 3.33; b = 7;
a = 4.7; a = 8;    a = -17;  b = 5;
a = 6;   a = 2;    a = 8;    b = -5.2;
//By calculator: x=0,5063353908 y=-1,215147148 z=-0,725964756

a = 12;       // to make matrix diadonal dominant
a = -1;

int n = 3;

//initialization of x array
for (int i = 0; i < n; i++) {
x[i] = b[i] / a[i][i];
}

//start Gauss-Seidel method
int iter = 0;
int n_iter_max = 20;
do {
for (int i = 0; i < n; i++) {
//xo is used to check accuracy
xo[i] = x[i];
x[i] = b[i];
for (int j = 0; j < n; j++) {
if (i != j) x[i] -= (a[i][j] * x[j]);
}
x[i] /= a[i][i];
}
iter++;
//check the accuracy
} while ((std::abs(x - xo) > epsilon ||
std::abs(x - xo) > epsilon ||
std::abs(x - xo) > epsilon) && (iter <= n_iter_max));

std::cout << "x = " << x << std::endl;
std::cout << "y = " << x << std::endl;
std::cout << "z = " << x << std::endl;

//  Final check

double b_cal;
for (int i = 0; i < n; ++i) {
b_cal[i] = 0.0;
for (int j = 0; j < n; j++) {
b_cal[i] += a[i][j] * x[j];
}
}
std::cout << "A*x = ";
for (int i = 0; i < n; ++i) {
std::cout << b_cal[i] << " (" << b_cal[i] - b[i] << ")  ";
}
std::cout << "n";

return 0;
}
```