function [b,iflag] = utriangsl( A, b);
%
% UTRIANGSL  solves an N by N linear system with
% upper triangular matrix.
%
% Usage
%       [b,iflag] = utriangsl( A, b)
%
% input:
%        A:     the upper triangular N by N matrix A
%
%        b:     the right hand side b. b can have more than one column. 
%               In this case a system with multiple right hand sides is
%               solved.
%
%
% output:
%        b:     the solution of the system A x = b, if iflag = 0
%
%        iflag: error flag
%               iflag = 0  solution could be computed and is
%                          stored in b
%               iflag = 1  dimension of A or of b is not correct
%               iflag > 1  zero diagonal element detected in row iflag+1
%
%


iflag = 0;

% get size of A and check dimensions
[m,n]    = size(A);
if ( m ~= n | n ~= size(b,1) )
   iflag = 1;
   return
end

% Start back substitution (column form)

for i = n:-1:1
   if( A(i,i) == 0 )
       iflag = i+1;
       return
   end 
   b(i,:) = b(i,:) / A(i,i);
   b(1:i-1,:) = b(1:i-1,:) - A(1:i-1,i) * b(i,:);
end