Implement sum

src/BandedMatrices.jl

@@ -7,7 +7,7 @@ import Base: axes, axes1, getproperty, getindex, setindex!, *, +, -, ==, <, <=,
                >=, /, \, adjoint, transpose, showerror, convert, size, view,
                unsafe_indices, first, last, size, length, unsafe_length, step, to_indices,
                to_index, show, fill!, similar, copy, promote_rule, real, imag,
-               copyto!, Array
+               copyto!, Array, sum, sum!
 
 using Base.Broadcast: AbstractArrayStyle, DefaultArrayStyle, Broadcasted
 import Base.Broadcast: BroadcastStyle, broadcasted
@@ -99,4 +99,5 @@ end
 include("precompile.jl")
 
 
+
 end #module

src/banded/BandedMatrix.jl

@@ -1007,3 +1007,4 @@ function resize(A::BandedSubBandedMatrix, n::Integer, m::Integer)
     l,u = bandwidths(A)
     _BandedMatrix(reshape(resize!(vec(copy(bandeddata(A))), (l+u+1)*m), l+u+1, m), n, l,u)
 end
+

src/generic/AbstractBandedMatrix.jl

@@ -336,3 +336,66 @@
     copy(A::Adjoint{T,<:AbstractBandedMatrix}) where T = copy(parent(A))'
     copy(A::Transpose{T,<:AbstractBandedMatrix}) where T = transpose(copy(parent(A)))
 end
+
+function sum!(ret::AbstractArray, A::AbstractBandedMatrix)
+    #Behaves similarly to Base.sum!
+    ret .= 0
+    n,m = size(A)
+    s = size(ret)
+    l = length(s)
+    #Check for singleton dimension and perform respective sum
+    if s[1] == 1 && (l == 1 || s[2]==1)
+        for j = 1:m, i = colrange(A, j)
+            ret .+= A[i, j]
+        end
+    elseif s[1] == n && (l == 1 || s[2]==1)
+        for i = 1:n, j = rowrange(A, i)
+            ret[i, 1] += A[i, j]
+        end
+    elseif s[1] == 1 && s[2] == m
+        for j = 1:m, i = colrange(A, j)
+            ret[1, j] += A[i, j]
+        end
+    elseif s[1] == n && s[2] == m
+        ret = A
+    else
+        throw(DimensionMismatch("reduction on matrix of size ($n, $m) with output size $s"))
+    end
+end
+
+function sum(A::AbstractBandedMatrix; dims=:)
+    if dims isa Colon
+        l, u = bandwidths(A)
+        ret = zero(eltype(A))
+        if l + u < 0
+            return ret
+        end
+        n, m = size(A)
+        for j = 1:m, i = colrange(A, j)
+            ret += A[i, j]
+        end
+        ret
+    elseif dims > 2
+        A
+    elseif dims == 2
+        l, u = bandwidths(A)
+        n, m = size(A)
+        ret = zeros(eltype(A), n, 1)
+        if l + u < 0
+            return ret
+        end
+        sum!(ret, A)
+        ret
+    elseif dims == 1
+        l, u = bandwidths(A)
+        n, m = size(A)
+        ret = zeros(eltype(A), 1, m)
+        if l + u < 0
+            return ret
+        end
+        sum!(ret, A)
+        ret
+    else
+        throw(ArgumentError("dimension must be ≥ 1, got $dims"))
+    end
+end

test/runtests.jl

@@ -26,3 +26,4 @@ include("test_symbanded.jl")
 include("test_tribanded.jl")
 include("test_interface.jl")
 include("test_miscs.jl")
+include("test_sum.jl")

test/test_sum.jl

@@ -0,0 +1,20 @@
+module TestSum
+
+using Test, BandedMatrices, Random
+
+r = brand(Float64,rand(1:10_000),rand(1:10_000),rand(-20:100),rand(-20:100))
+empty_r = brand(Float64,rand(1:1_000),rand(1:1_000),rand(1:100),rand(-200:-101))
+n,m = size(empty_r)
+matr = Matrix(r)
+@testset "sum" begin
+    @test sum(empty_r) == 0
+    @test sum(empty_r; dims = 2) == zeros(n,1)
+    @test sum(empty_r; dims = 1) == zeros(1,m)
+    @test sum(r) == sum(matr)
+    @test sum(r; dims=2) == sum(matr; dims=2)
+    @test sum(r; dims=1) == sum(matr; dims=1)
+    @test sum(r; dims=3) == r
+    @test_throws ArgumentError sum(r; dims=0)
+end
+
+end