! Intrinsic acosh violates 2008 Standard rule 13.7.5 line 5 ! "If the result is complex the imaginary part is expressed in radians and lies ! in the range 0<= AIMAG(ACOSH(X)) <= pi" pi is approx 3.14 ! Imaginary part must be >=0 and <= pi so it may not be negative complex :: cc=(-1.0,-1.0) print *,acosh(cc) ! Imaginary part negative print *,log(cc+sqrt(cc**2-1)) ! this is correct acosh(cc) end (1.06127512,-2.23703575) (-1.06127524,2.23703575) Fortran acosh provides negative of correct answer.
I think there are two issues here, one is the glibc also puts the result in the wrong quadrant. The other issue is the GCC's constant folding does too. What is interesting is they both put in the same quadrant though. I cannot comment if this is a bug because I don't have a copy of the IEEE spec and/or C11 spec too.
(In reply to Andrew Pinski from comment #1) > I think there are two issues here, one is the glibc also puts the result in > the wrong quadrant. The other issue is the GCC's constant folding does too. > What is interesting is they both put in the same quadrant though. > > I cannot comment if this is a bug because I don't have a copy of the IEEE > spec and/or C11 spec too. It's a bug in the Fortran Standard. OP should send an interpretation request to J3. n1256.pdf has 7.3.6 Hyperbolic functions 7.3.6.1 The cacosh functions The cacosh functions compute the complex arc hyperbolic cosine of z, with a branch cut at values less than 1 along the real axis. The cacosh functions return the complex arc hyperbolic cosine value, in the range of a half-strip of non-negative values along the real axis and in the interval [-i pi, +i pi] along the imaginary axis.
(In reply to kargl from comment #2) > (In reply to Andrew Pinski from comment #1) > > I think there are two issues here, one is the glibc also puts the result in > > the wrong quadrant. The other issue is the GCC's constant folding does too. > > What is interesting is they both put in the same quadrant though. > > > > I cannot comment if this is a bug because I don't have a copy of the IEEE > > spec and/or C11 spec too. > > It's a bug in the Fortran Standard. OP should send an interpretation > request to J3. n1256.pdf has > > 7.3.6 Hyperbolic functions > > 7.3.6.1 The cacosh functions > > The cacosh functions compute the complex arc hyperbolic cosine > of z, with a branch cut at values less than 1 along the real axis. > > The cacosh functions return the complex arc hyperbolic cosine value, > in the range of a half-strip of non-negative values along the real > axis and in the interval [-i pi, +i pi] along the imaginary axis. Checking FreeBSD libm sources, which differ from glibc, one finds the comment /* * cacosh(z) = I*cacos(z) or -I*cacos(z) * where the sign is chosen so Re(cacosh(z)) >= 0. */ which means libm chooses the Riemann sheet with the REAL part always positive. Fortran must be using a convention for choosing a different Riemann sheet such that AIMAG part is always positive.
The standard branch cut for acosh (not just a C standard, but as at https://dlmf.nist.gov/4.37 for example) follows from the principles that (a) acosh(conj(x)) = conj(acosh(x)) and (b) complex acosh should take the same value as real acosh for those real arguments for which acosh has a real value. Fixing the sign of the imaginary part of the result is inconsistent with (a) (which is a general principle for complex libm functions, not just acosh).