Deprecate cumsum, cumprod, cumsum_kbn, and accumulate when dim isn't specified (#24684)

* Deprecate cumsum, cumprod, cumsum_kbn, and accumulate when dim isn't
specified and dimension is larger than one

* Add docstrings for accumulation methods on vectors

Fix bug in accumulate when initial value is provided since the implementation
doesn't support multidimentional arrays.

Move the accumulate methods with initial values to the end

* Adjust docstrings to match methods

* axis -> dim

* Update NEWS.md

Author: andreasnoack

NEWS.md

@@ -663,6 +663,9 @@ Deprecated or removed
   * `a:b` is deprecated for constructing a `StepRange` when `a` and `b` have physical units
     (Dates and Times). Use `a:s:b`, where `s = Dates.Day(1)` or `s = Dates.Second(1)`.
 
+  * `cumsum`, `cumprod`, `accumulate`, and their mutating versions now require a `dim`
+    argument instead of defaulting to using the first dimension ([#24684]).
+
 Command-line option changes
 ---------------------------
 

base/abstractarraymath.jl

@@ -265,12 +265,12 @@ end
 # TODO: Needs a separate IndexCartesian method, this is only fast for IndexLinear
 
 """
-    cumsum_kbn(A, [dim::Integer=1])
+    cumsum_kbn(A, dim::Integer)
 
 Cumulative sum along a dimension, using the Kahan-Babuska-Neumaier compensated summation
-algorithm for additional accuracy. The dimension defaults to 1.
+algorithm for additional accuracy.
 """
-function cumsum_kbn(A::AbstractArray{T}, axis::Integer=1) where T<:AbstractFloat
+function cumsum_kbn(A::AbstractArray{T}, axis::Integer) where T<:AbstractFloat
     dimsA = size(A)
     ndimsA = ndims(A)
     axis_size = dimsA[axis]

base/deprecated.jl

@@ -2105,6 +2105,9 @@ end
     @deprecate chol!(x::Number, uplo) chol(x) false
 end
 
+@deprecate cumsum(A::AbstractArray)     cumsum(A, 1)
+@deprecate cumsum_kbn(A::AbstractArray) cumsum_kbn(A, 1)
+@deprecate cumprod(A::AbstractArray)    cumprod(A, 1)
 
 # issue #16307
 @deprecate finalizer(o, f::Function) finalizer(f, o)

base/multidimensional.jl

@@ -675,25 +675,25 @@ function accumulate_pairwise(op, v::AbstractVector{T}) where T
     return accumulate_pairwise!(op, out, v)
 end
 
-function cumsum!(out, v::AbstractVector, axis::Integer=1)
+function cumsum!(out, v::AbstractVector, dim::Integer)
     # we dispatch on the possibility of numerical stability issues
-    _cumsum!(out, v, axis, TypeArithmetic(eltype(out)))
+    _cumsum!(out, v, dim, TypeArithmetic(eltype(out)))
 end
 
-function _cumsum!(out, v, axis, ::ArithmeticRounds)
-    axis == 1 ? accumulate_pairwise!(+, out, v) : copy!(out, v)
+function _cumsum!(out, v, dim, ::ArithmeticRounds)
+    dim == 1 ? accumulate_pairwise!(+, out, v) : copy!(out, v)
 end
-function _cumsum!(out, v, axis, ::ArithmeticUnknown)
-    _cumsum!(out, v, axis, ArithmeticRounds())
+function _cumsum!(out, v, dim, ::ArithmeticUnknown)
+    _cumsum!(out, v, dim, ArithmeticRounds())
 end
-function _cumsum!(out, v, axis, ::TypeArithmetic)
-    axis == 1 ? accumulate!(+, out, v) : copy!(out, v)
+function _cumsum!(out, v, dim, ::TypeArithmetic)
+    dim == 1 ? accumulate!(+, out, v) : copy!(out, v)
 end
 
 """
-    cumsum(A, dim=1)
+    cumsum(A, dim::Integer)
 
-Cumulative sum along a dimension `dim`. See also [`cumsum!`](@ref)
+Cumulative sum along the dimension `dim`. See also [`cumsum!`](@ref)
 to use a preallocated output array, both for performance and to control the precision of the
 output (e.g. to avoid overflow).
 
@@ -714,22 +714,52 @@ julia> cumsum(a,2)
  4  9  15
 ```
 """
-function cumsum(A::AbstractArray{T}, axis::Integer=1) where T
+function cumsum(A::AbstractArray{T}, dim::Integer) where T
     out = similar(A, rcum_promote_type(+, T))
-    cumsum!(out, A, axis)
+    cumsum!(out, A, dim)
 end
 
 """
-    cumsum!(B, A, dim::Integer=1)
+    cumsum(x::AbstractVector)
 
-Cumulative sum of `A` along a dimension, storing the result in `B`. See also [`cumsum`](@ref).
+Cumulative sum a vector. See also [`cumsum!`](@ref)
+to use a preallocated output array, both for performance and to control the precision of the
+output (e.g. to avoid overflow).
+
+```jldoctest
+julia> cumsum([1, 1, 1])
+3-element Array{Int64,1}:
+ 1
+ 2
+ 3
+
+julia> cumsum([fill(1, 2) for i in 1:3])
+3-element Array{Array{Int64,1},1}:
+ [1, 1]
+ [2, 2]
+ [3, 3]
+```
+"""
+cumsum(x::AbstractVector) = cumsum(x, 1)
+
+"""
+    cumsum!(B, A, dim::Integer)
+
+Cumulative sum of `A` along the dimension `dim`, storing the result in `B`. See also [`cumsum`](@ref).
+"""
+cumsum!(B, A, dim::Integer) = accumulate!(+, B, A, dim)
+
+"""
+    cumsum!(y::AbstractVector, x::AbstractVector)
+
+Cumulative sum of a vector `x`, storing the result in `y`. See also [`cumsum`](@ref).
 """
-cumsum!(B, A, axis::Integer=1) = accumulate!(+, B, A, axis)
+cumsum!(y::AbstractVector, x::AbstractVector) = cumsum!(y, x, 1)
 
 """
-    cumprod(A, dim=1)
+    cumprod(A, dim::Integer)
 
-Cumulative product along a dimension `dim`. See also
+Cumulative product along the dimension `dim`. See also
 [`cumprod!`](@ref) to use a preallocated output array, both for performance and
 to control the precision of the output (e.g. to avoid overflow).
 
@@ -750,20 +780,79 @@ julia> cumprod(a,2)
  4  20  120
 ```
 """
-cumprod(A::AbstractArray, axis::Integer=1) = accumulate(*, A, axis)
+cumprod(A::AbstractArray, dim::Integer) = accumulate(*, A, dim)
+
+"""
+    cumprod(x::AbstractVector)
+
+Cumulative product of a vector. See also
+[`cumprod!`](@ref) to use a preallocated output array, both for performance and
+to control the precision of the output (e.g. to avoid overflow).
+
+```jldoctest
+julia> cumprod(fill(1//2, 3))
+3-element Array{Rational{Int64},1}:
+ 1//2
+ 1//4
+ 1//8
+
+julia> cumprod([fill(1//3, 2, 2) for i in 1:3])
+3-element Array{Array{Rational{Int64},2},1}:
+ Rational{Int64}[1//3 1//3; 1//3 1//3]
+ Rational{Int64}[2//9 2//9; 2//9 2//9]
+ Rational{Int64}[4//27 4//27; 4//27 4//27]
+```
+"""
+cumprod(x::AbstractVector) = cumprod(x, 1)
 
 """
-    cumprod!(B, A, dim::Integer=1)
+    cumprod!(B, A, dim::Integer)
 
-Cumulative product of `A` along a dimension, storing the result in `B`.
+Cumulative product of `A` along the dimension `dim`, storing the result in `B`.
 See also [`cumprod`](@ref).
 """
-cumprod!(B, A, axis::Integer=1) = accumulate!(*, B, A, axis)
+cumprod!(B, A, dim::Integer) = accumulate!(*, B, A, dim)
 
 """
-    accumulate(op, A, dim=1)
+    cumprod!(y::AbstractVector, x::AbstractVector)
 
-Cumulative operation `op` along a dimension `dim`. See also
+Cumulative product of a vector `x`, storing the result in `y`.
+See also [`cumprod`](@ref).
+"""
+cumprod!(y::AbstractVector, x::AbstractVector) = cumprod!(y, x, 1)
+
+"""
+    accumulate(op, A, dim::Integer)
+
+Cumulative operation `op` along the dimension `dim`. See also
+[`accumulate!`](@ref) to use a preallocated output array, both for performance and
+to control the precision of the output (e.g. to avoid overflow). For common operations
+there are specialized variants of `accumulate`, see:
+[`cumsum`](@ref), [`cumprod`](@ref)
+
+```jldoctest
+julia> accumulate(+, fill(1, 3, 3), 1)
+3×3 Array{Int64,2}:
+ 1  1  1
+ 2  2  2
+ 3  3  3
+
+julia> accumulate(+, fill(1, 3, 3), 2)
+3×3 Array{Int64,2}:
+ 1  2  3
+ 1  2  3
+ 1  2  3
+```
+"""
+function accumulate(op, A, dim::Integer)
+    out = similar(A, rcum_promote_type(op, eltype(A)))
+    accumulate!(op, out, A, dim)
+end
+
+"""
+    accumulate(op, x::AbstractVector)
+
+Cumulative operation `op` on a vector. See also
 [`accumulate!`](@ref) to use a preallocated output array, both for performance and
 to control the precision of the output (e.g. to avoid overflow). For common operations
 there are specialized variants of `accumulate`, see:
@@ -783,14 +872,67 @@ julia> accumulate(*, [1,2,3])
  6
 ```
 """
-function accumulate(op, A, axis::Integer=1)
-    out = similar(A, rcum_promote_type(op, eltype(A)))
-    accumulate!(op, out, A, axis)
+accumulate(op, x::AbstractVector) = accumulate(op, x, 1)
+
+"""
+    accumulate!(op, B, A, dim::Integer)
+
+Cumulative operation `op` on `A` along the dimension `dim`, storing the result in `B`.
+See also [`accumulate`](@ref).
+"""
+function accumulate!(op, B, A, dim::Integer)
+    dim > 0 || throw(ArgumentError("dim must be a positive integer"))
+    inds_t = indices(A)
+    indices(B) == inds_t || throw(DimensionMismatch("shape of B must match A"))
+    dim > ndims(A) && return copy!(B, A)
+    isempty(inds_t[dim]) && return B
+    if dim == 1
+        # We can accumulate to a temporary variable, which allows
+        # register usage and will be slightly faster
+        ind1 = inds_t[1]
+        @inbounds for I in CartesianRange(tail(inds_t))
+            tmp = convert(eltype(B), A[first(ind1), I])
+            B[first(ind1), I] = tmp
+            for i_1 = first(ind1)+1:last(ind1)
+                tmp = op(tmp, A[i_1, I])
+                B[i_1, I] = tmp
+            end
+        end
+    else
+        R1 = CartesianRange(indices(A)[1:dim-1])   # not type-stable
+        R2 = CartesianRange(indices(A)[dim+1:end])
+        _accumulate!(op, B, A, R1, inds_t[dim], R2) # use function barrier
+    end
+    return B
 end
 
+"""
+    accumulate!(op, y, x::AbstractVector)
 
+Cumulative operation `op` on a vector `x`, storing the result in `y`.
+See also [`accumulate`](@ref).
 """
-    accumulate(op, v0, A)
+function accumulate!(op::Op, y, x::AbstractVector) where Op
+    isempty(x) && return y
+    v1 = first(x)
+    _accumulate1!(op, y, v1, x, 1)
+end
+
+@noinline function _accumulate!(op, B, A, R1, ind, R2)
+    # Copy the initial element in each 1d vector along dimension `dim`
+    ii = first(ind)
+    @inbounds for J in R2, I in R1
+        B[I, ii, J] = A[I, ii, J]
+    end
+    # Accumulate
+    @inbounds for J in R2, i in first(ind)+1:last(ind), I in R1
+        B[I, i, J] = op(B[I, i-1, J], A[I, i, J])
+    end
+    B
+end
+
+"""
+    accumulate(op, v0, x::AbstractVector)
 
 Like `accumulate`, but using a starting element `v0`. The first entry of the result will be
 `op(v0, first(A))`.
@@ -810,30 +952,23 @@ julia> accumulate(min, 0, [1,2,-1])
  -1
 ```
 """
-function accumulate(op, v0, A, axis::Integer=1)
-    T = rcum_promote_type(op, typeof(v0), eltype(A))
-    out = similar(A, T)
-    accumulate!(op, out, v0, A, 1)
-end
-
-function accumulate!(op::Op, B, A::AbstractVector, axis::Integer=1) where Op
-    isempty(A) && return B
-    v1 = first(A)
-    _accumulate1!(op, B, v1, A, axis)
+function accumulate(op, v0, x::AbstractVector)
+    T = rcum_promote_type(op, typeof(v0), eltype(x))
+    out = similar(x, T)
+    accumulate!(op, out, v0, x)
 end
 
-function accumulate!(op, B, v0, A::AbstractVector, axis::Integer=1)
-    isempty(A) && return B
-    v1 = op(v0, first(A))
-    _accumulate1!(op, B, v1, A, axis)
+function accumulate!(op, y, v0, x::AbstractVector)
+    isempty(x) && return y
+    v1 = op(v0, first(x))
+    _accumulate1!(op, y, v1, x, 1)
 end
 
-
-function _accumulate1!(op, B, v1, A::AbstractVector, axis::Integer=1)
-    axis > 0 || throw(ArgumentError("axis must be a positive integer"))
+function _accumulate1!(op, B, v1, A::AbstractVector, dim::Integer)
+    dim > 0 || throw(ArgumentError("dim must be a positive integer"))
     inds = linearindices(A)
     inds == linearindices(B) || throw(DimensionMismatch("linearindices of A and B don't match"))
-    axis > 1 && return copy!(B, A)
+    dim > 1 && return copy!(B, A)
     i1 = inds[1]
     cur_val = v1
     B[i1] = cur_val
@@ -844,51 +979,6 @@ function _accumulate1!(op, B, v1, A::AbstractVector, axis::Integer=1)
     return B
 end
 
-"""
-    accumulate!(op, B, A, dim=1)
-
-Cumulative operation `op` on `A` along a dimension, storing the result in `B`.
-See also [`accumulate`](@ref).
-"""
-function accumulate!(op, B, A, axis::Integer=1)
-    axis > 0 || throw(ArgumentError("axis must be a positive integer"))
-    inds_t = indices(A)
-    indices(B) == inds_t || throw(DimensionMismatch("shape of B must match A"))
-    axis > ndims(A) && return copy!(B, A)
-    isempty(inds_t[axis]) && return B
-    if axis == 1
-        # We can accumulate to a temporary variable, which allows
-        # register usage and will be slightly faster
-        ind1 = inds_t[1]
-        @inbounds for I in CartesianRange(tail(inds_t))
-            tmp = convert(eltype(B), A[first(ind1), I])
-            B[first(ind1), I] = tmp
-            for i_1 = first(ind1)+1:last(ind1)
-                tmp = op(tmp, A[i_1, I])
-                B[i_1, I] = tmp
-            end
-        end
-    else
-        R1 = CartesianRange(indices(A)[1:axis-1])   # not type-stable
-        R2 = CartesianRange(indices(A)[axis+1:end])
-        _accumulate!(op, B, A, R1, inds_t[axis], R2) # use function barrier
-    end
-    return B
-end
-
-@noinline function _accumulate!(op, B, A, R1, ind, R2)
-    # Copy the initial element in each 1d vector along dimension `axis`
-    ii = first(ind)
-    @inbounds for J in R2, I in R1
-        B[I, ii, J] = A[I, ii, J]
-    end
-    # Accumulate
-    @inbounds for J in R2, i in first(ind)+1:last(ind), I in R1
-        B[I, i, J] = op(B[I, i-1, J], A[I, i, J])
-    end
-    B
-end
-
 ### from abstractarray.jl
 
 """