# Why is LinearIndices hardcoded to use Int?

**URL:** https://discourse.julialang.org/t/why-is-linearindices-hardcoded-to-use-int/17235
**Category:** General Usage
**Created:** [November 7, 2018, 1:25am UTC](https://discourse.julialang.org/t/why-is-linearindices-hardcoded-to-use-int/17235 "2018-11-07T01:25:17Z")
**Posts on this page:** 9
**Page:** 1

<div class="post-metadata">

### Author: ![Pbellive](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/pbellive/32/3604_2.png) [@Pbellive](https://discourse.julialang.org/u/Pbellive)
#### Post date: [November 7, 2018, 1:25am UTC](https://discourse.julialang.org/t/why-is-linearindices-hardcoded-to-use-int/17235/1 "2018-11-07T01:25:17Z")

</div>

Hi,  
The `LinearIndices` type is defined in Base at [line 382 of indices.jl](https://github.com/JuliaLang/julia/blob/0d4edb38b0b624421a77e5b89cfe22087345e9cd/base/indices.jl#L382) as

```julia
struct LinearIndices{N,R<:NTuple{N,AbstractUnitRange{Int}}} <: AbstractArray{Int,N}
    indices::R
end

```

I’m wondering why it’s defined to hold an `NTuple{N,AbstractUnitRange{Int}}}` rather than being defined along the lines of:

```julia
struct MyLinearIndices{N,Tn,R<:NTuple{N,AbstractUnitRange{Tn}}} <: AbstractArray{Tn,N}
    indices::R
end

```

The issue with the current implementation is that you can’t convert a Cartesian index to a linear index if the resulting linear index is too big to represented by an `Int`. For example (with `Int`=`Int64` on my system, running on current master branch julia):

```julia
julia> l = Int128(10_000_000);

julia> sz = (l,l,l)
(10000000, 10000000, 10000000)

julia> LI = LinearIndices(sz);

julia> LI[l,l,l]
3875820019684212736

julia> l^3
1000000000000000000000

```

I do need to convert cartesian indices to such large linear indices in practice. It’s part of a sparse 3D array implementation. I was able to get what I was looking for by extending a few methods for the `MyLinearIndices` type defined above:

```julia
MyLinearIndices(sz::NTuple{N,Tn}) where {N,Tn<:Integer} = MyLinearIndices(map(Base.OneTo, sz))

IndexStyle(::Type{<:MyLinearIndices}) = IndexLinear()
axes(iter::MyLinearIndices) = map(axes1, iter.indices)
size(iter::MyLinearIndices) = map(unsafe_length, iter.indices)
function getindex(iter::MyLinearIndices, i::Integer)
    @_inline_meta
    @boundscheck checkbounds(iter, i)
    i
end

function getindex(iter::MyLinearIndices, i::AbstractRange{<:Integer})
    @_inline_meta
    @boundscheck checkbounds(iter, i)
    @inbounds isa(iter, LinearIndices{1}) ? iter.indices[1][i] : (first(iter):last(iter))[i]
end

```

Using this new type gives me

```julia
julia> myLI = MyLinearIndices(sz);

julia> myLI[l,l,l]
1000000000000000000000

```

I don’t know the array indexing code very well so I don’t know whether or not there’s a good reason to force `LinearIndices` to use `Int`s rather than working with arbitrary integer types? Is there a good reason or was this just an oversight? If it was just an oversight then I can prepare a PR with a fix based on the `MyLinearIndices` stuff in this post. I believe it could count as a bugfix and not as a breaking change. On the other hand, if the current behaviour is intended I would be grateful if someone could explain the reasoning behind it a little.

Thanks very much, Patrick

---

<div class="post-metadata">

### Author: ![StefanKarpinski](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stefankarpinski/32/24_2.png) [@StefanKarpinski](https://discourse.julialang.org/u/StefanKarpinski)
#### Post date: [November 7, 2018, 4:36am UTC](https://discourse.julialang.org/t/why-is-linearindices-hardcoded-to-use-int/17235/2 "2018-11-07T04:36:25Z")

</div>

Linear indices generally represent offsets into regions of memory, which cannot be larger than the system integer type. Using linear indices is not very efficient otherwise since the linear index must be decoded before being used (involving shockingly slow integer division operations). Are you sure you need linear indices here? Can’t you use cartesian indices? Especially for sparse matrices, linear indexing is really inefficient.

---

<div class="post-metadata">

### Author: ![Pbellive](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/pbellive/32/3604_2.png) [@Pbellive](https://discourse.julialang.org/u/Pbellive)
#### Post date: [November 7, 2018, 6:03am UTC](https://discourse.julialang.org/t/why-is-linearindices-hardcoded-to-use-int/17235/3 "2018-11-07T06:03:20Z")

</div>

Thanks for the reply Stefan. That explanation makes sense. Let me give a bit more context on the application. As I said, this relates to a three dimensional sparse array implementation. I say array but it’s not technically a subtype of AbstractArray. These 3D “arrays”, called `SparseArray3D`, are part of a private project I’m involved with. In building and manipulating `SparseArray3D`s we work with Cartesian indices but under the hood the non-zero entries are stored in a dictionary keyed by the linear index. So say we have a sparse 3D array `N`. We’ll get a value from the array using

```julia
val = N[i,j,k]

```

where i,j,k will only be single integers, not ranges or anything else like that. `getindex` for `SparseArray3D`s will convert the i,j,k cartesian index to a linear index and grab the value from the dictionary. So the `LinearIndices` object is only ever used to convert a single Cartesian index to a linear index. Just a drop in replacement for `sub2ind`. Clearly we can do this manually. I was wondering though why it wasn’t possible using `LinearIndices[i,j,k]`. I see from your reply that there’s a good reason for this. Thanks again.

Regarding performance, we’ve found that for our typical access patterns our implementation is faster in practice than using a dictionary keyed by an (i,j,k) tuple or using a 3D analogue of the compressed sparse column/row format used for sparse matrices.

---

<div class="post-metadata">

### Author: ![mauro3](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/mauro3/32/292_2.png) [@mauro3](https://discourse.julialang.org/u/mauro3)
#### Post date: [November 7, 2018, 8:20am UTC](https://discourse.julialang.org/t/why-is-linearindices-hardcoded-to-use-int/17235/4 "2018-11-07T08:20:46Z")

</div>

Do you know why performance is worse with `(i,j,k)` keys? I thought that the hashing might be it, but no. Comparison is a bit slower but not enough to make up for the hashing:

```julia
julia> l = Int128(10_000_000); ll = Int(l); t = (ll,ll,ll);

julia> @btime hash($l);
  19.996 ns (0 allocations: 0 bytes)

julia> @btime hash($ll);
  2.965 ns (0 allocations: 0 bytes)

julia> @btime hash($t);
  11.252 ns (0 allocations: 0 bytes)

julia> @btime isequal($l,$l);
  0.018 ns (0 allocations: 0 bytes)

julia> @btime isequal($ll,$ll);
  0.018 ns (0 allocations: 0 bytes)

julia> @btime isequal($t,$t);
  1.968 ns (0 allocations: 0 bytes)

```

---

<div class="post-metadata">

### Author: ![Pbellive](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/pbellive/32/3604_2.png) [@Pbellive](https://discourse.julialang.org/u/Pbellive)
#### Post date: [November 7, 2018, 4:43pm UTC](https://discourse.julialang.org/t/why-is-linearindices-hardcoded-to-use-int/17235/5 "2018-11-07T16:43:43Z")

</div>

Henh. I’m not sure @mauro3. I’ll have to check with my colleague who did the benchmarking. I wasn’t actually involved in that. I had thought it was a matter of hashing performance but that doesn’t jive with what you show here.

---

<div class="post-metadata">

### Author: ![kristoffer.carlsson](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/kristoffer.carlsson/32/22_2.png) [@kristoffer.carlsson](https://discourse.julialang.org/u/kristoffer.carlsson)
#### Post date: [November 7, 2018, 4:46pm UTC](https://discourse.julialang.org/t/why-is-linearindices-hardcoded-to-use-int/17235/6 "2018-11-07T16:46:21Z")

</div>

> [@mauro3](#):
>
> ```julia
> julia> @btime isequal($ll,$ll); 0.018 ns (0 allocations: 0 bytes)
> 
> ```

FWIW, this is a benchmarking artifact (function call gets optimized away or constant folder). Also, benchmarking something like `isequal` between two Ints is almost useless because the timing of such a code will heavily depend on the context it is used in.

---

<div class="post-metadata">

### Author: ![StefanKarpinski](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stefankarpinski/32/24_2.png) [@StefanKarpinski](https://discourse.julialang.org/u/StefanKarpinski)
#### Post date: [November 7, 2018, 4:56pm UTC](https://discourse.julialang.org/t/why-is-linearindices-hardcoded-to-use-int/17235/7 "2018-11-07T16:56:57Z")

</div>

> [@Pbellive](#):
>
> under the hood the non-zero entries are stored in a dictionary keyed by the linear index.

Would it be possible to just skip the linear indexing operation and use a dictionary keyed by the `(i,j,k)` index tuple directly?

> [@Pbellive](#):
>
> Regarding performance, we’ve found that for our typical access patterns our implementation is faster in practice than using a dictionary keyed by an (i,j,k) tuple or using a 3D analogue of the compressed sparse column/row format used for sparse matrices.

I guess that you tried that already and found that it was slow… which is weird because it seems like what you’re currently doing is equivalent to using a hash function that consists of linearizing first and then hashing, whereas I would think that the built-in hash function which just mixes the hashes of the indices would be faster.

---

<div class="post-metadata">

### Author: ![Pbellive](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/pbellive/32/3604_2.png) [@Pbellive](https://discourse.julialang.org/u/Pbellive)
#### Post date: [November 7, 2018, 5:55pm UTC](https://discourse.julialang.org/t/why-is-linearindices-hardcoded-to-use-int/17235/8 "2018-11-07T17:55:46Z")

</div>

Thanks for the replies all. It looks like we should dig a bit deeper on the benchmarking and see if some greater insight can be gained on what’s going on here.

Cheers

---

<div class="post-metadata">

### Author: ![StefanKarpinski](https://sea2.discourse-cdn.com/julialang/user_avatar/discourse.julialang.org/stefankarpinski/32/24_2.png) [@StefanKarpinski](https://discourse.julialang.org/u/StefanKarpinski)
#### Post date: [November 7, 2018, 6:10pm UTC](https://discourse.julialang.org/t/why-is-linearindices-hardcoded-to-use-int/17235/9 "2018-11-07T18:10:14Z")

</div>

One thought is that instead of using linear indices, you could just use large random odd numbers as if they were strides. I.e. hash `(i, j, k)` as `c₁*i + c₂*j + c₃*k`—this should be faster than the default hashing but is a pretty decent hash function. You could, for example use these constants:

```julia
const c₁ = -5318449269026300911
const c₂ = 4566921428329651793
const c₃ = -6967723246364373735

```

In order to override the default hashing without monkey-patching, you’d need to wrap the `(i, j, k)` indices in a custom type like this:

```julia
struct Indices
    i::Int
    j::Int
    k::Int
end

Base.hash(x::Indices, h::UInt) = hash(c₁*x.i + c₂*x.j + c₃*x.k, h)

```

Untested but something like that.
