Add special methods for Cartesian indexing of CartesianIndices for certain permissible indices

[jishnub]

Closes #40314

Now a CartesianIndices, when indexed with ranges, Colon and compatible CartesianIndices returns another CartesianIndices instead of an Array of CartesianIndexes. The behavior is replicated while constructing views. This might have performance benefits.

julia> a = ones(2, 2)
2×2 Matrix{Float64}:
 1.0  1.0
 1.0  1.0

julia> c = CartesianIndices(a)
2×2 CartesianIndices{2, Tuple{Base.OneTo{Int64}, Base.OneTo{Int64}}}:
 CartesianIndex(1, 1)  CartesianIndex(1, 2)
 CartesianIndex(2, 1)  CartesianIndex(2, 2)

julia> c[1:2, :]
2×2 CartesianIndices{2, Tuple{UnitRange{Int64}, Base.OneTo{Int64}}}:
 CartesianIndex(1, 1)  CartesianIndex(1, 2)
 CartesianIndex(2, 1)  CartesianIndex(2, 2)

julia> c[c]
2×2 CartesianIndices{2, Tuple{Base.OneTo{Int64}, Base.OneTo{Int64}}}:
 CartesianIndex(1, 1)  CartesianIndex(1, 2)
 CartesianIndex(2, 1)  CartesianIndex(2, 2)

julia> c[c] === c
true

julia> c[:, :] === c
true

julia> @view c[:, :]
2×2 CartesianIndices{2, Tuple{Base.OneTo{Int64}, Base.OneTo{Int64}}}:
 CartesianIndex(1, 1)  CartesianIndex(1, 2)
 CartesianIndex(2, 1)  CartesianIndex(2, 2)

julia> (@view c[:, :]) === c
true

This doesn't work for combinations of range/colon and CartesianIndices though, as it dispatches on the ndims.

[johnnychen94]

A clean implementation with carefully written tests. Most importantly, this sets up a solid ground for tile operations performance optimization that can be found in image processing and related fields. (Previously, view is the workaround for it).

R = CartesianIndices((2, 3, 4));
indices = (1:2, 1:2, 1:2);
@btime getindex($R, $indices...);
# master  58.120 ns (1 allocation: 272 bytes)
# PR 8.454 ns (0 allocations: 0 bytes)

[jishnub]

Test failure seems unrelated

[jishnub]

gentle bump

[stevengj]

Yes, the error in testset errorshow has cropped up a few times now (e.g. #30233, @ararslan), and seems unrelated.