dexprl.f
1.13 KB
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
double precision function dexprl (x)
c july 1980 edition. w. fullerton, c3, los alamos scientific lab.
c
c evaluate exprel(x) = (exp(x) - 1.0) / x. for small abs(x) the
c taylor series is used. if x is negative the reflection formula
c exprel(x) = exp(x) * exprel(abs(x))
c may be used this reflection formula will be of use when the
c evaluation for small abs(x) is done by chebyshev series rather than
c taylor series.
c
double precision x, absx, alneps, xbnd, xln, xn, d1mach,
1 dexp, dlog
external d1mach
data nterms, xbnd / 0, 0.0d0 /
c
if (nterms.ne.0) go to 10
alneps = dlog(d1mach(3))
xn = 3.72d0 - 0.3d0*alneps
xln = dlog((xn+1.0d0)/1.36d0)
nterms = xn - (xn*xln+alneps)/(xln+1.36d0) + 1.5d0
xbnd = d1mach(3)
c
10 absx = dabs(x)
if (absx.le.0.5d0) go to 20
call entsrc (irold, 1)
dexprl = (dexp(x)-1.0d0)/x
call entsrc (ir, irold)
return
c
20 dexprl = 1.0d0
if (absx.lt.xbnd) return
c
dexprl = 0.0d0
do 30 i=1,nterms
dexprl = 1.0d0 + dexprl*x/dble(float(nterms+2-i))
30 continue
c
return
end