Today we are the 2nd of April. No more joke, the Boss said "I want the culprit, now!".
We asked the Experts, they said: "The HP BCD math [in Spice/Voyager/Coconut machines] implements
arithmetic algorithms described as '
infinite precision truncated to 13 digits', then rounded
once to the user 10-digit result. All is fully deterministic, no approximation, no random rounding error, every machine is giving
exactly the same result for a given operation, and the result is the best representation of the result with 10 digits. No doubt, no any piece of randomness, no place for different answers of any kind. Period".
So we are taking our list again, and assuming the Experts are right, we take as an hypothesis, in
apparent contradiction with the OP's claim, that the arithmetic is right and without flaw. What can we still get out of our list of suspects?
Let's skip over +, -, *, /, 1/x, x².
Sqrt is basically (according to the same Experts) similarly implementing a pseudo-division (as the oldest of us may have learn at school loooong ago) and thus out of suspicion.
Let's forget %, INT and FRAC too.
PI is a constant, and we can easily check that its numeric value 3.141592654 is the best 10-digit representation of pi.
What is remaining now? The >DEG and >RAD function that are converting between radians and degrees using a constant ratio.
Can we really suspect them ? Well, "
when you have eliminated the impossible, whatever remains, however improbable, must be the truth".
Let's see, the converting ratio is (pi/180) or (180/pi) of course, depending on the direction the conversion is done.
How can we check it?
Well, 1 >DEG gives 57.29577951, the best 10-digit representation of (180/pi), and
1 >RAD gives 0.01745329252 , the best 10-digit representation of (pi/180).
All this is perfectly correct but let's try now to do 7 >DEG. We should get 7*180/pi and we indeed get 401.0704566, the best representation of this value with 10 digits. Now notice that this is different from the user-level calculation (RPN style) 7 180 * PI / that gives 401.070456
5 with an error of 1 ULP.
This proves that the >DEG and >RAD functions are using a high accuracy ratio value, presumably an internal 13-digit number, to guarantee the result.
Is it possible that some machines have a ratio value slightly wrong in the last places, that shows up only with certain numbers? How can we check it?
If this hypothesis is correct, then the problem would not be an arithmetic issue as such, but a simple programming flaw with no mystery.
Now that we have a plausible suspect, it's time to go back to the crime scene to look for evidences.
However, I need to stop my investigations for now, because I have a
sextuple integral to compute for tomorrow morning :-)
J-F