HP41: accuracy of 13-digit routines

[Dieter]

(09-02-2015 11:19 AM)Michael Fehlhammer Wrote:  ... Therefore you should not assume that results of log, sin etc. yield more than 10 correct digits.

Wait a minute. There is a chapter on accuracy in the HP15C Advanced Functions Handbook. The internal algorithms of the 41C and 15C arithmetic and transcendental functions should be very similar at least, if not identical, so let's assume the following also applies to the 41 series. The appendix (chapter "A Hierarchy of Errors") distinguishes several error levels, among them these three:

Quote:Operations that deliver "correctly rounded" results whose error cannot exceed ½ unit in their last (10th) significant digit include the real algebraic operations +, –, *, ÷, x², √x, 1/x, and %, (...) and →H.MS. These results are the best that 10 significant digits can represent (...).

Operations that can suffer a slightly larger error, but still significantly smaller than one unit in the 10th significant digit of the result, include ∆%, →H, →RAD, →DEG, Px,y and Cx,y, LN, LOG, 10^x, ASIN, ACOS, ATAN (...) and finally SIN, COS, and TAN for real arguments in Degrees and Grads modes (...).

A function that grows to ∞ or decays to 0 exponentially fast as its argument approaches ±∞ may suffer an error larger than one unit in its 10th significant digit, but only if its magnitude is smaller than 10^–20 or larger than 10^20; and though the relative error gets worse as the result gets more extreme (small or large), the error stays below three units in the last (10th) significant digit. (...) Functions so affected, are e^x, y^x, x! (for non-integer x) (...).

So the accuracy level depends on the function, and sometimes even on its argument. Since the four basic arithmetic functions as well as x² and √x deliver correctly rounded results they must be virtually exact internally, i.e. in all 13 internal digits, maybe the last one truncated. This would yield a max. error of 0,501 ULP.

The logs, 10^x and the trig functions (in DEG/GRAD) and their inverses may exhibit a slightly larger error, yet "significantly less" than 1 ULP. Maybe we can translate this to 0,6 ULP or 1 unit in the 11th place, or a relative error within 1E–11, maybe 0,7 ULP (rel. error < 2E–11). So there are 11 exact or nearly-exact digits.

And finally the error may spill into the last of the 10 digits presented to the user, where the 10th digit may be off by 2 or even 3 units.

What does this mean for a Gamma implementation in 13-digit MCode? The functions used in the Lanczos algorithm are mostly +, –, x and / which are (nearly) exact. On the other hand the logs may only have 11 valid digits. So this is what can be expected: about 11 exact digits, and a final error similar to that of the internal log or trig functions.

(09-02-2015 12:03 PM)Ángel Martin Wrote:  Yep, and even if they were how would you output it using just the standard X-register approach?

That's a different story. The question is what accuracy level can be expected from the internal 13-digit routines. Since some functions may not return more than 11 valid digits, this seems to be the limit for MCode using these limited-accuracy functions.

What do you think?

Dieter