Skip to content

proper errors in libSparse instead of asserts #93

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.

Sign up for GitHub

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

Merged
ViralBShah merged 2 commits into from
Dec 8, 2025
Merged

proper errors in libSparse instead of asserts #93

Show file tree
Hide file tree
Changes from all commits
Commits
Show all changes
2 commits
Select commit Hold shift + click to select a range
d7b9f3c
proper errors in libSparse instead of asserts
luke-kiernan Nov 27, 2025
adc3654
fix failing tests
luke-kiernan Nov 27, 2025
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
23 changes: 15 additions & 8 deletions src/libSparse/AAFactorizations.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -57,17 +57,22 @@ function factor!(aa_fact::AAFactorization{T},
end

function solve(aa_fact::AAFactorization{T}, b::StridedVecOrMat{T}) where T<:vTypes
@assert size(aa_fact.matrixObj)[2] == size(b, 1)
size(aa_fact.matrixObj)[2] != size(b, 1) && throw(DimensionMismatch(
"Matrix and right-hand side size mismatch: got "
* "$(size(aa_fact.matrixObj)[2]) and $(size(b, 1))"))
factor!(aa_fact)
x = Array{T}(undef, size(aa_fact.matrixObj)[2], size(b)[2:end]...)
SparseSolve(aa_fact._factorization, b, x)
return x
end

function solve!(aa_fact::AAFactorization{T}, xb::StridedVecOrMat{T}) where T<:vTypes
@assert (xb isa StridedVector) ||
(size(aa_fact.matrixObj)[1] == size(aa_fact.matrixObj)[2]) "Can't in-place " *
"solve: x and b are different sizes and Julia cannot resize a matrix."
((xb isa StridedMatrix) &&
(size(aa_fact.matrixObj)[1] != size(aa_fact.matrixObj)[2])) &&
throw(ArgumentError("Can't in-place " *
"solve: x and b are different sizes and Julia cannot resize a matrix."
)
)
factor!(aa_fact)
SparseSolve(aa_fact._factorization, xb)
return xb # written in imitation of KLU.jl, which also returns
Expand All @@ -79,11 +84,13 @@ LinearAlgebra.ldiv!(aa_fact::AAFactorization{T}, xb::StridedVecOrMat{T}) where T
function LinearAlgebra.ldiv!(x::StridedVecOrMat{T},
aa_fact::AAFactorization{T},
b::StridedVecOrMat{T}) where T<:vTypes
@assert size(aa_fact.matrixObj)[2] == size(b, 1)
size(aa_fact.matrixObj)[2] != size(b, 1) && throw(DimensionMismatch(
"Matrix and right-hand side size mismatch: got "
* "$(size(aa_fact.matrixObj)[2]) and $(size(b, 1))"))
size(aa_fact.matrixObj)[2] != size(x, 1) && throw(DimensionMismatch(
"Matrix and output size mismatch: got "
* "$(size(aa_fact.matrixObj)[2]) and $(size(x, 1))"))
factor!(aa_fact)
SparseSolve(aa_fact._factorization, b, x)
return x
end



25 changes: 17 additions & 8 deletions src/libSparse/AASparseMatrices.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -23,7 +23,8 @@ function AASparseMatrix(n::Int, m::Int,
col::StridedVector{Clong}, row::StridedVector{Cint},
data::StridedVector{T},
attributes::att_type = ATT_ORDINARY) where T<:vTypes
@assert stride(col, 1) == 1 && stride(row, 1) == 1 && stride(data, 1) == 1
(stride(col, 1) != 1 || stride(row, 1) != 1 || stride(data, 1) != 1) &&
throw(ArgumentError("col, row, and data must have stride 1"))
# I'm assuming here that pointer(row) == pointer(_row_inds),
# ie that col, row, and data are passed by reference, not by value.
s = SparseMatrixStructure(n, m, pointer(col),
Expand Down Expand Up @@ -62,7 +63,7 @@ LinearAlgebra.istril(M::AASparseMatrix) = istri(M) && (MM.matrix.structure.attri
& ATT_TRIANGLE_MASK == ATT_UPPER_TRIANGLE)

function Base.getindex(M::AASparseMatrix, i::Int, j::Int)
@assert all((1, 1) .<= (i,j) .<= size(M))
((size(M)[1] >= i >= 1) && (size(M)[2] >= j >= 1)) || throw(BoundsError(M, (i, j)))
(startCol, endCol) = (M._colptr[j], M._colptr[j+1]-1) .+ 1
rowsInCol = @view M._rowval[startCol:endCol]
ind = searchsortedfirst(rowsInCol, i-1)
Expand All @@ -73,7 +74,7 @@ function Base.getindex(M::AASparseMatrix, i::Int, j::Int)
end

function Base.getindex(M::AASparseMatrix, i::Int)
@assert 1 <= i <= size(M)[1]*size(M)[2]
1 <= i <= size(M)[1]*size(M)[2] || throw(BoundsError(M, i))
return M[(i-1) % size(M)[1] + 1, div(i-1, size(M)[1]) + 1]
end
# Creates a new structure, referring to the same data,
Expand All @@ -83,15 +84,19 @@ Base.transpose(M::AASparseMatrix) = AASparseMatrix(SparseGetTranspose(M.matrix),
M._colptr, M._rowval, M._nzval)

function Base.:(*)(A::AASparseMatrix{T}, x::StridedVecOrMat{T}) where T<:vTypes
@assert size(x)[1] == size(A)[2]
size(x)[1] == size(A)[2] || throw(DimensionMismatch(
"Matrix and right-hand side size mismatch: got "
* "$(size(A)[2]) and $(size(x, 1))"))
y = Array{T}(undef, size(A)[1], size(x)[2:end]...)
SparseMultiply(A.matrix, x, y)
return y
end

function Base.:(*)(alpha::T, A::AASparseMatrix{T},
x::StridedVecOrMat{T}) where T<:vTypes
@assert size(x)[1] == size(A)[2]
size(x)[1] == size(A)[2] || throw(DimensionMismatch(
"Matrix and right-hand side size mismatch: got "
* "$(size(A)[2]) and $(size(x, 1))"))
y = Array{T}(undef, size(A)[1], size(x)[2:end]...)
SparseMultiply(alpha, A.matrix, x, y)
return y
Expand All @@ -102,7 +107,9 @@ Computes y += A*x in place. Note that this modifies its LAST argument.
"""
function muladd!(A::AASparseMatrix{T}, x::StridedVecOrMat{T},
y::StridedVecOrMat{T}) where T<:vTypes
@assert size(x) == size(y) && size(x)[1] == size(A)[2]
(size(x) == size(y) && size(x)[1] == size(A)[2]) || throw(DimensionMismatch(
"Matrix and right-hand side size mismatch: got "
* "$(size(A)[2]) and $(size(x, 1))"))
SparseMultiplyAdd(A.matrix, x, y)
end

Expand All @@ -111,6 +118,8 @@ Computes y += alpha*A*x in place. Note that this modifies its LAST argument.
"""
function muladd!(alpha::T, A::AASparseMatrix{T},
x::StridedVecOrMat{T}, y::StridedVecOrMat{T}) where T<:vTypes
@assert size(x) == size(y) && size(x)[1] == size(A)[2]
(size(x) == size(y) && size(x)[1] == size(A)[2]) || throw(DimensionMismatch(
"Matrix and right-hand side size mismatch: got "
* "$(size(A)[2]) and $(size(x, 1))"))
SparseMultiplyAdd(alpha, A.matrix, x, y)
end
end
8 changes: 4 additions & 4 deletions src/libSparse/wrappers.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -195,15 +195,15 @@ const LIBSPARSE = "/System/Library/Frameworks/Accelerate.framework/Versions/A/Fr
Base.cconvert(::Type{DenseMatrix{T}}, m::StridedMatrix{T}) where T<:vTypes = m

function Base.unsafe_convert(::Type{DenseMatrix{T}}, m::StridedMatrix{T}) where T<:vTypes
@assert stride(m, 1) == 1
stride(m, 1) == 1 || throw(ArgumentError("Stride of first dimension must be 1"))
# hard-coded ATT_ORDINARY
return DenseMatrix{T}(size(m)..., stride(m, 2), ATT_ORDINARY, pointer(m))
end

Base.cconvert(::Type{DenseVector{T}}, v::StridedVector{T}) where T<:vTypes = v

function Base.unsafe_convert(::Type{DenseVector{T}}, v::StridedVector{T}) where T<:vTypes
@assert stride(v, 1) == 1
stride(v, 1) == 1 || throw(ArgumentError("Stride of first dimension must be 1"))
return DenseVector{T}(size(v)[1], pointer(v))
end

Expand Down Expand Up @@ -356,7 +356,7 @@ for T in (Cfloat, Cdouble)
:_Z12SparseFactorh19SparseMatrix_Double
@eval function SparseFactorNoErrors(arg1::SparseFactorization_t,
arg2::SparseMatrix{$T})::SparseOpaqueFactorization{$T}
@assert arg1 != SparseFactorizationTBD
arg1 != SparseFactorizationTBD || throw(ArgumentError("Factorization type must be specified"))
@ccall(LIBSPARSE.$sparseFactorMatrix(arg1::Cuint,
arg2::SparseMatrix{$T})::SparseOpaqueFactorization{$T})
end
Expand All @@ -367,7 +367,7 @@ for T in (Cfloat, Cdouble)
arg2::SparseMatrix{$T},
arg3::SparseSymbolicFactorOptions,
arg4::SparseNumericFactorOptions)
@assert arg1 != SparseFactorizationTBD
arg1 != SparseFactorizationTBD || throw(ArgumentError("Factorization type must be specified"))
@ccall(LIBSPARSE.$sparseFactorMatrixOpts(
arg1::Cuint,
arg2::SparseMatrix{$T},
Expand Down
47 changes: 44 additions & 3 deletions test/libSparseTests.jl
View file
Open in desktop
Original file line number Diff line number Diff line change
Expand Up @@ -172,6 +172,7 @@ using LinearAlgebra
ordinary = AASparseMatrix(sparseM)
@test !issymmetric(ordinary) && !istril(ordinary) && !istriu(ordinary)
symmetricM = sparse(sparseM + sparseM')
# verify that the julia matrix is has the desired properties.
@assert issymmetric(symmetricM) && !istril(symmetricM) && !istriu(symmetricM)
symmetric = AASparseMatrix(symmetricM)
@test issymmetric(symmetric) && !istril(symmetric) && !istriu(symmetric)
Expand Down Expand Up @@ -215,7 +216,7 @@ using LinearAlgebra
factor!(f)
catch err1
end
@test sprint(showerror, err1) == "The matrix is singular."
@test occursin("singular", sprint(showerror, err1))

err2 = nothing
temp = sprand(N, N, 0.5)
Expand All @@ -233,11 +234,51 @@ using LinearAlgebra
nonSymmetric = sparse(temp*temp' + diagm(rand(N)) + singular)
try
f = AAFactorization(nonSymmetric)
@assert size(f.matrixObj)[1] == size(f.matrixObj)[2]
# Verify the matrix is square
@test size(f.matrixObj)[1] == size(f.matrixObj)[2]
factor!(f, AppleAccelerate.SparseFactorizationCholesky)
catch err3
end
@test sprint(showerror, err3) == "Cannot perform symmetric matrix factorization of non-symmetric matrix.\n"
@test startswith(sprint(showerror, err3),
"Cannot perform symmetric matrix factorization"
)
end

@testset "DimensionMismatch errors" begin
N = 4
M = 5
jlA = sprand(Float64, N, M, 0.5)
while rank(Array(jlA)) < min(N, M)
jlA = sprand(Float64, N, M, 0.5)
end
f = AAFactorization(jlA)

# Test solve with wrong dimension vector
wrong_b = rand(Float64, N) # should be M
@test_throws DimensionMismatch solve(f, wrong_b)

# Test solve with wrong dimension matrix
wrong_B = rand(Float64, N, 3) # should be M x 3
@test_throws DimensionMismatch solve(f, wrong_B)

# Test ldiv! with wrong dimensions
x = rand(Float64, N)
b = rand(Float64, M)
@test_throws DimensionMismatch LinearAlgebra.ldiv!(x, f, b)
end

@testset "ArgumentError for in-place solve" begin
N = 4
M = 5
# Test that in-place solve with matrix throws ArgumentError for non-square factorization
jlA = sprand(Float64, N, M, 0.9)
while rank(Array(jlA)) < min(N, M)
jlA = sprand(Float64, N, M, 0.9)
end
f = AAFactorization(jlA)

xb_matrix = rand(Float64, M, 3)
@test_throws ArgumentError solve!(f, xb_matrix)
end

@testset "Cholesky" begin
Expand Down