Date: Thu, 17 May 2018 18:11:53 +0200 From: Szabolcs Nagy <nsz@...t70.net> To: musl@...ts.openwall.com Subject: Re: What should be the result of CACOSH(F)(CCOSH(F)(-2 + 1j))? * Paolo Mantegazza <paolo.mantegazza@...imi.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.
Powered by blists - more mailing lists
Powered by Openwall GNU/*/Linux - Powered by OpenVZ