function x=mnhf_cramer(A,b)
%MNHF_CRAMER implements Cramer's Rule solution of Ax=b
%
% A - the matrix in question
% b - right-hand side vector
%
% Usage  mnhf_cramer(A,b);

tolerance = 1e-8;

[r,c] = size(A);       % Compute the size of A.
[rr,cc] = size(b);     % Compute the size of b.
if cc==1
 if isequal(r,c)
  if isequal(c,rr)
   
   % Initialize array, x.
   x = zeros(r,1);
   
   for ii=1:c
    % Define the matrix Aii as A with the appropriate column replaced 
    % by b. Use : to facilitate matrix index referencing.
    if ii==1
     Aii = [b A(:,2:c)];
    elseif (ii>1 && ii<length(b))
     Aii = [A(:,1:ii-1) b A(:,ii+1:c)];
    else
     Aii = [A(:,1:c-1) b];
    end    
    % Solve for x(ii) using the Cramer's rule formula.
    x(ii) = mnhf_determinant(Aii)/mnhf_determinant(A);
   end
   
   % Check residual.
   resid = norm(A*x-b,2);
   if resid>tolerance
    error('Residual exceeds tolerance.')
   end
   
  else
   fprintf('Size mismatch between matrix and vector.') 
  end
 else
  fprintf('Matrix is not square.')
 end
else
 fprintf('b is not a column vector.')
end