C99/POSIX define remquo:
https://pubs.opengroup.org/onlinepubs/9799919799/functions/remquo.html
The basic idea is that return value r is a remainder, and in q is stored
the low-order k bits of the quotient, but with q having the same sign as
the quotient.
The point is to reduce possibly large angles before using sin/cos.
I'm doing all this on netbsd-10, amd64, system gcc 10, intel 9th
generation i7. The following test program calls remquo with (-90, 90).
This is from a test in geodesiclib in proj of a high-accuracy sincos
implementation, which reduces angles to [45,-45] and then adjusts sign
by quadrant. But, it can be understood from the remquo specification.
The remainder is 0, and that's ok (-0, but not worrying about that).
The true quotient value is -1, and the low order bits should be
not only congruent to -1, but < 0.
The test program calls remquo twice, one with the -90 in a double, and
once with it as a literal. The first time gets a positive q -- which
messes up the sin calculation, wwhich is how I got here. The second is
correct.
Adding -ffloat-store didn't change the output.
I got the same behavior on netbsd-9 earmv7hf-el.
src/lib/libm/src/s_remquo.c is somewhat hard to understand....
Does anyone consider themselves the guardian of numerical correctness?
You are hereby nerdsniped!
----------------------------------------
#include <math.h>
#include <stdio.h>
int
main()
{
double x = -90.0;
double r; int q = 0;
r = remquo(x, 90, &q);
printf("remquo(x == -90.0, 90) r %f q %d", r, q);
if (q > 0)
printf("\tERROR: q positive");
printf("\n");
r = remquo(-90.0, 90, &q);
printf("remquo(-90.0, 90) r %f q %d", r, q);
if (q > 0)
printf("\tERROR: q positive");
printf("\n");
}
----------------------------------------
remquo(x == -90.0, 90) r -0.000000 q 1 ERROR: q positive
remquo(-90.0, 90) r -0.000000 q -1
----------------------------------------
