sparse2.jl

# Insert a value into a sorted list
# Uses linear search
# Returns 1 if a new value was inserted else 0
function insertsorted!(listrow, listval, row, val)
    if isempty(listrow)
        insert!(listrow, 1, row)
        insert!(listval, 1, val)
        return 1
    end
            
    it = 1
    @inbounds while it <= length(listrow) && listrow[it] < row
        it += 1
    end

    if !(it > length(listrow)) && listrow[it] == row
        listval[it] += val
        return 0
    end

    insert!(listrow, it, row)
    insert!(listval, it, val)
    return 1
end


dimlub(I) = length(I)==0 ? 0 : Int(maximum(I))
sparse2(I, J, V, nzperrow = 0) = sparse2(I, J, V, dimlub(I), dimlub(J), nzperrow)
function sparse2(I, J, V, m, n, nzperrow = 20)
    # Create buckets
    rows = [zeros(Int, 0) for i = 1:n]
    vals = [zeros(0) for i = 1:n]
    
    for r in rows
        sizehint!(r, nzperrow)
    end

    # Insert into bucket
    nzvals = 0
    for i in 1:length(I)
        nzvals += insertsorted!(rows[J[i]], vals[J[i]], I[i], V[i])
    end

    rowval = zeros(Int, nzvals)
    nzval = zeros(nzvals)
    colptr = zeros(Int, length(rows) +1)

    colptr[1] = 1
    offset = 0
    @inbounds for i in 1:length(rows)
        row = rows[i]
        val = vals[i]
        ns = length(rows[i])
        colptr[i + 1] = colptr[i] + ns
        for j in 1:ns
            rowval[offset + j] = row[j]
            nzval[offset + j] = val[j]
        end
        offset += ns
    end
    return SparseMatrixCSC(m, n, colptr, rowval, nzval);
end

n = 10^5
sparsity = 20 / n
accums = 5
nzs_1 = repmat(rand(1:n, Int(sparsity * n^2)), accums)
nzs_2 = repmat(rand(1:n, Int(sparsity * n^2)), accums)
vals = rand(Float64, length(nzs_1))

A2 = @time sparse2(nzs_1, nzs_2, vals)
A = @time sparse(nzs_1, nzs_2, vals)
A == A2