math.exp underflow (was: exp behaviour: is it a bug ?)

[Tim Peters]

[George N. White]
> There are a couple issues in this thread:
>
> 1)  cross-platform consistency
>
> 2)  treatment of exceptions: underflow, overflow, invalid, ...
>
> For 1), Sun's freely distributable math library should eliminate the
> effects of differences in vendor-supplied math libraries

Except that, as you say later, "speed wins", and fdlibm is almost certainly
slower than the native libm.  I happened to be writing a math library for a
HW vendor at the time K.C.Ng was writing fdlibm, and by exploiting platform
quirks I was able to get better accuracy at 10x the speed.  There's really
no good portable & safe way to deal with the mass of 754 special cases
without casting tricks and mounds of conditionals.

> For 2), support for reliable reporting of f.p. exceptions when they
> occur conflicts with techniques (e.g., combined multiply+add: y=a*x+b)
> used to get better f.p. performance in new processor designs.  Since
> the marketplace favors speed over all else, the conflict will likely
> be resolved in favor of performance.  This is reflected in the lack of
> support for f.p. exceptions in language standards.  It is much faster
> and easier to simply store a 0, +/-Inf, or NaN and keep going.

Stangely enough, on a Pentium Inf and NaN inputs cause + - * / to take up to
10x longer than normal operands.  I've never seen that documented, BTW -- it
was an empirical observation after tracking down why an FFT took about 5x
longer than I thought it should (turned out that the test harness generated
Infs and NaN inputs by mistake).

> Until we have benchmarks that explicity test underflow and overflow
> behaviour, vendors are going to focus on performance for the simplest
> cases.

Yup.  The saving grace for Python is that it's so slow on fp arithmetic
anyway that nobody would notice if it took the time to do it right <0.6
wink>.

> In practical terms, users need to a) write code that doesn't rely
> on exception handlers

You can always get at them, albeit in wildly differing platform-specific
ways.  I'd like Python to abstract that and so make it usable via one
spelling.

> and b) push for support in language stds. for features needed to work with
> signed 0, NaN's and Inf's in (e.g.: don't optimize away expressions like
> 'x != x').

Note that the Numerical C Extension Group's 754 recommendations were adopted
pretty much in whole by C9X ("the next" iteration of the C std, and already
ISO's current official version).

> As it happens, just this week I encountered a subroutine that was
> returning NaN's from a calculation of the form: y*log(y) where y=exp(-x)
> and x was large enough to cause underflow, giving 0*(-Inf)=NaN.  I tried
> trapping the underflow, but there were _thousands_ of log(y) expressions
> used in a way that did not cause any difficulty before I reached the one
> that was returning a NaN.  In this case, trapping underflows gave so many
> false positives that it didn't help identify the real problem area in the
> code.

Why didn't you trap invalid operation instead?  Every 754 operation that
produces a NaN from non-NaN inputs (like 0 and -Inf) signals inv. op.

> Perhaps a better approach for Johann's specific example would be to use a
> special non-underflow function, exp_p(-x), that calls an error routine
> when x exceeds some threshold, or add a test before exp(-x) is called.

W. Kahan's pithy definition of an exceptional operation is any behavior that
someone may take exception to <wink>.

math-libraries-are-both-the-best-defined-and-least-agreed-upon-ly y'rs  -
tim