Calculating e^x-1 on classic HPs

[Dieter]

(01-16-2016 12:40 PM)Gerson W. Barbosa Wrote:  This appeared to work, but some cancellation occurs at a critical region:
(...)
What's a few ULPs among friends anyway? :-)

Or a few hundred, or a few thousand... ;-)

(01-16-2016 12:40 PM)Gerson W. Barbosa Wrote:  (...)
That should be expected by just looking at the identity I have used:

\[e^{x}-1=\frac{e^{x}-e^{-x}}{e^{-x}+1}\]

It works better after replacing e^x – e^–x with 2 sinh(x) – provided an accurate sinh function is available. I tried this on the 34s and in 100.000 random samples there were no errors beyond ±7 ULP. About 75% of the results were within ±1 ULP.

Finally I did a test run with one million random values in [ln 0.9, ln 2] with the method I suggested, and the results are essentially the same as before:

Code:

seed = 47,11   n = 1000000

ULP      #cases
--------------------------
 -5          12 ─── 0,001%
 -4        1092 ─┐
 -3         674  ├─ 0,24%
 -2         676 ─┘
 -1      167703 ─┐
 ±0      659710  ├─ 99,5%
 +1      167784 ─┘
 +2         666 ─┐
 +3         677  ├─ 0,23%
 +4         992 ─┘
 +5          14 ─── 0,001%

Again, there are two peaks at +4 and –4 ULP that also showed up in earlier tests. I wonder where these come from. Any idea?

As already mentioned, for arguments in [ln10, ln11[ or [ln100, ln101[ etc. the last digit is lost, which leads to an interesting error distribution. Here are some results for x in [ln10, ln11]:

Code:

seed = 47,11    n = 100000

ULP      #cases
--------------------------
 -5       4990 ─── 5%
 -4      10025 ─┐
 -3       9952  ├─ 30%
 -2       9975 ─┘
 -1      10002 ─┐
 ±0      10024  ├─ 30%
 +1      10170 ─┘
 +2       9888 ─┐
 +3       9941  ├─ 30%
 +4      10042 ─┘
 +5       4991 ─── 5%

As expected, the errors are evenly distributed in this interval. Again, they do not exceed ±5 ULP.

Dieter