tridiagldl.m

Interpolation routines in matlab

%
% tridiagLDL.m
%
% Compute the LDL^T decomposition without a ymmetric
% tridiagonal matrix A. The subdiagonal, the diagonal,
% and the superdiagonal of A are stored as follows
%
%            d(1)   e(1)
%            e(1)   d(2)   e(2)
%                   e(2)   d(3)   e(3)      
%     
%                        :      :      :
%
%                           e(n-2)   d(n-1)   e(n-1)
%                                    e(n-1)   d(n)   
%
%
%
%
%
% function [d, e, iflag] = tridiagLDL( d, e )
%
% input:
%        d:     diagonal
%
%        e:     superdiagonal
%
% output:
%        e:     subdiagonal of L
%
%        d:     diagonal od D
%
%        iflag: error flag
%               iflag = 0  LU decomposition could be computed
%               iflag = 1  dimension of c or e do not match dimension of d
%               iflag > 1  zero pivot element detected in row iflag+1
%

function [d, e, iflag] = tridiagLDL( d, e );

iflag = 0;

% get size of A and check dimensions
n  = size(d(:),1);
if ( n-1 > size(e(:),1) )
   iflag = 1;
   return
end


for i = 1:n-1
   if( d(i) == 0 ) 
       iflag = i+1;
       return
   end 
   tmp    = e(i);
   e(i)   = e(i) / d(i);
   d(i+1) = d(i+1) - tmp * e(i);

end