-
git-svn-id: https://lxsd.femto-st.fr/svn/fvn@70 b657c933-2333-4658-acf2-d3c7c2708721
zlnrel.f
1.04 KB
complex(kind(1.d0)) function zlnrel (z)
implicit none
c april 1977 version. w. fullerton, c3, los alamos scientific lab.
c
c clnrel(z) = clog(1+z) with relative error accuracy near z = 0.
c let rho = cabs(z) and
c r**2 = cabs(1+z)**2 = (1+x)**2 + y**2 = 1 + 2*x + rho**2 .
c now if rho is small we may evaluate clnrel(z) accurately by
c clog(1+z) = cmplx (alog(r), carg(1+z))
c = cmplx (0.5*alog(r**2), carg(1+z))
c = cmplx (0.5*alnrel(2*x+rho**2), carg(1+z))
c
complex(kind(1.d0)) z
real(kind(1.d0)) dlnrel,d1mach,zarg,sqeps,x,rho
external dlnrel, zarg, d1mach
data sqeps /0.0d0/
c
if (sqeps.eq.0.) sqeps = sqrt (d1mach(4))
c
if ( abs(1.+z).lt.sqeps) call seteru (
1 54hzlnrel answer lt half precision because z too near -1, 54,
2 1, 1)
c
rho = abs(z)
if (rho.gt.0.375) zlnrel = log (1.0+z)
if (rho.gt.0.375) return
c
x = real(z)
zlnrel = cmplx (0.5*dlnrel(2.*x+rho**2), zarg(1.0+z), kind(1.d0))
c
return
end