From 31c67c2f615e72b4375e794d24f1b18afdead0b6 Mon Sep 17 00:00:00 2001 From: Andreas Noack Date: Thu, 28 Jun 2018 14:45:42 +0200 Subject: [PATCH] Adjust initialization in maximum and minimum Fixes #27836 --- base/reducedim.jl | 48 ++++++++++++---------- stdlib/SparseArrays/src/sparsematrix.jl | 2 - stdlib/SparseArrays/test/higherorderfns.jl | 3 ++ test/reducedim.jl | 8 ++-- 4 files changed, 34 insertions(+), 27 deletions(-) diff --git a/base/reducedim.jl b/base/reducedim.jl index c8d5e720b6f899..90b0b06dec54b6 100644 --- a/base/reducedim.jl +++ b/base/reducedim.jl @@ -79,25 +79,6 @@ end reducedim_initarray(A::AbstractArray, region, v0, ::Type{R}) where {R} = fill!(similar(A,R,reduced_indices(A,region)), v0) reducedim_initarray(A::AbstractArray, region, v0::T) where {T} = reducedim_initarray(A, region, v0, T) -function reducedim_initarray0(A::AbstractArray{T}, region, f, ops) where T - ri = reduced_indices0(A, region) - if isempty(A) - if prod(length, reduced_indices(A, region)) != 0 - reducedim_initarray0_empty(A, region, f, ops) # ops over empty slice of A - else - R = f == identity ? T : Core.Compiler.return_type(f, (T,)) - similar(A, R, ri) - end - else - R = f == identity ? T : typeof(f(first(A))) - si = similar(A, R, ri) - mapfirst!(f, si, A) - end -end - -reducedim_initarray0_empty(A::AbstractArray, region, f, ops) = mapslices(x->ops(f.(x)), A, region) -reducedim_initarray0_empty(A::AbstractArray, region,::typeof(identity), ops) = mapslices(ops, A, region) - # TODO: better way to handle reducedim initialization # # The current scheme is basically following Steven G. Johnson's original implementation @@ -124,8 +105,33 @@ function _reducedim_init(f, op, fv, fop, A, region) return reducedim_initarray(A, region, z, Tr) end -reducedim_init(f, op::typeof(max), A::AbstractArray{T}, region) where {T} = reducedim_initarray0(A, region, f, maximum) -reducedim_init(f, op::typeof(min), A::AbstractArray{T}, region) where {T} = reducedim_initarray0(A, region, f, minimum) +# initialization when computing minima and maxima requires a little care +for (f1, f2, supval) in ((min, max, Inf), (max, min, -Inf)) + function reducedim_init(f, op::typeof(f1), A::AbstractArray, region) + # First compute the reduce indices. This will throw an ArgumentError + # if any region is invalid + ri = Base.reduced_indices(A, region) + + # Next, throw if reduction is over a region with length zero + any(i -> iszero(size(A, i)), region) && _empty_reduce_error() + + # Make a view of the first slice of the region + A1 = view(A, ri...) + + if isempty(A1) + # If the slice is empty just return non-view version as the initial array + return copy(A1) + else + # otherwise use the min/max of the first slice as initial value + v0 = mapreduce(f, f2, A1) + + # but NaNs need to be avoided as intial values + v0 = v0 != v0 ? typeof(v0)(supval) : v0 + + return reducedim_initarray(A, region, v0) + end + end +end reducedim_init(f::Union{typeof(abs),typeof(abs2)}, op::typeof(max), A::AbstractArray{T}, region) where {T} = reducedim_initarray(A, region, zero(f(zero(T)))) diff --git a/stdlib/SparseArrays/src/sparsematrix.jl b/stdlib/SparseArrays/src/sparsematrix.jl index 47fe84a7cd194f..cf0cf6b6ea77d3 100644 --- a/stdlib/SparseArrays/src/sparsematrix.jl +++ b/stdlib/SparseArrays/src/sparsematrix.jl @@ -1608,8 +1608,6 @@ end # non-trivial, so use Arrays for output Base.reducedim_initarray(A::SparseMatrixCSC, region, v0, ::Type{R}) where {R} = fill(v0, Base.reduced_indices(A,region)) -Base.reducedim_initarray0(A::SparseMatrixCSC, region, v0, ::Type{R}) where {R} = - fill(v0, Base.reduced_indices0(A,region)) # General mapreduce function _mapreducezeros(f, op, ::Type{T}, nzeros::Int, v0) where T diff --git a/stdlib/SparseArrays/test/higherorderfns.jl b/stdlib/SparseArrays/test/higherorderfns.jl index f1bd3f5d25c9ca..5871e2d96e1e23 100644 --- a/stdlib/SparseArrays/test/higherorderfns.jl +++ b/stdlib/SparseArrays/test/higherorderfns.jl @@ -628,5 +628,8 @@ end @test Diagonal(spzeros(5)) \ view(rand(10), 1:5) == [Inf,Inf,Inf,Inf,Inf] end +@testset "Issue #27836" begin + @test minimum(sparse([1, 2], [1, 2], ones(Int32, 2)), dims = 1) isa Matrix +end end # module diff --git a/test/reducedim.jl b/test/reducedim.jl index 2c9f622cb49058..59dc1d202e823f 100644 --- a/test/reducedim.jl +++ b/test/reducedim.jl @@ -18,7 +18,7 @@ safe_maxabs(A::Array{T}, region) where {T} = safe_mapslices(maximum, abs.(A), re safe_minabs(A::Array{T}, region) where {T} = safe_mapslices(minimum, abs.(A), region) Areduc = rand(3, 4, 5, 6) -for region in Any[ +@testset "test reductions over region: $region" 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") @@ -318,10 +318,10 @@ end @test sum(Any[1 2;3 4], dims=1) == [4 6] @test sum(Vector{Int}[[1,2],[4,3]], dims=1)[1] == [5,5] -# issue #10461 -Areduc = rand(3, 4, 5, 6) -for region in Any[-1, 0, (-1, 2), [0, 1], (1,-2,3), [0 1; +@testset "Issue #10461. region=$region" for region in Any[-1, 0, (-1, 2), [0, 1], (1,-2,3), [0 1; 2 3], "hello"] + Areduc = rand(3, 4, 5, 6) + @test_throws ArgumentError sum(Areduc, dims=region) @test_throws ArgumentError prod(Areduc, dims=region) @test_throws ArgumentError maximum(Areduc, dims=region)