Tim Peters <[email protected]> added the comment:
> ...
> one at a time subtract a**2 (an argument squared) in descending
> order of magnitude
> ...
But that doesn't really help unless the sum of squares was computed without
care to begin with. Can do as well by skipping that but instead computing the
original sum of squares in increasing order of magnitude.
Cheapest way I know of that "seemingly always" reproduces the Decimal result
(when that's rounded back to float) combines fsum(), Veltkamp splitting, and
the correction trick from the paper:
def split(x, T27=ldexp(1.0, 27)+1.0):
t = x * T27
hi = t - (t - x)
lo = x - hi
assert hi + lo == x
return hi, lo
# result := hypot(*xs)
parts = []
for x in xs:
a, b = split(x)
parts.append(a*a)
parts.append(2.0*a*b)
parts.append(b*b)
result = sqrt(fsum(parts))
a, b = split(result)
parts.append(-a*a)
parts.append(-2.0*a*b)
parts.append(-b*b)
x = fsum(parts) # best float approx to sum(x_i ** 2) - result**2
result += x / (2.0 * result)
----------
_______________________________________
Python tracker <[email protected]>
<https://bugs.python.org/issue41513>
_______________________________________
_______________________________________________
Python-bugs-list mailing list
Unsubscribe:
https://mail.python.org/mailman/options/python-bugs-list/archive%40mail-archive.com