What should be the result of CACOSH(F)(CCOSH(F)(-2 + 1j))?

What should be the result of CACOSH(F)(CCOSH(F)(-2 + 1j))?

[Paolo Mantegazza]

Hi,

calling: j = sqrt(-1),

MUSL        answer is: -0.200000e+1 + 0.100000e+1j;

GLIBC       answer is: +0.20000e+1  -  0.100000e+1j;

SCILAB     answer is: +0.20000e+1  -  0.100000e+1j;

MATLAB   answer is: +0.20000e+1  -  0.100000e+1j;

UCLIB-NG answer is: +0.20000e+1  -  0.100000e+1j;

Math is not democracy so maybe MUSL's answer is the right one. In fact, 
with infinite precision at least, one should expect that, by applying 
the inverse of a function to a function, the result should be the used 
function argument.

So, does it either show a missed correct principal value or that MUSL is 
the smartest one?

In any case, following: 
http://mathworld.wolfram.com/InverseHyperbolicCosine.html, a way to have 
MUSL match GLIBC:SCILAB:MATLAB:UCLIB-NG is to change the two code lines 
in MUSL ./src/complex/cacosh.c

      z = cacos(z);
      return CMPLX(-cimag(z), creal(z));  // AKA j*cacos(z)

into

     return = clog(z + csqrt(z + 1) * csqrt(z - 1));  // AKA a 
definition of cacosh

As a further info, NEWLIB cacosh.c (not tested here) recites:

#if 0 /* does not give the principal value */
         w = I * cacos(z);
#else
         w = clog(z + csqrt(z + 1) * csqrt(z - 1));
#endif
         return w;

it should be remarked that, to provide the correct principal value using 
just cacos(z), the above mentioned link addresses the need of testing 
the cacosh argument in order to appropriately use either  j*cacos(z) or  
-j*cacos(z). It is therefore likely that the fix to be chosen, if any, 
should be based on computational efficiency.

Once more math is not democracy, so the answer must be left to the math 
savviest.

Regards, Paolo Mantegazza.

Re: What should be the result of CACOSH(F)(CCOSH(F)(-2 + 1j))?

[Szabolcs Nagy]

* Paolo Mantegazza <paolo.mantegazza@polimi.it> [2018-05-17 14:32:40 +0200]:
> Hi,
> 
> calling: j = sqrt(-1),
> 
> MUSL        answer is: -0.200000e+1 + 0.100000e+1j;
> 
> GLIBC       answer is: +0.20000e+1  -  0.100000e+1j;
> 
> SCILAB     answer is: +0.20000e+1  -  0.100000e+1j;
> 
> MATLAB   answer is: +0.20000e+1  -  0.100000e+1j;
> 
> UCLIB-NG answer is: +0.20000e+1  -  0.100000e+1j;
> 
> Math is not democracy so maybe MUSL's answer is the right one. In fact, with
> infinite precision at least, one should expect that, by applying the inverse
> of a function to a function, the result should be the used function
> argument.
> 
> So, does it either show a missed correct principal value or that MUSL is the
> smartest one?
> 

musl does not try very hard to provide correct complex/
so it's easily possible that the answer is wrong.

we can try to get these edge cases right, but there are
lot of large ulp error cases and intermediate value
overflows cases that are not trivial to fix.
(that's why we use the simplest expression in most cases
which is wrong but in an obvious way)

> In any case, following:
> http://mathworld.wolfram.com/InverseHyperbolicCosine.html, a way to have
> MUSL match GLIBC:SCILAB:MATLAB:UCLIB-NG is to change the two code lines in
> MUSL ./src/complex/cacosh.c
> 
>      z = cacos(z);
>      return CMPLX(-cimag(z), creal(z));  // AKA j*cacos(z)
> 
> into
> 
>     return = clog(z + csqrt(z + 1) * csqrt(z - 1));  // AKA a definition of
> cacosh
> 
> As a further info, NEWLIB cacosh.c (not tested here) recites:
> 
> #if 0 /* does not give the principal value */
>         w = I * cacos(z);
> #else
>         w = clog(z + csqrt(z + 1) * csqrt(z - 1));
> #endif
>         return w;
> 
> it should be remarked that, to provide the correct principal value using
> just cacos(z), the above mentioned link addresses the need of testing the
> cacosh argument in order to appropriately use either  j*cacos(z) or 
> -j*cacos(z). It is therefore likely that the fix to be chosen, if any,
> should be based on computational efficiency.
> 
> Once more math is not democracy, so the answer must be left to the math
> savviest.
> 
> Regards, Paolo Mantegazza.