Moritz Buhl wrote: > ... I noticed that some floating point operations cause failures of other > tests. > ... > Many edge cases for complex floating point operations are not covered at all.
Hi, https://marc.info/?l=openbsd-tech&m=150737856618497&w=2 is another example of an edge case for complex floating point operations. https://github.com/wch/r-source/blob/trunk/src/main/complex.c#L452-L455 gives a solution by checking if the imaginary part of the input complex number is too large (as otherwise sinh() is called which grows exponentially (see e.g. https://www.wolframalpha.com/input/?i=sinh(x) ) resulting in an overflow.) Note that the ctan() implementation in R is under GPL, so I am unsure if the check can be taken as is and committed to OpenBSD. s_ctanf.c probably needs a similar treatment. Best regards, Ingo Index: s_ctan.c =================================================================== RCS file: /cvs/src/lib/libm/src/s_ctan.c,v retrieving revision 1.7 diff -u -p -r1.7 s_ctan.c --- s_ctan.c 12 Sep 2016 19:47:02 -0000 1.7 +++ s_ctan.c 11 Jul 2019 12:31:41 -0000 @@ -135,9 +135,11 @@ double complex ctan(double complex z) { double complex w; - double d; + double d, wy, x, y; - d = cos (2.0 * creal (z)) + cosh (2.0 * cimag (z)); + x = 2.0 * creal(z); + y = 2.0 * cimag(z); + d = cos(x) + cosh(y); if (fabs(d) < 0.25) d = _ctans (z); @@ -148,7 +150,12 @@ ctan(double complex z) return (w); } - w = sin (2.0 * creal(z)) / d + (sinh (2.0 * cimag(z)) / d) * I; + if (isnan(y) || fabs(y) < 50.0) + wy = sinh(y) / d; + else + wy = (y < 0 ? -1.0 : 1.0); + + w = sin(x) / d + wy * I; return (w); } DEF_STD(ctan);
