deprecate (then remove) generalized linear indexing

[IainNZ]

I'm sure this has been discussed elsewhere, but maybe not recently or in a focussed/isolated fashion. If this shouldn't be an issue report, I'll move it to julia-users.

Essentially, I would say this behavior is surprising:

julia> x = reshape(1:3^3, (3,3,3));

julia> size(x)
(3,3,3)

julia> x[3,4]
12

This was a cause of a hard-to-find bug for me recently, where I accidentally forgot an index.
I suppose this is in some way a generalization of the fact that x[12] works, but I'm not really sure why or where that behavior is useful for 2D-or-higher indices. Is there some logic for the current behavior?

EDIT: most similar discussion: #5396

[timholy]

The last index is effectively interpreted as a linear index, if the dimensionality is insufficient. (More properly, for LinearFast AbstractArrays, everything gets converted to an linear index in the end, although there is intermediate bounds-checking.) Extra 1s are also dropped. It's been this way "forever," but I won't close this issue if you're hoping for reconsideration. Suffice it to say that I've sometimes found this behavior useful, but I also understand how it could be a source of bugs.

Relevant code:
LinearFast (converts all to a linear index):

julia/base/abstractarray.jl

Line 508 in 0c20e64

unsafe_getindex(A, sub2ind(size(A), J...))

LinearSlow (expands to full dimensionality):

julia/base/abstractarray.jl

Lines 525 to 551 in 0c20e64

	elseif N > AN
	# Drop trailing ones
	Isplat = Expr[:(I[$d]) for d = 1:AN]
	Osplat = Expr[:(to_index(I[$d]) == 1) for d = AN+1:N]
	quote
	$(Expr(:meta, :inline, :propagate_inbounds))
	(&)($(Osplat...)) || throw_boundserror(A, I)
	getindex(A, $(Isplat...))
	end
	else
	# Expand the last index into the appropriate number of indices
	Isplat = Expr[:(I[$d]) for d = 1:N-1]
	i = 0
	for d=N:AN
	push!(Isplat, :(s[$(i+=1)]))
	end
	sz = Expr(:tuple)
	sz.args = Expr[:(size(A, $d)) for d=N:AN]
	szcheck = Expr[:(size(A, $d) > 0) for d=N:AN]
	quote
	$(Expr(:meta, :inline, :propagate_inbounds))
	# ind2sub requires all dimensions to be > 0:
	(&)($(szcheck...)) || throw_boundserror(A, I)
	s = ind2sub($sz, to_index(I[$N]))
	getindex(A, $(Isplat...))
	end
	end

[eschnett]

Whether the functionality is useful is not quite the same question as whether this should be the default behaviour. There could be a function linearindexing with this behaviour, or one could create a subarray...

[johnmyleswhite]

FWIW, I tend to feel like this behavior is too odd to be the default.

[IainNZ]

I would say that I'm looking for reconsideration, yes. I think its got a high surprising-ness to it, relative to its (minimal-but-nonzero) utility. I guess my main initial goal was understanding how much utility it does have, and to whether "fixing it" would have nontrivial performance applications.

In terms of options for "fixing it", after I admit relatively little consideration, it doesn't seem like anything but requiring all dimensions would be consistent. e.g. consider

julia> x = reshape(1:3^3, (3,3,3));
julia> x[3,3]
9

In some ways I find that even more surprising than my other example (but would have found it returning x[3,3,:] less so).

[BobPortmann]

I agree and feel strongly that it should be an error to index an array with fewer indices than the rank of the array. It just makes it too easy to introduce hard to find errors into code when using higher rank arrays (say 4-10 dimensional). It is easy enough to use reshape to make a copy free reference if one wants to use the trailing linear indexing feature (and this leaves the intent clear in the code). And one could easily write some helper functions to make trailing linear indexing easier to use. I would also note that other languages (e.g., IDL and Fortran) make this an error so the present behavior will be surprising to many.

The case of 1-D linear indexing is perhaps OK since indexing as a one dimensional array is visually distinct (compared with, e.g., indexing a 8-D array as 7-D). However, since it is so easy to use vec I am not sure even this is necessary.

This was discussed a bit in JuliaLang/LinearAlgebra.jl#42 but I don't know how to link to specific comments and that is a long issue.

[rfourquet]

@BobPortmann you can link to specific comments by right clicking on the timestamp (next to the author name in the header of the comment), and selecting "copy link location" or equivalent.

[toivoh]

+1 to reconsidering this, it just seems so much more likely that this would be used by mistake than on purpose. Even for linear indexing I think that we should consider if it should have its own method instead of being overloaded onto regular indexing.

[johnmyleswhite]

I kind of love the idea of giving linear indexing a special method.

[timholy]

What about A[i, LinearIndex(j)]?

[johnmyleswhite]

I'd be ok with that.

[JeffBezanson]

+1 to dropping this behavior and having separate syntax of some kind for linear indexing. This behavior was just taken as the default choice since it's what some other environments do. I would describe it as an over-generalization of a premature optimization.

[StefanKarpinski]

I would describe it as an over-generalization of a premature optimization.

Precisely. Let's ditch this for sure.

[tknopp]

-1 for dropping the 1D version. Its very natural to loop over all elements in a multidimensional array using

for n=1:length(x)
  x[n] = ...
end

I know this can be achieve in other ways (e.g. enumerate) but the above is quite easy to remember.

[lobingera]

@JeffBezanson, i disagree with the statement about over-generalization of a premature optimization.

Linear indexing is the natural indexing scheme and anything providing multi dimensional access is actually some infrastructure to make it easier for the person writing code or design an algorithm.
Taking linear indexing away completely just feels wrong to me, introducing additional syntax (which means to me, i have to hit more keys on the REPL before getting output) seems to me additional burden; a[i] is just linear, a[i,j] is dimensional; where is the problem?

Now to the interesting problem that started this: Having a a[i,j] access when a is actually defined 3D. I was not aware that this existed, but actually this is something that i'd like to have in a language.

The initial problem for me (and i think the title of the issue is missleading) is rather in the area of a compiler warning or some code checking (like lint).