Define reducedim_initarray(A, region, ::UndefInitializer, ::Type) for easier initialization

also simplify min/max init

Replace `UndefInitializer` with `_InitialValue`

Author: N5N3

base/reducedim.jl

@@ -91,7 +91,7 @@ end
 # reducedim_initarray is called by
 reducedim_initarray(A::AbstractArrayOrBroadcasted, region, init, ::Type{R}) where {R} = fill!(similar(A,R,reduced_indices(A,region)), init)
 reducedim_initarray(A::AbstractArrayOrBroadcasted, region, init::T) where {T} = reducedim_initarray(A, region, init, T)
-
+reducedim_initarray(A::AbstractArrayOrBroadcasted, region, ::_InitialValue, ::Type{R}) where {R} = similar(A,R,reduced_indices(A,region))
 # TODO: better way to handle reducedim initialization
 #
 # The current scheme is basically following Steven G. Johnson's original implementation
@@ -125,46 +125,38 @@ function _reducedim_init(f, op, fv, fop, A, region)
 end
 
 # initialization when computing minima and maxima requires a little care
-for (f1, f2, initval, typeextreme) in ((:min, :max, :Inf, :typemax), (:max, :min, :(-Inf), :typemin))
-    @eval 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 = reduced_indices(A, region)
+function reducedim_init(f::F, ::Union{typeof(min),typeof(max)}, A::AbstractArray, region) where {F}
+    # First compute the reduce indices. This will throw an ArgumentError
+    # if any region is invalid
+    ri = reduced_indices(A, region)
 
-        # Next, throw if reduction is over a region with length zero
-        any(i -> isempty(axes(A, i)), region) && _empty_reduce_error()
+    # Next, throw if reduction is over a region with length zero
+    any(i -> isempty(axes(A, i)), region) && _empty_reduce_error()
 
-        # Make a view of the first slice of the region
-        A1 = view(A, ri...)
+    # 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 map(f, A1)
-        else
-            # otherwise use the min/max of the first slice as initial value
-            v0 = mapreduce(f, $f2, A1)
-
-            T = _realtype(f, promote_union(eltype(A)))
-            Tr = v0 isa T ? T : typeof(v0)
-
-            # but NaNs and missing need to be avoided as initial values
-            if v0 isa Number && isnan(v0)
-                # v0 is NaN
-                v0 = oftype(v0, $initval)
-            elseif isunordered(v0)
-                # v0 is missing or a third-party unordered value
-                Tnm = nonmissingtype(Tr)
-                # TODO: Some types, like BigInt, don't support typemin/typemax.
-                # So a Matrix{Union{BigInt, Missing}} can still error here.
-                v0 = $typeextreme(Tnm)
-            end
-            # v0 may have changed type.
-            Tr = v0 isa T ? T : typeof(v0)
+    # calculate the output type
+    T = promote_typejoin_union(_return_type(f, Tuple{eltype(A)}))
 
-            return reducedim_initarray(A, region, v0, Tr)
-        end
+    map!(f, reducedim_initarray(A,region,_InitialValue(),T), A1)
+end
+
+function reducedim_init(f::ExtremaMap, ::typeof(_extrema_rf), A::AbstractArray, region)
+    ri = reduced_indices(A, region)
+    any(i -> isempty(axes(A, i)), region) && _empty_reduce_error()
+    A1 = view(A, ri...)
+    IT = eltype(A)
+    if missing isa IT
+        RT = promote_typejoin_union(_return_type(f.f, Tuple{nonmissingtype(IT)}))
+        T = Union{Tuple{RT,RT},Tuple{Missing,Missing}}
+    else
+        RT = promote_typejoin_union(_return_type(f.f, Tuple{IT}))
+        T = Union{Tuple{RT,RT}}
     end
+    map!(f, reducedim_initarray(A,region,_InitialValue(),T), A1)
 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))), _realtype(f, T))
 

stdlib/LinearAlgebra/src/special.jl

@@ -7,7 +7,8 @@
 # Diagonal/Bidiagonal/Tridiagonal/SymTridiagonal matrix. However, reducedim should
 # yield a dense vector to increase performance.
 Base.reducedim_initarray(A::Union{Diagonal,Bidiagonal,Tridiagonal,SymTridiagonal}, region, init, ::Type{R}) where {R} = fill(convert(R, init), Base.reduced_indices(A,region))
-
+Base.reducedim_initarray(A::Union{Diagonal,Bidiagonal,Tridiagonal,SymTridiagonal}, region, ::Base._InitialValue, ::Type{R}) where {R} =
+    Array{R}(undef, Base.to_shape(Base.reduced_indices(A, region)))
 
 # Interconversion between special matrix types