[Python-ideas] Float range class

Nathaniel Smith njs at pobox.com
Sat Jan 10 16:48:14 CET 2015


On Sat, Jan 10, 2015 at 8:37 AM, Andrew Barnert
<abarnert at yahoo.com.dmarc.invalid> wrote:
>
> Actually, here's what numpy does:
>
>     step = (stop-start)/float((num-1))
>     y = arange(0, num) * step + start
>     y[-1] = stop
>
> De-vectorize that, and it's exactly your too-naive algorithm (except for the
> off-by-1), with a simple fix for inexact stop.

Numpy in general doesn't spend a lot of time worrying about 1 ulp
errors versus 2 ulp errors or whatever. We try to avoid using
particularly buggy platform transcendental functions, but otherwise,
well, real computations accumulate small rounding errors in dozens of
places, so worrying about 1 ulp here or there is usually a waste of
time. We even use the naive summation algorithm in many cases (though
we recently started upgrading to pairwise summation whenever it's
cheap to do so), and don't even have an implementation of Kahan
summation, never mind something as fancy as math.fsum. We'll probably
get around to it eventually, but it's not a high priority.

Basically what I'm saying is that numpy is not necessarily an
authority on the perfect high-precision implementation of little
numeric cores :-).

I can say that AFAIK we've never gotten any bug reports about rounding
errors in linspace.

> From a quick check, the timing is almost identical (6.89us vs. 6.49us to
> enumerate linspace(1, 2, 5000).
>
> It doesn't _normally_ accumulate rounding errors. For example, whether I use
> the same test, or much smaller numbers, every 5th element, and sometimes the
> one immediately after, is off by one digit, but the intervening elements are
> identical. However, if you get down to denormal values, then the underflow
> _does_ accumulate. For example, linspace(0, 2e-323, 10) should give (0, 0,
> 5e-324, 5e-324. ...), and the smart version does that, but the naive
> version--and numpy--gives all 0s. And if you build a simple 4.2 fixed-point
> type and ask for linspace(0, .1, 100), again you get all 0s. On the other
> hand, the naive version avoids overflow problems at the other end of the
> scale, so maybe that's a wash?

Okay, I guess I spoke too soon :-) https://github.com/numpy/numpy/issues/5437

-n

-- 
Nathaniel J. Smith
Postdoctoral researcher - Informatics - University of Edinburgh
http://vorpus.org


More information about the Python-ideas mailing list