function [x]=substs(a,b)
n=length(b);
x=zeros(n,1);
x(n)=b(n)/a(n,n);
for i=n-1:-1:1
   soma=a(i,i+1:n)*x(i+1:n);
   x(i)=(b(i)-soma)/a(i,i);
end
return