diff --git a/base/multidimensional.jl b/base/multidimensional.jl index aa3c3f3028884..f37d999ce2725 100644 --- a/base/multidimensional.jl +++ b/base/multidimensional.jl @@ -345,33 +345,35 @@ end # the corresponding cartesian index into the parent, and then uses # dims to convert back to a linear index into the parent array. # -# However, a common case is linindex::UnitRange. -# Since div is slow and in(j::Int, linindex::UnitRange) is fast, +# However, a common case is linindex::Range. +# Since div is slow and in(j::Int, linindex::Range) is fast, # it can be much faster to generate all possibilities and # then test whether the corresponding linear index is in linindex. # One exception occurs when only a small subset of the total # is desired, in which case we fall back to the div-based algorithm. -stagedfunction merge_indexes(V, indexes::NTuple, dims::Dims, linindex::UnitRange{Int}) - N = length(indexes) +#stagedfunction merge_indexes{T<:Integer}(V, parentindexes::NTuple, parentdims::Dims, linindex::Union(Colon,Range{T}), lindim) +stagedfunction merge_indexes_in{TT}(V, parentindexes::TT, parentdims::Dims, linindex, lindim) + N = length(parentindexes) # number of parent axes we're merging N > 0 || throw(ArgumentError("cannot merge empty indexes")) + lengthexpr = linindex == Colon ? (:(prod(size(V)[lindim:end]))) : (:(length(linindex))) + L = symbol(string("Istride_", N+1)) # length of V's trailing dimensions quote - n = length(linindex) - Base.Cartesian.@nexprs $N d->(I_d = indexes[d]) - L = 1 - dimoffset = ndims(V.parent) - length(dims) - Base.Cartesian.@nexprs $N d->(L *= dimsize(V.parent, d+dimoffset, I_d)) - if n < 0.1L # this has not been carefully tuned - return merge_indexes_div(V, indexes, dims, linindex) + n = $lengthexpr + Base.Cartesian.@nexprs $N d->(I_d = parentindexes[d]) + pdimoffset = ndims(V.parent) - length(parentdims) + Istride_1 = 1 # parentindexes strides + Base.Cartesian.@nexprs $N d->(Istride_{d+1} = Istride_d*dimsize(V.parent, d+pdimoffset, I_d)) + Istridet = Base.Cartesian.@ntuple $(N+1) d->Istride_d + if n < 0.1*$L # this has not been carefully tuned + return merge_indexes_div(V, parentindexes, parentdims, linindex, lindim) end Pstride_1 = 1 # parent strides - Base.Cartesian.@nexprs $(N-1) d->(Pstride_{d+1} = Pstride_d*dims[d]) - Istride_1 = 1 # indexes strides - Base.Cartesian.@nexprs $(N-1) d->(Istride_{d+1} = Istride_d*dimsize(V, d+dimoffset, I_d)) - Base.Cartesian.@nexprs $N d->(counter_d = 1) # counter_0 is a linear index into indexes + Base.Cartesian.@nexprs $(N-1) d->(Pstride_{d+1} = Pstride_d*parentdims[d]) + Base.Cartesian.@nexprs $N d->(counter_d = 1) # counter_0 is a linear index into parentindexes Base.Cartesian.@nexprs $N d->(offset_d = 1) # offset_0 is a linear index into parent k = 0 index = Array(Int, n) - Base.Cartesian.@nloops $N i d->(1:dimsize(V, d+dimoffset, I_d)) d->(offset_{d-1} = offset_d + (I_d[i_d]-1)*Pstride_d; counter_{d-1} = counter_d + (i_d-1)*Istride_d) begin + Base.Cartesian.@nloops $N i d->(1:dimsize(V.parent, d+pdimoffset, I_d)) d->(offset_{d-1} = offset_d + (I_d[i_d]-1)*Pstride_d; counter_{d-1} = counter_d + (i_d-1)*Istride_d) begin if in(counter_0, linindex) index[k+=1] = offset_0 end @@ -379,27 +381,35 @@ stagedfunction merge_indexes(V, indexes::NTuple, dims::Dims, linindex::UnitRange index end end -merge_indexes(V, indexes::NTuple, dims::Dims, linindex) = merge_indexes_div(V, indexes, dims, linindex) + +# HACK: dispatch seemingly wasn't working properly +function merge_indexes(V, parentindexes::NTuple, parentdims::Dims, linindex, lindim) + if isa(linindex, Colon) || isa(linindex, Range) + return merge_indexes_in(V, parentindexes, parentdims, linindex, lindim) + end + merge_indexes_div(V, parentindexes, parentdims, linindex, lindim) +end # This could be written as a regular function, but performance # will be better using Cartesian macros to avoid the heap and # an extra loop. -stagedfunction merge_indexes_div(V, indexes::NTuple, dims::Dims, linindex) - N = length(indexes) +stagedfunction merge_indexes_div{TT}(V, parentindexes::TT, parentdims::Dims, linindex, lindim) + N = length(parentindexes) N > 0 || throw(ArgumentError("cannot merge empty indexes")) Istride_N = symbol("Istride_$N") + lengthexpr = :(length(linindex)) quote - Base.Cartesian.@nexprs $N d->(I_d = indexes[d]) + Base.Cartesian.@nexprs $N d->(I_d = parentindexes[d]) Pstride_1 = 1 # parent strides - Base.Cartesian.@nexprs $(N-1) d->(Pstride_{d+1} = Pstride_d*dims[d]) - Istride_1 = 1 # indexes strides - dimoffset = ndims(V.parent) - length(dims) - Base.Cartesian.@nexprs $(N-1) d->(Istride_{d+1} = Istride_d*dimsize(V.parent, d+dimoffset, I_d)) - n = length(linindex) - L = $(Istride_N) * dimsize(V.parent, $N+dimoffset, indexes[end]) + Base.Cartesian.@nexprs $(N-1) d->(Pstride_{d+1} = Pstride_d*parentdims[d]) + Istride_1 = 1 # parentindexes strides + pdimoffset = ndims(V.parent) - length(parentdims) + Base.Cartesian.@nexprs $(N-1) d->(Istride_{d+1} = Istride_d*dimsize(V.parent, d+pdimoffset, I_d)) + n = $lengthexpr + L = $(Istride_N) * dimsize(V.parent, $N+pdimoffset, parentindexes[end]) index = Array(Int, n) for i = 1:n - k = linindex[i] # k is the indexes-centered linear index + k = linindex[i] # k is the parentindexes-centered linear index 1 <= k <= L || throw(BoundsError()) k -= 1 j = 0 # j will be the new parent-centered linear index diff --git a/base/subarray.jl b/base/subarray.jl index b55144b45387f..fbcfa111f0aac 100644 --- a/base/subarray.jl +++ b/base/subarray.jl @@ -3,7 +3,7 @@ typealias ViewIndex Union(Int, NonSliceIndex) typealias RangeIndex Union(Int, Range{Int}, UnitRange{Int}, Colon) # LD is the last dimension up through which this object has efficient -# linear indexing. If LD==N, then the object itself has efficient +# linear indexing. If LD==length(I), then the object itself has efficient # linear indexing. immutable SubArray{T,N,P<:AbstractArray,I<:(ViewIndex...),LD} <: AbstractArray{T,N} parent::P @@ -40,28 +40,45 @@ function _slice(A, I) slice_unsafe(A, I) end -stagedfunction slice_unsafe{T,NP,IndTypes}(A::AbstractArray{T,NP}, I::IndTypes) +# The most complicated part of this is matching the axes between the +# input index tuples (denoted by J), the index tuples that get stored +# in the view (denoted by I), and the overall dimensionality of the +# view. +# The complexities increase when you create a view-of-a-view, because +# then there is also the index tuple of the parent view (denoted IV) +# to consider. +# +# Examples: +# S1 = sub(A::Matrix, 2, 3:5) ndims(S1) == length(I) == length(J) == 2 +# S2 = slice(A::Matrix, 2, 3:5) ndims(S2) == 1, length(I) == length(J) == 2 +# S3 = sub(A::Matrix, 4:17) ndims(S3) == length(I) == length(J) == 1 +# S4 = sub(S2, 1:2) ndims(S4) == length(J) == 1, length(I) == 2 +# S3 addresses the trailing dimensions of the parent by linear indexing. +# For S4, J[1] corresponds to I[2], because of the slice along +# dimension 1 in S2 + +stagedfunction slice_unsafe{T,NP,IndTypes}(A::AbstractArray{T,NP}, J::IndTypes) N = 0 sizeexprs = Array(Any, 0) - for k = 1:length(I) - i = I[k] - if !(i <: Real) + for Jindex = 1:length(J) + j = J[Jindex] + if !(j <: Real) N += 1 - push!(sizeexprs, dimsizeexpr(I[k], k, length(I), :A, :I)) + push!(sizeexprs, dimsizeexpr(j, Jindex, length(J), :A, :J)) end end dims = :(tuple($(sizeexprs...))) - LD = subarray_linearindexing_dim(A, I) - strideexpr = stride1expr(A, I, :A, :I, LD) - exfirst = first_index_expr(:A, :I, length(I)) + LD = subarray_linearindexing_dim(A, J) + strideexpr = stride1expr(A, J, :A, :J, LD) + exfirst = first_index_expr(:A, :J, length(J)) quote $exfirst - SubArray{$T,$N,$A,$I,$LD}(A, I, $dims, f, $strideexpr) + SubArray{$T,$N,$A,$J,$LD}(A, J, $dims, f, $strideexpr) end end # Conventional style (drop trailing singleton dimensions, keep any -# other singletons) +# other singletons by converting them to ranges, e.g., 3:3) sub(A::AbstractArray, I::ViewIndex...) = _sub(A, I) sub(A::AbstractArray, I::(ViewIndex...)) = _sub(A, I) function _sub(A, I) @@ -69,32 +86,33 @@ function _sub(A, I) sub_unsafe(A, I) end -stagedfunction sub_unsafe{T,NP,IndTypes}(A::AbstractArray{T,NP}, I::IndTypes) +stagedfunction sub_unsafe{T,NP,IndTypes}(A::AbstractArray{T,NP}, J::IndTypes) sizeexprs = Array(Any, 0) Itypes = Array(Any, 0) Iexprs = Array(Any, 0) - N = length(I) - while N > 0 && I[N] <: Real + N = length(J) + while N > 0 && J[N] <: Real N -= 1 end - for k = 1:length(I) - if k <= N - push!(sizeexprs, dimsizeexpr(I[k], k, length(I), :A, :I)) + for Jindex = 1:length(J) + j = J[Jindex] + if Jindex <= N + push!(sizeexprs, dimsizeexpr(j, Jindex, length(J), :A, :J)) end - if k < N && I[k] <: Real + if Jindex < N && j <: Real push!(Itypes, UnitRange{Int}) - push!(Iexprs, :(Int(I[$k]):Int(I[$k]))) + push!(Iexprs, :(Int(J[$Jindex]):Int(J[$Jindex]))) else - push!(Itypes, I[k]) - push!(Iexprs, :(I[$k])) + push!(Itypes, j) + push!(Iexprs, :(J[$Jindex])) end end dims = :(tuple($(sizeexprs...))) Iext = :(tuple($(Iexprs...))) It = tuple(Itypes...) - LD = subarray_linearindexing_dim(A, I) - strideexpr = stride1expr(A, I, :A, :I, LD) - exfirst = first_index_expr(:A, :I, length(It)) + LD = subarray_linearindexing_dim(A, J) + strideexpr = stride1expr(A, J, :A, :J, LD) + exfirst = first_index_expr(:A, :J, length(It)) quote $exfirst SubArray{$T,$N,$A,$It,$LD}(A, $Iext, $dims, f, $strideexpr) @@ -103,68 +121,87 @@ end # Constructing from another SubArray # This "pops" the old SubArray and creates a more compact one -stagedfunction slice_unsafe{T,NV,PV,IV,PLD,IndTypes}(V::SubArray{T,NV,PV,IV,PLD}, I::IndTypes) +stagedfunction slice_unsafe{T,NV,PV,IV,PLD,IndTypes}(V::SubArray{T,NV,PV,IV,PLD}, J::IndTypes) N = 0 sizeexprs = Array(Any, 0) indexexprs = Array(Any, 0) Itypes = Array(Any, 0) - k = 0 # index into I - LD, die_next_vector, Ilast, isLDdone = 0, false, Void, false # for linear indexing inference - for j = 1:length(IV) - if IV[j] <: Real - push!(indexexprs, :(V.indexes[$j])) - push!(Itypes, IV[j]) + # The next two Ints, if nonzero, record information about the place + # in the index tuple at which trailing dimensions got packed into a + # single Vector{Int}. For stride1 computation, we need to keep track + # of whether the index that triggered this had uniform stride. + # Iindex_lin is the spot in the resulting index tuple + # Jindex_lin is the corresponding spot in the input index tuple + Iindex_lin = Jindex_lin = 0 + # Linear indexing inference makes use of the following variables: + # LD: the last dimension up through which linear indexing is efficient + # isLDdone: true if we've quit incrementing LD + # die_next_vector: if true, stop incrementing LD on the next + # "extended" input index + # jprev: holds the previous input index type + LD, die_next_vector, jprev, isLDdone = 0, false, Void, false # for linear indexing inference + Jindex = 0 + for IVindex = 1:length(IV) + iv = IV[IVindex] + if iv <: Real + push!(indexexprs, :(V.indexes[$IVindex])) + push!(Itypes, iv) if !isLDdone LD += 1 end else - k += 1 + Jindex += 1 + j = J[Jindex] + if Jindex < length(J) || Jindex == NV || IVindex == length(IV) + if !(j <: Real) + N += 1 + push!(sizeexprs, dimsizeexpr(j, Jindex, length(J), :V, :J)) + end + push!(indexexprs, :(reindex(V.indexes[$IVindex], J[$Jindex]))) + push!(Itypes, rangetype(iv, j)) + else + # We have a linear index that spans more than one + # dimension of the parent + N += 1 + push!(sizeexprs, dimsizeexpr(j, Jindex, length(J), :V, :J)) + push!(indexexprs, :(merge_indexes(V, V.indexes[$IVindex:end], size(V.parent)[$IVindex:end], J[$Jindex], $Jindex))) + push!(Itypes, Array{Int, 1}) + Iindex_lin = length(Itypes) + Jindex_lin = Jindex + break + end if !isLDdone if LD < PLD LD += 1 - Ilast, LD, die_next_vector, isdone = nextLD(Ilast, I[k], LD, die_next_vector) + jprev, LD, die_next_vector, isdone = nextLD(jprev, j, LD, die_next_vector) isLDdone |= isdone else - if I[k] <: Real + if j <: Real LD += 1 else isLDdone = true end end end - if k < length(I) || k == NV || j == length(IV) - if !(I[k] <: Real) - N += 1 - push!(sizeexprs, dimsizeexpr(I[k], k, length(I), :V, :I)) - end - push!(indexexprs, :(reindex(V.indexes[$j], I[$k]))) - push!(Itypes, rangetype(IV[j], I[k])) - else - # We have a linear index that spans more than one dimension of the parent - N += 1 - push!(sizeexprs, dimsizeexpr(I[k], k, length(I), :V, :I)) - push!(indexexprs, :(merge_indexes(V, V.indexes[$j:end], size(V.parent)[$j:end], I[$k]))) - push!(Itypes, Array{Int, 1}) - break - end end end - for i = k+1:length(I) - if !(I[i] <: Real) + for Jind = Jindex+1:length(J) + j = J[Jind] + if !(j <: Real) N += 1 - push!(sizeexprs, dimsizeexpr(I[i], i, length(I), :V, :I)) + push!(sizeexprs, dimsizeexpr(j, Jind, length(J), :V, :J)) isLDdone = true elseif !isLDdone LD += 1 end - push!(indexexprs, :(I[$i])) - push!(Itypes, I[i]) + push!(indexexprs, :(J[$Jind])) + push!(Itypes, j) end Inew = :(tuple($(indexexprs...))) dims = :(tuple($(sizeexprs...))) It = tuple(Itypes...) LD = max(LD, subarray_linearindexing_dim(PV, It)) - strideexpr = stride1expr(PV, It, :(V.parent), :Inew, LD) + strideexpr = stride1expr(PV, It, :(V.parent), :Inew, LD, :J, Iindex_lin, Jindex_lin) exfirst = first_index_expr(:(V.parent), :Inew, length(It)) quote Inew = $Inew @@ -173,52 +210,59 @@ stagedfunction slice_unsafe{T,NV,PV,IV,PLD,IndTypes}(V::SubArray{T,NV,PV,IV,PLD} end end -stagedfunction sub_unsafe{T,NV,PV,IV,PLD,IndTypes}(V::SubArray{T,NV,PV,IV,PLD}, I::IndTypes) - N = length(I) - while N > 0 && I[N] <: Real +stagedfunction sub_unsafe{T,NV,PV,IV,PLD,IndTypes}(V::SubArray{T,NV,PV,IV,PLD}, J::IndTypes) + N = length(J) + while N > 0 && J[N] <: Real N -= 1 end sizeexprs = Array(Any, 0) indexexprs = Array(Any, 0) Itypes = Array(Any, 0) + ItypesLD = Array(Any, 0) preexprs = Array(Any, 0) - k = 0 - LD, die_next_vector, Ilast, isLDdone = 0, false, Void, false # for linear indexing inference - for j = 1:length(IV) - if IV[j] <: Real - push!(indexexprs, :(V.indexes[$j])) - push!(Itypes, IV[j]) + LD, die_next_vector, jprev, isLDdone = 0, false, Void, false + Jindex = 0 + for IVindex = 1:length(IV) + iv = IV[IVindex] + if iv <: Real + push!(indexexprs, :(V.indexes[$IVindex])) + push!(Itypes, iv) + push!(ItypesLD, iv) if !isLDdone LD += 1 end else - k += 1 - if k <= N - push!(sizeexprs, dimsizeexpr(I[k], k, length(I), :V, :I)) + Jindex += 1 + j = J[Jindex] + if Jindex <= N + push!(sizeexprs, dimsizeexpr(j, Jindex, length(J), :V, :J)) end - if k < N && I[k] <: Real + if Jindex < N && j <: Real # convert scalar to a range sym = gensym() - push!(preexprs, :($sym = reindex(V.indexes[$j], Int(I[$k])))) + push!(preexprs, :($sym = reindex(V.indexes[$IVindex], Int(J[$Jindex])))) push!(indexexprs, :($sym:$sym)) push!(Itypes, UnitRange{Int}) - elseif k < length(I) || k == NV || j == length(IV) + push!(ItypesLD, j) + elseif Jindex < length(J) || Jindex == NV || IVindex == length(IV) # simple indexing - push!(indexexprs, :(reindex(V.indexes[$j], I[$k]))) - push!(Itypes, rangetype(IV[j], I[k])) + push!(indexexprs, :(reindex(V.indexes[$IVindex], J[$Jindex]))) + push!(Itypes, rangetype(iv, j)) + push!(ItypesLD, Itypes[end]) else # We have a linear index that spans more than one dimension of the parent - push!(indexexprs, :(merge_indexes(V, V.indexes[$j:end], size(V.parent)[$j:end], I[$k]))) + push!(indexexprs, :(merge_indexes(V, V.indexes[$IVindex:end], size(V.parent)[$IVindex:end], J[$Jindex], $Jindex))) push!(Itypes, Array{Int, 1}) + push!(ItypesLD, Itypes[end]) break end if !isLDdone if LD < PLD LD += 1 - Ilast, LD, die_next_vector, isdone = nextLD(Ilast, I[k], LD, die_next_vector) + jprev, LD, die_next_vector, isdone = nextLD(jprev, j, LD, die_next_vector) isLDdone |= isdone else - if I[k] <: Real + if j <: Real LD += 1 else isLDdone = true @@ -227,18 +271,21 @@ stagedfunction sub_unsafe{T,NV,PV,IV,PLD,IndTypes}(V::SubArray{T,NV,PV,IV,PLD}, end end end - for i = k+1:length(I) - if i <= N - push!(sizeexprs, dimsizeexpr(I[i], i, length(I), :V, :I)) + for Jind = Jindex+1:length(J) + j = J[Jind] + if Jind <= N + push!(sizeexprs, dimsizeexpr(j, Jind, length(J), :V, :J)) end - push!(indexexprs, :(I[$i])) - push!(Itypes, I[i]) + push!(indexexprs, :(J[$Jind])) + push!(Itypes, j) + push!(ItypesLD, Itypes[end]) end Inew = :(tuple($(indexexprs...))) dims = :(tuple($(sizeexprs...))) It = tuple(Itypes...) + ItLD = tuple(ItypesLD...) LD = max(LD, subarray_linearindexing_dim(PV, It)) - strideexpr = stride1expr(PV, It, :(V.parent), :Inew, LD) + strideexpr = stride1expr(PV, ItLD, :(V.parent), :Inew, LD) preex = isempty(preexprs) ? nothing : Expr(:block, preexprs...) exfirst = first_index_expr(:(V.parent), :Inew, length(It)) quote @@ -294,6 +341,7 @@ getindex(::Colon, i) = i step(::Colon) = 1 first(::Colon) = 1 +in(::Int, ::Colon) = true ## Strides stagedfunction strides{T,N,P,I}(V::SubArray{T,N,P,I}) @@ -315,7 +363,7 @@ end stride(V::SubArray, d::Integer) = d <= ndims(V) ? strides(V)[d] : strides(V)[end] * size(V)[end] -function stride1expr(Atype::Type, Itypes::Tuple, Aexpr, Inewsym, LD) +function stride1expr(Atype::Type, Itypes::Tuple, Aexpr, Isym, LD, Jsym=Isym, Iindex_lin=0, Jindex_lin=0) if LD == 0 return 0 end @@ -325,13 +373,20 @@ function stride1expr(Atype::Type, Itypes::Tuple, Aexpr, Inewsym, LD) if I <: Real ex = :($ex * size($Aexpr, $d)) else - ex = :($ex * step($Inewsym[$d])) + if d == Iindex_lin + ex = :($ex * step_sa($Jsym[$Jindex_lin])) + else + ex = :($ex * step($Isym[$d])) + end break end end ex end +step_sa(arg) = step(arg) +step_sa(::Integer) = 1 + # This might be conservative, but better safe than sorry function iscontiguous{T,N,P,I,LD}(::Type{SubArray{T,N,P,I,LD}}) LD == length(I) || return false @@ -373,20 +428,20 @@ function first_index_expr(Asym, Isym::Symbol, n::Int) end # Detecting whether one can support fast linear indexing -function nextLD(Ilast, I, LD, die_next_vector) +function nextLD(jprev, j, LD, die_next_vector) isdone = false - if I <: Real - if Ilast != Void && !(Ilast <: Real) + if j <: Real + if jprev != Void && !(jprev <: Real) die_next_vector = true end - elseif die_next_vector || I <: Vector + elseif die_next_vector || j <: Vector LD -= 1 isdone = true - elseif I == Colon - elseif I <: UnitRange + elseif j == Colon + elseif j <: UnitRange die_next_vector = true - elseif I <: Range - if !(Ilast == Void || Ilast <: Real) + elseif j <: Range + if !(jprev == Void || jprev <: Real) LD -= 1 isdone = true end @@ -394,20 +449,20 @@ function nextLD(Ilast, I, LD, die_next_vector) else error("This shouldn't happen (linear indexing inference)") end - Ilast = I - return Ilast, LD, die_next_vector, isdone + jprev = j + return jprev, LD, die_next_vector, isdone end function subarray_linearindexing_dim{A<:AbstractArray}(::Type{A}, It::Tuple) isa(Base.linearindexing(A), Base.LinearSlow) && return 0 isempty(It) && return 0 - Ilast = Void + jprev = Void LD = 0 die_next_vector = false while LD < length(It) LD += 1 I = It[LD] - Ilast, LD, die_next_vector, isdone = nextLD(Ilast, I, LD, die_next_vector) + jprev, LD, die_next_vector, isdone = nextLD(jprev, I, LD, die_next_vector) if isdone break end diff --git a/test/subarray.jl b/test/subarray.jl index b3780b95afa75..b7be037b24745 100644 --- a/test/subarray.jl +++ b/test/subarray.jl @@ -167,8 +167,8 @@ function test_linear(A, B) end # "mixed" means 2 indexes even for N-dimensional arrays -test_mixed{T}(A::AbstractArray{T,1}, B::Array) = nothing -test_mixed{T}(A::AbstractArray{T,2}, B::Array) = nothing +test_mixed{T}(::AbstractArray{T,1}, ::Array) = nothing +test_mixed{T}(::AbstractArray{T,2}, ::Array) = nothing test_mixed(A, B::Array) = _test_mixed(A, reshape(B, size(A))) function _test_mixed(A, B) L = length(A) @@ -194,10 +194,11 @@ function err_li(I::Tuple, ld::Int, ldc::Int) error("Linear indexing inference mismatch") end -function err_li(S::SubArray, ldc::Int) +function err_li(S::SubArray, ld::Int, szC) println(summary(S)) @show S.indexes - @show ldc + @show ld + @show szC error("Linear indexing inference mismatch") end @@ -228,18 +229,42 @@ function runtests(A::SubArray, I...) AA = copy_to_array(A) # Direct test of linear indexing inference C = Agen_nodrop(AA, I...) - ld = single_stride_dim(C) + Cld = ld = single_stride_dim(C) + Cdim = AIindex = 0 + while Cdim <= Cld && AIindex < length(A.indexes) + AIindex += 1 + if isa(A.indexes[AIindex], Real) + ld += 1 + else + Cdim += 1 + end + end # sub - S = sub(A, I...) + local S + try + S = sub(A, I...) + catch err + @show typeof(A) + @show A.indexes + @show I + rethrow(err) + end ldc = getLD(S) - ldc <= ld || err_li(S, ld) + ldc <= ld || err_li(S, ld, size(C)) test_linear(S, C) test_cartesian(S, C) test_mixed(S, C) # slice - S = slice(A, I...) + try + S = slice(A, I...) + catch err + @show typeof(A) + @show A.indexes + @show I + rethrow(err) + end ldc = getLD(S) - ldc <= ld || err_li(S, ld) + ldc <= ld || err_li(S, ld, size(C)) test_linear(S, C) test_cartesian(S, C) test_mixed(S, C) @@ -248,28 +273,28 @@ end # indexN is a cartesian index, indexNN is a linear index for 2 dimensions, and indexNNN is a linear index for 3 dimensions function runviews{T}(SB::AbstractArray{T,3}, indexN, indexNN, indexNNN) for i3 in indexN, i2 in indexN, i1 in indexN - runtests(A, i1, i2, i3) + runtests(SB, i1, i2, i3) end for i2 in indexNN, i1 in indexN - runtests(A, i1, i2) + runtests(SB, i1, i2) end for i1 in indexNNN - runtests(A, i1) + runtests(SB, i1) end end function runviews{T}(SB::AbstractArray{T,2}, indexN, indexNN, indexNNN) for i2 in indexN, i1 in indexN - runtests(A, i1, i2) + runtests(SB, i1, i2) end for i1 in indexNN - runtests(A, i1) + runtests(SB, i1) end end function runviews{T}(SB::AbstractArray{T,1}, indexN, indexNN, indexNNN) for i1 in indexN - runtests(A, i1) + runtests(SB, i1) end end @@ -277,25 +302,55 @@ runviews{T}(SB::AbstractArray{T,0}, indexN, indexNN, indexNNN) = nothing ######### Tests ######### +testfull = Bool(parseint(get(ENV, "JULIA_TESTFULL", "0"))) + ### Views from Arrays ### -A = reshape(1:5*7*11, 11, 7, 5) index5 = (2, :, 2:5, 1:2:5, [4,1,5]) # all work with at least size 5 index25 = (8, :, 2:11, 12:3:22, [4,1,5,9]) index125 = (113, :, 85:121, 2:15:92, [99,14,103]) -runviews(A, index5, index25, index125) + +if testfull + let A = reshape(1:5*7*11, 11, 7, 5) + runviews(A, index5, index25, index125) + end +end ### Views from views ### -B = reshape(1:13^3, 13, 13, 13) -# "outer" indexes create snips that have at least size 5 along each dimension, with the exception of Int-slicing +# "outer" indexes create snips that have at least size 5 along each dimension, +# with the exception of Int-slicing oindex = (:, 6, 3:7, 13:-2:1, [8,4,6,12,5,7]) -for o3 in oindex, o2 in oindex, o1 in oindex - sliceB = slice(B, o1, o2, o3) - runviews(sliceB, index5, index25, index125) +if testfull + let B = reshape(1:13^3, 13, 13, 13) + for o3 in oindex, o2 in oindex, o1 in oindex + sliceB = slice(B, o1, o2, o3) + runviews(sliceB, index5, index25, index125) + end + end end +if !testfull + let B = reshape(1:13^3, 13, 13, 13) + for oind in ((:,:,:), + (:,:,6), + (:,6,:), + (6,:,:), + (:,3:7,:), + (3:7,:,:), + (3:7,6,:), + (3:7,6,6), + (6,3:7,3:7), + (13:-2:1,:,:), + ([8,4,6,12,5,7],:,3:7), + (6,6,[8,4,6,12,5,7])) + runtests(B, oind...) + sliceB = slice(B, oind) + runviews(sliceB, index5, index25, index125) + end + end +end ####### "Classical" tests #######