function [x,L,U,er]=gauss(A,b)
N=length(A);
x=zeros(N,1);
L=zeros(size(A));
U=zeros(size(A));
er='';
for k=1:N-1
    if A(k,k)==0 
      er='Matrice singolare';
      return
    end
  for i=k+1:N    
    A(i,k)=A(i,k)/A(k,k);
    for j=k+1:N
      A(i,j)=A(i,j)-A(i,k)*A(k,j);
    end
    b(i)=b(i)-A(i,k)*b(k);
  end
end
x=zeros(N,1);
U=triu(A);
L=eye(N)+tril(A,-1);
for i=N:-1:1
  x(i)=(b(i)-U(i,:)*x)/U(i,i);
end