| [abd,lda,n,ml,mu,ipvt,det]=cgbdi(abd,lda,n,ml,mu,ipvt,det); |
function [abd,lda,n,ml,mu,ipvt,det]=cgbdi(abd,lda,n,ml,mu,ipvt,det);
%***BEGIN PROLOGUE CGBDI
%***PURPOSE Compute the determinant of a complex band matrix using the
% factors from CGBCO or CGBFA.
%***LIBRARY SLATEC (LINPACK)
%***CATEGORY D3C2
%***TYPE COMPLEX (SGBDI-S, DGBDI-D, CGBDI-C)
%***KEYWORDS BANDED, DETERMINANT, INVERSE, LINEAR ALGEBRA, LINPACK,
% MATRIX
%***AUTHOR Moler, C. B., (U. of New Mexico)
%***DESCRIPTION
%
% CGBDI computes the determinant of a band matrix
% using the factors computed by CGBCO or CGBFA.
% If the inverse is needed, use CGBSL N times.
%
% On Entry
%
% ABD COMPLEX(LDA, N)
% the output from CGBCO or CGBFA.
%
% LDA INTEGER
% the leading dimension of the array ABD .
%
% N INTEGER
% the order of the original matrix.
%
% ML INTEGER
% number of diagonals below the main diagonal.
%
% MU INTEGER
% number of diagonals above the main diagonal.
%
% IPVT INTEGER(N)
% the pivot vector from CGBCO or CGBFA.
%
% On Return
%
% DET COMPLEX(2)
% determinant of original matrix.
% Determinant = DET(1) * 10.0**DET(2)
% with 1.0 .LE. CABS1(DET(1)) .LT. 10.0
% or DET(1) = 0.0 .
%
%***REFERENCES J. J. Dongarra, J. R. Bunch, C. B. Moler, and G. W.
% Stewart, LINPACK Users' Guide, SIAM, 1979.
%***ROUTINES CALLED (NONE)
%***REVISION HISTORY (YYMMDD)
% 780814 DATE WRITTEN
% 890831 Modified array declarations. (WRB)
% 890831 REVISION DATE from Version 3.2
% 891214 Prologue converted to Version 4.0 format. (BAB)
% 900326 Removed duplicate information from DESCRIPTION section.
% (WRB)
% 920501 Reformatted the REFERENCES section. (WRB)
%***end PROLOGUE CGBDI
persistent i m ten zdum ;
ipvt_shape=size(ipvt);ipvt=reshape(ipvt,1,[]);
abd_shape=size(abd);abd=reshape([abd(:).',zeros(1,ceil(numel(abd)./prod([lda])).*prod([lda])-numel(abd))],lda,[]);
%
if isempty(ten), ten=0; end;
if isempty(i), i=0; end;
if isempty(m), m=0; end;
if isempty(zdum), zdum=0; end;
% cabs1= @(zdum) abs(real(zdum)) + abs(aimag(zdum));real :: cabs1;
%
cabs1= @(zdum) abs(real(zdum)) + abs(imag(zdum));
%***FIRST EXECUTABLE STATEMENT CGBDI
m = fix(ml + mu + 1);
det(1) =complex(1.0e0,0.0e0);
det(2) =complex(0.0e0,0.0e0);
ten = 10.0e0;
for i = 1 : n;
if( ipvt(i)~=i )
det(1) = -det(1);
end;
det(1) = abd(m,i).*det(1);
cabs1= @(zdum) abs(real(zdum)) + abs(imag(zdum));
if( cabs1(det(1))==0.0e0 )
break;
end;
cabs1= @(zdum) abs(real(zdum)) + abs(imag(zdum));
while( cabs1(det(1))<1.0e0 );
det(1) = complex(ten,0.0e0).*det(1);
det(2) = det(2) -complex(1.0e0,0.0e0);
end;
cabs1= @(zdum) abs(real(zdum)) + abs(imag(zdum));
while( cabs1(det(1))>=ten );
det(1) = det(1)./complex(ten,0.0e0);
det(2) = det(2) +complex(1.0e0,0.0e0);
end;
end;
ipvt_shape=zeros(ipvt_shape);ipvt_shape(:)=ipvt(1:numel(ipvt_shape));ipvt=ipvt_shape;
abd_shape=zeros(abd_shape);abd_shape(:)=abd(1:numel(abd_shape));abd=abd_shape;
end %subroutine cgbdi
%DECK CGBFA
|
|