On 02/05/2012 11:08 AM, James Courtier-Dutton wrote:
Hi,
I looked at this a bit closer.
sin(1.0e22) is outside the +-2^63 range, so FPREM1 is used to bring it
inside the range.
So, I looked at FPREM1 a bit closer.
#include<stdio.h>
#include<math.h>
int main (void)
{
long double x, r, m;
x = 1.0e22;
// x = 5.26300791462049950360708478127784;<- This is what the answer
should be give or take 2PI.
m = M_PIl * 2.0;
r = remainderl(x, m); // Utilizes FPREM1
printf ("x = %.17Lf\n", x);
printf ("m = %.17Lf\n", m);
printf ("r = %.17Lf\n", r);
return 1;
}
This outputs:
x = 10000000000000000000000.00000000000000000
m = 6.28318530717958648
r = 2.66065232182161996
But, r should be
5.26300791462049950360708478127784... or
-1.020177392559086973318201985281...
according to wolfram alpha and most arbitrary maths libs I tried.
I need to do a bit more digging, but this might point to a bug in the
cpu instruction FPREM1
Kind Regards
James
As I recall, the remaindering instruction was documented as using a
66-bit rounded approximation fo PI, in case that is what you refer to.
--
Tim Prince
