float ranges: mismatch between length and bounds error

[StefanKarpinski]

julia> r = 2.0^53-2:2.0^53+3
9.00719925474099000e+15:9.00719925474099600e+15

julia> r[1]
9.00719925474099000e+15

julia> r[2]
9.00719925474099100e+15

julia> r[3]
9.00719925474099200e+15

julia> r[4]
9.00719925474099200e+15

julia> r[5]
9.00719925474099400e+15

julia> r[6]
9.00719925474099600e+15

julia> r[7]
9.00719925474099600e+15

julia> r[8]
9.00719925474099600e+15

julia> r[9]
BoundsError()
in ref, /Users/stefan/projects/julia/j/range.j:76
in run_repl, /Users/stefan/projects/julia/j/client.j:22
in _start, /Users/stefan/projects/julia/j/client.j:134

julia> length(r)
7

[JeffBezanson]

There are also corner cases at the other end of the dynamic range:

julia> r = 1.0:1e-19:2.0
1.0:9.99999999999999975e-20:2.0

julia> r[1]
1.0

julia> r[2]
1.0

julia> length(r)
0

[JeffBezanson]

Maybe the fix is to generalize the step == 0 check to step/start < eps(start) to disallow such ranges.

[StefanKarpinski]

So the solution to this is just to disallow the range 2.0^53-2:2.0^53+3? That seems like an incredible copout to me. It seems to me that the correct answer is to have a range object that produces the following sequence of values:

julia> for x = float(2^53-2:2^53+3)
         println(x)
       end
9.00719925474099000e+15
9.00719925474099100e+15
9.00719925474099200e+15
9.00719925474099200e+15
9.00719925474099400e+15
9.00719925474099600e+15

[StefanKarpinski]

I find it highly questionable that range lengths are mapped to Int. This means that on a 32-bit machine, we have length(-2.0^53:2.0^53) == 1. That's very broken.

[JeffBezanson]

OK, Float64 ranges (or all ranges) can always have Int64 lengths if you want.

I'm worried about the change to done, especially for Range1. It seems like it might be a lot slower. Maybe it's time to cache the lengths of ranges.

Ranges where the step is less than eps(start) are questionable: their lengths can overflow 64-bit counters, and they can have adjacent elements that are equal. But if your fix performs well and you're happy with it I'll go along with it.

[StefanKarpinski]

Computing the length upon construction might make a lot of sense. I agree that ranges with adjacent elements that are equal are weird, but it seems like the best solution. Making the range 2.0^53-2:2.0^53+3 illegal seems kind of bogus. Another argument for storing length explicitly is something like this:

julia> r = -2:3
-2:3

julia> length(r)
6

julia> r += 2.0^53
9.00719925474099000e+15:9.00719925474099600e+15

julia> length(r)
7

If the length of the former was part of its state, then this could work sanely but copying the length of the range. Maybe we should store the start, step and length or a range instead. I think this came up once before and we felt that keeping the stop value around was a good idea to avoid confusion. But maybe we shouldn't bother and just compute the last value for display and make sure that if a range is input as it is printed, you'll get back the same length range.

[JeffBezanson]

Things are definitely better now. But we still have stuff like

julia> length(0:1e-20:1)
0

Even with 64-bit lengths, you can't handle everything.

I'm not sure it's so bad to use Int32 for length on 32-bit systems. You don't generally expect to have more than 2^32-1 of anything on a 32-bit system.

[StefanKarpinski]

See 44b1aaf and e35cc22.

[JeffBezanson]

Can we close this now?