Merge pull request #19100 from JuliaLang/teh/backport_19015

[release-0.5] Backport #19015

base/abstractarray.jl

@@ -202,6 +202,14 @@ function trailingsize(A, n)
  end
  return s
 end
+function trailingsize(inds::Indices, n)
+ s = 1
+ for i=n:length(inds)
+ s *= unsafe_length(inds[i])
+ end
+ return s
+end
+# This version is type-stable even if inds is heterogeneous
 function trailingsize(inds::Indices)
  @_inline_meta
  prod(map(unsafe_length, inds))

base/multidimensional.jl

@@ -207,10 +207,10 @@ index_ndims() = ()
 
 # Recursively compute the lengths of a list of indices, without dropping scalars
 # These need to be inlined for more than 3 indexes
-index_lengths(A::AbstractArray, I::Colon) = (length(A),)
+index_lengths(A::AbstractArray, I::Colon) = (_length(A),)
 @inline index_lengths(A::AbstractArray, I...) = index_lengths_dim(A, 1, I...)
 index_lengths_dim(A, dim) = ()
-index_lengths_dim(A, dim, ::Colon) = (trailingsize(A, dim),)
+index_lengths_dim(A, dim, ::Colon) = (trailingsize(indices(A), dim),)
 @inline index_lengths_dim(A, dim, ::Colon, i, I...) = (_length(indices(A, dim)), index_lengths_dim(A, dim+1, i, I...)...)
 @inline index_lengths_dim(A, dim, ::Real, I...) = (1, index_lengths_dim(A, dim+1, I...)...)
 @inline index_lengths_dim{N}(A, dim, ::CartesianIndex{N}, I...) = (1, index_lengths_dim(A, dim+N, I...)...)

base/reducedim.jl

@@ -2,71 +2,86 @@
 
 ## Functions to compute the reduced shape
 
+function reduced_dims(a::AbstractArray, region)
+ to_shape(reduced_indices(a, region)) # to_shape keeps the return-type consistent, when it's possible to do so
+end
+
+function reduced_dims0(a::AbstractArray, region)
+ to_shape(reduced_indices0(a, region))
+end
+
 # for reductions that expand 0 dims to 1
-reduced_dims(a::AbstractArray, region) = reduced_dims(size(a), region)
+reduced_indices(a::AbstractArray, region) = reduced_indices(indices(a), region)
 
 # for reductions that keep 0 dims as 0
-reduced_dims0(a::AbstractArray, region) = reduced_dims0(size(a), region)
-
-function reduced_dims{N}(siz::NTuple{N,Int}, d::Int, rd::Int)
- if d < 1
- throw(ArgumentError("dimension must be ≥ 1, got $d"))
- elseif d == 1
- return tuple(rd, siz[d+1:N]...)::typeof(siz)
- elseif 1 < d < N
- return tuple(siz[1:d-1]..., rd, siz[d+1:N]...)::typeof(siz)
- elseif d == N
- return tuple(siz[1:N-1]..., rd)::typeof(siz)
+reduced_indices0(a::AbstractArray, region) = reduced_indices0(indices(a), region)
+
+function reduced_indices{N}(inds::Indices{N}, d::Int, rd::AbstractUnitRange)
+ d < 1 && throw(ArgumentError("dimension must be ≥ 1, got $d"))
+ if d == 1
+ return (oftype(inds[1], rd), tail(inds)...)
+ elseif 1 < d <= N
+ return tuple(inds[1:d-1]..., oftype(inds[d], rd), inds[d+1:N]...)::typeof(inds)
  else
- return siz
+ return inds
  end
 end
-reduced_dims{N}(siz::NTuple{N,Int}, d::Int) = reduced_dims(siz, d, 1)
+reduced_indices{N}(inds::Indices{N}, d::Int) = reduced_indices(inds, d, OneTo(1))
 
-function reduced_dims0{N}(siz::NTuple{N,Int}, d::Int)
- if d < 1
- throw(ArgumentError("dimension must be ≥ 1, got $d"))
- elseif d <= N
- return reduced_dims(siz, d, (siz[d] == 0 ? 0 : 1))
+function reduced_indices0{N}(inds::Indices{N}, d::Int)
+ d < 1 && throw(ArgumentError("dimension must be ≥ 1, got $d"))
+ if d <= N
+ return reduced_indices(inds, d, (inds[d] == OneTo(0) ? OneTo(0) : OneTo(1)))
  else
- return siz
+ return inds
  end
 end
 
-function reduced_dims{N}(siz::NTuple{N,Int}, region)
- rsiz = [siz...]
+function reduced_indices{N}(inds::Indices{N}, region)
+ rinds = [inds...]
  for i in region
  isa(i, Integer) || throw(ArgumentError("reduced dimension(s) must be integers"))
- d = convert(Int, i)::Int
+ d = Int(i)
  if d < 1
  throw(ArgumentError("region dimension(s) must be ≥ 1, got $d"))
  elseif d <= N
- rsiz[d] = 1
+ rinds[d] = oftype(rinds[d], OneTo(1))
  end
  end
- tuple(rsiz...)::typeof(siz)
+ tuple(rinds...)::typeof(inds)
 end
 
-function reduced_dims0{N}(siz::NTuple{N,Int}, region)
- rsiz = [siz...]
+function reduced_indices0{N}(inds::Indices{N}, region)
+ rinds = [inds...]
  for i in region
  isa(i, Integer) || throw(ArgumentError("reduced dimension(s) must be integers"))
- d = convert(Int, i)::Int
+ d = Int(i)
  if d < 1
  throw(ArgumentError("region dimension(s) must be ≥ 1, got $d"))
  elseif d <= N
- rsiz[d] = (rsiz[d] == 0 ? 0 : 1)
+ rinds[d] = oftype(rinds[d], (rinds[d] == OneTo(0) ? OneTo(0) : OneTo(1)))
  end
  end
- tuple(rsiz...)::typeof(siz)
+ tuple(rinds...)::typeof(inds)
 end
 
-function regionsize(a, region)
- s = 1
- for d in region
- s *= size(a,d)
- end
- s
+# The following are deprecated in julia-0.6
+function reduced_dims(::Tuple{}, d::Int)
+ d < 1 && throw(ArgumentError("dimension must be ≥ 1, got $d"))
+ ()
+end
+reduced_dims(::Tuple{}, region) = ()
+function reduced_dims(dims::Dims, region)
+ map(last, reduced_dims(map(n->OneTo(n), dims), region))
+end
+
+function reduced_dims0(::Tuple{}, d::Int)
+ d < 1 && throw(ArgumentError("dimension must be ≥ 1, got $d"))
+ ()
+end
+reduced_dims0(::Tuple{}, region) = ()
+function reduced_dims0(dims::Dims, region)
+ map(last, reduced_dims0(map(n->OneTo(n), dims), region))
 end
 
 
@@ -82,10 +97,10 @@ for (Op, initval) in ((:(typeof(&)), true), (:(typeof(|)), false))
  @eval initarray!(a::AbstractArray, ::$(Op), init::Bool) = (init && fill!(a, $initval); a)
 end
 
-reducedim_initarray{R}(A::AbstractArray, region, v0, ::Type{R}) = fill!(similar(A,R,reduced_dims(A,region)), v0)
+reducedim_initarray{R}(A::AbstractArray, region, v0, ::Type{R}) = fill!(similar(A,R,reduced_indices(A,region)), v0)
 reducedim_initarray{T}(A::AbstractArray, region, v0::T) = reducedim_initarray(A, region, v0, T)
 
-reducedim_initarray0{R}(A::AbstractArray, region, v0, ::Type{R}) = fill!(similar(A,R,reduced_dims0(A,region)), v0)
+reducedim_initarray0{R}(A::AbstractArray, region, v0, ::Type{R}) = fill!(similar(A,R,reduced_indices0(A,region)), v0)
 reducedim_initarray0{T}(A::AbstractArray, region, v0::T) = reducedim_initarray0(A, region, v0, T)
 
 # TODO: better way to handle reducedim initialization
@@ -387,11 +402,11 @@ For an array input, returns the value and index of the minimum over the given re
 """
 function findmin{T}(A::AbstractArray{T}, region)
  if isempty(A)
- return (similar(A, reduced_dims0(A, region)),
- zeros(Int, reduced_dims0(A, region)))
+ return (similar(A, reduced_indices0(A, region)),
+ similar(dims->zeros(Int, dims), reduced_indices0(A, region)))
  end
  return findminmax!(<, reducedim_initarray0(A, region, typemax(T)),
- zeros(Int, reduced_dims0(A, region)), A)
+ similar(dims->zeros(Int, dims), reduced_indices0(A, region)), A)
 end
 
 """
@@ -414,11 +429,11 @@ For an array input, returns the value and index of the maximum over the given re
 """
 function findmax{T}(A::AbstractArray{T}, region)
  if isempty(A)
- return (similar(A, reduced_dims0(A,region)),
- zeros(Int, reduced_dims0(A,region)))
+ return (similar(A, reduced_indices0(A,region)),
+ similar(dims->zeros(Int, dims), reduced_indices0(A,region)))
  end
  return findminmax!(>, reducedim_initarray0(A, region, typemin(T)),
- zeros(Int, reduced_dims0(A, region)), A)
+ similar(dims->zeros(Int, dims), reduced_indices0(A, region)), A)
 end
 
 reducedim1(R, A) = length(indices1(R)) == 1

base/sort.jl

@@ -2,7 +2,7 @@
 
 module Sort
 
-using Base: Order, copymutable, linearindices, linearindexing, viewindexing, LinearFast
+using Base: Order, copymutable, linearindices, linearindexing, viewindexing, LinearFast, _length
 
 import
  Base.sort,
@@ -59,7 +59,8 @@ issorted(itr;
  issorted(itr, ord(lt,by,rev,order))
 
 function select!(v::AbstractVector, k::Union{Int,OrdinalRange}, o::Ordering)
- sort!(v, 1, length(v), PartialQuickSort(k), o)
+ inds = indices(v, 1)
+ sort!(v, first(inds), last(inds), PartialQuickSort(k), o)
  v[k]
 end
 select!(v::AbstractVector, k::Union{Int,OrdinalRange};
@@ -180,7 +181,7 @@ searchsorted{T<:Real}(a::Range{T}, x::Real, o::DirectOrdering) =
 
 for s in [:searchsortedfirst, :searchsortedlast, :searchsorted]
  @eval begin
- $s(v::AbstractVector, x, o::Ordering) = $s(v,x,1,length(v),o)
+ $s(v::AbstractVector, x, o::Ordering) = (inds = indices(v, 1); $s(v,x,first(inds),last(inds),o))
  $s(v::AbstractVector, x;
  lt=isless, by=identity, rev::Bool=false, order::Ordering=Forward) =
  $s(v,x,ord(lt,by,rev,order))
@@ -420,7 +421,7 @@ sort(v::AbstractVector; kws...) = sort!(copymutable(v); kws...)
 ## selectperm: the permutation to sort the first k elements of an array ##
 
 selectperm(v::AbstractVector, k::Union{Integer,OrdinalRange}; kwargs...) =
- selectperm!(Vector{eltype(k)}(length(v)), v, k; kwargs..., initialized=false)
+ selectperm!(similar(Vector{eltype(k)}, indices(v,1)), v, k; kwargs..., initialized=false)
 
 function selectperm!{I<:Integer}(ix::AbstractVector{I}, v::AbstractVector,
  k::Union{Int, OrdinalRange};
@@ -463,7 +464,7 @@ function sortperm!{I<:Integer}(x::AbstractVector{I}, v::AbstractVector;
  order::Ordering=Forward,
  initialized::Bool=false)
  if indices(x,1) != indices(v,1)
- throw(ArgumentError("index vector must be the same length as the source vector, $(indices(x,1)) != $(indices(v,1))"))
+ throw(ArgumentError("index vector must have the same indices as the source vector, $(indices(x,1)) != $(indices(v,1))"))
  end
  if !initialized
  @inbounds for i = indices(v,1)

base/sparse/sparsematrix.jl

@@ -1793,9 +1793,9 @@ end
 # and computing reductions along columns into SparseMatrixCSC is
 # non-trivial, so use Arrays for output
 Base.reducedim_initarray{R}(A::SparseMatrixCSC, region, v0, ::Type{R}) =
- fill!(Array{R}(Base.reduced_dims(A,region)), v0)
+ fill!(similar(dims->Array{R}(dims), Base.reduced_indices(A,region)), v0)
 Base.reducedim_initarray0{R}(A::SparseMatrixCSC, region, v0, ::Type{R}) =
- fill!(Array{R}(Base.reduced_dims0(A,region)), v0)
+ fill!(similar(dims->Array{R}(dims), Base.reduced_indices0(A,region)), v0)
 
 # General mapreduce
 function _mapreducezeros(f, op, T::Type, nzeros::Int, v0)

base/statistics.jl

@@ -24,12 +24,12 @@ function mean(f::Callable, iterable)
  return total/count
 end
 mean(iterable) = mean(identity, iterable)
-mean(f::Callable, A::AbstractArray) = sum(f, A) / length(A)
-mean(A::AbstractArray) = sum(A) / length(A)
+mean(f::Callable, A::AbstractArray) = sum(f, A) / _length(A)
+mean(A::AbstractArray) = sum(A) / _length(A)
 
 function mean!{T}(R::AbstractArray{T}, A::AbstractArray)
  sum!(R, A; init=true)
- scale!(R, length(R) / length(A))
+ scale!(R, _length(R) / _length(A))
  return R
 end
 
@@ -140,7 +140,7 @@ function centralize_sumabs2!{S,T,N}(R::AbstractArray{S}, A::AbstractArray{T,N},
 end
 
 function varm{T}(A::AbstractArray{T}, m::Number; corrected::Bool=true)
- n = length(A)
+ n = _length(A)
  n == 0 && return convert(real(momenttype(T)), NaN)
  n == 1 && return convert(real(momenttype(T)), abs2(A[1] - m)/(1 - Int(corrected)))
  return centralize_sumabs2(A, m) / (n - Int(corrected))
@@ -150,7 +150,7 @@ function varm!{S}(R::AbstractArray{S}, A::AbstractArray, m::AbstractArray; corre
  if isempty(A)
  fill!(R, convert(S, NaN))
  else
- rn = div(length(A), length(R)) - Int(corrected)
+ rn = div(_length(A), _length(R)) - Int(corrected)
  scale!(centralize_sumabs2!(R, A, m), convert(S, 1/rn))
  end
  return R
@@ -266,7 +266,7 @@ stdm(iterable, m::Number; corrected::Bool=true) =
 _conj{T<:Real}(x::AbstractArray{T}) = x
 _conj(x::AbstractArray) = conj(x)
 
-_getnobs(x::AbstractVector, vardim::Int) = length(x)
+_getnobs(x::AbstractVector, vardim::Int) = _length(x)
 _getnobs(x::AbstractMatrix, vardim::Int) = size(x, vardim)
 
 function _getnobs(x::AbstractVecOrMat, y::AbstractVecOrMat, vardim::Int)
@@ -293,11 +293,11 @@ unscaled_covzm(x::AbstractMatrix, y::AbstractMatrix, vardim::Int) =
 
 # covzm (with centered data)
 
-covzm(x::AbstractVector, corrected::Bool=true) = unscaled_covzm(x) / (length(x) - Int(corrected))
+covzm(x::AbstractVector, corrected::Bool=true) = unscaled_covzm(x) / (_length(x) - Int(corrected))
 covzm(x::AbstractMatrix, vardim::Int=1, corrected::Bool=true) =
  scale!(unscaled_covzm(x, vardim), inv(size(x,vardim) - Int(corrected)))
 covzm(x::AbstractVector, y::AbstractVector, corrected::Bool=true) =
- unscaled_covzm(x, y) / (length(x) - Int(corrected))
+ unscaled_covzm(x, y) / (_length(x) - Int(corrected))
 covzm(x::AbstractVecOrMat, y::AbstractVecOrMat, vardim::Int=1, corrected::Bool=true) =
  scale!(unscaled_covzm(x, y, vardim), inv(_getnobs(x, y, vardim) - Int(corrected)))
 
@@ -553,16 +553,18 @@ function median!{T}(v::AbstractVector{T})
  isnan(x) && return x
  end
  end
- n = length(v)
+ inds = indices(v, 1)
+ n = length(inds)
+ mid = div(first(inds)+last(inds),2)
  if isodd(n)
- return middle(select!(v,div(n+1,2)))
+ return middle(select!(v,mid))
  else
- m = select!(v, div(n,2):div(n,2)+1)
+ m = select!(v, mid:mid+1)
  return middle(m[1], m[2])
  end
 end
 median!{T}(v::AbstractArray{T}) = median!(vec(v))
-median{T}(v::AbstractArray{T}) = median!(copy!(Array{T,1}(length(v)), v))
+median{T}(v::AbstractArray{T}) = median!(copy!(Array{T,1}(_length(v)), v))
 
 """
  median(v[, region])

test/offsetarray.jl

@@ -42,6 +42,20 @@ S = OffsetArray(view(A0, 1:2, 1:2), (-1,2)) # LinearSlow
 @test A[1, [4,3]] == S[1, [4,3]] == [4,2]
 @test A[:, :] == S[:, :] == A
 
+A_3_3 = OffsetArray(Array{Int}(3,3), (-2,-1))
+A_3_3[:, :] = reshape(1:9, 3, 3)
+for i = 1:9 @test A_3_3[i] == i end
+A_3_3[-1:1, 0:2] = reshape(1:9, 3, 3)
+for i = 1:9 @test A_3_3[i] == i end
+A_3_3[:, :] = 1:9
+for i = 1:9 @test A_3_3[i] == i end
+A_3_3[-1:1, 0:2] = 1:9
+for i = 1:9 @test A_3_3[i] == i end
+A_3_3[:] = 1:9
+for i = 1:9 @test A_3_3[i] == i end
+A_3_3[1:9] = 1:9
+for i = 1:9 @test A_3_3[i] == i end
+
 # CartesianIndexing
 @test A[CartesianIndex((0,3))] == S[CartesianIndex((0,3))] == 1
 @test_throws BoundsError A[CartesianIndex(1,1)]
@@ -259,6 +273,9 @@ A = OffsetArray(rand(4,4), (-3,5))
 @test maximum(A) == maximum(parent(A))
 @test minimum(A) == minimum(parent(A))
 @test extrema(A) == extrema(parent(A))
+@test maximum(A, 1) == OffsetArray(maximum(parent(A), 1), (0,A.offsets[2]))
+@test maximum(A, 2) == OffsetArray(maximum(parent(A), 2), (A.offsets[1],0))
+@test maximum(A, 1:2) == maximum(parent(A), 1:2)
 C = similar(A)
 cumsum!(C, A, 1)
 @test parent(C) == cumsum(parent(A), 1)
@@ -293,6 +310,13 @@ I,J,N = findnz(z)
 @test find(x->x>0, h) == [-1,1]
 @test find(x->x<0, h) == [-2,0]
 @test find(x->x==0, h) == [2]
+@test mean(A_3_3) == median(A_3_3) == 5
+@test mean(x->2x, A_3_3) == 10
+@test mean(A_3_3, 1) == median(A_3_3, 1) == OffsetArray([2 5 8], (0,A_3_3.offsets[2]))
+@test mean(A_3_3, 2) == median(A_3_3, 2) == OffsetArray([4,5,6]'', (A_3_3.offsets[1],0))
+@test var(A_3_3) == 7.5
+@test std(A_3_3, 1) == OffsetArray([1 1 1], (0,A_3_3.offsets[2]))
+@test std(A_3_3, 2) == OffsetArray([3,3,3]'', (A_3_3.offsets[1],0))
 
 @test_approx_eq vecnorm(v) vecnorm(parent(v))
 @test_approx_eq vecnorm(A) vecnorm(parent(A))

test/reducedim.jl

@@ -20,7 +20,7 @@ for region in Any[
  1, 2, 3, 4, 5, (1, 2), (1, 3), (1, 4), (2, 3), (2, 4), (3, 4),
  (1, 2, 3), (1, 3, 4), (2, 3, 4), (1, 2, 3, 4)]
  # println("region = $region")
- r = fill(NaN, Base.reduced_dims(size(Areduc), region))
+ r = fill(NaN, map(length, Base.reduced_indices(indices(Areduc), region)))
  @test sum!(r, Areduc) ≈ safe_sum(Areduc, region)
  @test prod!(r, Areduc) ≈ safe_prod(Areduc, region)
  @test maximum!(r, Areduc) ≈ safe_maximum(Areduc, region)
@@ -62,7 +62,7 @@ end
 # Test reduction along first dimension; this is special-cased for
 # size(A, 1) >= 16
 Breduc = rand(64, 3)
-r = fill(NaN, Base.reduced_dims(size(Breduc), 1))
+r = fill(NaN, map(length, Base.reduced_indices(indices(Breduc), 1)))
 @test sum!(r, Breduc) ≈ safe_sum(Breduc, 1)
 @test sumabs!(r, Breduc) ≈ safe_sumabs(Breduc, 1)
 @test sumabs2!(r, Breduc) ≈ safe_sumabs2(Breduc, 1)