function [dbesi1result,x]=dbesi1(x);
dbesi1result=[];
persistent bi1cs first firstCall nti1 xmax xmin xsml y ; if isempty(firstCall),firstCall=1;end;
;
if isempty(nti1), nti1=0; end;
%***BEGIN PROLOGUE DBESI1
%***PURPOSE Compute the modified (hyperbolic) Bessel function of the
% first kind of order one.
%***LIBRARY SLATEC (FNLIB)
%***CATEGORY C10B1
%***TYPE doubleprecision (BESI1-S, DBESI1-D)
%***KEYWORDS FIRST KIND, FNLIB, HYPERBOLIC BESSEL FUNCTION,
% MODIFIED BESSEL FUNCTION, ORDER ONE, SPECIAL FUNCTIONS
%***AUTHOR Fullerton, W., (LANL)
%***DESCRIPTION
%
% DBESI1(X) calculates the doubleprecision modified (hyperbolic)
% Bessel function of the first kind of order one and doubleprecision
% argument X.
%
% Series for BI1 on the interval 0. to 9.00000E+00
% with weighted error 1.44E-32
% log weighted error 31.84
% significant figures required 31.45
% decimal places required 32.46
%
%***REFERENCES (NONE)
%***ROUTINES CALLED D1MACH, DBSI1E, DCSEVL, INITDS, XERMSG
%***REVISION HISTORY (YYMMDD)
% 770701 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)
%***end PROLOGUE DBESI1
if isempty(bi1cs), bi1cs=zeros(1,17); end;
if isempty(xmax), xmax=0; end;
if isempty(xmin), xmin=0; end;
if isempty(xsml), xsml=0; end;
if isempty(y), y=0; end;
if isempty(first), first=false; end;
if firstCall, bi1cs(1)=[-.19717132610998597316138503218149d-2]; end;
if firstCall, bi1cs(2)=[+.40734887667546480608155393652014d+0]; end;
if firstCall, bi1cs(3)=[+.34838994299959455866245037783787d-1]; end;
if firstCall, bi1cs(4)=[+.15453945563001236038598401058489d-2]; end;
if firstCall, bi1cs(5)=[+.41888521098377784129458832004120d-4]; end;
if firstCall, bi1cs(6)=[+.76490267648362114741959703966069d-6]; end;
if firstCall, bi1cs(7)=[+.10042493924741178689179808037238d-7]; end;
if firstCall, bi1cs(8)=[+.99322077919238106481371298054863d-10]; end;
if firstCall, bi1cs(9)=[+.76638017918447637275200171681349d-12]; end;
if firstCall, bi1cs(10)=[+.47414189238167394980388091948160d-14]; end;
if firstCall, bi1cs(11)=[+.24041144040745181799863172032000d-16]; end;
if firstCall, bi1cs(12)=[+.10171505007093713649121100799999d-18]; end;
if firstCall, bi1cs(13)=[+.36450935657866949458491733333333d-21]; end;
if firstCall, bi1cs(14)=[+.11205749502562039344810666666666d-23]; end;
if firstCall, bi1cs(15)=[+.29875441934468088832000000000000d-26]; end;
if firstCall, bi1cs(16)=[+.69732310939194709333333333333333d-29]; end;
if firstCall, bi1cs(17)=[+.14367948220620800000000000000000d-31]; end;
if firstCall, first=[true]; end;
firstCall=0;
%***FIRST EXECUTABLE STATEMENT DBESI1
if( first )
[nti1 ,bi1cs]=initds(bi1cs,17,0.1.*real(d1mach(3)));
xmin = 2.0d0.*d1mach(1);
xsml = sqrt(4.5d0.*d1mach(3));
xmax = log(d1mach(2));
end;
first = false;
%
y = abs(x);
if( y>3.0d0 )
%
if( y>xmax )
xermsg('SLATEC','DBESI1','ABS(X) SO BIG I1 OVERFLOWS',2,2);
end;
%
dbesi1result = exp(y).*dbsi1e(x);
else;
%
dbesi1result = 0.0d0;
if( y==0.0d0 )
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;
%
if( y<=xmin )
xermsg('SLATEC','DBESI1','ABS(X) SO SMALL I1 UNDERFLOWS',1,1);
end;
if( y>xmin )
dbesi1result = 0.5d0.*x;
end;
if( y>xsml )
dbesi1result = x.*(0.875d0+dcsevl(y.*y./4.5d0-1.0d0,bi1cs,nti1));
end;
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 DBESI
