Mail Archives: djgpp-workers/1999/06/27/10:26:28
Sorry, it seems we are not there yet: I ran the Cygnus test suite with
the new math functions, and found a few additional problems.
- The NaN is loaded by the functions from a long whose bit pattern
is 0xFFC00000. This produces a negative NaN, with the sign bit
set. I know that the distinction is silly, but "%+f" will print
it as "-NaN", which isn't nice.
I think the bit pattern should be 0x7ff80000. Do you agree?
- A more serious problem is with asin and acos for arguments whose
absolute value is between 1 and 1.5 (inclusive). These manage to
pass the tests for badarg, wind up in fpatan and return -NaN,
which has the wrong sign, and errno is not set.
- atan2 returns NaN when both arguments are +/-Inf. The docs keeps
silent about this case (atan2 from libm only returns NaN when y is
Inf). It seems to me that if y is finite, x can be Inf and the
x87 will produce a finite result, so perhaps we need to modify the
checks of x?
- atanh(+/-1.25) returns -/+10 instead of a NaN. This seems like
another case of problems with checking arguments whose absolute
value is between 1 and 2, like in the case of acos and asin.
- cos suffers a large loss of significant digits near Pi/2: for
example cos(-1.571) loses about 40 bits of the mantissa (they come
out as zero). The same happens to sin near Pi. Is this some flaw
of the x87? If so, perhaps we should instead compute sin(Pi/2 - x)
and sin(Pi - x) for such arguments?
- cosh doesn't set errno to ERANGE for an argument that is +/-Inf,
unlike the docs says. I think the problem is that the flow jumps
from abarg: to infarg:, but the latter doesn't set errno to 2.
- ldexp sets errno when its first argument is -0.0. I guess we have
another case of +0 vs -0 controversy to handle. Seems like
clearing the high bit when testing x for zero should do the trick.
- modf doesn't set the sign bit of the integral part if it is zero.
For example, -0.81 sets the integral part to +0, not -0. I'm not
sure if we should bother.
- when the first argument of modf is NaN or Inf, perhaps we should
set the integral portion to NaN and Inf, respectively, not to zero
like the code now does, and return zero. At least for the Inf
argument, it makes more sense to say that its integral part is
infinite, not its fractional part.
- what should pow(-0.,+1) return? your last version returns +0,
which might be okay, but in that case we should document it, since
it contradicts the rule that x^1 == x for any x. The same goes
for pow(-0.,-1): should it be +Inf or -Inf?
Btw, libm functions consistently do NOT set errno for NaN arguments
(except for pow(NaN,0)). Our behavior is well-documented, so this
is only FYI. But I wonder whether some widely-accepted standard
rules this.
- Raw text -