Skip to content
New issue

Have a question about this project? Sign up for a free GitHub account to open an issue and contact its maintainers and the community.

By clicking “Sign up for GitHub”, you agree to our terms of service and privacy statement. We’ll occasionally send you account related emails.

Already on GitHub? Sign in to your account

Resize buffers in sparse! to satisfy buffer checks in constructor #314

Merged
merged 1 commit into from
Dec 12, 2022
Merged
Show file tree
Hide file tree
Changes from all commits
Commits
File filter

Filter by extension

Filter by extension

Conversations
Failed to load comments.
Loading
Jump to
Jump to file
Failed to load files.
Loading
Diff view
Diff view
18 changes: 12 additions & 6 deletions src/sparsematrix.jl
Original file line number Diff line number Diff line change
Expand Up @@ -1068,9 +1068,8 @@ intermediate CSR forms and require `length(csrrowptr) >= m + 1`,
`length(csrcolval) >= length(I)`, and `length(csrnzval >= length(I))`. Input
array `klasttouch`, workspace for the second stage, requires `length(klasttouch) >= n`.
Optional input arrays `csccolptr`, `cscrowval`, and `cscnzval` constitute storage for the
returned CSC form `S`. `csccolptr` requires `length(csccolptr) >= n + 1`. If necessary,
`cscrowval` and `cscnzval` are automatically resized to satisfy
`length(cscrowval) >= nnz(S)` and `length(cscnzval) >= nnz(S)`; hence, if `nnz(S)` is
returned CSC form `S`. If necessary, these are resized automatically to satisfy
`length(csccolptr) = n + 1`, `length(cscrowval) = nnz(S)` and `length(cscnzval) = nnz(S)`; hence, if `nnz(S)` is
unknown at the outset, passing in empty vectors of the appropriate type (`Vector{Ti}()`
and `Vector{Tv}()` respectively) suffices, or calling the `sparse!` method
neglecting `cscrowval` and `cscnzval`.
Expand All @@ -1081,6 +1080,7 @@ representation of the result's transpose.
You may reuse the input arrays' storage (`I`, `J`, `V`) for the output arrays
(`csccolptr`, `cscrowval`, `cscnzval`). For example, you may call
`sparse!(I, J, V, csrrowptr, csrcolval, csrnzval, I, J, V)`.
Note that they will be resized to satisfy the conditions above.

For the sake of efficiency, this method performs no argument checking beyond
`1 <= I[k] <= m` and `1 <= J[k] <= n`. Use with care. Testing with `--check-bounds=yes`
Expand Down Expand Up @@ -1140,6 +1140,9 @@ function sparse!(I::AbstractVector{Ti}, J::AbstractVector{Ti}, V::Union{Tv,Abstr
end
# This completes the unsorted-row, has-repeats CSR form's construction

# The output array csccolptr can now be resized safely even if aliased with I
resize!(csccolptr, n + 1)

# Sweep through the CSR form, simultaneously (1) calculating the CSC form's column
# counts and storing them shifted forward by one in csccolptr; (2) detecting repeated
# entries; and (3) repacking the CSR form with the repeated entries combined.
Expand Down Expand Up @@ -1188,10 +1191,13 @@ function sparse!(I::AbstractVector{Ti}, J::AbstractVector{Ti}, V::Union{Tv,Abstr
Base.hastypemax(Ti) && (countsum <= typemax(Ti) || throw(ArgumentError("more than typemax(Ti)-1 == $(typemax(Ti)-1) entries")))
end

# Now knowing the CSC form's entry count, resize cscrowval and cscnzval if necessary
# Now knowing the CSC form's entry count, resize cscrowval and cscnzval
# Note: This is done unconditionally to appease the buffer checks in the SparseMatrixCSC
# constructor. If these checks are lifted this resizing is only needed if the
# buffers are too short. csccolptr is resized above.
cscnnz = countsum - Tj(1)
length(cscrowval) < cscnnz && resize!(cscrowval, cscnnz)
length(cscnzval) < cscnnz && resize!(cscnzval, cscnnz)
resize!(cscrowval, cscnnz)
resize!(cscnzval, cscnnz)

# Finally counting-sort the row and nonzero values from the CSR form into cscrowval and
# cscnzval. Tracking write positions in csccolptr corrects the column pointers.
Expand Down
78 changes: 75 additions & 3 deletions test/sparsematrix_constructors_indexing.jl
Original file line number Diff line number Diff line change
Expand Up @@ -12,6 +12,10 @@ using Dates
include("forbidproperties.jl")
include("simplesmatrix.jl")

function same_structure(A, B)
return all(getfield(A, f) == getfield(B, f) for f in (:m, :n, :colptr, :rowval))
end

@testset "uniform scaling should not change type #103" begin
A = spzeros(Float32, Int8, 5, 5)
B = I - A
Expand Down Expand Up @@ -59,9 +63,6 @@ end
end

@testset "spzeros for pattern creation (structural zeros)" begin
function same_structure(A, B)
return all(getfield(A, f) == getfield(B, f) for f in (:m, :n, :colptr, :rowval))
end
I = [1, 2, 3]
J = [1, 3, 4]
V = zeros(length(I))
Expand Down Expand Up @@ -1625,4 +1626,75 @@ end
@test_throws ArgumentError SparseArrays.expandptr([2; 3])
end

@testset "sparse!" begin
using SparseArrays: sparse!, getcolptr, getrowval, nonzeros

function allocate_arrays(m, n)
N = round(Int, 0.5 * m * n)
Tv, Ti = Float64, Int
I = Ti[rand(1:m) for _ in 1:N]; I = Ti[I; I]
J = Ti[rand(1:n) for _ in 1:N]; J = Ti[J; J]
V = Tv.(I)
csrrowptr = Vector{Ti}(undef, m + 1)
csrcolval = Vector{Ti}(undef, length(I))
csrnzval = Vector{Tv}(undef, length(I))
klasttouch = Vector{Ti}(undef, n)
csccolptr = Vector{Ti}(undef, n + 1)
cscrowval = Vector{Ti}()
cscnzval = Vector{Tv}()
return I, J, V, klasttouch, csrrowptr, csrcolval, csrnzval, csccolptr, cscrowval, cscnzval
end

for (m, n) in ((10, 5), (5, 10), (10, 10))
# Passing csr vectors
I, J, V, klasttouch, csrrowptr, csrcolval, csrnzval = allocate_arrays(m, n)
S = sparse(I, J, V, m, n)
S! = sparse!(I, J, V, m, n, +, klasttouch, csrrowptr, csrcolval, csrnzval)
@test S == S!
@test same_structure(S, S!)

# Passing csr vectors + csccolptr
I, J, V, klasttouch, csrrowptr, csrcolval, csrnzval, csccolptr = allocate_arrays(m, n)
S = sparse(I, J, V, m, n)
S! = sparse!(I, J, V, m, n, +, klasttouch, csrrowptr, csrcolval, csrnzval, csccolptr)
@test S == S!
@test same_structure(S, S!)
@test getcolptr(S!) === csccolptr

# Passing csr vectors, and csc vectors
I, J, V, klasttouch, csrrowptr, csrcolval, csrnzval, csccolptr, cscrowval, cscnzval =
allocate_arrays(m, n)
S = sparse(I, J, V, m, n)
S! = sparse!(I, J, V, m, n, +, klasttouch, csrrowptr, csrcolval, csrnzval,
csccolptr, cscrowval, cscnzval)
@test S == S!
@test same_structure(S, S!)
@test getcolptr(S!) === csccolptr
@test getrowval(S!) === cscrowval
@test nonzeros(S!) === cscnzval

# Passing csr vectors, and csc vectors of insufficient lengths
I, J, V, klasttouch, csrrowptr, csrcolval, csrnzval, csccolptr, cscrowval, cscnzval =
allocate_arrays(m, n)
S = sparse(I, J, V, m, n)
S! = sparse!(I, J, V, m, n, +, klasttouch, csrrowptr, csrcolval, csrnzval,
resize!(csccolptr, 0), resize!(cscrowval, 0), resize!(cscnzval, 0))
@test S == S!
@test same_structure(S, S!)
@test getcolptr(S!) === csccolptr
@test getrowval(S!) === cscrowval
@test nonzeros(S!) === cscnzval

# Passing csr vectors, and csc vectors aliased with I, J, V
I, J, V, klasttouch, csrrowptr, csrcolval, csrnzval = allocate_arrays(m, n)
S = sparse(I, J, V, m, n)
S! = sparse!(I, J, V, m, n, +, klasttouch, csrrowptr, csrcolval, csrnzval, I, J, V)
@test S == S!
@test same_structure(S, S!)
@test getcolptr(S!) === I
@test getrowval(S!) === J
@test nonzeros(S!) === V
end
end

end