Code covered by the BSD License  

Highlights from
slatec

slatec

by

 

06 Apr 2007 (Updated )

The slatec library converted into matlab functions.

[besk0eresult,x]=besk0e(x);
function [besk0eresult,x]=besk0e(x);
besk0eresult=[];
persistent ak02cs ak0cs besk0e bk0cs first firstCall ntak0 ntak02 ntk0 xsml y ; if isempty(firstCall),firstCall=1;end; 

if isempty(ak02cs), ak02cs=zeros(1,14); end;
if isempty(ak0cs), ak0cs=zeros(1,17); end;
if isempty(besk0eresult), besk0eresult=0; end;
if isempty(bk0cs), bk0cs=zeros(1,11); end;
if isempty(xsml), xsml=0; end;
if isempty(y), y=0; end;
if isempty(ntak0), ntak0=0; end;
if isempty(ntak02), ntak02=0; end;
if isempty(ntk0), ntk0=0; end;
%***BEGIN PROLOGUE  BESK0E
%***PURPOSE  Compute the exponentially scaled modified (hyperbolic)
%            Bessel function of the third kind of order zero.
%***LIBRARY   SLATEC (FNLIB)
%***CATEGORY  C10B1
%***TYPE      SINGLE PRECISION (BESK0E-S, DBSK0E-D)
%***KEYWORDS  EXPONENTIALLY SCALED, FNLIB, HYPERBOLIC BESSEL FUNCTION,
%             MODIFIED BESSEL FUNCTION, ORDER ZERO, SPECIAL FUNCTIONS,
%             THIRD KIND
%***AUTHOR  Fullerton, W., (LANL)
%***DESCRIPTION
%
% BESK0E(X) computes the exponentially scaled modified (hyperbolic)
% Bessel function of third kind of order zero for real argument
% X .GT. 0.0, i.e., EXP(X)*K0(X).
%
% Series for BK0        on the interval  0.          to  4.00000D+00
%                                        with weighted error   3.57E-19
%                                         log weighted error  18.45
%                               significant figures required  17.99
%                                    decimal places required  18.97
%
% Series for AK0        on the interval  1.25000D-01 to  5.00000D-01
%                                        with weighted error   5.34E-17
%                                         log weighted error  16.27
%                               significant figures required  14.92
%                                    decimal places required  16.89
%
% Series for AK02       on the interval  0.          to  1.25000D-01
%                                        with weighted error   2.34E-17
%                                         log weighted error  16.63
%                               significant figures required  14.67
%                                    decimal places required  17.20
%
%***REFERENCES  (NONE)
%***ROUTINES CALLED  BESI0, CSEVL, INITS, R1MACH, XERMSG
%***REVISION HISTORY  (YYMMDD)
%   770401  DATE WRITTEN
%   890531  Changed all specific intrinsics to generic.  (WRB)
%   890531  REVISION DATE from Version 3.2
%   891214  Prologue converted to Version 4.0 format.  (BAB)
%   900315  CALLs to XERROR changed to CALLs to XERMSG.  (THJ)
%   900326  Removed duplicate information from DESCRIPTION section.
%           (WRB)
%***end PROLOGUE  BESK0E
if isempty(first), first=false; end;
if firstCall,   bk0cs(1)=[-.03532739323390276872e0];  end;
if firstCall,   bk0cs(2)=[.3442898999246284869e0];  end;
if firstCall,   bk0cs(3)=[.03597993651536150163e0];  end;
if firstCall,   bk0cs(4)=[.00126461541144692592e0];  end;
if firstCall,   bk0cs(5)=[.00002286212103119451e0];  end;
if firstCall,   bk0cs(6)=[.00000025347910790261e0];  end;
if firstCall,   bk0cs(7)=[.00000000190451637722e0];  end;
if firstCall,   bk0cs(8)=[.00000000001034969525e0];  end;
if firstCall,   bk0cs(9)=[.00000000000004259816e0];  end;
if firstCall,   bk0cs(10)=[.00000000000000013744e0];  end;
if firstCall,   bk0cs(11)=[.00000000000000000035e0];  end;
if firstCall,   ak0cs(1)=[-.07643947903327941e0];  end;
if firstCall,   ak0cs(2)=[-.02235652605699819e0];  end;
if firstCall,   ak0cs(3)=[.00077341811546938e0];  end;
if firstCall,   ak0cs(4)=[-.00004281006688886e0];  end;
if firstCall,   ak0cs(5)=[.00000308170017386e0];  end;
if firstCall,   ak0cs(6)=[-.00000026393672220e0];  end;
if firstCall,   ak0cs(7)=[.00000002563713036e0];  end;
if firstCall,   ak0cs(8)=[-.00000000274270554e0];  end;
if firstCall,   ak0cs(9)=[.00000000031694296e0];  end;
if firstCall,   ak0cs(10)=[-.00000000003902353e0];  end;
if firstCall,   ak0cs(11)=[.00000000000506804e0];  end;
if firstCall,   ak0cs(12)=[-.00000000000068895e0];  end;
if firstCall,   ak0cs(13)=[.00000000000009744e0];  end;
if firstCall,   ak0cs(14)=[-.00000000000001427e0];  end;
if firstCall,   ak0cs(15)=[.00000000000000215e0];  end;
if firstCall,   ak0cs(16)=[-.00000000000000033e0];  end;
if firstCall,   ak0cs(17)=[.00000000000000005e0];  end;
if firstCall,   ak02cs(1)=[-.01201869826307592e0];  end;
if firstCall,   ak02cs(2)=[-.00917485269102569e0];  end;
if firstCall,   ak02cs(3)=[.00014445509317750e0];  end;
if firstCall,   ak02cs(4)=[-.00000401361417543e0];  end;
if firstCall,   ak02cs(5)=[.00000015678318108e0];  end;
if firstCall,   ak02cs(6)=[-.00000000777011043e0];  end;
if firstCall,   ak02cs(7)=[.00000000046111825e0];  end;
if firstCall,   ak02cs(8)=[-.00000000003158592e0];  end;
if firstCall,   ak02cs(9)=[.00000000000243501e0];  end;
if firstCall,   ak02cs(10)=[-.00000000000020743e0];  end;
if firstCall,   ak02cs(11)=[.00000000000001925e0];  end;
if firstCall,   ak02cs(12)=[-.00000000000000192e0];  end;
if firstCall,   ak02cs(13)=[.00000000000000020e0];  end;
if firstCall,   ak02cs(14)=[-.00000000000000002e0];  end;
if firstCall,   first=[true];  end;
firstCall=0;
%***FIRST EXECUTABLE STATEMENT  BESK0E
if( first )
[ntk0 ,bk0cs]=inits(bk0cs,11,0.1.*r1mach(3));
[ntak0 ,ak0cs]=inits(ak0cs,17,0.1.*r1mach(3));
[ntak02 ,ak02cs]=inits(ak02cs,14,0.1.*r1mach(3));
xsml = sqrt(4.0.*r1mach(3));
end;
first = false;
%
if( x<=0. )
xermsg('SLATEC','BESK0E','X IS ZERO OR NEGATIVE',2,2);
end;
if( x>2. )
%
if( x<=8. )
besk0eresult =(1.25+csevl((16../x-5.)./3.,ak0cs,ntak0))./sqrt(x);
end;
if( x>8. )
besk0eresult =(1.25+csevl(16../x-1.,ak02cs,ntak02))./sqrt(x);
end;
else;
%
y = 0.;
if( x>xsml )
y = x.*x;
end;
besk0eresult = exp(x).*(-log(0.5.*x).*besi0(x)-.25+csevl(.5.*y-1.,bk0cs,ntk0));
csnil=dbstack(1); csnil=csnil(1).name(1)~='@';
if csnil&&~isempty(inputname(1)), assignin('caller','FUntemp',x); evalin('caller',[inputname(1),'=FUntemp;']); end
return;
end;
%
csnil=dbstack(1); csnil=csnil(1).name(1)~='@';
if csnil&&~isempty(inputname(1)), assignin('caller','FUntemp',x); evalin('caller',[inputname(1),'=FUntemp;']); end
end
%DECK BESK0

Contact us