Corrects offset1 computing of SubArray for IdOffsetRange

[johnnychen94]

I believe this patch somehow fixes issue #100 but it needs careful review since I might misunderstand the design of SubArray and IdOffsetRange...

issue #100 seems only happens when it's a linear index

P.S. Tests labeled with issue #100 were broken without this patch

cc: @mcabbott

[johnnychen94]

Hmmm. the failures in the nightly test belong to another issue.

[johnnychen94]

JuliaLang/julia#34049 break the test here, looks like it's another issue wrt linear index.

[timholy]

Really amazing work! I'm super-impressed.

My inclination would be to mimic https://github.com/JuliaLang/julia/blob/4d9a3349170d9beea8584cc9a730d37609280ac3/base/subarray.jl#L354-L355, with the thinking that IdOffsetRange is another index-preserving case. Care to play with that?

[johnnychen94]

Before I dig deeper and confuse myself, is there any rule that says linear index for OffsetArray should start from 1?

Oh yes. https://docs.julialang.org/en/v1/devdocs/offset-arrays/#Linear-indexing-(LinearIndices)-1 It turns out it's time for me to step into the developer docs.

Regardless of the array's native indices, linear indices always range from 1:length(A)

A = zeros(-3:1, -2:0)
I = axes(A, 1)
A[I] # BoundsError

[timholy]

Those docs immediately go on to say

By this definition, 1-dimensional arrays always use Cartesian indexing with the array's native indices.

If it's not sufficiently clear we should reword it (suggestions welcome).

[johnnychen94]

"Obviously" this PR doesn't fix all the things.

julia> out = zeros(-2:1, -3:-1)
4×3 OffsetArray(::Array{Float64,2}, -2:1, -3:-1) with eltype Float64 with indices -2:1×-3:-1:
 0.0  0.0  0.0
 0.0  0.0  0.0
 0.0  0.0  0.0
 0.0  0.0  0.0

julia> out[1] = 1
1

julia> ind = axes(out, 1)
-2:1

julia> SubArray(out, (ind, 1)) # it should raise a BoundsError but it doesn't
4-element view(OffsetArray(::Array{Float64,2}, -2:1, -3:-1), -2:1, 1) with eltype Float64 with indices -2:1:
 2.3851344346e-314
 2.385134482e-314
 2.3851346243e-314
 2.3851333437e-314

Well, seems like it's a SubArray issue:

julia> SubArray(zeros(4, 4), (-4:-1, 1))
4-element view(::Array{Float64,2}, -4:-1, 1) with eltype Float64:
 2.0e-323
 2.0e-323
 2.33858836e-314
 2.5e-323

[fredrikekre]

I don't think the SubArray constructor is meant to be used like that, it completely skips argument checking (presumably for efficiency reasons). You should probably use view (or @view), which indeed throws a BoundsError:

julia> view(out, ind, 1)
ERROR: BoundsError: attempt to access 4×3 OffsetArray(::Array{Float64,2}, -2:1, -3:-1) with eltype Float64 with indices -2:1×-3:-1 at index [-2:1, 1]

[timholy]

Test failure on master is not your fault, so we can merge this and fix it in another PR. (EDIT: oh, I see you're way ahead of me, I saw #102.)

However...with apologies I'm now wondering if this should be fixed in Base (for older Julia versions we will need this, so this PR should be merged). The issue is this: if we fix it here, then other array types that define ranges with axes starting at values other than 1 will have to make the same fix.

Looking at https://github.com/JuliaLang/julia/blob/4d9a3349170d9beea8584cc9a730d37609280ac3/base/subarray.jl#L354-L355, now I'm noting that the method below it is identical in cases where first(inds) == 1. Since that's the default for ranges anyway, should we just delete the method below it?

The hard part will be writing a test. If you're interested in tackling this, you could copy much of OffsetArrays/src/axes.jl to test/testhelpers/OffsetArrays.jl if you wish. Or define something like

# Use a range type not defined in Base
struct StrangeRange <: AbstractUnitRange{Int}
    stop::Int
end
Base.first(r::StrangeRange) = -r.stop
Base.last(r::StrangeRange)  =  r.stop
Base.axes(r::StrangeRange)  = (0:2*r.stop,)
...

Note this is a case where the values and the axes are deliberately misaligned AND the first element of neither is 1. If that doesn't shake out bugs...

If you don't want to tackle Base, just let me know and I can do it.

[codecov-io]

Codecov Report

Merging #101 into master will decrease coverage by 0.46%.
The diff coverage is 0%.

@@            Coverage Diff             @@
##           master     #101      +/-   ##
==========================================
- Coverage   84.44%   83.97%   -0.47%     
==========================================
  Files           2        2              
  Lines         180      181       +1     
==========================================
  Hits          152      152              
- Misses         28       29       +1

Impacted Files	Coverage Δ
src/axes.jl	72.72% <0%> (-2.28%)	⬇️

---

Continue to review full report at Codecov.

Legend - Click here to learn more
Δ = absolute <relative> (impact), ø = not affected, ? = missing data
Powered by Codecov. Last update ff1acb3...f683dc4. Read the comment docs.

[johnnychen94]

Looking at JuliaLang/julia:base/subarray.jl@4d9a334#L354-L355, now I'm noting that the method below it is identical in cases where first(inds) == 1. Since that's the default for ranges anyway, should we just delete the method below it?

The benchmark says we'd better not do so...

julia> A = rand(5, 5);

julia> I = Base.IdentityUnitRange(3:5)
Base.IdentityUnitRange(3:5)

julia> @btime Base.compute_offset1(A, 1, (I, 1));
  120.822 ns (1 allocation: 32 bytes)

julia> @btime Base.compute_offset1(A, 1, (3:5, 1));
  30.544 ns (0 allocations: 0 bytes)

Oops, looks like I should do it with $:

julia> @btime Base.compute_offset1($A, 1, ($I, 1));
  1.420 ns (0 allocations: 0 bytes)

julia> @btime Base.compute_offset1($A, 1, (3:5, 1));
  1.421 ns (0 allocations: 0 bytes)