% Generate test matrices for Gaussian elimination
% Inputs:
%   ncase = integer from 1 to 5, indicating type of matrix.
%          See code below for documentation.
%   n    = dimension of matrix A to generate
% Outputs:
%   A    = n-by-n matrix
%   x    = n-by-1 vector
%   b    = A*x.  
% function [A,x,b] = matgen(ncase,n)
function [A,x,b] = matgen(ncase,n)
rand('state',0);       % reset random number generate, to always generate
                       % the same random matrices
if     (ncase == 1),    % random Hessenberg orthogonal matrix; very well-conditioned
   A = gallery('randhess',n);
elseif (ncase == 2),    % random matrix; usually well-conditioned
   A = rand(n,n)-.5;
elseif (ncase == 3),    % Hilbert matrix; very ill-conditioned 
   A = hilb(n);
elseif (ncase == 4)     % well conditioned, but bad pivot growth with GEPP or GENP, not GECP
   A = -tril(ones(n,n),-1)+eye(n); A(:,n) = ones(n,1);
elseif (ncase == 5)     % well conditioned, bad pivot growth with GENP 
   n2 = fix(n/2);
   A = [[1e-14*randn(n2,n2),randn(n2,n-n2)];[randn(n-n2,n2),randn(n-n2,n-n2)]];
end
x = randn(n,1);
b = A*x;