|Re: What is unique about 39?|
Message #12 Posted by Rodger Rosenbaum on 12 May 2007, 1:15 a.m.,
in response to message #11 by Palmer O. Hanson, Jr.
When we were doing comparative tests on machines back in the 1980's it was said that some machines (I think it was the HP's, but I'm not sure, and the results which follow suggest that the supposition was incorrect) had "smarter" y^x algorithms in the sense that if the exponent was an integer, then the machine did repeated multiplications rather than, it was said, taking the logarithm of y, multiplying by x and taking the antilog. The result quoted above was for a program with either the square root multiplied by itself or for the square root squared with the x^2 function. When I did the problem with the square root raised to the second power using the y^x function the result was 3.2137E-05 with 205 zeroes. I didn't mistype. The thousandths digit is really a three.
I wrote a program which compared the results for each integer calculated using the two methods. When the integer was 39 the result with the x^2 method was zero but was 1E-8 with the y^x method. The results were identical for the other 499 integer values. I extended the search for differences. The differences I found before I stopped were at 965, 2813, 1902, 3441 and 3446.
If the the HP41 and the HP-12C were using the "smarter" y^x algorithm then I would have thought that there would have been no differences. I haven't identified anything unusual about 39, 965, etc.
What is happening here is the same phenomenon that causes the HP-71 not to get the "expected" result in the Savage benchmark.
To take the case of 39. The Y^X function takes the natural log of Y, multiplies it by X and then takes the EXP of that result. This all takes place with 13 digit arithmetic on the HP41, HP15, etc.
Do this same calculation with exact arithmetic using a PC program like Derive, Maple, Mathematica, or such, so that you can be sure you're getting correct results. If you start with the 10 digit SQRT(39), take the 13 digit LN, multiply by 2, take the 13 digit EXP you should get exactly 38.99999999502. Notice the last 3 digits. Those digits will be dropped when the result is rounded to 10 digits, but they determine which way the rounding will go. If we round 38.99999999502 to 10 digits we get 40.00000000. Why didn't the HP41 get 40.00000000? When the calculation was done exactly, the result was 38.99999999502, but the HP41 must have gotten a 13 digit internal result a little less than 38.99999999500, and therefore it didn't round the result up. If the HP41's internal result was in error by only -3 ULP's, that would account for the error. This is perfectly excusable. It's why the calculator uses 3 more digits internally than are returned to the user. This way the internal errors are thrown away. It's only when the exactly correct 13 digit value is right on the borderline where it will round up or down that this kind of error occurs.
I wrote a Mathematica program to search for all the values resulting from the n SQRT 2 y^x sequence where the last 3 digits of the exactly correct 13 digit value are 495, 496, 497, 498, 499, 500, 501, 502, 503, 504, or 505, with n ranging from 1 to 100. The values are:
n Last 3 of 13 digits
Any of these could have exhibited the error, but the HP41 often gets it right. But if the error is going to happen, it will be one of these cases where the rounding boundary is very close.
I searched up to n=7000 and found a bunch more possible candidates for the error. The only way to be sure is to test them all on the HP41 to see if the error occurs. The additional ones I found that make the error are:
And, note that the error the HP41 makes isn't always negative; it can be in either direction. For example, the exact 13 digit result for n=3446 is 3446.000000497 so the HP41's internal result must have had the last 3 digits greater than 500, causing it to round up when it shouldn't have.
This can also happen on a Saturn machine, of course. It will happen when the internal 15 digit result has the last 3 digits very close to 500. I searched for these and got a list of candidates which I then tested on the HP50. The ones that made the error are:
The HP-33S should make the error with these values of n.
I can't observe what the HP41 is doing internally, but it's pretty easy to see what the internal long real results are on a Saturn machine. So, to show that my explanation is correct, here is what I found on an HP48G using Jim Donnelly's program to do 15 digit arithmetic. The exact result for n=415 should be 415.000000001495, but the HP48G gets 415.000000001502, and in fact when you execute 415 SQRT 2 y^x on the keyboard, you get 415.000000002 instead of 415.000000001. The HP48G got 502 for the last 3 digits of its internal result instead of the 495 it should have gotten, so it rounded up when it shouldn't have.