Merge pull request #13612 from JuliaLang/mb/drop-bear

RFC: Drop dimensions indexed by scalars

Author: andreasnoack

NEWS.md

@@ -47,6 +47,9 @@ Library improvements
 
   * Linear algebra:
 
+    * All dimensions indexed by scalars are now dropped, whereas previously only
+      trailing scalar dimensions would be omitted from the result.
+
     * New `normalize` and `normalize!` convenience functions for normalizing
       vectors ([#13681]).
 

base/linalg/diagonal.jl

@@ -121,7 +121,10 @@ function A_ldiv_B!{T}(D::Diagonal{T}, V::AbstractMatrix{T})
         if d == zero(T)
             throw(SingularException(i))
         end
-        V[i,:] *= inv(d)
+        d⁻¹ = inv(d)
+        for j=1:size(V,2)
+            @inbounds V[i,j] *= d⁻¹
+        end
     end
     V
 end

base/multidimensional.jl

@@ -166,15 +166,15 @@ index_lengths_dim(A, dim, ::Colon) = (trailingsize(A, dim),)
 @inline index_lengths_dim(A, dim, i::AbstractArray{Bool}, I...) = (sum(i), index_lengths_dim(A, dim+1, I...)...)
 @inline index_lengths_dim(A, dim, i::AbstractArray, I...) = (length(i), index_lengths_dim(A, dim+1, I...)...)
 
-# shape of array to create for getindex() with indexes I, dropping trailing scalars
+# shape of array to create for getindex() with indexes I, dropping scalars
 index_shape(A::AbstractArray, I::AbstractArray) = size(I) # Linear index reshape
 index_shape(A::AbstractArray, I::AbstractArray{Bool}) = (sum(I),) # Logical index
 index_shape(A::AbstractArray, I::Colon) = (length(A),)
 @inline index_shape(A::AbstractArray, I...) = index_shape_dim(A, 1, I...)
 index_shape_dim(A, dim, I::Real...) = ()
 index_shape_dim(A, dim, ::Colon) = (trailingsize(A, dim),)
 @inline index_shape_dim(A, dim, ::Colon, i, I...) = (size(A, dim), index_shape_dim(A, dim+1, i, I...)...)
-@inline index_shape_dim(A, dim, ::Real, I...) = (1, index_shape_dim(A, dim+1, I...)...)
+@inline index_shape_dim(A, dim, ::Real, I...) = (index_shape_dim(A, dim+1, I...)...)
 @inline index_shape_dim(A, dim, i::AbstractVector{Bool}, I...) = (sum(i), index_shape_dim(A, dim+1, I...)...)
 @inline index_shape_dim(A, dim, i::AbstractVector, I...) = (length(i), index_shape_dim(A, dim+1, I...)...)
 
@@ -238,7 +238,8 @@ end
         $(Expr(:meta, :inline))
         D = eachindex(dest)
         Ds = start(D)
-        @nloops $N i dest d->(j_d = unsafe_getindex(I[d], i_d)) begin
+        idxlens = index_lengths(src, I...) # TODO: unsplat?
+        @nloops $N i d->(1:idxlens[d]) d->(j_d = unsafe_getindex(I[d], i_d)) begin
             d, Ds = next(D, Ds)
             v = @ncall $N unsafe_getindex src j
             unsafe_setindex!(dest, v, d)
@@ -248,18 +249,18 @@ end
 end
 
 # checksize ensures the output array A is the correct size for the given indices
-checksize(A::AbstractArray, I::AbstractArray) = size(A) == size(I) || throw(DimensionMismatch("index 1 has size $(size(I)), but size(A) = $(size(A))"))
-checksize(A::AbstractArray, I::AbstractArray{Bool}) = length(A) == sum(I) || throw(DimensionMismatch("index 1 selects $(sum(I)) elements, but length(A) = $(length(A))"))
-@generated function checksize(A::AbstractArray, I...)
-    N = length(I)
-    quote
-        @nexprs $N d->(_checksize(A, d, I[d]) || throw(DimensionMismatch("index $d selects $(length(I[d])) elements, but size(A, $d) = $(size(A,d))")))
-    end
-end
-_checksize(A::AbstractArray, dim, I) = size(A, dim) == length(I)
-_checksize(A::AbstractArray, dim, I::AbstractVector{Bool}) = size(A, dim) == sum(I)
-_checksize(A::AbstractArray, dim, ::Colon) = true
-_checksize(A::AbstractArray, dim, ::Real) = size(A, dim) == 1
+@noinline throw_checksize_error(arr, dim, idx) = throw(DimensionMismatch("index $d selects $(length(I[d])) elements, but size(A, $d) = $(size(A,d))"))
+
+checksize(A::AbstractArray, I::AbstractArray) = size(A) == size(I) || throw_checksize_error(A, 1, I)
+checksize(A::AbstractArray, I::AbstractArray{Bool}) = length(A) == sum(I) || throw_checksize_error(A, 1, I)
+
+checksize(A::AbstractArray, I...) = _checksize(A, 1, I...)
+_checksize(A::AbstractArray, dim) = true
+# Skip scalars
+_checksize(A::AbstractArray, dim, ::Real, J...) = _checksize(A, dim, J...)
+_checksize(A::AbstractArray, dim, I, J...) = (size(A, dim) == length(I) || throw_checksize_error(A, dim, I); _checksize(A, dim+1, J...))
+_checksize(A::AbstractArray, dim, I::AbstractVector{Bool}, J...) = (size(A, dim) == sum(I) || throw_checksize_error(A, dim, I); _checksize(A, dim+1, J...))
+_checksize(A::AbstractArray, dim, ::Colon, J...) = _checksize(A, dim+1, J...)
 
 @inline unsafe_setindex!(v::BitArray, x::Bool, ind::Int) = (Base.unsafe_bitsetindex!(v.chunks, x, ind); v)
 @inline unsafe_setindex!(v::BitArray, x, ind::Real) = (Base.unsafe_bitsetindex!(v.chunks, convert(Bool, x), to_index(ind)); v)
@@ -585,7 +586,8 @@ end
 
         storeind = 1
         Xc, Bc = X.chunks, B.chunks
-        @nloops($N, i, d->1:size(X, d+1),
+        idxlens = index_lengths(B, I0, I...) # TODO: unsplat?
+        @nloops($N, i, d->(1:idxlens[d+1]),
                 d->nothing, # PRE
                 d->(ind += stride_lst_d - gap_lst_d), # POST
                 begin # BODY
@@ -608,7 +610,8 @@ end
         $(symbol(:offset_, N)) = 1
         ind = 0
         Xc, Bc = X.chunks, B.chunks
-        @nloops $N i X d->(offset_{d-1} = offset_d + (unsafe_getindex(I[d], i_d)-1)*stride_d) begin
+        idxlens = index_lengths(B, I...) # TODO: unsplat?
+        @nloops $N i d->(1:idxlens[d]) d->(offset_{d-1} = offset_d + (unsafe_getindex(I[d], i_d)-1)*stride_d) begin
             ind += 1
             unsafe_bitsetindex!(Xc, unsafe_bitgetindex(Bc, offset_0), ind)
         end

base/sparse/sparsematrix.jl

@@ -1352,8 +1352,6 @@ function getindex{T}(A::SparseMatrixCSC{T}, i0::Integer, i1::Integer)
     ((r1 > r2) || (A.rowval[r1] != i0)) ? zero(T) : A.nzval[r1]
 end
 
-getindex{T<:Integer}(A::SparseMatrixCSC, i::Integer, J::AbstractVector{T}) = getindex(A,[i],J)
-
 # Colon translation
 getindex(A::SparseMatrixCSC, ::Colon, ::Colon) = copy(A)
 getindex(A::SparseMatrixCSC, i, ::Colon)       = getindex(A, i, 1:size(A, 2))

base/sparse/sparsevector.jl

@@ -345,6 +345,7 @@ sprandbool(r::AbstractRNG, n::Integer, p::AbstractFloat) = sprand(r, n, p, trueb
 
 ## Indexing into Matrices can return SparseVectors
 
+# Column slices
 function getindex(x::SparseMatrixCSC, ::Colon, j::Integer)
     checkbounds(x, :, j)
     r1 = convert(Int, x.colptr[j])
@@ -369,6 +370,38 @@ end
     SparseVector(M.m, M.rowval, M.nzval)
 end
 
+# Row slices
+getindex(A::SparseMatrixCSC, i::Integer, ::Colon) = A[i, 1:end]
+getindex{Tv,Ti}(A::SparseMatrixCSC{Tv,Ti}, i::Integer, J::AbstractVector{Bool}) = A[i, find(J)]
+function Base.getindex{Tv,Ti}(A::SparseMatrixCSC{Tv,Ti}, i::Integer, J::AbstractVector)
+    checkbounds(A, i, J)
+    nJ = length(J)
+    colptrA = A.colptr; rowvalA = A.rowval; nzvalA = A.nzval
+
+    nzinds = Array(Ti, 0)
+    nzvals = Array(Tv, 0)
+
+    # adapted from SparseMatrixCSC's sorted_bsearch_A
+    ptrI = 1
+    @inbounds for j = 1:nJ
+        col = J[j]
+        rowI = i
+        ptrA = colptrA[col]
+        stopA = colptrA[col+1]-1
+        if ptrA <= stopA
+            if rowvalA[ptrA] <= rowI
+                ptrA = searchsortedfirst(rowvalA, rowI, ptrA, stopA, Base.Order.Forward)
+                if ptrA <= stopA && rowvalA[ptrA] == rowI
+                    push!(nzinds, ptrI)
+                    push!(nzvals, nzvalA[ptrA])
+                end
+            end
+            ptrI += 1
+        end
+    end
+    return SparseVector(nJ, nzinds, nzvals)
+end
+
 
 # Logical and linear indexing into SparseMatrices
 getindex{Tv}(A::SparseMatrixCSC{Tv}, I::AbstractVector{Bool}) = _logical_index(A, I) # Ambiguities

base/subarray.jl

@@ -604,8 +604,8 @@ end
 
 # Indexing with non-scalars. For now, this returns a copy, but changing that
 # is just a matter of deleting the explicit call to copy.
-getindex{T,N,P,IV}(V::SubArray{T,N,P,IV}, I::ViewIndex...) = copy(sub(V, I...))
-unsafe_getindex{T,N,P,IV}(V::SubArray{T,N,P,IV}, I::ViewIndex...) = copy(sub_unsafe(V, I))
+getindex{T,N,P,IV}(V::SubArray{T,N,P,IV}, I::ViewIndex...) = copy(slice(V, I...))
+unsafe_getindex{T,N,P,IV}(V::SubArray{T,N,P,IV}, I::ViewIndex...) = copy(slice_unsafe(V, I))
 
 # Nonscalar setindex! falls back to the AbstractArray versions
 

doc/manual/arrays.rst

@@ -233,8 +233,8 @@ where each ``I_k`` may be:
 The result ``X`` generally has dimensions
 ``(length(I_1), length(I_2), ..., length(I_n))``, with location
 ``(i_1, i_2, ..., i_n)`` of ``X`` containing the value
-``A[I_1[i_1], I_2[i_2], ..., I_n[i_n]]``. Trailing dimensions
-indexed with scalars are dropped. For example, the dimensions of ``A[I, 1]`` will be
+``A[I_1[i_1], I_2[i_2], ..., I_n[i_n]]``. All dimensions indexed with scalars are
+dropped. For example, the result of ``A[2, I, 3]`` will be a vector with size
 ``(length(I),)``. Boolean vectors are first transformed with ``find``; the size of
 a dimension indexed by a boolean vector will be the number of true values in the vector.
 As a special part of this syntax, the ``end`` keyword may be used to represent the last

test/arrayops.jl

@@ -93,7 +93,7 @@ a = rand(1, 1, 8, 8, 1)
 sz = (5,8,7)
 A = reshape(1:prod(sz),sz...)
 @test A[2:6] == [2:6;]
-@test A[1:3,2,2:4] == cat(3,46:48,86:88,126:128)
+@test A[1:3,2,2:4] == cat(2,46:48,86:88,126:128)
 @test A[:,7:-3:1,5] == [191 176 161; 192 177 162; 193 178 163; 194 179 164; 195 180 165]
 @test A[:,3:9] == reshape(11:45,5,7)
 rng = (2,2:3,2:2:5)
@@ -111,7 +111,7 @@ tmp = cat([1,3],blk,blk)
 
 x = rand(2,2)
 b = x[1,:]
-@test isequal(size(b), (1, 2))
+@test isequal(size(b), (2,))
 b = x[:,1]
 @test isequal(size(b), (2,))
 
@@ -129,7 +129,7 @@ B[[3,1],[2,4]] = [21 22; 23 24]
 B[4,[2,3]] = 7
 @test B == [0 23 1 24 0; 11 12 13 14 15; 0 21 3 22 0; 0 7 7 0 0]
 
-@test isequal(reshape(1:27, 3, 3, 3)[1,:], [1  4  7  10  13  16  19  22  25])
+@test isequal(reshape(1:27, 3, 3, 3)[1,:], [1,  4,  7,  10,  13,  16,  19,  22,  25])
 
 a = [3, 5, -7, 6]
 b = [4, 6, 2, -7, 1]
@@ -575,12 +575,12 @@ let
 
     @test isequal(c[:,1], cv)
     @test isequal(c[:,3], cv2)
-    @test isequal(c[4,:], [2.0 2.0 2.0 2.0]*1000)
+    @test isequal(c[4,:], [2.0, 2.0, 2.0, 2.0]*1000)
 
     c = cumsum_kbn(A, 2)
 
-    @test isequal(c[1,:], cv2')
-    @test isequal(c[3,:], cv')
+    @test isequal(c[1,:], cv2)
+    @test isequal(c[3,:], cv)
     @test isequal(c[:,4], [2.0,2.0,2.0,2.0]*1000)
 
 end

test/bitarray.jl

@@ -171,10 +171,10 @@ function gen_getindex_data()
     m1, m2 = rand_m1m2()
     produce((m1, m2, Bool))
     m1, m2 = rand_m1m2()
-    produce((m1, 1:m2, BitMatrix))
-    produce((m1, :, BitMatrix))
+    produce((m1, 1:m2, BitVector))
+    produce((m1, :, BitVector))
     m1, m2 = rand_m1m2()
-    produce((m1, randperm(m2), BitMatrix))
+    produce((m1, randperm(m2), BitVector))
     m1, m2 = rand_m1m2()
     produce((1:m1, m2, BitVector))
     produce((:, m2, BitVector))

test/parallel_exec.jl

@@ -242,7 +242,7 @@ map!(x->1, d)
 @test fill!(d, 2.) == fill(2, 10, 10)
 @test d[:] == fill(2, 100)
 @test d[:,1] == fill(2, 10)
-@test d[1,:] == fill(2, 1, 10)
+@test d[1,:] == fill(2, 10)
 
 # Boundary cases where length(S) <= length(pids)
 @test 2.0 == remotecall_fetch(D->D[2], id_other, Base.shmem_fill(2.0, 2; pids=[id_me, id_other]))