Commit 671a9a2c authored by ulrich_y's avatar ulrich_y

Made branch cuts in plog1 more precise

parent ef6e73d9
......@@ -459,8 +459,19 @@ CONTAINS
implicit none
type(inum) :: a,b
complex(kind=prec) plog1
!TODO!!
plog1 = log(1.-a%c/b%c)
if ( (abs(aimag(a)).lt.zero).and.(abs(aimag(b)).lt.zero) ) then
! Both arguments are real, only here does the ieps matter
plog1 = log(abs(1.-a%c/b%c))
! this does not depend on the sign of a
if (real(a).gt.real(b) .and. real(b).gt. 0) then
plog1 = plog1 + b%i0*i_*pi
elseif(real(a).lt.real(b) .and. real(b).lt. 0) then
plog1 = plog1 - b%i0*i_*pi
endif
else
plog1 = log(1.-a%c/b%c)
endif
END FUNCTION
#ifndef NOCACHE
......
Markdown is supported
0% or
You are about to add 0 people to the discussion. Proceed with caution.
Finish editing this message first!
Please register or to comment