Floating point exception handling #38
Description
Problem
Currently, the DACE does not check for floating point exceptions in any of its calculations, nor are non-normal values taken into account when DA calculations are performed. This has some strange consequences, for example:
- Anything multiplied by a zero DA becomes zero, even if the other DA contains infinities (overflows)
- Tanh becomes zero due to this for large arguments causing overflow during sinh/cosh (see Wrong hyperbolic tangent for argument >> 0 #37)
- most certainly other functions may suffer from similar effects
Open Questions
This is a very difficult problem for which we first have to answer several very basic question before a solution can be implemented.
What constitutes an error?
To really solve this, we have to decide what FP exceptions actually mean in DA context. Is an overflow in order 20 reason to make the entire DA invalid? Do we care about underflows?
How to signal errors
Should these throw visible exceptions/errors the user must handle right there (like sqrt(-1)) or should DAs have a silent error flag? Should there be several error flags indicating different errors (inf, nan, ...)? How is this exposed to the different DACE interfaces?
How to propagate errors
How does each and every built-in operation propagate errors / non-normal DA input? E.g. 1/inf = 0 is perfectly valid in FP, but inf/inf yields a nan.
How do we detect these problems in the first place?
- Range checks for each function to make sure it can't overflow? How would this make sure the high order coefficients don't overflow either?
- Checking each and every coefficient in a result for non-normal? Surely this will be incredibly slow...
- Checking the FP environment flags after each DA operation for errors? Is this portable/reliable/fast? How does this interact with other libraries manipulating the FP environment (see here)