Bug 3866 - pow(double, double) does not fulfil several of C99's requirements
Summary: pow(double, double) does not fulfil several of C99's requirements
Status: RESOLVED FIXED
Alias: None
Product: glibc
Classification: Unclassified
Component: math (show other bugs)
Version: 2.4
: P2 normal
Target Milestone: ---
Assignee: Not yet assigned to anyone
URL:
Keywords:
Depends on:
Blocks:
 
Reported: 2007-01-12 23:30 IST by Richard B. Kreckel
Modified: 2014-05-28 19:46 IST (History)
3 users (show)

See Also:
Host: i386-linux-gnu
Target:
Build:
Last reconfirmed:


Attachments
test program (2.75 KB, text/plain)
2008-08-05 14:05 IST, Michael Kerrisk
Details

Note You need to log in before you can comment on or make changes to this bug.
Description Richard B. Kreckel 2007-01-12 23:30:54 IST
Here is a list of bugs in GlibCs pow implementation.


1) pow(-1, 1e100):

Returns NaN and raises the "invalid" floating-point exception.

It should return 1 and raise no exception. (This is because all large
floating-point numbers are even integers if the radix is even. This
motivates pow(-1, inf) = +1, which is treated correctly by GlibC).


2) pow(-1, -1e100):

Same as above.


3) pow(-3.141592, 1e100):

Returns NaN and raises the "invalid" floating-point exception.

It should return +inf and raise the "overflow" floating-point exception flag.


4) pow(-3.141592, -1e100):

Returns NaN and raises the "invalid" floating-point exception.

It should return +0 and raise the "underflow" floating-point exception flag.



4) pow(-0, -1e100):

Raises the "invalid" floating point exception.

It shouldn't do that: this is just a case of exact division by
zero. Division by zero is correctly flagged and +inf is returned, but
the "invalid" flag is wrong.



5) pow(-inf, -1e100):

Raises the "invalid" floating point exception.

It shouldn't do that. After all, it correctly returns 0.


6) pow(-inf, 1e100):

Raises the "invalid" floating point exception.

It shouldn't do that. After all, it correctly returns +inf.


7) pow(-0, -9.0071991544401990e+15)

This returns +inf.

It should return -inf, since the double precision exponent is an odd number
(fabs(fmod(-9.0071991544401990e+15,2.0))==1.0 exactly).
Comment 1 Michael Kerrisk 2008-08-05 14:05:57 IST
Created attachment 2892 [details]
test program

Just to add a little detail here regarding these four points:

1) pow(-1, 1e100):
2) pow(-1, -1e100):
3) pow(-3.141592, 1e100):
4) pow(-3.141592, -1e100):

If x is less than -0, then (on x86, glibc 2.8) pow() fails with EDOM/FE_INVALID
for any y whose absolute value is greater than about 9.223373e18. 
(EDOM/FE_INVALID == pole error, which is the error that occurs if x == 0, and y
is negative)
Comment 2 Michael Kerrisk 2008-08-05 14:19:15 IST
Regarding:
4) pow(-0, -1e100)

this problem seems to be absent on glibc 2.8.  There, pow(-0, -1e100) does the
right thing, returning -Inf and giving EDOM/FE_DIVBYZERO.
Comment 3 Michael Kerrisk 2008-08-05 14:21:10 IST
Regarding:
5) pow(-inf, -1e100):
6) pow(-inf, 1e100):

The problem persists in glibc 2.8, and again the cut-off y value that
demonstrates this behavior is |y| greater than about  9.223373e18.
Comment 4 Richard B. Kreckel 2008-08-05 22:04:52 IST
(In reply to comment #2)
> 4) pow(-0, -1e100)
> 
> this problem seems to be absent on glibc 2.8.  There, pow(-0, -1e100) does the
> right thing, returning -Inf and giving EDOM/FE_DIVBYZERO.

But the right thing to do would return +Inf and giving EDOM/FE_DIVBYZERO, as
-1e100 is < 0 and not an odd integer (specified in. F.9.4.4, second dash.)
Comment 5 Michael Kerrisk 2008-08-08 05:48:08 IST
(In reply to comment #4)
> (In reply to comment #2)
> > 4) pow(-0, -1e100)
> > 
> > this problem seems to be absent on glibc 2.8.  There, pow(-0, -1e100) does 
the
> > right thing, returning -Inf and giving EDOM/FE_DIVBYZERO.
> But the right thing to do would return +Inf and giving EDOM/FE_DIVBYZERO, as
> -1e100 is < 0 and not an odd integer (specified in. F.9.4.4, second dash.)

Oops.  Yes, of course you are right.  And glibc 2.8 is doing the right thing, 
I just wrote the wrong thing. pow(-0, -1e100) gives +Inf with 
EDOM/FE_DIVBYZERO.
Comment 6 Joseph Myers 2012-02-28 17:00:42 IST
Confirmed still present on x86: (1), (2), (3), (4a), (4b), (5), (6), (7).

On x86_64: (1), (2) pass, (3) has the right result but no OVERFLOW exception, (4a) has the right result but no UNDERFLOW exception, (4b) has the right result but no DIVBYZERO exception, (5), (6), (7) pass.

There are two tests numbered (4) which I have marked (4a) and (4b).  I do not reproduce with current sources the report of (4b) pow(-0, -1e100) passing.
Comment 7 Joseph Myers 2012-03-21 12:20:27 IST
(4b), (5), (6), (7) fixed by:

commit 2460d3aa21f04cdf28497683bd3e29183189f779
Author: Joseph Myers <joseph@codesourcery.com>
Date:   Wed Mar 21 12:14:57 2012 +0000

    Fix pow of zero and infinity to large powers.

(1), (2), (3), (4a) remain to be fixed.
Comment 8 Joseph Myers 2012-03-28 15:03:59 IST
Remaining cases fixed by:

commit d6270972f79fe89a96fa7a3909991dad2e317033
Author: Joseph Myers <joseph@codesourcery.com>
Date:   Wed Mar 28 14:57:58 2012 +0000

    Fix pow of negative numbers to integer exponents (bugs 369, 2678, 3866).
Comment 9 Jackie Rosen 2014-02-16 18:27:12 IST Comment hidden (spam)