[Numpy-discussion] arange(start, stop, step) and floating point (Ticket #8)

Tim Hochberg tim.hochberg at cox.net
Wed Feb 8 15:41:31 CST 2006

Sasha wrote:

>On 2/8/06, Robert Kern <robert.kern at gmail.com> wrote:
>> ...
>>arange() does allow for fractional steps unlinke range(). You may fix the
>>docstring if you like. However, I don't think it is possible to ensure that
>>invariant in the face of floating point. That's why we have linspace().
>There is certainly a way to ensure that arange(..., stop, ...)[-1] <
>stop in the face of floating point -- just repeat start += step with
>start in a volatile double variable until it exceeds stop to get the
>length of the result.   There might be an O(1) solution as well, but
>it may require some assumptions about the floating point unit.
Isn't that bad numerically? That is, isn't (n*step) much more accurate 
than (step + step + ....)? It also seems needlessly inefficient; you 
should be able to to it in at most a few steps:

length = (stop - start)/step
while length * step < stop:
    length += 1
while length * step >= stop:
   length -= 1

Fix errors, convert to C and enjoy. It should normally take only a few 
tries to get the right N.

I see that the internals of range use repeated adding to make the range. 
I imagine that is why you proposed the repeated adding. I think that 
results in error that's on the order of length ULP, while multiplying 
would result in error on the order of 1 ULP. So perhaps we should fix 
XXX_fill to be more accurate if nothing else.

>In any case, I can do one of the following depending on a vote:
>1 (default). Document length=ceil((stop - start)/step) in the arange docstring
That has the virtue of being easy to explain.

>2. Change arange to be a fast equivalent of array(range(start, stop,
>step), dtype).
No thank you.

>3. Change arange to ensure that arange(..., stop, ...)[-1] < stop.

I see that Travis has vetoed this in any event, but perhaps we should 
fix up the fill functions to be more accurate and maybe most of the 
problem would just magically go away.


More information about the Numpy-discussion mailing list