kusano 7d535a
function info = try2(A,pin);
kusano 7d535a
% TRY2 : test SUPERLU with 2 outputs
kusano 7d535a
%
kusano 7d535a
% info = try2(A,pin);
kusano 7d535a
% normally info is the residual norm; 
kusano 7d535a
% but info is at least 10^6 if U is not triangular
kusano 7d535a
% or if the factors contain explicit zeros.
kusano 7d535a
% Copyright (c) 1995 by Xerox Corporation.  All rights reserved.
kusano 7d535a
% HELP COPYRIGHT for complete copyright and licensing notice.
kusano 7d535a
kusano 7d535a
kusano 7d535a
n = max(size(A));
kusano 7d535a
info = 0;
kusano 7d535a
kusano 7d535a
if nargin == 1
kusano 7d535a
    [l,u] = superlu(A);
kusano 7d535a
elseif nargin == 2
kusano 7d535a
    [l,u] = superlu(A,pin);
kusano 7d535a
else 
kusano 7d535a
    error('requires 1 or 2 inputs');
kusano 7d535a
end;
kusano 7d535a
kusano 7d535a
figure(2); spy(l); title('L');
kusano 7d535a
figure(1); spy(u); title('U');
kusano 7d535a
drawnow;
kusano 7d535a
kusano 7d535a
format compact
kusano 7d535a
disp('SUPERLU with 2 outputs:');
kusano 7d535a
if any(any(triu(l,1)))
kusano 7d535a
    disp('L is *NOT* lower triangular.');
kusano 7d535a
else
kusano 7d535a
    disp('L is lower triangular.');
kusano 7d535a
end;
kusano 7d535a
if nnz(l) == nnz(l+l)
kusano 7d535a
    disp('L has no explicit zeros.');
kusano 7d535a
else
kusano 7d535a
    disp('L contains explicit zeros.');
kusano 7d535a
    info = info+10^6;
kusano 7d535a
end;
kusano 7d535a
if any(any(tril(u,-1)))
kusano 7d535a
    disp('U is *NOT* upper triangular.');
kusano 7d535a
    info = info + 10^6;
kusano 7d535a
else
kusano 7d535a
    disp('U is upper triangular.');
kusano 7d535a
end;
kusano 7d535a
if nnz(u) == nnz(u+u)
kusano 7d535a
    disp('U has no explicit zeros.');
kusano 7d535a
else
kusano 7d535a
    disp('U contains explicit zeros.');
kusano 7d535a
    info = info+10^6;
kusano 7d535a
end;
kusano 7d535a
kusano 7d535a
if nargin == 1
kusano 7d535a
    rnorm = norm(A - l*u,inf);
kusano 7d535a
    fprintf(1,'||A -  L*U|| = %d\n', rnorm);
kusano 7d535a
elseif size(pin) == [n n]
kusano 7d535a
    rnorm = norm(A*pin' - l*u,inf);
kusano 7d535a
    fprintf(1,'||A*PIN'' -  L*U|| = %d\n', rnorm);
kusano 7d535a
elseif max(size(pin)) == n
kusano 7d535a
    rnorm = norm(A(:,pin) - l*u,inf);
kusano 7d535a
    fprintf(1,'||A(:,PIN) -  L*U|| = %d\n', rnorm);
kusano 7d535a
else
kusano 7d535a
    rnorm = norm(A - l*u,inf);
kusano 7d535a
    fprintf(1,'||A -  L*U|| = %d\n', rnorm);
kusano 7d535a
end;
kusano 7d535a
kusano 7d535a
info = info + rnorm;
kusano 7d535a
disp(' ');