Sum of squares for sparse matrix

sumsq_sparse.jl

# These two functions are defined only in Julia 0.4
rowvals(S::SparseMatrixCSC) = S.rowval
nzrange(S::SparseMatrixCSC, col::Integer) = S.colptr[col]:(S.colptr[col+1]-1)

# Based on the example from http://julia.readthedocs.org/en/latest/stdlib/sparse/#Base.nzrange
function sumsq(A::SparseMatrixCSC)
    vals = nonzeros(A)
    m, n = size(A)
    count = zero(eltype(A))
    for i = 1:n
       for j in nzrange(A, i)
          count += vals[j]^2
       end
    end
    count
end

function sumsq_vectorized(A::SparseMatrixCSC)
    sum(A .* A)
end


A = sprand(2000, 15037, .1)

@assert abs(sumsq(A) - sumsq_vectorized(A)) < 1e-7

@time for i in 1:100
    sumsq(A)
end
# elapsed time: 0.615286572 seconds (1600 bytes allocated)

@time for i in 1:100
    sumsq_vectorized(A)
end
# elapsed time: 7.953261595 seconds (9650144000 bytes allocated, 46.79% gc time)