D issues are now tracked on GitHub. This Bugzilla instance remains as a read-only archive.
Issue 5628 - std.math unittest disabled - roundoff error in pow() on SSE2
Summary: std.math unittest disabled - roundoff error in pow() on SSE2
Status: RESOLVED FIXED
Alias: None
Product: D
Classification: Unclassified
Component: phobos (show other issues)
Version: D2
Hardware: x86_64 All
: P2 normal
Assignee: No Owner
URL:
Keywords: pull
Depends on:
Blocks:
 
Reported: 2011-02-20 14:25 UTC by Brad Roberts
Modified: 2019-12-25 23:37 UTC (History)
2 users (show)

See Also:


Attachments

Note You need to log in before you can comment on or make changes to this issue.
Description Brad Roberts 2011-02-20 14:25:05 UTC
The test either takes an enormous amount of time or it goes into an infinite loop somewhere.
Comment 1 Brad Roberts 2011-04-30 22:21:19 UTC
Reduced test case:

Must be built with debugging turned on, otherwise it doesn't loop:
  dmd -m64 -gc  bug-pow.d

module bug;

real pow(real x, ubyte n) @trusted pure nothrow
{
    real p = 1.0;
    ubyte m = n;

    switch (n)
    {
    default:
    }

    while (1)
    {
        if (n & 1)
            p *= x;
        n >>= 1;
        if (!n)
            break;
        x *= x;
    }
    return p;
}

int main()
{
    immutable real x = 46;
    immutable ubyte three = 3;
    assert(pow(x,three) == x * x * x);

    return 0;
}

Extracted from std/math.d, function:

typeof(Unqual!(F).init * Unqual!(G).init) pow(F, G)(F x, G n) @trusted pure nothrow
if (isIntegral!(F) && isIntegral!(G))

And it's following unittest.
Comment 2 Brad Roberts 2012-01-01 21:24:28 UTC
The bug in comment 1 is fixed.  There are 3 asserts left that fail w/in std.math:

   assert(pow(xd, neg2) == 1 / (x * x));
   assert(pow(xf, neg8) == 1 / ((x * x) * (x * x) * (x * x) * (x * x)));

   assert(feqrel(real.min_normal/8,real.min_normal/17)==3);

They're currently versioned out for x86_64
Comment 3 Don 2013-04-09 02:28:41 UTC
The remaining bug in comment 2 is just a rounding error.
The last bit of 1/ x*x  is different when the intermediate values are 80 bit reals, vs when they are 64 bit doubles.

It is a bug, but it's not a compiler bug, just a fairly minor Phobos one.

Dropping severity to normal, and changing to Phobos.
Comment 4 Martin Nowak 2013-11-14 12:17:23 UTC
The test is wrong because of excess precision.
To fix this bug we need to replace == with approxEqual and determine the allowed error, right?
Comment 5 Iain Buclaw 2013-11-29 04:10:40 UTC
(In reply to comment #4)
> The test is wrong because of excess precision.
> To fix this bug we need to replace == with approxEqual and determine the
> allowed error, right?

Right, I've come some other failing tests in gdc that I've had to tweak to allow an extra bit of rounding accuracy.

feqrel() should do the trick with this.
Comment 6 Vladimir Panteleev 2017-07-02 18:12:31 UTC
Here is the extracted test case from std.math:

void main()
{
    import std.math;

    immutable real x = 46;
    immutable float xf = x;
    immutable double xd = x;

    immutable short neg2 = -2;
    immutable long neg8 = -8;

    assert(pow(xd, neg2) == 1 / (x * x));
    assert(pow(xf, neg8) == 1 / ((x * x) * (x * x) * (x * x) * (x * x)));
}
Comment 7 Dlang Bot 2019-12-15 10:41:04 UTC
@berni44 created dlang/phobos pull request #7321 "Fix Issue 5628 - std.math unittest disabled - roundoff error in pow()" fixing this issue:

- Fix Issue 5628 - std.math unittest disabled - roundoff error in pow() on SSE2

https://github.com/dlang/phobos/pull/7321
Comment 8 Dlang Bot 2019-12-25 23:37:33 UTC
dlang/phobos pull request #7321 "Fix Issue 5628 - std.math unittest disabled - roundoff error in pow()" was merged into master:

- 7b9c47452152217dccc83ad8a35945ee472dda1a by Bernhard Seckinger:
  Fix Issue 5628 - std.math unittest disabled - roundoff error in pow() on SSE2

https://github.com/dlang/phobos/pull/7321