function [x, nit] = gausssidel(A, b, x0, err, maxnit)
  n = length(b);
  D = diag(diag(A));
  L = tril(A, -1);
  U = A - D - L;
  M = D + L;
  N = -U;
  T = M^(-1)*N;
  c = M^(-1)*b;
  x = x0;
  for k = 1:maxnit
    old_x = x;
    x = T*x+c;
    if norm(x-old_x,inf) <= ((1-norm(T,inf))/norm(T,inf))*err
      nit = k;
      return;
    endif
  endfor
  error('Too dificult');
end;