[issue2819] Full precision summation

[Mark Dickinson]

Mark Dickinson <dickinsm at gmail.com> added the comment:

Some comments/questions:

(1) It seems wasteful to wrap every addition in PyFPE_START/END_PROTECT, 
and to check for NaNs and infinities after every addition.  I'd wrap the 
whole thing in a single PyFPE_START/END_PROTECT, replace _math_sum_add 
with an inline addition, and just let NaNs and Infs sort themselves out.

If the result comes out finite (as it surely will in almost all 
applications), then all the summands were necessarily finite and there's 
nothing more to do.  If the result comes out as an infinity or NaN, you 
need to decide whether it's appropriate to return a NaN, an infinity, or 
to raise OverflowError or ValueError.  I'm not sure it's worth trying to 
do the right thing for all special value cases, but if you do want to 
follow 'spirit of IEEE 754' rules for special values, they should look 
something like this:

  (1) if the summands include a NaN, return a NaN
  (2) else if the summands include infinities of both signs, raise
      ValueError,
  (3) else if the summands include infinities of only one sign, return
      infinity with that sign,
  (4) else (all summands are finite) if the result is infinite, raise
      OverflowError.  (The result can never be a NaN if all summands
      are finite.)

Note that some sums involving overflow won't be computed correctly:  
e.g. [1e308, 1e308, -1e308] will likely sum to infinity instead of
returning 1e308.  I don't see any easy way around this, and it's 
probably not worth worrying about.

(2) The algorithm is only guaranteed to work correctly assuming IEEE 754 
semantics.  Python currently doesn't insist on IEEE 754 floating point, 
so what should happen on non IEEE-754 machines?

(3) Rather than duplicating the math module code in cmathmodule.c, why 
not have the complex version simply sum real parts and imaginary parts 
separately, using a version of the code that's already in mathmodule.c?

----------
nosy: +marketdickinson

__________________________________
Tracker <report at bugs.python.org>
<http://bugs.python.org/issue2819>
__________________________________