Tim Peters <t...@python.org> added the comment:
Well, that can't work: the most likely result for a long input is 0.0 (try
it!). frexp() forces the mantissa into range [0.5, 1.0). Multiply N of those,
and the result _can_ be as small as 2**-N. So, as in Mark's code, every
thousand times (2**-1000 is nearing the subnormal boundary) frexp() is used
again to force the product-so-far back into range. That's straightforward when
going "left to right".
With fancier reduction schemes, "it varies". Aiming for "obviously correct"
rather than for maximal cleverness ;-) , here I'll associate each partial
product with an integer e such that it's guaranteed (in the absence of
infinities, NaNs, zeroes), abs(partial_product) >= 2^^-e. Then quick integer
arithmetic can be used in advance to guard against a partial product
underflowing:
def frpair(seq):
from math import frexp, ldexp
stack = []
exp = 0
for i, x in enumerate(seq, start=1):
m, e = frexp(x)
exp += e
stack += [(m, 1)]
while not i&1:
i >>= 1
(x, e1), (y, e2) = stack[-2:]
esum = e1 + e2
if esum >= 1000:
x, e = frexp(x)
exp += e
y, e = frexp(y)
exp += e
esum = 2
stack[-2:] = [(x * y, esum)]
total = 1.0
totale = 0
while stack:
x, e = stack.pop()
totale += e
if totale >= 1000:
total, e = frexp(total)
exp += e
x, e = frexp(x)
exp += e
totale = 2
total *= x
return ldexp(total, exp)
But I see no obvious improvement in accuracy over "left to right" for the
uniformly random test cases I've tried.
----------
_______________________________________
Python tracker <rep...@bugs.python.org>
<https://bugs.python.org/issue41458>
_______________________________________
_______________________________________________
Python-bugs-list mailing list
Unsubscribe:
https://mail.python.org/mailman/options/python-bugs-list/archive%40mail-archive.com
