% Matlab code testge
%
% Test Gaussian elimination with partial and complete pivoting
%
% function result = testge(A,x,b)
%
%    inputs: 
%      A = matrix
%      x = true solution
%      b = A*x
%    output: (all norms are 1-norms)
%      result(1) = dimension n of A
%      result(2) = log_10(2-norm condition number(A))
%      result(3) = backward error of factorization from gepp 
%      result(4) = pivot growth from gepp = norm(Upp,max)/norm(A,max)
%      result(5) = backward error of factorization from gecp 
%      result(6) = pivot growth from gecp = norm(Ucp,max)/norm(A,max)
%      result(7) = backward error of solution from gepp 
%      result(8) = backward error of solution from gecp 
%      result(9) = condition estimate from gepp
%      result(10) = condition estimate from gecp
%      result(11) = explicit computation of condition number
%      result(12) = true error norm(xgepp - x,1)/norm(x,1)
%      result(13) = error bound on norm(xgepp - x,1)/norm(x,1)
%      result(14) = true error norm(xgecp - x,1)/norm(x,1)
%      result(15) = error bound on norm(xgecp - x,1)/norm(x,1)
%
function result = testge(A,x,b)
% 
result(1) = length(x);
result(2) = log10(cond(A));
% Call GEPP, solve Ax=b
[Lpp,Upp,Prpp] = lu(A);   %  A = Prpp' * Lpp * Upp
xgepp = A\b;
%
% Call GECPa, solve Ax=b
[Prcp,Lcp,Ucp,Pccp] = gecp(A); % A = Prcp * Lcp * Ucp * Pccp'
xgecp = Pccp*(Ucp\(Lcp\(Prcp'*b)));
%
result(3) = norm(A - Prpp'*Lpp*Upp,1)/norm(A,1);
result(4) = max(max(abs(Upp)))/max(max(abs(A)));
%
result(5) = norm(A - Prcp*Lcp*Ucp*Pccp',1)/norm(A,1);
result(6) = max(max(abs(Ucp)))/max(max(abs(A)));
%
result(7) = norm(A*xgepp-b,1)/(norm(A,1)*norm(x,1));
result(8) = norm(A*xgecp-b,1)/(norm(A,1)*norm(x,1));
%
result(9) = 1/rcond(A);
result(10) = 1/rcond(A*Pccp);
[U,S,V]=svd(A);
result(11) = norm(A,1)*norm(V*inv(S)*U,1);
% 
result(12)  = norm(xgepp - x, 1)/norm(x,1);
result(13)  = result(7)*result(9);   % backward error * cond estimate
result(14) = norm(xgecp - x, 1)/norm(x,1);
result(15) = result(8)*result(10);  % backward error * cond estimate
%
