function [x] = gaussSolve(A,r,b)
% gaussSolve solve system using gauss
% Syntax: x = gaussSolve(A,r,b)
% A is the n by n matrix produced by gauss
% b is the right hand side vector (not changed by gaussSolve)
% x is the solution (column vector) 
   
   [m,n] = size(A);
   x = zeros(n,1);

   % since we didn't hsve a right hand side vector in gauss
   % then we must ake the modifications to this vector here
   % using the pivot vector r and matrix A returned from gauss

   for k = 1 : n-1
      for i = k+1 : n
         b(r(i)) = b(r(i)) - A(r(i),k)*b(r(k));
      end
   end

   % back substitution using pivot vector

   x(n) = b(r(n)) / A(r(n),n);
   for i = n-1 : -1 : 1;
      s = b(r(i));
      for j = i+1 : n
         s = s - A(r(i),j)*x(j);
      end
      x(i) = s / A(r(i),i);
   end
end % gaussSolve